Télécharger nlin_burger.dgibi

Retour à la liste

* fichier :  nlin_burger.dgibi
*************************************************************************
*
* NOM         : nlin_burger.dgibi
*
* DESCRIPTION : We compute the 2D Burger problem.
*               Scalar numerical viscosity is added.
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Alberto BECCANTINI (CEA/DEN/DM2S/SFME/LTMF)
*               mél : beccantini@cea.fr
*************************************************************************
* Prière de PRENDRE LE TEMPS de compléter les commentaires
* en cas de modification de ce sous-programme afin de faciliter
* la maintenance !
*************************************************************************
*************************************************************************
*************************************************************************
******* PERSONAL PROCEDURES *********************************************
*************************************************************************
*************************************************************************
*
*BEGINPROCEDUR matmas
*
*************************************************************************
*
* NOM         : MATMAS
*
* DESCRIPTION : We compute 
*
*               (Np , Nd )
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Alberto BECCANTINI (CEA/DEN/DM2S/SFME/LTMF)
*               mél : beccantini@cea.fr
*************************************************************************
* Prière de PRENDRE LE TEMPS de compléter les commentaires
* en cas de modification de ce sous-programme afin de faciliter
* la maintenance !
*************************************************************************
*
*
'DEBPROC' MATMAS ;
*
'ARGUMENT'  _mt*'MAILLAGE'    ;
'ARGUMENT'  nomp*'MOT     '  ;
'ARGUMENT'  nomd*'MOT     '  ;
'ARGUMENT'  discp*'MOT     '  ;
'ARGUMENT'  discd*'MOT     ' ;
*
*  Arguments
*  _mt   = MAILLAGE, quaf mesh
*  nomp  = MOT, name of the primal variable
*  nomd  = MOT, name of the dual variable
*  discp = MOT, name of the discretization of the primal variable
*  discd = MOT, name of the discretization of the dual variable
*
idim   = 'VALEUR' 'DIME' ;
*
*   Matrix A for the operator NLIN
*
numop  = 1 ;
numder = idim ;
numvar = 1 ;
* The variable involved in A is the primal variable T
numdat = 0 ;
numcof = 0 ;
*
discg = 'LINE' ;
* Geometrical discretisation
methgau = 'GAU7' ;
* Gauss points involved
*
A = ININLIN numop numvar numdat numcof numder ;   
A . 'VAR' . 1 . 'NOMDDL' = 'MOTS' nomp ;
A . 'VAR' . 1 . 'DISC'   = discp ;
*
*   Matrix B for the operator NLIN
*
B = ININLIN numop numvar 0 0 numder ;
* For B the coef no coefficients are involved
* but 1
B . 'VAR' . 1 . 'NOMDDL' = 'MOTS' nomd ;
B . 'VAR' . 1 . 'DISC'   = discd ;
*
A . 1 . 1 . 0 = 'LECT' ;
B . 1 . 1 . 0 = 'LECT'  ;
*
mat = 'NLIN' discg _mt A B methgau ;
'RESPRO' mat ;
'FINPROC' ;
*
*ENDPROCEDUR matmas
*BEGINPROCEDUR matkv
*
*************************************************************************
*
* NOM         : MATKV
*
* DESCRIPTION : We compute the matrix
*
*               ((cx \dep{u_x}{x} '+' cy \dep{u_x}{y}), Nvx)
*               ((cx \dep{u_y}{x} '+' cy \dep{u_y}{y}), Nvy)
*
*               where Nv = shape function for (u_x,u_y)^T
*
*               Names of primal variables = names of the dual variables
*                = 'UX' 'UY'
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Alberto BECCANTINI (CEA/DEN/DM2S/SFME/LTMF)
*               mél : beccantini@cea.fr
*************************************************************************
* Prière de PRENDRE LE TEMPS de compléter les commentaires
* en cas de modification de ce sous-programme afin de faciliter
* la maintenance !
*************************************************************************
*
'DEBPROC' MATKV ;
'ARGUMENT'   _mt*'MAILLAGE' ;
'ARGUMENT' discv*'MOT     ' ;
'ARGUMENT' discc*'MOT     ' ;
'ARGUMENT'     cx*'CHPOINT' ;
'ARGUMENT'     cy*'CHPOINT' ;
*
*  _mt   = solid quaf mesh
*  discv = discretization type ux, uy
*  discc = discretization type for cx, cy
*  cx    = CHPOINT ('SCAL')
*  cy    = CHPOINT ('SCAL')
*
idim   = 'VALEUR' 'DIME' ;
*
numop  = 2 ;
numder = idim ;
numvar = 2 ;
* Two variables, ux, uy
numdat = 2 ;
* cx, cy
numcof = 2 ;
* Coef: cx  cy (two identity functions)
discg = 'LINE' ;
* Geometrical discretisation
methgau = 'GAU7' ;
* Gauss points involved
*
*   Matrix A for the operator NLIN
*
A = ININLIN numop numvar numdat numcof numder ;   
A . 'VAR' . 1 . 'NOMDDL' = 'MOTS' 'UX' ;
A . 'VAR' . 1 . 'DISC'   = discv  ;
A . 'VAR' . 2 .  'NOMDDL' = 'MOTS' 'UY' ;
A . 'VAR' . 2 .  'DISC'   = discv  ;
* cx
A . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
A . 'DAT' . 1 . 'DISC' =  discc ;
A . 'DAT' . 1 . 'VALEUR' = cx ;
A . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
A . 'COF' . 1 . 'LDAT' = 'LECT' 1 ;
* cy 
A . 'DAT' . 2 . 'NOMDDL' = 'MOTS' 'SCAL' ;
A . 'DAT' . 2 . 'DISC' =  discc ;
A . 'DAT' . 2 . 'VALEUR' = cy ;
A . 'COF' . 2 . 'COMPOR' = 'IDEN' ;
A . 'COF' . 2 . 'LDAT' = 'LECT' 2 ;
* 
*   Matrix B for the operator NLIN
*
numvar = 2 ;
B = ININLIN numop numvar 0 0 numder ;
*
* For B, no coefficients are involved
* but 1
*
B . 'VAR' . 1 . 'NOMDDL' = 'MOTS' 'UX' ;
B . 'VAR' . 1 . 'DISC'   = discv ;
B . 'VAR' . 2 . 'NOMDDL' = 'MOTS' 'UY' ;
B . 'VAR' . 2 . 'DISC'   = discv ;
*
* Contribution
*
* cx \dep{ux}{x} Nvx
A . 1 . 1 . 1 = 'LECT' 1 ;
B . 1 . 1 . 0 = 'LECT'   ;
* cy \dep{ux}{y}  Nvx
A . 1 . 1 . 2 = 'LECT' 2 ;
* cx \dep{uy}{x} Nvy
A . 2 . 2 . 1 = 'LECT' 1 ;
B . 2 . 2 . 0 = 'LECT'   ;
* cy \dep{uy}{y}  Nd
A . 2 . 2 . 2 = 'LECT' 2 ;
*
mat = 'NLIN' discg _mt A B  methgau ;
'RESPRO' mat ;
'FINPROC' ;
*
*ENDPROCEDUR matkv
*BEGINPROCEDUR matlap
*
*************************************************************************
*
* NOM         : MATLAP
*
* DESCRIPTION : We compute the matrix arising from the computation of
*               the scalar product
*
*               ( -(coef \div( \alpha \grad T)) , Nd )
*
*               where Nd = test function for the dual variable.
*
*               -\int_{\Omega} coef \div( \alpha \grad T) Nd dV =
*                  -\int_{\delta \Omega}
*                     coef \alpha (\grad T \cdot n) Nd dS '+'
*                  +\int_{\Omega}
*                     coef \alpha (\grad T \cdot \grad Nd) dV '+'
*                  +\int_{\Omega}
*                     Nd \alpha (\grad T \cdot \grad coef) dV
*
*               The matrix is issue from the 2nd and 3rd contribution
*               (volume integrals).
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Alberto BECCANTINI (CEA/DEN/DM2S/SFME/LTMF)
*               mél : beccantini@cea.fr
*************************************************************************
* Prière de PRENDRE LE TEMPS de compléter les commentaires
* en cas de modification de ce sous-programme afin de faciliter
* la maintenance !
*************************************************************************
*
*
'DEBPROC' MATLAP ;
*
'ARGUMENT'  _mt*'MAILLAGE'    ;
'ARGUMENT'   nomt*'MOT     '  ;
'ARGUMENT'  disct*'MOT     '  ;
'ARGUMENT'   nomd*'MOT     '  ;
'ARGUMENT'  discd*'MOT     '  ;
'ARGUMENT'  discal*'MOT     ' ;
'ARGUMENT'  alpha*'CHPOINT'   ;
'ARGUMENT'  discof*'MOT     ' ;
'ARGUMENT'  coef*'CHPOINT'    ;
*
*    Arguments
*
*    _mt    = surface QUAF mesh
*    nomt   = name of the primal variable T
*    disct  = discretization type of T
*    nomd   = name of the dual variable
*    discd  = discretization type of the dual variable
*    discal = discretization type of alpha
*    discof = discretization type of coef
*    alpha  = CHPOINT ('SCAL')
*    coef   = CHPOINT ('SCAL')
*
idim   = 'VALEUR' 'DIME' ;
*
*    Four contributions
*
*    (coef alpha \dep{T}{x} \dep{Nd}{x}
*   + coef alpha \dep{T}{y} \dep{Nd}{y}
*   + Nd alpha \dep{T}{x} \dep{coef}{x}
*   + Nd alpha \dep{T}{y} \dep{coef}{y})
*
*   Matrix A for the operator NLIN
*
numop  = 4 ;
numder = idim ;
numvar = 1 ;
* The variable involved in A is the primal variable T
numdat = 2 ;
* Two data: alpha, coef
numcof = 4 ;
* Four functions: f(alpha) = alpha
*                 f(coef)  = coef
*                 f(coef) = \dep{coef}{x}
*                 f(coef) = \dep{coef}{y}
discg = 'LINE' ;
* Geometrical discretisation
methgau = 'GAU7' ;
* Gauss points involved
*
A = ININLIN numop numvar numdat numcof numder ;   
A . 'VAR' . 1 . 'NOMDDL' = 'MOTS' nomt ;
A . 'VAR' . 1 . 'DISC'   = disct ;
* alpha
A . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
A . 'DAT' . 1 . 'DISC' = discal ;
A . 'DAT' . 1 . 'VALEUR' = alpha ;
* coef
A . 'DAT' . 2 . 'NOMDDL' = 'MOTS' 'SCAL' ;
A . 'DAT' . 2 . 'DISC' = discof ;
A . 'DAT' . 2 . 'VALEUR' = coef ;
* Function alpha
A . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
A . 'COF' . 1 . 'LDAT' = 'LECT' 1 ;
* Function coef
A . 'COF' . 2 . 'COMPOR' = 'IDEN' ;
A . 'COF' . 2 . 'LDAT' = 'LECT' 2 ;
* Function \dep{coef}{x}
A . 'COF' . 3 . 'COMPOR' = 'D/DX1' ;
A . 'COF' . 3 . 'LDAT' = 'LECT' 2  ;
* Function \dep{coef}{y}
A . 'COF' . 4 . 'COMPOR' = 'D/DX2' ;
A . 'COF' . 4 . 'LDAT' = 'LECT' 2 ;
*
*   Matrix B for the operator NLIN
*
B = ININLIN numop numvar 0 0 numder ;
* For B the coef no coefficients are involved
* but 1
B . 'VAR' . 1 . 'NOMDDL' = 'MOTS' nomd ;
B . 'VAR' . 1 . 'DISC'   = discd ;
*
* Contribution
*
* coef alpha \dep{T}{x} \dep{Nd}{x}
A . 1 . 1 . 1 = 'LECT' 1 2 ;
B . 1 . 1 . 1 = 'LECT'  ;
* coef alpha \dep{T}{y} \dep{Nd}{y}
A . 2 . 1 . 2 = 'LECT' 1 2 ;
B . 2 . 1 . 2 = 'LECT'  ;
* Nd alpha \dep{T}{x} \dep{coef}{x}
A . 3 . 1 . 1 = 'LECT' 1 3 ;
B . 3 . 1 . 0 = 'LECT'  ;
* Nd alpha \dep{T}{y} \dep{coef}{y}
A . 4 . 1 . 2 = 'LECT' 1 4 ;
B . 4 . 1 . 0 = 'LECT'  ;
mlapn = 'NLIN' discg _mt A B  methgau ;
'RESPRO' mlapn ;
'FINPROC' ;
*ENDPROCEDUR matlap
*BEGINPROCEDUR resoup
*************************************************************************
**** Resolution of a linear system **************************************
*************************************************************************
*
'DEBPROC' RESOUP ;
   'ARGUMENT' LOGTPS*'LOGIQUE ' ;
   'ARGUMENT' mat*'MATRIK  ' ;
   'ARGUMENT' matpre/'MATRIK  ' ;
   'ARGUMENT' smb*'CHPOINT ' ;
   'ARGUMENT' ccl*'CHPOINT ' ;
   'ARGUMENT' res*'FLOTTANT' ;
   'ARGUMENT' nit*'ENTIER'   ;
   'ARGUMENT' pre/'ENTIER'   ;
   'ARGUMENT' gggtre/'TABLE   ' ;
   'ARGUMENT' gggtcv/'LOGIQUE ' ;
   'ARGUMENT' typslv/'MOT     ' ;
*
*  Det. tt s.t.
*
*  mat . tt = smb
*  tt(boundary) = ccl
*
*  
*  mat = matrix to "inverse"
*  matpre = matrix for which preconditioner has been computed (optional)
*  smb = right hand side
*  ccl = imposed boundary conditions
*  res = rvk . 'RESID'
*  nit = max number of linear iteration (for an iterative solver)
*  pre = type of preconditioning matrix (for an iterative solver)
*  gggtre = table of resolution to store 'XINIT'
*  gggtcv = graphic of the history for the iterative solver
*           (optional FAUX)
*  typslv = 'DIRECT'/'ITER' (optional, ITER)
*
   'SI' ('NON' ('EXISTE' typslv)) ;
      tslv = VRAI ;
   'SINON' ;      
      'SI' ('EGA' typslv 'DIRECT') ;
         tslv = FAUX ;
      'SINON' ;
         'SI' ('EGA' typslv 'ITER') ;
            tslv = VRAI ;
         'SINON' ;
            'ERREUR' 'Solveur DIRECT ou ITERatif ?' ;
         'FINSI' ;
      'FINSI' ;
   'FINSI' ;
    gggniv = 1 ;
   'SI' ('NON' ('EXISTE' pre)) ;
      pre = 5 ;
   'FINSI' ;
   'SI' ('OU' ('NON' ('EXISTE' gggtcv)) ('NON' tslv)) ;
      gggtcv = FAUX ;
   'FINSI' ;      
   rv = 'EQEX' ;
   rvk = rv . 'METHINV' ;
   'SI' ('EXISTE' matpre) ;
      rvk . 'MATASS' = matpre ;
      rvk . 'MAPREC' = matpre ;
   'SINON' ;      
      rvk . 'MATASS' = mat ;
      rvk . 'MAPREC' = mat ;
   'FINSI' ;      
   'SI' tslv ;
      rvk . 'TYPINV' = 3 ;
   'SINON' ;      
      rvk . 'TYPINV' = 1 ;
   'FINSI' ;      
   rvk . 'SCALING' = 1 ;
*  Scaling factor   
   rvk . 'OUBMAT' = 1 ;
*  oubmat = 1 -> We eliminate the elementary matrix
*  'ILUTPPIV' -> 1 We always search for a new pivot
*             -> 0 We do not search for a new pivot
   rvk . 'ILUTPPIV' = 0.1D0 ;
*   rvk . 'ILUTPPIV' = 0.01D0 ;
*   rvk . 'ILUTPPIV' = 0.00001D0 ;
*   rvk . 'ILUTPPIV' = 0.D0 ;
   rvk . 'IMPINV' = gggniv ;
   rvk . 'TYRENU' = 'SLOA' ;
   rvk . 'PCMLAG' = 'APR2' ;
   rvk . 'NITMAX' = nit ;
   rvk . 'RESID' = res ;
   rvk . 'PRECOND' = pre ;
   rvk . 'BCGSBTOL' = 1.D-200 ;
   rvk . 'ILUTLFIL' = 1.5 ;
*
   'SI' ('EXISTE' gggtre) ;
      'SI' ('EXISTE' gggtre 'XINIT') ;
         rvk . 'XINIT' = gggtre . 'XINIT' ;
      'FINSI' ;
   'FINSI' ;      
   'TEMPS' 'ZERO' ;
   tt = 'KRES' mat 'TYPI' rvk 'SMBR' smb 'CLIM' ccl ;
   TABTPS = TEMP 'NOEC';
   tcpu = TABTPS.'TEMPS_CPU'.'INITIAL';
   tcpus = '/' ('FLOTTANT' tcpu) 100.D0 ;
   'SI' LOGTPS ;
      'MESSAGE' ('CHAINE' 'tcpu (s) = ' tcpus) ;
   'FINSI' ;
   'SI' ('ET' gggtcv ('EXISTE' rvk 'CONVINV')) ;
      lcvg = rvk . 'CONVINV' ;
      nit = 'DIME' lcvg ;
      'SI' ('>' nit 1) ;
         lit = 'PROG' 1.D0 PAS 1.D0 NPAS ('-' nit 1) ;
         evit = 'EVOL' 'MANU' 'iter' lit 'log residu'
          ('/' ('LOG' lcvg) ('LOG' 10.D0)) ;
         'DESSIN' evit 'TITR' titglob ;
      'FINSI' ;         
   'FINSI' ;
   'SI' ('EXISTE' gggtre) ;
      gggtre . 'XINIT' = tt ;
   'FINSI' ;      
   'RESPRO' tt ;
*
* End of procedure file RESOU
*
'FINPROC' ;
*
*************************************************************************
**** End of the Resolution of a linear system ***************************
*************************************************************************
*ENDPROCEDUR resoup
*************************************************************************
*************************************************************************
******* END OF PERSONAL PROCEDURES **************************************
*************************************************************************
*************************************************************************
**
*
 'OPTION' 'DIME' 2 'ELEM' 'QUA4' ;
 GRAPH = FAUX ;
*
**** Mesh
*
*
* A4         A3
*
*
*
*
* A1         A2
*
 L  = 1.0 ;
 A1 = 0.0 0.0 ;
 A2 = L   0.0 ;
 A3 = L   L ;
 A4 = 0.0 L ;
*
 DX = L '/' 20 ;
 DY = DX ;
*
 A1A2 = A1 'DROIT' A2 'DINI' DX 'DFIN' DX ;
 A2A3 = A2 'DROIT' A3 'DINI' DY 'DFIN' DY ;
 A3A4 = A3 'DROIT' A4 'DINI' DX 'DFIN' DX ;
 A4A1 = A4 'DROIT' A1 'DINI' DY 'DFIN' DY ;
*
 DOM1 = 'SURFACE' (A1A2 'ET' A2A3 'ET' A3A4 'ET' A4A1) 'PLAN' ;
 PDOM1 = 'CHANGER' DOM1 'POI1' ;
 'SI' GRAPH ;
    'TRACER' DOM1 'TITRE' 'Mesh' ;
 'FINSI' ;
 QDOM1 = 'CHANGER' DOM1 'QUAF' ;
 QA4A1 = 'CHANGER' A4A1 'QUAF' ;
 QA1A2 = 'CHANGER' A1A2 'QUAF' ;
 QA2A3 = 'CHANGER' A2A3 'QUAF' ;
 QA3A4 = 'CHANGER' A3A4 'QUAF' ;
 'ELIMINATION' (QA4A1 'ET' QDOM1) 1.0D-4 ;
 'ELIMINATION' (QA1A2 'ET' QDOM1) 1.0D-4 ;
 'ELIMINATION' (QA2A3 'ET' QDOM1) 1.0D-4 ;
 'ELIMINATION' (QA3A4 'ET' QDOM1) 1.0D-4 ;
*
* Boundary conditions
*
 LIMDIR1  = 'DIFF' ('CHANGER' 'POI1' A4A1) ('MANUEL' 'POI1' A1) ;
 LIMDIR2  = 'DIFF' ('CHANGER' 'POI1' A1A2) ('MANUEL' 'POI1' A1) ;
 LIMDIR3  = ('MANUEL' 'POI1' A1) ;
 UN1 = 'MANUEL' 'CHPO' LIMDIR1 2 'UX' 1.0 'UY' 0.0 ;
 UN2 = 'MANUEL' 'CHPO' LIMDIR2 2 'UX' 0.0 'UY' 1.0 ;
 UN3 = 'MANUEL' 'CHPO' LIMDIR3 2 'UX' 1.0 'UY' 1.0 ;
 UNLIM = UN1 '+' UN2 '+' UN3 ;
*
* IC
*
 UN = ('MANUEL' 'CHPO' DOM1 2 'UX' 0.0 'UY' 0.0) '+' UN1 '+'
       UN2 '+' UN3 ;
*
 VECN = 'VECTEUR' UN ;
 'SI' GRAPH ;
    'TRACER' DOM1 VECN ;
 'FINSI' ;
*
* Time step
*
 DT = 1.0D16 '*' DX '/' 1.0 ;
*
 NU = 'MANUEL' 'CHPO' DOM1 1 'SCAL' (0.25 '*' DX) ;
 CHONE = 'MANUEL' 'CHPO' DOM1 1 'SCAL' 1.0 ;
*
* Matrices
*
 MMAS1 = MATMAS QDOM1 'UX' 'UX' 'LINE' 'LINE' ;
 MMAS2 = MATMAS QDOM1 'UY' 'UY' 'LINE' 'LINE' ;
 MMAS = MMAS1 'ET' MMAS2 ;
 MKON = MATKV QDOM1 'LINE' 'LINE' ('EXCO' 'UX' UN) ('EXCO' 'UY' UN) ;
 MLAP1 = MATLAP  QDOM1 'UX' 'LINE' 'UX' 'LINE' 'LINE' NU 'LINE' CHONE ;
 MLAP2 = MATLAP  QDOM1 'UY' 'LINE' 'UY' 'LINE' 'LINE' NU 'LINE' CHONE ;
 MLAP = MLAP1 'ET' MLAP2 ;
*
*
* Total matrix
*
 MMASSDT = MMAS '/' DT ;
 mconst = MMASSDT 'ET' MLAP ;
 mtot = mconst 'ET' MKON ;
 mtotm = 'KOPS' 'RIMA' mtot  ;
 smbv pipi = 'KOPS' 'MATRIK' ;
*
* Resolution
*
* TAB1=TABLE;
* TAB1.'TITRE'= TABLE ;
* TAB1.1='TIRR ';
* TAB1.'TITRE' . 1  = 'Exact' ;
* TAB1.2='MARQ CROI';
* TAB1.'TITRE' . 2  = 'Numerical' ;
*
 TPS = 0.0 ;
 UN = 'COPIER' UN ;
 UNM = 'COPIER' UN ;
 LIT = 'PROG' ;
 LER = 'PROG' ;
 ERRO = 1.0D10 ;
*
 TABRESUN = 'TABLE' ;
 TABRESUN . 'XINIT' = UN ;
 NITERU = 100 ;
 PREC = 5 ;
 METINV = MOT 'DIRECT' ;
* 
 'REPETER' BL1 100 ;
    TPS = TPS '+' DT ;
    qv  = MMASSDT '*' UN ;
    UN  = RESOUP FAUX mtotm mtotm qv UNLIM 1.0D-60 NITERU
                PREC TABRESUN METINV ;
    'SI' (ERRO < 1.0D-10) ;
       'QUITTER' BL1 ;
    'FINSI' ;
    'SI' (((&BL1 '/' 5) '*' 5) 'EGA' &BL1) ;
        ERRO = 'MAXIMUM' (UN '-' UNM) 'ABS' ;
        'MESSAGE' ('CHAINE' 'ITER = ' &BL1 ',  TPS = ' TPS  ',  ERRO = '
           ERRO) ;
        UNM = 'COPIER' UN ;
        LIT = LIT 'ET' ('PROG' &BL1) ;
        L10 = 'LOG' 10. ;
        LER = LER 'ET' ('PROG' (('LOG' (ERRO '+' 1.0D-20)) '/' L10))  ;
        EVER = 'EVOL' 'MANU' 'it' LIT 'Log(Err)' LER ;
        'SI' GRAPH ;
           'DESSIN' EVER 'TITR' 'Convergence history' 'NCLK' ;
        'FINSI' ;
    'FINSI' ;
*
****** We update the convective matrix
*
    MKON = MATKV QDOM1 'LINE' 'LINE' ('EXCO' 'UX' UN) ('EXCO' 'UY' UN) ;
    MTOT = MCONST 'ET' MKON ;
    mtotm = 'KOPS' 'RIMA' mtot  ;
 'FIN' BL1 ;
*
 UX = 'EXCO' UN 'UX' ;
 UY = 'EXCO' UN 'UY' ;
*
 VECN = 'VECTEUR' UN ;
 'SI' GRAPH ;
    'TRACER' UX DOM1 VECN 'TITR' 'ux' ;  
    'TRACER' UY DOM1 VECN 'TITR' 'uy' ;
 'FINSI' ;
*
* Test
*
 NITER = &BL1 ;
*
 'SI' (NITER > 40) ;
     'MESSAGE' 'Convergence not reached' ;
     'ERREUR' 5 ;
 'FINSI' ;
*
 XX YY = 'COORDONNEE' DOM1 ;
 DOMH = (YY '-' XX) 'POIN' 'SUPERIEUR' 0.1 ;
 erro = 'MAXIMUM' (('EXCO' 'UX' ('REDU' UN DOMH)) '-' 1.0) 'ABS' ;
 'SI' (erro > 1.0D-2) ;
     'ERREUR' 5 ;
 'FINSI' ;
 erro = 'MAXIMUM' (('EXCO' 'UY' ('REDU' UN DOMH))) 'ABS' ;
 'SI' (erro > 1.0D-2) ;
     'ERREUR' 5 ;
 'FINSI' ;
*
 'FIN' ;
 
 
 
 
 

© Cast3M 2003 - Tous droits réservés.
Mentions légales