Télécharger defila2.dgibi

Retour à la liste

*
************************************************************************
************************************************************************
* fichier : defila2.dgibi
*
'OPTI' 'ECHO' 0 ;
*
'SAUTER' 2 'LIGNE' ;
'MESSAGE' ' Execution de defila2.dgibi' ;
'SAUTER' 2 'LIGNE' ;
*
graph    = faux ;
complet  = faux ;
interact = faux ;
*
************************************************************************
* NOM         : DEFILA2
* DESCRIPTION : Ecoulement sous une surface libre soumise à une pression
*               On a tenté de simplifier par rapport à defila :
*               - plus de MATRIK
*               - procédures de post-traitements (bilans de force)
*                 plus courtes et plus générales.
*               - on teste différents solveurs (itératifs, direct,
*                 méthode de projection ou couplée)
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 30/09/2011, version initiale
* HISTORIQUE : v1, 30/09/2011, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* Prière de PRENDRE LE TEMPS de compléter les commentaires
* en cas de modification de ce sous-programme afin de faciliter
* la maintenance !
************************************************************************
*
 'OPTION' 'DIME' 2 'MODE' 'PLAN' 'ISOV' 'SULI' 'ELEM' 'QUA4' ;
* 'OPTION' 'DIME' 3 'MODE' 'TRID' 'ISOV' 'SURF' 'ELEM' 'CUB8' ;
vdim = 'VALEUR' 'DIME' ;
vmod = 'VALEUR' 'MODE' ;
'SI' ('NON' interact) ;
   'OPTI' 'TRACER' 'PSC' ;
'SINON' ;
   'OPTI' 'TRACER' 'X' ;
'FINSI' ;
lcmp  = 'MOTS' 'LX1' 'LX2' 'LX3' 'LX4' 'LXP' 'LXPD' ;
opti inco lcmp lcmp ;
 
debug = faux ;
graphqdm = faux ;
graphsur = faux ;
lsauv = faux ;
*
*fic = 'CHAINE' '/test4/gounand/kong/'
*   'defila22dRe150Eu1FE1We1r4it100.sauv' ;
*'MESSAGE' ('CHAINE' 'Loading ' fic '...') ;
*'OPTI' 'RESTITUER' fic ;
*'RESTITUER' ;
*
************************************************************************
*
*
*   PROCEDURES
*
*
************************************************************************
************************************************************************
* NOM         : NLINP
* DESCRIPTION :   Gestion du parallélisme avec les assistants avant
*  l'appel à NLIN
*
'DEBPROC' NLINP ;
'ARGUMENT' gdisc*'MOT' ;
'ARGUMENT' _mt*'MAILLAGE' ;
'ARGUMENT' A*'TABLE' ;
'ARGUMENT' B*'TABLE' ;
'ARGUMENT' mot1*'MOT' ;
nmot = 1 ;
'ARGUMENT' mot2/'MOT' ;
'SI' ('EXISTE' mot2) ;
   nmot = '+' nmot 1 ;
   'ARGUMENT' mot3/'MOT' ;
   'SI' ('EXISTE' mot3) ;
      nmot = '+' nmot 1 ;
   'FINSI' ;
'FINSI' ;
*
nbpart = 'VALEUR' 'ASSI' ;
crit1 = '>' nbpart 1 ;
crit2 = '>' ('/' ('NBEL' _mt) nbpart) 1000 ;
crit3 =  vrai ;
* Appel systématique pour test du parallelisme avec ASSI
*'SI' (crit1 'ET' crit2 'ET' crit3) ;
'SI' vrai ;
   mapart = 'PART' 'ARLE' _mt nbpart ;
   'SI' ('EGA' nmot 1) ;
      tmat    = 'ASSI' 'TOUS' 'NLIN' gdisc mapart A B mot1 ;
   'FINSI' ;
   'SI' ('EGA' nmot 2) ;
      tmat    = 'ASSI' 'TOUS' 'NLIN' gdisc mapart A B mot1 mot2 ;
   'FINSI' ;
   'SI' ('EGA' nmot 3) ;
      tmat    = 'ASSI' 'TOUS' 'NLIN' gdisc mapart A B mot1
                mot2 mot3 ;
   'FINSI' ;
   mat = 'ETG' tmat ;
'SINON' ;
   'SI' ('EGA' nmot 1) ;
      mat = 'NLIN' gdisc _mt A B mot1 ;
   'FINSI' ;
   'SI' ('EGA' nmot 2) ;
      mat = 'NLIN' gdisc _mt A B mot1 mot2 ;
   'FINSI' ;
   'SI' ('EGA' nmot 3) ;
      mat = 'NLIN' gdisc _mt A B mot1 mot2 mot3 ;
   'FINSI' ;
'FINSI' ;
'RESPRO' mat ;
*
* End of procedure file NLINP
*
'FINPROC' ;
*
* Procédure pour résoudre un système linéaire
*
'DEBPROC' MONRES ;
'ARGUMENT' mat*'RIGIDITE' ;
'ARGUMENT' smb*'CHPOINT' ;
'ARGUMENT' tres/'TABLE' ;
'SI' ('EXISTE' tres) ;
   'SI' ('EGA' (tres . 'TYPINV') 0) ;
      'OPTI' impi 0 ;
      sol = 'RESOUD' mat smb 'NOID' ;
      'OPTI' impi 0 ;
   'SINON' ;
      'SI' ('EXISTE' tres 'LTIME') ;
         ltime = tres . 'LTIME' ;
      'SINON' ;
         ltime = FAUX ;
      'FINSI' ;
*
      'SI' ('EGA' ltime vrai) ;
         sol tt  = 'KRES' mat smb 'TYPI' tres ;
         'LISTE' tt ;
      'SINON' ;
         sol  = 'KRES' mat smb 'TYPI' tres ;
      'FINSI' ;
   'FINSI' ;
'SINON' ;
   sol = 'RESOUD' mat smb ;
'FINSI' ;
'RESPRO' sol ;
'FINPROC' ;
*
* Procédure pour afficher la résultante d'un chpoint et de sa valeur
* absolue (utile pour debugger sans lister le CHPOINT en question)
*
'DEBPROC' HCHPO ;
'ARGUMENT' chp*'CHPOINT' ;
'ARGUMENT' tit*'MOT' ;
*
achp = 'ABS' chp ;
nachp = 'MAXIMUM' achp ;
*'SI' ('<EG' nachp 1.D-60) ; nachp= 1.D0 ; 'FINSI' ;
nachp = 1.D0 ;
cc = '/' achp nachp ;
*hh = '+' cc 1. ;
hh = 'ABS' cc ;
rcc = 'RESULT' cc ;
rhh = 'RESULT' hh ;
'MESSAGE' ('CHAINE' 'Resultante de ' tit) ;
'LISTE' rcc ;
*'LISTE' rhh ;
'FINPROC' ;
*
* Procédure pour attendre un moment
*
'DEBPROC' attente ;
'ARGUMENT' s*'FLOTTANT' ;
'SI' interact ;
   'REPETER' i ('+' ('ENTIER' ('*' 100000 s)) 1) ;
   'FIN' i ;
'FINSI' ;
'FINPROC' ;
*BEGINPROCEDUR affvar
************************************************************************
* NOM         : AFFVAR
* DESCRIPTION : Affiche des variables
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
*
*
'DEBPROC' AFFVAR ;
'REPETER' bcl ;
   'ARGUMENT' x/'FLOTTANT' ;
   'SI' ('EXISTE' x) ;
      'ARGUMENT' lx*'MOT' ;
      'MESSAGE' ('CHAINE' lx '=' x) ;
   'SINON' ;
      'QUITTER' bcl ;
   'FINSI' ;
'FIN' bcl ;
'FINPROC' ;
*
* End of procedure file AFFVAR
*
*ENDPROCEDUR affvar
*BEGINPROCEDUR append
************************************************************************
* NOM         : APPEND
* DESCRIPTION : Rajoute :
*               - un entier à un listentier
*               - un réel   à un listreel
*               - un objet (liste, evolution, matrice ou chpoint)
*                 à un indice de table ('MOT' ou 'ENTIER')
*                    * si l'indice n'existe pas
*                    * 'ET' si l'indice existe
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 10/09/2004, version initiale
* HISTORIQUE : v1, 10/09/2004, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' APPEND ;
'ARGUMENT' tab/'TABLE' ;
'SI' ('EXISTE' tab) ;
   'ARGUMENT' itab/'MOT' ;
   'SI' ('NON' ('EXISTE' itab)) ;
      'ARGUMENT' itab*'ENTIER' ;
   'FINSI' ;
   lobj = FAUX ;
   'SI' ('NON' lobj) ;
      'ARGUMENT' lr/'LISTREEL' ;
      'SI' ('EXISTE' lr) ;
         obj = lr ; lobj = VRAI ;
      'FINSI' ;
   'FINSI' ;
   'SI' ('NON' lobj) ;
      'ARGUMENT' le/'LISTENTI' ;
      'SI' ('EXISTE' le) ;
         obj = le ; lobj = VRAI ;
      'FINSI' ;
   'FINSI' ;
   'SI' ('NON' lobj) ;
      'ARGUMENT' lev/'EVOLUTION' ;
      'SI' ('EXISTE' lev) ;
         obj = lev ; lobj = VRAI ;
      'FINSI' ;
   'FINSI' ;
   'SI' ('NON' lobj) ;
      'ARGUMENT' lm/'MAILLAGE' ;
      'SI' ('EXISTE' lm) ;
         obj = lm ; lobj = VRAI ;
      'FINSI' ;
   'FINSI' ;
   'SI' ('NON' lobj) ;
      'ARGUMENT' chpo/'CHPOINT' ;
      'SI' ('EXISTE' chpo) ;
         obj = chpo ; lobj = VRAI ;
      'FINSI' ;
   'FINSI' ;
   'SI' ('NON' lobj) ;
      'ARGUMENT' rig/'RIGIDITE' ;
      'SI' ('EXISTE' rig) ;
         obj = rig ; lobj = VRAI ;
      'FINSI' ;
   'FINSI' ;
   'SI' ('NON' lobj) ;
      'ARGUMENT' matk/'MATRIK' ;
      'SI' ('EXISTE' matk) ;
         obj = matk ; lobj = VRAI ;
      'FINSI' ;
   'FINSI' ;
   'SI' ('NON' lobj) ;
      cherr = 'CHAINE'
       'Il faut fournir un objet liste, evolution, matrice ou chpoint.'
        ;
      'ERREUR' cherr ;
   'FINSI' ;
   'SI' ('EXISTE' tab itab) ;
      'SI' ('EGA' ('TYPE' obj) 'CHPOINT') ;
         tab . itab = '+' (tab . itab) obj ;
      'SINON' ;
         tab . itab = 'ET' (tab . itab) obj ;
      'FINSI' ;
   'SINON' ;
      tab . itab = obj ;
   'FINSI' ;
   'RESPRO' tab ;
'FINSI' ;
'ARGUMENT' lenti/'LISTENTI' ;
'ARGUMENT' lreel/'LISTREEL' ;
'SI' ('EXISTE' lenti) ;
   'ARGUMENT' enti*'ENTIER' ;
   lenti = 'ET' lenti ('LECT' enti) ;
   'RESPRO' lenti ;
'FINSI' ;
'SI' ('EXISTE' lreel) ;
   'ARGUMENT' reel*'FLOTTANT' ;
   lreel = 'ET' lreel ('PROG' reel) ;
   'RESPRO' lreel ;
'FINSI' ;
*
* End of procedure file APPEND
*
'FINPROC' ;
*ENDPROCEDUR append
*BEGINPROCEDUR calimet
************************************************************************
* NOM         : CALIMET
* DESCRIPTION :
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, ??/??/2007, version initiale
* HISTORIQUE : v1, ??/??/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' CALIMET ;
'ARGUMENT' _mt*'MAILLAGE' ;
'ARGUMENT' gdisc*'MOT' ;
'ARGUMENT' methgau*'MOT' ;
*methgau = 'GAU7' ;
tmot = 'TABLE' ;
tnom = 'TABLE' ;
idim = 'VALEUR' 'DIME' ;
vdim  = DEADUTIL 'DIMM' _mt ;
idx = 0 ;
cim = 'CHAINE' 'IMET' ;
cg  = 'CHAINE' 'G' ;
'REPETER' iidim idim ;
   'REPETER' jidim idim ;
      'SI' ('>EG' &jidim &iidim) ;
         idx = '+' idx 1 ;
         tmot . idx = 'CHAINE' cim &iidim &jidim ;
* Convention Castem opposee à convention NLIN
         tnom . idx = 'CHAINE' cg &jidim &iidim ;
      'FINSI' ;
   'FIN' jidim ;
'FIN' iidim ;
*
lvid = 'LECT' ;
dtm = 'DIME' tmot ;
tchpo = 'TABLE' 'ESCLAVE' ;
idx = 0 ;
'REPETER' itm dtm ;
   mcm = tmot . &itm ;
   numop  = 1 ;
   numder = vdim ;
   numvar = 1 ;
   numdat = 0 ;
   numcof = 1 ;
   A = ININLIN numop numvar numdat numcof numder ;
   A . 'VAR' . 1 . 'NOMDDL' = 'MOTS' 'DUMM' ;
   A . 'VAR' . 1 . 'DISC'   = 'CSTE' ;
   A . 'VAR' . 1 . 'VALEUR' = 1. ;
   A . 'COF' . 1 . 'COMPOR' = mcm ;
   A . 'COF' . 1 . 'LDAT' = lvid ;
   A . 1 . 1 . 0 = 'LECT' 1 ;
   numdat = 0 ;
   numcof = 0 ;
   B = ININLIN numop numvar numdat numcof numder ;
   B . 'VAR' . 1 . 'NOMDDL' = 'MOTS' 'DUMM' ;
   B . 'VAR' . 1 . 'DISC'   = 'CSTE' ;
   B . 'VAR' . 1 . 'VALEUR' = 1. ;
   B . 1 . 1 . 0 = lvid ;
   cpo = NLINP gdisc _mt A B 'ERF1' methgau ;
   cpo = 'NOMC' (tnom . &itm) cpo ;
   idx = '+' idx 1 ;
   tchpo . idx = cpo ;
'FIN' itm ;
imet = 'ETG' tchpo ;
'RESPRO' imet ;
*
* End of procedure file CALIMET
*
'FINPROC' ;
*ENDPROCEDUR calimet
*BEGINPROCEDUR defdd
************************************************************************
* NOM         : DEFDD
* DESCRIPTION :
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, ??/??/2007, version initiale
* HISTORIQUE : v1, ??/??/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' DEFDD ;
*'ARGUMENT' _cmt*'MAILLAGE' ;
*'ARGUMENT' cmt*'MAILLAGE' ;
*'ARGUMENT' sur*'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
'ARGUMENT' idir/'ENTIER' ;
'SI' ('NON' ('EXISTE' idir)) ;
   idir = 0 ;
'FINSI' ;
*
NOMDEP = @STBL (TDISC . 'XN' . 'NOMPRI') ;
NOMFOR = @STBL (TDISC . 'XN' . 'NOMDUA') ;
*
vdim = 'VALEUR' 'DIME' ;
DISCG  = TDISC . 'GEOM' . 'DISC' ;
_hau = tdisc . 'hau' . 'QUAF' ;
hau = tdisc . 'hau' . discg ;
   vnor  = GNOR _hau tdisc 'NPRI' discg 'FPRI' 1. 'NDUA' 'XN' ;
   vnor = 'NOMC' nomfor nomdep vnor ;
*   'SI' ('EGA' vdim 3) ;
*      vnor = '*' vnor -1. ;
*   'FINSI' ;
   vnorn = '/' vnor ('**' ('PSCAL' vnor vnor nomdep nomdep) 0.5 ) ;
*   trvec hau vnorn 'Vnorn' ;
*  Correction de vnorn aux extrémités
   phau = 'CHANGER' 'POI1' hau ;
   'SI' ('EGA' idir 0) ;
      'SI' ('EGA' vdim 2) ;
         mcorr = ('POIN' hau 'INITIAL')
             'ET' ('POIN' hau 'FINAL') ;
      'SINON' ;
         bhau = tdisc . 'bhau' . discg ;
*         mcorr = 'CONTOUR' sur ;
*---------------pour 3D
         mcorr = bhau ;
*         'TRACER' (hau 'ET' ('COULEUR' bhau roug)) ;
      'FINSI' ;
   'SINON' ;
      mcorr = hau ;
   'FINSI' ;
   pmcorr = 'CHANGER' 'POI1' mcorr ;
   phaur = 'DIFF' phau pmcorr ;
   vnorn1 = 'REDU' vnorn phaur ;
   vvn = 'PROG' vdim * 0. ;
   'REMPLACER' vvn vdim -1. ;
**   vnorn2 = 'MANUEL' 'CHPO' mcorr 2 'UX' 0. 'UY' 1. ;
**   'LISTE' nomdep ;
**   'LISTE' vvn ;
   vnorn2 = 'MANUEL' 'CHPO' mcorr nomdep vvn ;
*   vnorn2 = 'MANUEL' 'CHPO' mcorr 1 ('EXTRAIRE' nomdep  vdim) 1. ;
   vnorn = vnorn1 '+' vnorn2 ;
*   trvec hau vnorn 'Vnorn2' ;
'RESPRO' vnorn ;
'FINPROC' ;
*
* End of procedure file DEFDD
*
*ENDPROCEDUR defdd
*BEGINPROCEDUR defmail
************************************************************************
* NOM         : DEFMAIL
* DESCRIPTION :
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 27/01/2011, version initiale
* HISTORIQUE : v1, 27/01/2011, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' DEFMAIL ;
'ARGUMENT' rsou*'FLOTTANT' ;
'ARGUMENT' lav*'FLOTTANT' ;
'ARGUMENT' lap*'FLOTTANT' ;
'ARGUMENT' prof*'FLOTTANT' ;
'ARGUMENT' raff*'ENTIER' ;
*
vdim = 'VALEUR' 'DIME' ;
vmod = 'VALEUR' 'MODE' ;
tdisc = MODENLIN 'GEOM' discg ;
'SI' ('EGA' vdim 2) ;
   'SI' ('EGA' vmod 'AXIS') ;
      tdisc = MODENLIN tdisc 'INCO' 'XN' discg
                                 'VECT' 'UR' 'UZ' 'FR' 'FZ'
                             'INCO' 'UN' discv
                                 'VECT' 'UR' 'UZ' 'FR' 'FZ' ;
   'SINON' ;
      tdisc = MODENLIN tdisc 'INCO' 'XN' discg
                                 'VECT' 'UX' 'UY' 'FX' 'FY'
                             'INCO' 'UN' discv
                                 'VECT' 'UX' 'UY' 'FX' 'FY' ;
   'FINSI' ;
'SINON' ;
      tdisc = MODENLIN tdisc 'INCO' 'XN' discg
                                 'VECT' 'UX' 'UY' 'UZ' 'FX' 'FY' 'FZ'
                             'INCO' 'UN' discv
                                 'VECT' 'UX' 'UY' 'UZ' 'FX' 'FY' 'FZ' ;
'FINSI' ;
'SI' ('EGA' discp 'LINM') ;
   'SI' ('EGA' vdim 2) ;
      nompre = 'MOTS' 'LX1' 'LX2' 'LX3' ;
   'FINSI' ;
   'SI' ('EGA' vdim 3) ;
      nompre = 'MOTS' 'LX1' 'LX2' 'LX3' 'LX4' ;
   'FINSI' ;
   tdisc = MODENLIN tdisc 'INCO' 'PN' discp 'SCAL' nompre nompre ;
'SINON' ;
   nompre = 'LXP' ;
   tdisc = MODENLIN tdisc 'INCO' 'PN' discp 'SCAL' nompre nompre ;
'FINSI' ;
*tdisc = JETMOD ;
raf2 = '*' 2 raff ;
raf4 = '*' 4 raff ;
raf8 = '*' 8 raff ;
raf16 = '*' 16 raff ;
raf32 = '*' 32 raff ;
*
dr = '/' rsou raff ;
dav = '/' lav raff ;
dap = '/' lap raff ;
drp = '/' prof raff ;
*
mlav = '*' lav -1. ;
mrsou = '*' rsou -1. ;
mprof = '*' prof -1. ;
'SI' ('EGA' vdim 2) ;
   pA = mlav mprof ; pB = mrsou mprof ; pC = 0. mprof ;
   pD = rsou mprof ; pE = lap mprof ;
*pJ = mlav 0. ; pII = mrsou 0. ; pH =  0. -0.5 ;
   pJ = mlav 0. ; pII = mrsou 0. ; pH =  0. 0. ;
   pG = rsou 0. ; pF = lap 0. ;
   'SI' ('EGA' raff -1) ; pF= lap ('*' mprof 0.5) ;  'FINSI' ;
*
   lAB = 'DROIT' pA pB 'DINI' dav 'DFIN' dr ;
   lBC = 'DROIT' pB pC 'DINI' dr 'DFIN' dr ;
   lCD = 'DROIT' pC pD 'DINI' dr 'DFIN' dr ;
   lDE = 'DROIT' pD pE 'DINI' dr 'DFIN' dap ;
   lEF = 'DROIT' pE pF 'DINI' drp 'DFIN' dr ;
   lFG = 'DROIT' pF pG 'DINI' dap 'DFIN' dr ;
   lGH = 'DROIT' pG pH 'DINI' dr 'DFIN' dr ;
   lHII = 'DROIT' pH pII 'DINI' dr 'DFIN' dr ;
   lIIJ = 'DROIT' pII pJ 'DINI' dr 'DFIN' dav ;
   lJA = 'DROIT' pJ pA 'DINI' dr 'DFIN' drp ;
*
   bas = lAB 'ET' lBC 'ET' lCD 'ET' lDE ;
   dro = lEF ;
   hau = LFG 'ET' lGH 'ET' lHII 'ET' lIIJ ;
   gau = lJA ;
*
   'SI' ('EGA' raff -1) ;
      lAE = 'DROIT' 1 pA pE ;
      lEF = 'DROIT' 1 pE pF ;
      lFJ = 'DROIT' 1 pF pJ ;
      lJA = 'DROIT' 1 pJ pA ;
      bas = lAE ;
      dro = lEF ;
      hau = LFJ ;
      gau = lJA ;
   'FINSI';
   mtw = 'DALLER' bas dro hau gau ;
   cmtw = 'CONTOUR' mtw ;
*'TRACER' mtw ;
'FINSI' ;
'SI' ('EGA' vdim 3) ;
   pA = mlav 0. mprof ; pB = mrsou 0. mprof ; pC = 0. 0. mprof ;
   pD = rsou 0. mprof ; pE = lap 0. mprof ;
*pJ = mlav 0. ; pII = mrsou 0. ; pH =  0. -0.5 ;
   pJ = mlav 0. 0. ; pII = mrsou 0. 0. ; pH =  0. 0. 0. ;
   pG = rsou 0. 0. ; pF = lap 0. 0. ;
*
   lAB = 'DROIT' pA pB 'DINI' dav 'DFIN' dr ;
   lBC = 'DROIT' pB pC 'DINI' dr 'DFIN' dr ;
   lCD = 'DROIT' pC pD 'DINI' dr 'DFIN' dr ;
   lDE = 'DROIT' pD pE 'DINI' dr 'DFIN' dap ;
   lEF = 'DROIT' pE pF 'DINI' drp 'DFIN' dr ;
   lFG = 'DROIT' pF pG 'DINI' dap 'DFIN' dr ;
   lGH = 'DROIT' pG pH 'DINI' dr 'DFIN' dr ;
   lHII = 'DROIT' pH pII 'DINI' dr 'DFIN' dr ;
   lIIJ = 'DROIT' pII pJ 'DINI' dr 'DFIN' dav ;
   lJA = 'DROIT' pJ pA 'DINI' dr 'DFIN' drp ;
*
   bas2 = lAB 'ET' lBC 'ET' lCD 'ET' lDE ;
   dro2 = lEF ;
   hau2 = LFG 'ET' lGH 'ET' lHII 'ET' lIIJ ;
   gau2 = lJA ;
   fro  = 'DALLER' bas2 dro2 hau2 gau2 ;
   dhau = 'INVERSE' hau2 ;
*   cmtw2 = 'CONTOUR' mtw2 ;
*
   vtran = 0. lav 0. ;
*   vtran = 0. ('*' rsou 2.) 0. ;
   pH2   = 'PLUS' pH (0. rsou 0.) ;
   pH3   = 'PLUS' pH vtran ;
   lgen1 = 'DROIT' pH  pH2 'DINI' dr 'DFIN' dr  ;
   lgen2 = 'DROIT' pH2 pH3 'DINI' dr 'DFIN' dav ;
   lgen  = lgen1 'ET' lgen2 ;
*   lgen  = 'DROIT' 1 pH pH3 ;
*
   rea = 'PLUS' fro vtran ;
   bas = 'GENERATRICE' bas2 lgen ;
   dro = 'GENERATRICE' dro2 lgen ;
   hau = 'GENERATRICE' hau2 lgen ;
   gau = 'GENERATRICE' gau2 lgen ;
   mtw = 'VOLUME' fro 'GENE' lgen ;
   cmtw = 'ENVELOPPE' mtw ;
   pint = 0. rsou ('/' mprof 2.) ;
   cmtw = 'ORIENTER' cmtw 'POIN' pint ;
   bhau1 = lgen 'PLUS' pJ ;
   bhau2 = lgen 'PLUS' pF ;
   bhau  = bhau1 'ET' bhau2 ;
*
*   'TRACER' mtw ;
*
'FINSI' ;
tdisc = QUAFME tdisc bas 'bas' dro 'dro' hau 'hau' gau 'gau'
    cmtw 'cmtw' mtw 'mtw' ('*' 1.D-5 rsou) ;
'SI' ('EGA' vdim  3) ;
   tdisc = QUAFME tdisc fro 'fro' rea 'rea' bhau 'bhau' dhau 'dhau'
       ('*' 1.D-5 rsou) ;
'FINSI' ;
*
*   vnor  = GNOR _cmtw tdisc 'NPRI' discg 'FPRI' 1. 'NDUA' 'XN' ;
*   NOMDEP = @STBL (TDISC . 'XN' . 'NOMINC') ;
*   vnorn = '/' vnor ('**' ('PSCAL' vnor vnor nomdep nomdep) 0.5 ) ;
*   trvec cmtw vnorn 'Vnorn' ;
*
'RESPRO' tdisc ;
*
* End of procedure file DEFMAIL
*
'FINPROC' ;
*ENDPROCEDUR defmail
*BEGINPROCEDUR defvit
************************************************************************
* NOM         : DEFVIT
* DESCRIPTION :
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, ??/??/2007, version initiale
* HISTORIQUE : v1, ??/??/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' DEFVIT ;
*'ARGUMENT' vit*'CHPOINT' ;
'ARGUMENT' lnclk/'LOGIQUE' ;
*
'SI' ('NON' ('EXISTE' lnclk)) ;
   lnclk = FAUX ;
'FINSI' ;
   vdim = 'VALEUR' 'DIME' ;
   'SI' ('EGA' vdim 2) ; rvit = vit ;
   'SINON' ; rvit = 'REDU' vit cmtw ; 'FINSI' ;
   TRVEC mtw nomvit rvit 'Vit' 'VECT' 'NORM' 'VNOR' ;
   'SI' ('EGA' vdim 2) ;
      mtvit = tdisc . 'mtw' . discv ;
      blcou = 'BLOQ' 'T' (tdisc . 'bas' . discv) ;
      fcou = FCOURANT mtvit vit blcou ;
      'TRACER' mtw ('*' fcou -1.) cmtw 'TITR' 'Fcourant' ;
   'FINSI' ;
   prep = GLNO _mtw tdisc 'PN' pre discv ;
   'TRACER' 'CACH' mtw prep cmtw 20 'TITR' ('CHAINE' 'Pression') ;
*      dvit = 'EXCO' nomdep dsol ;
*      TRVEC2 mtw nomvit  rvit 'Vit normee' ;
*      TRVEC cmtw nomvit rvit 'Vit' ;
*      TRVEC mtrr ('REDU' vit mtrr) 'Vit' ;
*      'TRACER' mtrr ('REDU' prep mtrr) cmtrr
*         'TITR' ('CHAINE' 'Pression') ;
*      JETPROF1 vit ;
*      JETPROF1 vit ('*' haut 0.4) ;
*
* End of procedure file DEFVIT
*
'FINPROC' ;
*ENDPROCEDUR defvit
*BEGINPROCEDUR defvsurf
************************************************************************
* NOM         : DEFVSURF
* DESCRIPTION :
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, ??/??/2007, version initiale
* HISTORIQUE : v1, ??/??/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' DEFVSURF ;
'ARGUMENT' rmax/'FLOTTANT' ;
'ARGUMENT' lnclk/'LOGIQUE' ;
*
'SI' ('NON' ('EXISTE' lnclk)) ;
   lnclk = FAUX ;
'FINSI' ;
vdim = 'VALEUR' 'DIME' ;
'SI' ('EGA' vdim 2) ;
      ihau = 'INVERSE' hau ;
'SINON' ;
      ihau = dhau ;
'FINSI' ;
      i=0 ;
      tabev = 'TABLE' ;
      tabt  = 'TABLE' ;
      i = '+' i 1 ;
      rhau = 'EXTRAIRE' ('EVOL' 'CHPO' ('COORDONNEE' 1 hau) 'SCAL' ihau)
         'ORDO' ;
      hv = 'REDU' vit ihau ;
NOMDEP = @STBL (TDISC . 'XN' . 'NOMPRI') ;
      nhv = '**' ('PSCAL' vit vit nomdep nomdep) 0.5  ;
      zhau = 'EXTRAIRE' ('EVOL' 'CHPO' nhv 'SCAL' ihau)   'ORDO' ;
      tabev . i = 'EVOL' 'MANU' rhau zhau ;
      tabt  . i = 'CHAINE' '|V| surf' ;
      tix = 's' ; tiy = '|V|' ; tit = 'CHAINE' tiy '(' tix ')' ;
      'SI' ('EXISTE' rmax) ;
         binf = '-' 0. rmax ;
         bsup = '+' 0. rmax ;
         dessevol (@STBL tabev) tabt tit tix tiy
            ('PROG' binf bsup) lnclk ;
      'SINON' ;
         dessevol (@STBL tabev) tabt tit tix tiy
            lnclk ;
      'FINSI' ;
*'SINON' ;
*   'SI' ('EXISTE' rmax) ;
*      xhau = 'COORDONNEE' 1 hau ;
*      yhau = 'COORDONNEE' 2 hau ;
*      rhau = '**' ('+' ('**' xhau 2) ('**' yhau 2)) 0.5 ;
*      phau = 'POIN' rhau 'INFERIEUR' rmax ;
*      redhau = 'ELEM' hau 'APPUYE' 'LARGEMENT' phau ;
*   'SINON' ;
*      redhau = hau ;
*   'FINSI' ;
*   'SI' lnclk ;
*      'TRACER' 'CACH' redhau 'TITR' 'Surface' 'NCLK' ;
*   'SINON' ;
*      'TRACER' 'CACH' redhau 'TITR' 'Surface' 'NCLK' ;
*   'FINSI' ;
*'FINSI' ;
*
* End of procedure file DEFVSURF
*
'FINPROC' ;
*ENDPROCEDUR defvsurf
*BEGINPROCEDUR dessevol
************************************************************************
* NOM         : DESSEVOL
* DESCRIPTION : Dessine des évolutions : choisit automatiquement
*               les options, marqueurs, couleurs...
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 16/11/2004, version initiale
* HISTORIQUE : v1, 16/11/2004, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' DESSEVOL ;
'ARGUMENT' evtot*'EVOLUTION' ;
'ARGUMENT' tabt*'TABLE' ;
'ARGUMENT' tit*'MOT' ;
'ARGUMENT' tix*'MOT' ;
'ARGUMENT' tiy*'MOT' ;
'ARGUMENT' lnclk/'LOGIQUE' ;
'ARGUMENT' nb/'ENTIER' ;
'ARGUMENT' lx/'LISTREEL' ;
*
'SI' ('NON' ('EXISTE' lnclk)) ;
   lnclk = FAUX ;
'FINSI' ;
*
* nb = 0 : noir et blanc
* nb = 1 : couleur
* nb = 2 : couleur + marqueurs
* nb = 3 : couleur + marqueurs + tirets
* nb = 4 : couleur + marqueurs regu
* nb = 5 : couleur + marqueurs regu + tirets
* nb = 6 : nb + marqueurs regu + tirets
* nb = 7 : nb + marqueurs
* nb = 8 : nb + marqueurs regu
*
'SI' ('NON' ('EXISTE' nb)) ;
   nb = 3 ;
'FINSI' ;
*
nt  = 'DIME' tabt ;
nev = 'DIME' evtot ;
*
* Attention, dans evtot, il y a une évolution avec des noms de points ?
*
*'SI' ('NEG' nev nt) ;
*   cherr = 'CHAINE' 'Evolution and title table : not same dim.' ;
*   'ERREUR' cherr ;
*'FINSI' ;
*
tev = 'TABLE' ;
tev . 'TITRE' = tabt ;
*
toto = 'TABLE' ;
*
*lcoul = 'MOTS' 'TURQ' 'VERT' 'JAUN' 'ROSE' 'ROUG' 'BLEU' ;
'SI' ('EGA' ('VALEUR' 'TRAC') 'PSC') ;
   lcoul = 'MOTS' 'BLEU' 'ROUG' 'VERT' 'VIOL' 'OLIV' 'ORAN'  ;
'SINON' ;
   lcoul = 'MOTS' 'TURQ' 'VERT' 'JAUN' 'ROSE' 'ROUG' 'BLEU' ;
'FINSI' ;
lmarq = 'MOTS' 'TRID' 'TRIU' 'LOSA' 'CARR' 'ETOI' 'PLUS' 'CROI'
               'TRIL' 'TRIR' ;
ltirr = 'MOTS' 'TIRR' 'TIRC' 'TIRL' 'TIRM' ;
*
lnb = 'LECT' 0 6 7 8 ;
*'SI' ('OU' ('EGA' nb 0) ('EGA' nb 6)) ;
'SI' (dans ('LECT' nb) lnb) ;
   ev2 = evtot ;
'SINON' ;
   icou = 0 ;
   'REPETER' iev nev ;
      ii = &iev ;
      evi = 'EXTRAIRE' evtot 'COUR' ii ;
      'SI' ('NEG' ('TYPE' ('EXTRAIRE' evi 'ORDO')) 'LISTMOTS') ;
         icou = '+' icou 1 ;
      'FINSI' ;
*      ii2 = '/' ('+' ii 1) 2 ;
*      ci  = EXMOMOD lcoul ii2 ;
*      ci  = EXMOMOD lcoul ii ;
      ci  = EXMOMOD lcoul icou ;
      APPEND toto 'EVOLUTION' ('COULEUR' evi ci) ;
   'FIN' iev ;
   ev2 = toto . 'EVOLUTION' ;
'FINSI' ;
*
'REPETER' iev nev ;
   ii = &iev ;
   mi = EXMOMOD lmarq ii ;
   ti = EXMOMOD ltirr ii ;
   'SI' ('OU' ('EGA' nb 2) ('EGA' nb 7)) ;
      tev . ii = 'CHAINE' 'MARQ' ' ' mi ;
   'FINSI' ;
   'SI' ('EGA' nb 3) ;
      tev . ii = 'CHAINE' 'MARQ' ' ' mi ' ' ti ;
   'FINSI' ;
   'SI' ('OU' ('EGA' nb 4) ('EGA' nb 8)) ;
      tev . ii = 'CHAINE' 'MARQ' ' ' mi ' REGU' ;
   'FINSI' ;
   'SI' ('OU' ('EGA' nb 5) ('EGA' nb 6)) ;
      tev . ii = 'CHAINE' 'MARQ' ' ' mi ' ' ti ' REGU' ;
   'FINSI' ;
'FIN' iev ;
*
'SI' ('EXISTE' lx) ;
   dim4 = 'EGA' ('DIME' lx) 4 ;
   xmin = 'EXTRAIRE' lx 1 ; xmax = 'EXTRAIRE' lx 2 ;
   'SI' dim4 ;
      ymin = 'EXTRAIRE' lx 3 ; ymax = 'EXTRAIRE' lx 4 ;
   'FINSI' ;
   'SI' dim4 ;
   'SI' lnclk ;
      'DESSIN' ev2 'TITR' tit 'TITX' tix 'TITY' tiy 'LEGE' 'MIMA' tev
      'XBOR' xmin xmax 'YBOR' ymin ymax 'NCLK' motopt ;
   'SINON' ;
      'DESSIN' ev2 'TITR' tit 'TITX' tix 'TITY' tiy 'LEGE' 'MIMA' tev
      'XBOR' xmin xmax 'YBOR' ymin ymax motopt ;
   'FINSI' ;
   'SINON' ;
   'SI' lnclk ;
      'DESSIN' ev2 'TITR' tit 'TITX' tix 'TITY' tiy 'LEGE' 'MIMA' tev
      'XBOR' xmin xmax 'NCLK' motopt ;
   'SINON' ;
      'DESSIN' ev2 'TITR' tit 'TITX' tix 'TITY' tiy 'LEGE' 'MIMA' tev
      'XBOR' xmin xmax motopt ;
   'FINSI' ;
   'FINSI' ;
'SINON' ;
   'SI' lnclk ;
      'DESSIN' ev2 'TITR' tit 'TITX' tix 'TITY' tiy 'LEGE' 'MIMA' tev
      'NCLK' motopt ;
   'SINON' ;
      'DESSIN' ev2 'TITR' tit 'TITX' tix 'TITY' tiy 'LEGE' 'MIMA' tev
      motopt ;
   'FINSI' ;
'FINSI' ;
*
* End of procedure file DESSEVOL
*
'FINPROC' ;
*ENDPROCEDUR dessevol
*BEGINPROCEDUR dessvec
************************************************************************
* NOM         : DESSVEC
* DESCRIPTION : Trace des bilans de forces projetés sur une direction
*               sur une ligne.
*               Utile pour tracer des bilans de forces
*
*               Display vector fields.
*               Useful for visualization of force balance.
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 26/09/2011, version initiale
* HISTORIQUE : v1, 26/09/2011, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
*
*
'DEBPROC' DESSVEC ;
'ARGUMENT' tdom*'MAILLAGE' ;
*'ARGUMENT' tdomq*'MAILLAGE' ;
'ARGUMENT' xdom*'CHPOINT'  ;
'ARGUMENT' tix*'MOT' ;
'ARGUMENT' tiy*'MOT' ;
'ARGUMENT' nocomp/'LISTMOTS' ;
'SI' ('EXISTE' nocomp) ;
   'ARGUMENT' cdir*'CHPOINT' ;
'FINSI' ;
*
'ARGUMENT' tvec/'TABLE' ;
'SI' ('EXISTE' tvec) ;
   'ARGUMENT' ttit*'TABLE' ;
   i = '-' ('DIME' tvec) 1 ;
'SINON' ;
   tvec = 'TABLE' ; ttit = 'TABLE' ;
   i = 0 ;
'FINSI' ;
*
'REPETER' livec ;
   'SI' ('EGA' i 0) ;
      'ARGUMENT' ccvec*'CHPOINT' ;
   'SINON' ;
      'ARGUMENT' ccvec/'CHPOINT' ;
   'FINSI' ;
   'SI' ('EXISTE' ccvec) ;
      'ARGUMENT' ttvec*'MOT' ;
   'SINON' ;
      'QUITTER' livec ;
   'FINSI' ;
   i = '+' i 1 ;
*   'MESSAGE' ('CHAINE' 'i=' i) ;
*   'LISTE' ccvec ;
*   'LISTE' tvec ;
   tvec . i = ccvec ;
   ttit . i = ttvec ;
'FIN' livec ;
*'ARGUMENT' echv/'FLOTTANT' ;
*mmas = gmass2 tdomq tdisc 'NPRI' discv 'NDUA' discv ;
'ARGUMENT' mmas/'RIGIDITE' ;
'ARGUMENT' lrmax/'LISTREEL' ;
'ARGUMENT' lnclk/'LOGIQUE' ;
'SI' ('NON' ('EXISTE' lnclk)) ;
   lnclk = faux ;
'FINSI' ;
*
tdim = 'DIME' tvec ;
'SI' ('EXISTE' tvec 'SOUSTYPE') ; tdim = '-' tdim 1 ; 'FINSI' ;
'SI' ('EXISTE' tvec 'CREATEUR') ; tdim = '-' tdim 1 ; 'FINSI' ;
tvecr = 'TABLE' ;
'REPETER' idim tdim ;
   tvecr . &idim = 'REDU' (tvec . &idim) tdom ;
'FIN' idim ;
*tvecr = 'ASSI' 'TOUS' 'REDU' tvec tdom ;
*tdim = '-' ('DIME' tvecr) 1 ;
'SI' ('EXISTE' nocomp) ;
   tvecp = 'TABLE' ;
   'REPETER' idim tdim ;
      tvecp . &idim = 'PSCAL' (tvec . &idim) cdir nocomp nomdep ;
   'FIN' idim ;
'SINON' ;
   tvecp = tvecr ;
'FINSI' ;
'SI' ('EXISTE' mmas) ;
   'REPETER' idim tdim ;
      tvecp . &idim = 'RESOUD' mmas (tvecp . &idim) ;
   'FIN' idim ;
'FINSI' ;
*
lx = 'EXTRAIRE' ('EVOL' 'CHPO' xdom tdom) 'ORDO' ;
'REPETER' idim tdim ;
   tvecp . &idim = 'EVOL' 'CHPO' (tvecp . &idim) tdom ;
'FIN' idim ;
'REPETER' idim tdim ;
   tvecp . &idim = 'EXTRAIRE' (tvecp . &idim) 'ORDO' ;
'FIN' idim ;
'REPETER' idim tdim ;
   tvecp . &idim = 'EVOL' 'MANUEL' lx (tvecp . &idim) ;
'FIN' idim ;
*tvecp = 'ASSI' 'TOUS' 'EVOL' 'CHPO' tvecp tdom ;
*tvecp = 'ASSI' 'TOUS' 'EXTRAIRE' tvecp 'ORDO' ;
*tvecp = 'ASSI' 'TOUS' 'EVOL' 'MANU' lx tvecp ;
tit = 'CHAINE' tiy '(' tix ')' ;
'SI' ('EXISTE' lrmax) ;
   dessevol (@STBL tvecp) ttit tit tix tiy
      lrmax lnclk ;
'SINON' ;
   dessevol (@STBL tvecp) ttit tit tix tiy
      lnclk ;
'FINSI' ;
*
* End of procedure file DESSVEC
*
'FINPROC' ;
*ENDPROCEDUR dessvec
*BEGINPROCEDUR errrel
************************************************************************
* NOM         : ERRREL
* DESCRIPTION : Calcul d'une erreur relative
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 23/04/2003, version initiale
* HISTORIQUE : v1, 23/04/2003, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' ERRREL ;
'ARGUMENT' val*'FLOTTANT' ;
'ARGUMENT' valref*'FLOTTANT' ;
*
'SI' ('<' ('ABS' valref) 1.D-10) ;
   echref = 1.D0 ;
'SINON' ;
   echref = valref ;
'FINSI' ;
*
errabs = 'ABS' ('/' ('-' val valref) echref);
*
'RESPRO' errabs ;
*
* End of procedure file ERRREL
*
'FINPROC' ;
*ENDPROCEDUR errrel
*BEGINPROCEDUR exmomod
************************************************************************
* NOM         : EXMOMOD
* DESCRIPTION : Extraction d'un mot d'un listmots
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 23/06/2003, version initiale
* HISTORIQUE : v1, 23/06/2003, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' EXMOMOD ;
'ARGUMENT' lm*'LISTMOTS' i*'ENTIER' ;
j = 'DIME' lm ;
k = '+' (MODULO ('-' i 1) j) 1 ;
lemot = 'EXTRAIRE' lm k ;
* Usage de l'opérateur text pour éviter que lemot
* ne soit interprété comme un opérateur
'RESPRO' 'TEXTE' lemot ;
*
* End of procedure file EXMOMOD
*
'FINPROC' ;
*ENDPROCEDUR exmomod
*BEGINPROCEDUR formar
************************************************************************
* NOM         : FORMAR
* DESCRIPTION : formate un réel de facon courte
*               pratique pour les noms de
*               sauvegarde
*               Exemples :
* 'MESSAGE' ('CHAINE' (formar 2.9e5   1)) ;
* 2.9E5
* 'MESSAGE' ('CHAINE' (formar -2.9e5  1)) ;
* -2.9E5
* 'MESSAGE' ('CHAINE' (formar 2.9e-5  1)) ;
* 2.9E-5
* 'MESSAGE' ('CHAINE' (formar -2.9e-5 1)) ;
* -2.9E-5
* 'MESSAGE' ('CHAINE' (formar 2.9     1)) ;
* 2.9
* 'MESSAGE' ('CHAINE' (formar -2.9    1)) ;
* -2.9
* 'MESSAGE' ('CHAINE' (formar 0       1)) ;
* 0
* 'MESSAGE' ('CHAINE' (formar 0       1)) ;
* 0
* 'MESSAGE' ('CHAINE' (formar 2.9e5   0)) ;
* 3E5
* 'MESSAGE' ('CHAINE' (formar -2.9e5  0)) ;
* -3E5
* 'MESSAGE' ('CHAINE' (formar 2.9e-5  0)) ;
* 3E-5
* 'MESSAGE' ('CHAINE' (formar -2.9e-5 0)) ;
* -3E-5
* 'MESSAGE' ('CHAINE' (formar 2.9     0)) ;
* 3
* 'MESSAGE' ('CHAINE' (formar -2.9    0)) ;
* -3
* 'MESSAGE' ('CHAINE' (formar 0       0)) ;
* 0
* 'MESSAGE' ('CHAINE' (formar 0       0)) ;
* 0
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 18/02/2003, version initiale
* HISTORIQUE : v1, 18/02/2003, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' FORMAR ;
'ARGUMENT' fl*'FLOTTANT' ;
'ARGUMENT' vir/'ENTIER  ' ;
'SI' ('NON' ('EXISTE' vir)) ;
   vir = 1 ;
'SINON' ;
   'SI' ('<' vir 0) ;
      'ERREUR' 'fournir un entier positif' ;
   'FINSI' ;
'FINSI' ;
'SI' ('<' ('ABS' fl) 10.D-100) ;
   chfl = 'CHAINE' '0' ;
'SINON' ;
*! sans le 1.D-10, ca ne fonctionne pas
*! qd on entre pile poil une puissance de 10
   lfl = LOG10 ('ABS' fl) ;
*   lfl = '+' (LOG10 ('ABS' fl)) 1.D-10  ;
   slfl = 'SIGNE' ('ENTIER' lfl) ;
   'SI' ('EGA' slfl 1) ;
      elfl = 'ENTIER' lfl ;
   'SINON' ;
      elfl = '-' ('ENTIER' lfl) 1 ;
   'FINSI' ;
   man = '/' fl ('**' 10.D0 elfl) ;
*
* Une verrue pour des histoires de précision...
*
   'SI' ('EGA' man 10.D0 ('**' 10.D0 ('*' vir -1.D0))) ;
      man  = '/' man 10.D0 ;
      elfl = '+' elfl 1 ;
   'FINSI' ;
*
   sman = 'SIGNE' man ;
   'SI' ('EGA' sman 1) ;
      fman = 'CHAINE' '(F' ('+' vir 2) '.0' vir ')' ;
   'SINON' ;
      fman = 'CHAINE' '(F' ('+' vir 3) '.0' vir ')' ;
   'FINSI' ;
   'SI' ('NEG' vir 0) ;
      'SI' ('NEG' elfl 0) ;
         chfl = 'CHAINE' 'FORMAT' fman man 'E' elfl ;
      'SINON' ;
         chfl = 'CHAINE' 'FORMAT' fman man ;
      'FINSI' ;
   'SINON' ;
      man2 = 'ENTIER' ('+' man ('*' 0.5D0 sman)) ;
      'SI' ('NEG' elfl 0) ;
         chfl = 'CHAINE' man2 'E' elfl ;
      'SINON' ;
         chfl = 'CHAINE' man2 ;
      'FINSI' ;
   'FINSI' ;
'FINSI' ;
'RESPRO' chfl ;
*
* End of procedure file FORMAR
*
'FINPROC' ;
*ENDPROCEDUR formar
*BEGINPROCEDUR gchpo
************************************************************************
* NOM         : GCHPO
* DESCRIPTION : Une matrice de masse
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v2, 14/03/2006, mise à jour NLIN évolué
* VERSION    : v1, 13/05/2004, version initiale
* HISTORIQUE : v1, 13/05/2004, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' GCHPO ;
'ARGUMENT'     _mt*'MAILLAGE' ;
'ARGUMENT'    nomq*'MOT     ' ;
'ARGUMENT'    coef*'FLOTTANT' ;
*
'SI' ('EGA' discq 'LINM') ;
   discq = 'CSTE' ;
'FINSI' ;
*
idim = DEADUTIL 'DIMM' _mt ;
vdim = 'VALEUR' 'DIME' ;
*
discg = TDISC . 'GEOM' . 'DISC' ;
tnomq = TDISC . nomq ;
*
ms = tnomq . 'NOMPRI' . 1 ;
numop = 1 ; numder = idim ; numvar = 1 ; numdat = 0 ; numcof = 0 ;
A = ININLIN numop numvar numdat numcof numder ;
A . 'VAR' . 1 . 'NOMDDL' = 'MOTS' 'DUMM' ;
A . 'VAR' . 1 . 'DISC'   = 'CSTE' ;
A . 'VAR' . 1 . 'VALEUR' = 0. ;
*
A . 1 . 1 . 0 = 'LECT' ;
*
B = ININLIN numop numvar numdat numcof numder ;
B . 'VAR' . 1 . 'NOMDDL' = ms ;
B . 'VAR' . 1 . 'DISC'   = tnomq . 'DISC' ;
*
B . 1 . 1 . 0 = 'LECT' ;
*
mgchpos = '+' (NLIN discg _mt A B 'ERF1' 'GAU7') coef ;
*
mgchpo = 'NOMC' ms (tnomq . 'NOMPRI' . 1) mgchpos ;
*
ninc = 'DIME' (tnomq . 'NOMPRI') ;
'SI' ('>' ninc 1) ;
   'REPETER' iinc ('-' ninc 1) ;
      mgchpo = '+' mgchpo
         ('NOMC' ms (tnomq . 'NOMPRI' . ('+' &iinc 1)) mgchpos) ;
   'FIN' iinc ;
'FINSI' ;
*
'RESPRO' mgchpo ;
'FINPROC' ;
*
* End of procedure file GCHPO
*
*ENDPROCEDUR gchpo
*BEGINPROCEDUR gdiv2
************************************************************************
* NOM         : GDIV2
* DESCRIPTION : Une matrice de masse
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v2, 14/03/2006, mise à jour NLIN évolué
* VERSION    : v1, 13/05/2004, version initiale
* HISTORIQUE : v1, 13/05/2004, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' GDIV2 ;
'ARGUMENT'     _mt*'MAILLAGE' ;
'ARGUMENT'     _smt/'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
*
* Lectures
*
debug = FAUX ;
lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
                 'NCOF' 'FCOF' 'CCOF' 'GBBT' 'GMBT' ;
* Il faut initialiser valt et valq, sinon on peut capturer ceux de
* la procédure appelante
valt = 'valt' ; valq = 'valq' ;
lbbt = 0 ;
*
'REPETER' imotcle ;
   'ARGUMENT' motcle/'MOT' ;
   'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
   'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
      cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
   'FINSI' ;
   'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'FPRI') ; 'ARGUMENT' valt*'LISTREEL' ; 'FINSI' ;
   'SI' ('EGA' motcle 'FDUA') ; 'ARGUMENT' valq*'FLOTTANT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'FCOF') ; 'ARGUMENT' valo*'FLOTTANT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'GBBT') ; lbbt = 1 ; 'FINSI' ;
   'SI' ('EGA' motcle 'GMBT') ; lbbt = 2 ; 'FINSI' ;
'FIN' imotcle ;
*
* Tests
*
discg = TDISC . 'GEOM' . 'DISC' ;
'SI' ('EXISTE' tdisc 'methgau') ;
   methgau = tdisc . 'methgau' . 'amor' ;
'SINON' ;
   methgau = 'GAU7' ;
'FINSI' ;
tnomt = TDISC . nomt ;
lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
tnomq = TDISC . nomq ;
lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
*
lcof = 'EXISTE' TDISC nomo ;
'SI' lcof ; ncof = 1 ; tcof = TDISC . nomo ;
'SINON' ; ncof = 0 ;
'FINSI' ;
*
'SI' debug ;
'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
'SINON' ; 'MESSAGE' 'pas de coef  detecte' ;
'FINSI' ;
'FINSI' ;
*
   vdim = 'VALEUR' 'DIME' ;
   vmod = 'VALEUR' 'MODE' ;
   idim = 0 ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
      idim = 2 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
      idim = 2 ;
      iaxi = VRAI ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
      idim = 3 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' vdim 1) ;
      idim = 1 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' idim 0) ;
      'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
   'FINSI' ;
   'SI' iaxi ;
      dp = ('*' PI 2.D0) ;
      rmt = 'COORDONNEE' 1 _mt ;
      ncof = ncof '+' 2 ;
   'FINSI' ;
* Scalaire ou vecteur
ninct = 'DIME' (tnomt . 'NOMPRI') ;
nincq = 'DIME' (tnomq . 'NOMPRI') ;
'SI' ('NEG' ninct idim) ;
   cherr = 'CHAINE'
   'la primale doit etre un vecteur' ;
   'ERREUR' cherr ;
'FINSI' ;
'SI' ('NEG' nincq 1) ;
   cherr = 'CHAINE'
   'la duale doit etre un scalaire' ;
   'ERREUR' cherr ;
'FINSI' ;
*
   numop = 1 ; numder = idim ; numvar = ninct ;
   numdat = ncof ; numcof = ncof ;
   A = ININLIN numop numvar numdat numcof numder ;
*
   lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
   'REPETER' iiinct ninct ;
      iinct = &iiinct ;
      A . 'VAR' . iinct . 'NOMDDL' = tnomt . 'NOMPRI' . iinct ;
      A . 'VAR' . iinct . 'DISC'   = tnomt . 'DISC'  ;
      'SI' lvalt ;
         'SI' lvt ;
            A . 'VAR' . iinct . 'VALEUR'   = 'EXTRAIRE' valt iinct ;
         'SINON' ;
            A . 'VAR' . iinct . 'VALEUR'   = valt ;
         'FINSI' ;
      'FINSI' ;
   'FIN' iiinct ;
*
   icof = 0 ;
   'SI' lcof ;
      icof = '+' icof 1 ;
      A . 'DAT' . icof . 'NOMDDL' = tcof . 'NOMPRI' . 1 ;
      A . 'DAT' . icof . 'DISC'   = tcof . 'DISC' ;
      A . 'DAT' . icof . 'VALEUR' = valo ;
      A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
      A . 'COF' . icof . 'LDAT'   = 'LECT' icof ;
      ll = 'LECT' 1 ;
   'SINON' ;
      ll = 'LECT' ;
   'FINSI' ;
*
   'SI' iaxi ;
      icof = '+' icof 1 ;
      A . 'DAT' . icof . 'NOMDDL' = 'MOTS' 'SCAL' ;
      A . 'DAT' . icof . 'DISC'   = 'CSTE' ;
      A . 'DAT' . icof . 'VALEUR' = dp ;
      A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
      A . 'COF' . icof . 'LDAT'   = 'LECT' icof ;
      icof = '+' icof 1 ;
      A . 'DAT' . icof . 'NOMDDL' = 'MOTS' 'SCAL' ;
      A . 'DAT' . icof . 'DISC'   = discg ;
      A . 'DAT' . icof . 'VALEUR' = rmt ;
      A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
      A . 'COF' . icof . 'LDAT'   = 'LECT' icof ;
      lldpr = ll 'ET' ('LECT' ('-' icof 1) icof) ;
      lldp = ll 'ET' ('LECT' ('-' icof 1)) ;
   'FINSI' ;
*
   'SI' iaxi ;
      'REPETER' iidim idim ;
         A . 1 . &iidim . &iidim  = lldpr ;
      'FIN' iidim ;
      A . 1 . 1 . 0 = lldp ;
   'SINON' ;
      'REPETER' iidim idim ;
         A . 1 . &iidim . &iidim  = ll ;
      'FIN' iidim ;
   'FINSI' ;
*
   numvar = 1 ;
   numdat = 0 ;
   numcof = 0 ;
*
   B = ININLIN numop numvar numdat numcof numder ;
   B . 'VAR' . 1 . 'NOMDDL' = tnomq . 'NOMDUA' . 1 ;
   B . 'VAR' . 1 . 'DISC'   = tnomq . 'DISC'  ;
   'SI' lvalq ;
      B . 'VAR' . 1 . 'VALEUR'   = valq ;
   'FINSI' ;
   B . 1 . 1 . 0  = 'LECT' ;
*
   'SI' ('OU' ('EGA' lbbt 0) ('EGA' lbbt 1)) ;
      'SI' ('EXISTE' _smt) ;
          mgdiv2 = 'NLIN' discg _mt _smt A B methgau ;
      'SINON' ;
          mgdiv2 = NLINP discg _mt A B methgau ;
      'FINSI' ;
   'FINSI' ;
   'SI' ('OU' ('EGA' lbbt 1) ('EGA' lbbt 2)) ;
      B . 'VAR' . 1 . 'NOMDDL' = tnomq . 'NOMPRI' . 1 ;
      'REPETER' iiinct ninct ;
         iinct = &iiinct ;
         A . 'VAR' . iinct . 'NOMDDL' = tnomt . 'NOMDUA' . iinct ;
      'FIN' iiinct ;
      'SI' ('EXISTE' _smt) ;
          mgdiv3 = 'NLIN' discg _mt _smt B A methgau ;
      'SINON' ;
          mgdiv3 = NLINP discg _mt B A methgau ;
      'FINSI' ;
   'FINSI' ;
   'SI' ('EGA' lbbt 0) ;
      mgdiv = mgdiv2 ;
   'FINSI' ;
   'SI' ('EGA' lbbt 1) ;
      mgdiv = mgdiv2 'ET' mgdiv3 ;
   'FINSI' ;
   'SI' ('EGA' lbbt 2) ;
      mgdiv = mgdiv3 ;
   'FINSI' ;
'RESPRO' mgdiv ;
'FINPROC' ;
*
* End of procedure file GDIV2
*
*ENDPROCEDUR gdiv2
*BEGINPROCEDUR getcoo
************************************************************************
* NOM         : GETCOO
* DESCRIPTION :
* Renvoie les coordonnées des points dans un champ type déplacement
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 22/04/2011, version initiale
* HISTORIQUE : v1, 22/04/2011, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
*
*
'DEBPROC' GETCOO ;
'ARGUMENT' mail*'MAILLAGE' ;
'ARGUMENT' incop*'LISTMOTS' ;
*
dim = 'VALEUR' 'DIME' ;
'REPETER' iidim dim ;
   idim= &iidim ;
   icoo = 'NOMC' ('EXTRAIRE' incop idim)
                 ('COORDONNEE' idim mail) ;
   'SI' ('EGA' idim 1) ;
      vcoo = icoo ;
   'SINON' ;
      vcoo = 'ET' vcoo icoo ;
   'FINSI' ;
'FIN' iidim ;
'RESPRO' vcoo ;
*
* End of procedure file GETCOO
*
'FINPROC' ;
*ENDPROCEDUR getcoo
*BEGINPROCEDUR gforc
************************************************************************
* NOM         : GFORC
* DESCRIPTION : Calcul de la force associée à une pression imposée
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, ??/??/2007, version initiale
* HISTORIQUE : v1, ??/??/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' GFORC ;
'ARGUMENT' _surf*'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
'ARGUMENT' pfor*'CHPOINT' ;
*
vdim = 'VALEUR' 'DIME' ;
DISCG  = TDISC . 'GEOM' . 'DISC' ;
fpfor = GNOR _surf tdisc 'NPRI' discg 'CPRI' pfor 'NDUA' 'XN' ;
'RESPRO' fpfor ;
*
* End of procedure file GFORC
*
'FINPROC' ;
*ENDPROCEDUR gforc
*BEGINPROCEDUR ggravi
************************************************************************
* NOM         : GGRAVI
* DESCRIPTION : Calcul de la force associée au potentiel gravitaire
*               (\rho g z si g vertical)
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 22/04/2011
* HISTORIQUE : v1, 22/04/2011, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
*
*
'DEBPROC' GGRAVI ;
'ARGUMENT' _surf*'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
'ARGUMENT' coef*'FLOTTANT' ;
'ARGUMENT' ang*'FLOTTANT' ;
*
vdim = 'VALEUR' 'DIME' ;
pgrax  = '*' ('COORDONNEE' 1 _surf) ('*' +1. ('SIN' ang)) ;
pgraz  = '*' ('COORDONNEE' vdim _surf) ('*' -1. ('COS' ang)) ;
DISCG  = TDISC . 'GEOM' . 'DISC' ;
fpgrax = GNOR _surf tdisc 'NPRI' discg 'CPRI' pgrax 'NDUA' 'XN' ;
fpgraz = GNOR _surf tdisc 'NPRI' discg 'CPRI' pgraz 'NDUA' 'XN' ;
fpgra = '+' fpgrax fpgraz ;
fpgra = '*' fpgra ('*' -1. coef) ;
'RESPRO' fpgra ;
*
* End of procedure file GGRAVI
*
'FINPROC' ;
*ENDPROCEDUR ggravi
*BEGINPROCEDUR gkforc
************************************************************************
* NOM         : GKFORC
* DESCRIPTION : Calcul de la matrice tangente associée à la force
*               de pression imposée. Cette matrice tangente est
*               partielle car elle ne prend en compte que le gradient
*               surfacique de pression. On peut difficilement
*               faire autrement si on ne considère que la surface.
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, ??/??/2007, version initiale
* HISTORIQUE : v1, ??/??/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' GKFORC ;
'ARGUMENT' _surf*'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
'ARGUMENT' pfor*'CHPOINT' ;
'ARGUMENT' ijaco/'ENTIER' ;
'SI' ('NON' ('EXISTE' ijaco)) ;
   ijaco = 0 ;
'FINSI' ;
vdim = 'VALEUR' 'DIME' ;
DISCG  = TDISC . 'GEOM' . 'DISC' ;
*fpfor = GNOR _surf tdisc 'NPRI' discg 'CPRI' pfor 'NDUA' 'XN' ;
   k1 = GNORGC _surf tdisc 'NPRI' 'XN'
                             'NCOF' discg 'CCOF' pfor
                             'NDUA' 'XN' ;
*   k1 = GNOR _surf tdisc 'NPRI' discg 'NDUA' 'XN' ;
*   k1 = '*' k1 -1. ;
*   k1 = 'CHANGER' 'INCO' k1 ('MOTS' 'SCAL')
*            ('MOTS' ('EXTRAIRE' NOMVIT vdim)) NOMVIT NOMVIT ;
   k2 = GNORKTAN _surf tdisc 'NPRI' 'XN'
        'NCOF' discg 'CCOF' pfor 'NDUA' 'XN' ;
   'SI' ('OU' ('EGA' ijaco 0) ('EGA' ijaco 3)) ;
      ktfor = k1 'ET' k2 ;
   'FINSI' ;
   'SI' ('EGA' ijaco 1) ;
      ktfor = k1 ;
   'FINSI' ;
   'SI' ('EGA' ijaco 2) ;
      ktfor = k2 ;
   'FINSI' ;
   ktfor = '*' ktfor -1. ;
'RESPRO' ktfor ;
*
* End of procedure file GKFORC
*
'FINPROC' ;
*ENDPROCEDUR gkforc
*BEGINPROCEDUR gkgravi
************************************************************************
* NOM         : GKGRAVI
* DESCRIPTION : Calcul de la matrice tangente de la force
*               associée au potentiel gravitaire (calculée par GGRAVI)
*               en fonction des déplacements des points de la surface.
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 22/04/2011
* HISTORIQUE : v1, 22/04/2011, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
*
*
'DEBPROC' GKGRAVI ;
'ARGUMENT' _surf*'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
'ARGUMENT' ijaco*'ENTIER' ;
*'SI' ('NON' ('EXISTE' ijaco)) ;
*   ijaco = 0 ;
*'FINSI' ;
'ARGUMENT' coef*'FLOTTANT' ;
'ARGUMENT' ang*'FLOTTANT' ;
*
vdim = 'VALEUR' 'DIME' ;
pgrax  = '*' ('COORDONNEE' 1 _surf) ('*' +1. ('SIN' ang)) ;
pgraz  = '*' ('COORDONNEE' vdim _surf) ('*' -1. ('COS' ang)) ;
*pgra  = '*' ('-' ('COORDONNEE' vdim _surf) H) -1. ;
*pgra  = '*' ('COORDONNEE' vdim _surf) -1. ;
DISCG  = TDISC . 'GEOM' . 'DISC' ;
NOMDEP = @STBL (TDISC . 'XN' . 'NOMPRI') ;
NOMFOR = @STBL (TDISC . 'XN' . 'NOMDUA') ;
*fpgra = GNOR _surf tdisc 'NPRI' discg 'CPRI' pgra 'NDUA' 'XN' ;
   k1x = GNOR _surf tdisc 'NPRI' discg 'NDUA' 'XN' ;
   k1x = '*' k1x ('*' +1. ('SIN' ang)) ;
   k1x = 'CHANGER' 'INCO' k1x ('MOTS' 'SCAL')
            ('MOTS' ('EXTRAIRE' NOMDEP 1)) NOMFOR NOMFOR ;
   k2x = GNORKTAN _surf tdisc 'NPRI' 'XN'
        'NCOF' discg 'CCOF' pgrax 'NDUA' 'XN' ;
   k1z = GNOR _surf tdisc 'NPRI' discg 'NDUA' 'XN' ;
   k1z = '*' k1z ('*' -1. ('COS' ang)) ;
   k1z = 'CHANGER' 'INCO' k1z ('MOTS' 'SCAL')
            ('MOTS' ('EXTRAIRE' NOMDEP vdim)) NOMFOR NOMFOR ;
   k2z = GNORKTAN _surf tdisc 'NPRI' 'XN'
        'NCOF' discg 'CCOF' pgraz 'NDUA' 'XN' ;
   'SI' ('EGA' ijaco 0) ;
      ktgra = k1x 'ET' k1z 'ET' k2x 'ET' k2z ;
   'FINSI' ;
   'SI' ('EGA' ijaco 1) ;
      ktgra = k1x 'ET' k1z ;
   'FINSI' ;
   'SI' ('EGA' ijaco 2) ;
      ktgra = k2x 'ET' k2z ;
   'FINSI' ;
   ktgra = '*' ktgra coef ;
'RESPRO' ktgra ;
*
* End of procedure file GKGRAVI
*
'FINPROC' ;
*ENDPROCEDUR gkgravi
*BEGINPROCEDUR gkvol
************************************************************************
* NOM         : GKVOL
* DESCRIPTION : Matrice tangente associée à la variation du volume
*               contenu dans une surface (calculé par GVOL)
*               en fonction des déplacements des points de la surface.
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 22/04/2011, version initiale
* HISTORIQUE : v1, 22/04/2011, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
*
*
'DEBPROC' GKVOL ;
'ARGUMENT' _surf*'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
'ARGUMENT' ijaco/'ENTIER' ;
'SI' ('NON' ('EXISTE' ijaco)) ;
   ijaco = 0 ;
'FINSI' ;
*  Vecteur position et calcul du volume
NOMVIT = @STBL (TDISC . 'XN' . 'NOMPRI') ;
DISCG  = TDISC . 'GEOM' . 'DISC' ;
vdim = 'VALEUR' 'DIME' ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
      fdim = 3 ;
   'SINON' ;
      fdim = vdim ;
   'FINSI' ;
   vpos = GETCOO _surf nomvit ;
   kvol1 = GNOR _surf tdisc 'NPRI' ('CHAINE' discg 'V')
                             'NDUA' 'XN' 'FDUA' ('PROG' vdim * 1.) ;
   kvol2 = GNORKTAN _surf tdisc 'NPRI' ('CHAINE' discg 'V')
            'NCOF' ('CHAINE' discg 'V') 'CCOF' vpos
            'NDUA' 'XN' 'FDUA' ('PROG' vdim * 1.) ;
   'SI' ('EGA' ijaco 0) ;
      kvol = '/' ('+' kvol1 kvol2) fdim ;
   'FINSI' ;
   'SI' ('EGA' ijaco 1) ;
      kvol = '/' kvol1 fdim ;
   'FINSI' ;
   'SI' ('EGA' ijaco 2) ;
      kvol = '/' kvol2 fdim ;
   'FINSI' ;
'RESPRO' kvol ;
*
* End of procedure file GKVOL
*
'FINPROC' ;
*ENDPROCEDUR gkvol
*BEGINPROCEDUR glno
************************************************************************
* NOM         : GLNO
* DESCRIPTION : Une matrice de masse
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v2, 14/03/2006, mise à jour NLIN évolué
* VERSION    : v1, 13/05/2004, version initiale
* HISTORIQUE : v1, 13/05/2004, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' GLNO ;
'ARGUMENT'     _mt*'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
*
'ARGUMENT' nomt*'MOT' ;
'ARGUMENT' valt*'CHPOINT' ;
'ARGUMENT' nomq/'MOT' ;
*
'SI' ('NON' ('EXISTE' nomq)) ;
   nomq = 'QUAF' ;
'FINSI' ;
*
* Second membre
*
smb = GMASS2 _mt tdisc 'NPRI' nomt 'CPRI' valt 'NDUA' nomq ;
mat = GMASS2 _mt tdisc 'NPRI' nomq 'NDUA' nomq ;
cpro = 'KRES' mat smb ;
*
'RESPRO' cpro ;
'FINPROC' ;
*
* End of procedure file GLNO
*
*ENDPROCEDUR glno
*BEGINPROCEDUR gmass2
************************************************************************
* NOM         : GMASS2
* DESCRIPTION : Une matrice de masse
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v2, 14/03/2006, mise à jour NLIN évolué
* VERSION    : v1, 13/05/2004, version initiale
* HISTORIQUE : v1, 13/05/2004, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' GMASS2 ;
'ARGUMENT'     _mt*'MAILLAGE' ;
'ARGUMENT'     _smt/'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
*
* Lectures
*
debug = FAUX ;
lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
                 'NCOF' 'FCOF' 'CCOF' ;
* Il faut initialiser valt et valq, sinon on peut capturer ceux de
* la procédure appelante
valt = 'valt' ; valq = 'valq' ;
'REPETER' imotcle ;
   'ARGUMENT' motcle/'MOT' ;
   'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
   'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
      cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
   'FINSI' ;
   'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
   tst1 = 'EGA' motcle 'FPRI' ; tst2 = 'EGA' motcle 'FDUA' ;
   tst = tst1 'OU' tst2 ;
   'SI' tst ;
      'SI' tst1 ; tt = TDISC . nomt ; 'FINSI' ;
      'SI' tst2 ; tt = TDISC . nomq ; 'FINSI' ;
      isvec = ('>' ('DIME' (tt. 'NOMPRI')) 1) ;
      'SI' isvec ; 'ARGUMENT' valv*'LISTREEL' ; 'SINON' ;
         'ARGUMENT' valv*'FLOTTANT' ;
      'FINSI' ;
      'SI' tst1 ; valt = valv ; 'FINSI' ;
      'SI' tst2 ; valq = valv ; 'FINSI' ;
   'FINSI' ;
   'SI' ('EGA' motcle 'FCOF') ; 'ARGUMENT' valo*'FLOTTANT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
'FIN' imotcle ;
*
* Tests
*
discg = TDISC . 'GEOM' . 'DISC' ;
'SI' ('EXISTE' tdisc 'methgau') ;
   methgau = tdisc . 'methgau' . 'mass' ;
'SINON' ;
   methgau = 'GAU7' ;
'FINSI' ;
tnomt = TDISC . nomt ;
lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
tnomq = TDISC . nomq ;
lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
* Scalaire ou vecteur
ninct = 'DIME' (tnomt . 'NOMPRI') ;
nincq = 'DIME' (tnomq . 'NOMPRI') ;
'SI' ('NEG' ninct nincq) ;
   cherr = 'CHAINE'
   'les primales et duales nont pas le meme nombre de composantes' ;
   'ERREUR' cherr ;
'FINSI' ;
ninc = ninct ;
*
lcof = 'EXISTE' TDISC nomo ;
'SI' lcof ; ncof = 1 ; tcof = TDISC . nomo ;
'SINON' ; ncof = 0 ;
'FINSI' ;
*
'SI' debug ;
'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
'SINON' ; 'MESSAGE' 'pas de coef  detecte' ;
'FINSI' ;
'FINSI' ;
*
   vdim = 'VALEUR' 'DIME' ;
   vmod = 'VALEUR' 'MODE' ;
   idim = 0 ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
      idim = 2 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
      idim = 2 ;
      iaxi = VRAI ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
      idim = 3 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' vdim 1) ;
      idim = 1 ;
      iaxi = FAUX ;
   'FINSI' ;
*   'MESSAGE' ('CHAINE' 'iaxi=' iaxi );
   'SI' ('EGA' idim 0) ;
      'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
   'FINSI' ;
   'SI' iaxi ;
      dprmt = '*' ('COORDONNEE' 1 _mt) ('*' PI 2.D0) ;
   'FINSI' ;
*
* Optimisation possible : construire la matrice par blocs
* qd valt et valq ne sont pas donnés
*
   numop = ninc ; numder = idim ; numvar = ninc ;
   numdat = ncof ; numcof = ncof ;
   A = ININLIN numop numvar numdat numcof numder ;
   'SI' lcof ;
      A . 'DAT' . 1 . 'NOMDDL' = tcof . 'NOMPRI' . 1 ;
      A . 'DAT' . 1 . 'DISC'   = tcof . 'DISC' ;
      A . 'DAT' . 1 . 'VALEUR' = valo ;
      A . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
      A . 'COF' . 1 . 'LDAT'   = 'LECT' 1 ;
      ll = 'LECT' 1 ;
   'SINON' ;
      ll = 'LECT' ;
   'FINSI' ;
   lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
   'REPETER' iiinc ninc ;
      iinc = &iiinc ;
      A . 'VAR' . iinc . 'NOMDDL' = tnomt . 'NOMPRI' . iinc ;
      A . 'VAR' . iinc . 'DISC'   = tnomt . 'DISC'  ;
      'SI' lvalt ;
         'SI' lvt ;
            A . 'VAR' . iinc . 'VALEUR'   = 'EXTRAIRE' valt iinc ;
         'SINON' ;
            A . 'VAR' . iinc . 'VALEUR'   = valt ;
         'FINSI' ;
      'FINSI' ;
      A . iinc . iinc . 0  = ll ;
   'FIN' iiinc ;
*
   'SI' iaxi ;
      numdat = 1 ;
      numcof = 1 ;
   'SINON' ;
      numdat = 0 ;
      numcof = 0 ;
   'FINSI' ;
   B = ININLIN numop numvar numdat numcof numder ;
   'SI' iaxi ;
      B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
      B . 'DAT' . 1 . 'DISC'   = discg ;
      B . 'DAT' . 1 . 'VALEUR' = dprmt ;
      B . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
      B . 'COF' . 1 . 'LDAT'   = 'LECT' 1 ;
      ll = 'LECT' 1 ;
   'SINON' ;
      ll = 'LECT' ;
   'FINSI' ;
   lvq = 'EGA' ('TYPE' valq) 'LISTREEL' ;
   'REPETER' iiinc ninc ;
      iinc = &iiinc ;
      B . 'VAR' . iinc . 'NOMDDL' = tnomq . 'NOMDUA' . iinc ;
      B . 'VAR' . iinc . 'DISC'   = tnomq . 'DISC'  ;
      'SI' lvalq ;
         'SI' lvq ;
            B . 'VAR' . iinc . 'VALEUR'   = 'EXTRAIRE' valq iinc ;
         'SINON' ;
            B . 'VAR' . iinc . 'VALEUR'   = valq ;
         'FINSI' ;
      'FINSI' ;
      B . iinc . iinc . 0  = ll ;
   'FIN' iiinc ;
*
   'SI' ('EXISTE' _smt) ;
       mgmass2 = 'NLIN' discg _mt _smt A B methgau ;
   'SINON' ;
       mgmass2 = NLINP discg _mt A B methgau ;
   'FINSI' ;
*
'RESPRO' mgmass2 ;
'FINPROC' ;
*
* End of procedure file GMASS2
*
*ENDPROCEDUR gmass2
*BEGINPROCEDUR gnorgc
************************************************************************
* NOM         : GNORGC
* DESCRIPTION : Une matrice de masse
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v2, 14/03/2006, mise à jour NLIN évolué
* VERSION    : v1, 13/05/2004, version initiale
* HISTORIQUE : v1, 13/05/2004, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' GNORGC ;
'ARGUMENT'     _mt*'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
*
* Lectures
*
vdim = 'VALEUR' 'DIME' ;
debug = FAUX ;
lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
                 'NCOF' 'FCOF' 'CCOF' ;
* Il faut initialiser valt et valq, sinon on peut capturer ceux de
* la procédure appelante
valt = 'valt' ; valq = 'valq' ;
'REPETER' imotcle ;
   'ARGUMENT' motcle/'MOT' ;
   'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
   'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
      cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
   'FINSI' ;
   'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
   tst1 = 'EGA' motcle 'FPRI' ; tst2 = 'EGA' motcle 'FDUA' ;
   tst3 = 'EGA' motcle 'FCOF' ;
   tst = tst1 'OU' tst2 'OU' tst3 ;
   'SI' tst ;
      'SI' tst1 ; tt = TDISC . nomt ; 'FINSI' ;
      'SI' tst2 ; tt = TDISC . nomq ; 'FINSI' ;
      'SI' tst3 ; tt = TDISC . nomo ; 'FINSI' ;
      isvec = ('>' ('DIME' (tt. 'NOMPRI')) 1) ;
      'SI' isvec ; 'ARGUMENT' valv*'LISTREEL' ; 'SINON' ;
         'ARGUMENT' valv*'FLOTTANT' ;
      'FINSI' ;
      'SI' tst1 ; valt = valv ; 'FINSI' ;
      'SI' tst2 ; valq = valv ; 'FINSI' ;
      'SI' tst3 ; valo = valv ; 'FINSI' ;
   'FINSI' ;
   'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
'FIN' imotcle ;
*
* Tests
*
discg = TDISC . 'GEOM' . 'DISC' ;
'SI' ('EXISTE' tdisc 'methgau') ;
   methgau = tdisc . 'methgau' . 'mass' ;
'SINON' ;
   methgau = 'GAU7' ;
'FINSI' ;
tnomt = TDISC . nomt ;
lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
tnomq = TDISC . nomq ;
lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
* Scalaire ou vecteur
ninct = 'DIME' (tnomt . 'NOMPRI') ;
nincq = 'DIME' (tnomq . 'NOMPRI') ;
'SI' ('NEG' ninct vdim) ;
   cherr = 'CHAINE'
   'la primale doit etre un vecteur' ;
   'ERREUR' cherr ;
'FINSI' ;
'SI' ('NEG' nincq vdim) ;
   cherr = 'CHAINE'
   'la duale doit etre un vecteur' ;
   'ERREUR' cherr ;
'FINSI' ;
*ninc = ninct ;
*
lcof = 'EXISTE' TDISC nomo ;
'SI' lcof ; tcof = TDISC . nomo ;
ncof = 'DIME' (tcof . 'NOMPRI') ;
'SINON' ; ncof = 0 ;
'FINSI' ;
'SI' ('NEG' ncof 1) ;
   cherr = 'CHAINE'
   'il faut un coefficient scalaire' ;
   'ERREUR' cherr ;
'FINSI' ;
*
'SI' debug ;
'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
'SINON' ; 'MESSAGE' 'pas de coef  detecte' ;
'FINSI' ;
'FINSI' ;
*
   vdim = 'VALEUR' 'DIME' ;
   vmod = 'VALEUR' 'MODE' ;
   idim = 0 ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
      idim = 2 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
      idim = 2 ;
      iaxi = VRAI ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
      idim = 3 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' vdim 1) ;
      idim = 1 ;
      iaxi = FAUX ;
   'FINSI' ;
*   'MESSAGE' ('CHAINE' 'iaxi=' iaxi );
   'SI' ('EGA' idim 0) ;
      'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
   'FINSI' ;
   'SI' iaxi ;
      dprmt = '*' ('COORDONNEE' 1 _mt) ('*' PI 2.D0) ;
   'FINSI' ;
*
* Optimisation possible : construire la matrice par blocs
* qd valt et valq ne sont pas donnés
*
   numop = ninct '*' nincq ; numder = idim ; numvar = ninct ;
   numdat = 1 ; numcof = idim ;
   A = ININLIN numop numvar numdat numcof numder ;
   A . 'DAT' . 1 . 'NOMDDL' = tcof . 'NOMPRI' . 1 ;
   A . 'DAT' . 1 . 'DISC'   = tcof . 'DISC' ;
   A . 'DAT' . 1 . 'VALEUR' = valo ;
   'REPETER' iicof numcof ;
      icof = &iicof ;
      A . 'COF' . icof . 'COMPOR' = 'CHAINE' 'D/DX' icof ;
      A . 'COF' . icof . 'LDAT'   = 'LECT' 1 ;
   'FIN' iicof ;
   lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
   iop = 0 ;
   'REPETER' iiinct ninct ;
      iinct = &iiinct ;
      A . 'VAR' . iinct . 'NOMDDL' = tnomt . 'NOMPRI' . iinct ;
      A . 'VAR' . iinct . 'DISC'   = tnomt . 'DISC'  ;
      'SI' lvalt ;
         'SI' lvt ;
            A . 'VAR' . iinct . 'VALEUR'   = 'EXTRAIRE' valt iinct ;
         'SINON' ;
            A . 'VAR' . iinct . 'VALEUR'   = valt ;
         'FINSI' ;
      'FINSI' ;
      ll = 'LECT' iinct ;
      'REPETER' iiincq nincq ;
         iop = '+' iop 1 ;
         A . iop . iinct . 0  = 'LECT' ;
*         A . iop . iinct . 0  = ll ;
      'FIN' iiincq ;
   'FIN' iiinct ;
*
   'SI' iaxi ;
      numdat = 1 ;
      numcof = idim '+' 1 ;
   'SINON' ;
      numdat = 0 ;
      numcof = idim ;
   'FINSI' ;
   numvar = nincq ;
   B = ININLIN numop numvar numdat numcof numder ;
   icof = 0 ;
   'REPETER' iiidim idim ;
      icof = '+' icof 1 ;
      B . 'COF' . icof . 'COMPOR' = 'CHAINE' 'VNOR' &iiidim ;
      B . 'COF' . icof . 'LDAT'   = 'LECT' ;
   'FIN' iiidim ;
*
   'SI' iaxi ;
      icof = '+' icof 1 ;
      B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
      B . 'DAT' . 1 . 'DISC'   = discg ;
      B . 'DAT' . 1 . 'VALEUR' = dprmt ;
      B . 'COF' . icof . 'COMPOR' = 'IDEN' ;
      B . 'COF' . icof . 'LDAT'   = 'LECT' 1 ;
      ll = 'LECT' icof ;
   'SINON' ;
      ll = 'LECT' ;
   'FINSI' ;
   lvq = 'EGA' ('TYPE' valq) 'LISTREEL' ;
   'REPETER' iiincq nincq ;
      iincq = &iiincq ;
      B . 'VAR' . iincq . 'NOMDDL' = tnomq . 'NOMDUA' . iincq ;
      B . 'VAR' . iincq . 'DISC'   = tnomq . 'DISC'  ;
      'SI' lvalq ;
         'SI' lvq ;
            B . 'VAR' . iincq . 'VALEUR'   = 'EXTRAIRE' valq iincq ;
         'SINON' ;
            B . 'VAR' . iincq . 'VALEUR'   = valq ;
         'FINSI' ;
      'FINSI' ;
   'FIN' iiincq ;
   iop = 0 ;
   'REPETER' iiinct ninct ;
*      iinct = &iiinct ;
      'REPETER' iiincq nincq ;
         iincq = &iiincq ;
         iop = '+' iop 1 ;
         B . iop . iincq . 0  = ('LECT' iincq) 'ET' ll ;
      'FIN' iiincq ;
   'FIN' iiinct ;
*
   mgnorgc = NLIN discg _mt A B methgau ;
*
'RESPRO' mgnorgc ;
'FINPROC' ;
*
* End of procedure file GNORGC
*
*ENDPROCEDUR gnorgc
*BEGINPROCEDUR gnorktan
************************************************************************
* NOM         : GNORKTAN
* DESCRIPTION : Matrice tangente associée à la variation de la normale
*               à une surface (calculée par GNOR)
*               en fonction des déplacements des points de la surface.
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 22/04/2011, version initiale
* HISTORIQUE : v1, 22/04/2011, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
*
*
'DEBPROC' GNORKTAN ;
'ARGUMENT'     _mt*'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
*
* Lectures
*
dim = 'VALEUR' 'DIME' ;
mdim = DEADUTIL 'DIMM' _mt ;
'SI' ('NEG' mdim ('-' dim 1)) ;
   'ERREUR' 'Dim. maillage .neq. dim. espace - 1' ;
'FINSI' ;
loi  = 'CHAINE' 'VNOJ' ;
debug = FAUX ;
lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
                 'NCOF' 'FCOF' 'CCOF' ;
* Il faut initialiser valt et valq, sinon on peut capturer ceux de
* la procédure appelante
valt = 'valt' ; valq = 'valq' ;
'REPETER' imotcle ;
   'ARGUMENT' motcle/'MOT' ;
   'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
   'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
      cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
   'FINSI' ;
   'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
   tst1 = 'EGA' motcle 'FPRI' ; tst2 = 'EGA' motcle 'FDUA' ;
   tst = tst1 'OU' tst2 ;
   'SI' tst ;
      'SI' tst1 ; tt = TDISC . nomt ; 'FINSI' ;
      'SI' tst2 ; tt = TDISC . nomq ; 'FINSI' ;
      isvec = ('>' ('DIME' (tt. 'NOMPRI')) 1) ;
      'SI' isvec ; 'ARGUMENT' valv*'LISTREEL' ; 'SINON' ;
         'ARGUMENT' valv*'FLOTTANT' ;
      'FINSI' ;
      'SI' tst1 ; valt = valv ; 'FINSI' ;
      'SI' tst2 ; valq = valv ; 'FINSI' ;
   'FINSI' ;
   'SI' ('EGA' motcle 'FCOF') ; 'ARGUMENT' valo*'FLOTTANT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
'FIN' imotcle ;
*
* Tests
*
discg = TDISC . 'GEOM' . 'DISC' ;
'SI' ('EXISTE' tdisc 'methgau') ;
   methgau = tdisc . 'methgau' . 'mass' ;
'SINON' ;
   methgau = 'GAU7' ;
'FINSI' ;
tnomt = TDISC . nomt ;
lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
tnomq = TDISC . nomq ;
lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
* Scalaire ou vecteur
ninct = 'DIME' (tnomt . 'NOMPRI') ;
nincq = 'DIME' (tnomq . 'NOMPRI') ;
'SI' ('NEG' ninct dim) ;
   cherr = 'CHAINE'
   'la primale doit etre un vecteur' ;
   'ERREUR' cherr ;
'FINSI' ;
'SI' ('NEG' nincq dim) ;
   cherr = 'CHAINE'
   'la duale doit etre un vecteur' ;
   'ERREUR' cherr ;
'FINSI' ;
ninc = dim ;
*
lcof = 'EXISTE' TDISC nomo ;
'SI' lcof ; tcof = TDISC . nomo ;
ncof = 'DIME' (tcof . 'NOMPRI') ;
'SINON' ; ncof = 0 ;
'FINSI' ;
*
'SI' debug ;
'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
'SINON' ; 'MESSAGE' 'pas de coef  detecte' ;
'FINSI' ;
'FINSI' ;
*
   vdim = 'VALEUR' 'DIME' ;
   vmod = 'VALEUR' 'MODE' ;
   idim = 0 ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
      idim = 2 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
      idim = 2 ;
      iaxi = VRAI ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
      idim = 3 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' vdim 1) ;
      idim = 1 ;
      iaxi = FAUX ;
   'FINSI' ;
*   'MESSAGE' ('CHAINE' 'iaxi=' iaxi );
   'SI' ('EGA' idim 0) ;
      'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
   'FINSI' ;
   'SI' iaxi ;
      deupi = '*' PI 2.D0 ;
      dprmt = '*' ('COORDONNEE' 1 _mt) deupi ;
   'FINSI' ;
*
* Optimisation possible : construire la matrice par blocs
* qd valt et valq ne sont pas donnés
*
   numop = idim '*' idim '*' idim ;
   'SI' iaxi ;
      numop = numop '+' idim ;
   'FINSI' ;
   numder = idim ; numvar = ninct ;
   numdat = ncof ; numcof = ncof ;
   A = ININLIN numop numvar numdat numcof numder ;
   'SI' lcof ;
      lvo = 'EGA' ('TYPE' valo) 'LISTREEL' ;
      'REPETER' iicof ncof ;
         icof = &iicof ;
         A . 'DAT' . icof . 'NOMDDL' = tcof . 'NOMPRI' . icof ;
         A . 'DAT' . icof . 'DISC'   = tcof . 'DISC' ;
         'SI' lvo ;
            A . 'DAT' . icof . 'VALEUR' = 'EXTRAIRE' valo icof ;
         'SINON' ;
            A . 'DAT' . icof . 'VALEUR' = valo ;
         'FINSI' ;
         A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
         A . 'COF' . icof . 'LDAT'   = 'LECT' icof ;
      'FIN' iicof ;
   'SINON' ;
      ll = 'LECT' ;
   'FINSI' ;
   lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
   iop = 0 ;
   'REPETER' iiinct ninct ;
      iinct = &iiinct ;
      A . 'VAR' . iinct . 'NOMDDL' = tnomt . 'NOMPRI' . iinct ;
      A . 'VAR' . iinct . 'DISC'   = tnomt . 'DISC'  ;
      'SI' lvalt ;
         'SI' lvt ;
            A . 'VAR' . iinct . 'VALEUR'   = 'EXTRAIRE' valt iinct ;
         'SINON' ;
            A . 'VAR' . iinct . 'VALEUR'   = valt ;
         'FINSI' ;
      'FINSI' ;
      'REPETER' iiincq nincq ;
         'REPETER' iiider numder ;
            iop = '+' iop 1 ;
            'SI' lcof ;
               icof = 'MINIMUM' ('LECT' &iiincq ncof) ;
               A . iop . iinct . &iiider = 'LECT' icof ;
            'SINON' ;
               A . iop . iinct . &iiider = ll ;
            'FINSI' ;
         'FIN' iiider ;
      'FIN' iiincq ;
   'FIN' iiinct ;
   'SI' iaxi ;
      'REPETER' iiincq nincq ;
         iop = '+' iop 1 ;
         'SI' lcof ;
            icof = 'MINIMUM' ('LECT' &iiincq ncof) ;
            A . iop . 1 . 0 = 'LECT' icof ;
         'SINON' ;
            A . iop . 1 . 0 = ll ;
         'FINSI' ;
      'FIN' iiincq ;
   'FINSI' ;
*
*   'SI' iaxi ;
*      numdat = 1 ;
*      numcof = dim '+' 1 ;
*   'SINON' ;
      numdat = 0 ;
      numcof = idim '*' idim '*' idim ;
*   'FINSI' ;
   'SI' iaxi ;
      numdat = '+' numdat 2 ;
      numcof = '+' numcof ('+' idim 2) ;
   'FINSI' ;
   numvar = nincq ;
   B = ININLIN numop numvar numdat numcof numder ;
*
   lvq = 'EGA' ('TYPE' valq) 'LISTREEL' ;
   'REPETER' iiinc nincq ;
      iinc = &iiinc ;
      B . 'VAR' . iinc . 'NOMDDL' = tnomq . 'NOMDUA' . iinc ;
      B . 'VAR' . iinc . 'DISC'   = tnomq . 'DISC'  ;
      'SI' lvalq ;
         'SI' lvq ;
            B . 'VAR' . iinc . 'VALEUR'   = 'EXTRAIRE' valq iinc ;
         'SINON' ;
            B . 'VAR' . iinc . 'VALEUR'   = valq ;
         'FINSI' ;
      'FINSI' ;
   'FIN' iiinc ;
   idat = 0 ;
   icof = 0 ;
   'SI' iaxi ;
      'REPETER' iiidim idim ;
         icof = '+' icof 1 ;
         B . 'COF' . icof . 'COMPOR' = 'CHAINE' 'VNOR' &iiidim ;
         B . 'COF' . icof . 'LDAT'   = 'LECT' ;
      'FIN' iiidim ;
      idat = '+' idat 1 ;
      icof = '+' icof 1 ;
      B . 'DAT' . idat . 'NOMDDL' = 'MOTS' 'SCAL' ;
      B . 'DAT' . idat . 'DISC'   = discg ;
      B . 'DAT' . idat . 'VALEUR' = dprmt ;
      B . 'COF' . icof . 'COMPOR' = 'IDEN' ;
      B . 'COF' . icof . 'LDAT'   = 'LECT' idat ;
      ll = 'LECT' icof ;
      idat = '+' idat 1 ;
      icof = '+' icof 1 ;
      B . 'DAT' . idat . 'NOMDDL' = 'MOTS' 'SCAL' ;
      B . 'DAT' . idat . 'DISC'   = 'CSTE' ;
      B . 'DAT' . idat . 'VALEUR' = deupi ;
      B . 'COF' . icof . 'COMPOR' = 'IDEN' ;
      B . 'COF' . icof . 'LDAT'   = 'LECT' idat ;
      ll2 = 'LECT' icof ;
   'SINON' ;
      ll = 'LECT' ;
   'FINSI' ;
*
  iop = 0 ;
  'REPETER' iiinct ninct ;
     'REPETER' iiincq nincq ;
         'REPETER' iiider numder ;
            iop  = '+' iop 1 ;
            icof = '+' icof 1 ;
            lcomp = 'CHAINE' loi &iiincq &iiinct &iiider ;
*            lcomp = 'CHAINE' loi &iiinct &iiincq &iiider ;
            B . 'COF' . icof . 'COMPOR' = lcomp ;
            B . 'COF' . icof . 'LDAT'   = 'LECT' ;
            B . iop . &iiincq . 0 = ('LECT' icof) 'ET' ll ;
         'FIN' iiider ;
      'FIN' iiincq ;
   'FIN' iiinct ;
   'SI' iaxi ;
      'REPETER' iiincq nincq ;
         iincq = &iiincq ;
         iop = '+' iop 1 ;
         B . iop . iincq . 0 = ('LECT' iincq) 'ET' ll2 ;
      'FIN' iiincq ;
   'FINSI' ;
*
*   mgnorkt = NLIN discg _mt A B 'CRES' methgau ;
   mgnorkt = NLIN discg _mt A B methgau ;
*
'RESPRO' mgnorkt ;
'FINPROC' ;
*
* End of procedure file GNORKTAN
*
*ENDPROCEDUR gnorktan
*BEGINPROCEDUR gnor
************************************************************************
* NOM         : GNOR
* DESCRIPTION : Calcule le champ de normales à une surface.
*               Peut servir à calculer une pression, un potentiel
*               lié à la gravité, un volume contenu dans une surface.
*               Attention à l'orientation de la surface !
*
*               Computes a field of normal to a surface.
*               Also useful to compute a pressure field,
*               a gravity potential field, a volume enclosed
*               by a surface.
*               WARNING : The orientation of the surface matters !
*
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 22/04/2011
* HISTORIQUE : v1, 22/04/2011, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
*
*
'DEBPROC' GNOR ;
'ARGUMENT'     _mt*'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
*
* Lectures
*
dim = 'VALEUR' 'DIME' ;
debug = FAUX ;
lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
                 'NCOF' 'FCOF' 'CCOF' ;
* Il faut initialiser valt et valq, sinon on peut capturer ceux de
* la procédure appelante
valt = 'valt' ; valq = 'valq' ;
'REPETER' imotcle ;
   'ARGUMENT' motcle/'MOT' ;
   'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
   'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
      cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
   'FINSI' ;
   'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
   tst1 = 'EGA' motcle 'FPRI' ; tst2 = 'EGA' motcle 'FDUA' ;
   tst3 = 'EGA' motcle 'FCOF' ;
   tst = tst1 'OU' tst2 'OU' tst3 ;
   'SI' tst ;
      'SI' tst1 ; tt = TDISC . nomt ; 'FINSI' ;
      'SI' tst2 ; tt = TDISC . nomq ; 'FINSI' ;
      'SI' tst3 ; tt = TDISC . nomo ; 'FINSI' ;
      isvec = ('>' ('DIME' (tt. 'NOMPRI')) 1) ;
      'SI' isvec ; 'ARGUMENT' valv*'LISTREEL' ; 'SINON' ;
         'ARGUMENT' valv*'FLOTTANT' ;
      'FINSI' ;
      'SI' tst1 ; valt = valv ; 'FINSI' ;
      'SI' tst2 ; valq = valv ; 'FINSI' ;
      'SI' tst3 ; valo = valv ; 'FINSI' ;
   'FINSI' ;
   'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
'FIN' imotcle ;
*
* Tests
*
discg = TDISC . 'GEOM' . 'DISC' ;
'SI' ('EXISTE' tdisc 'methgau') ;
   methgau = tdisc . 'methgau' . 'mass' ;
'SINON' ;
   methgau = 'GAU7' ;
'FINSI' ;
tnomt = TDISC . nomt ;
lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
tnomq = TDISC . nomq ;
lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
* Scalaire ou vecteur
ninct = 'DIME' (tnomt . 'NOMPRI') ;
nincq = 'DIME' (tnomq . 'NOMPRI') ;
'SI' ('ET' ('NEG' ninct 1) ('NEG' ninct dim)) ;
   cherr = 'CHAINE'
   'la primale doit etre un scalaire ou un vecteur' ;
   'ERREUR' cherr ;
'FINSI' ;
'SI' ('NEG' nincq dim) ;
   cherr = 'CHAINE'
   'la duale doit etre un vecteur' ;
   'ERREUR' cherr ;
'FINSI' ;
*ninc = ninct ;
*
lcof = 'EXISTE' TDISC nomo ;
'SI' lcof ; tcof = TDISC . nomo ;
ncof = 'DIME' (tcof . 'NOMPRI') ;
'SINON' ; ncof = 0 ;
'FINSI' ;
*
'SI' debug ;
'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
'SINON' ; 'MESSAGE' 'pas de coef  detecte' ;
'FINSI' ;
'FINSI' ;
*
   vdim = 'VALEUR' 'DIME' ;
   vmod = 'VALEUR' 'MODE' ;
   idim = 0 ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
      idim = 2 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
      idim = 2 ;
      iaxi = VRAI ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
      idim = 3 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' vdim 1) ;
      idim = 1 ;
      iaxi = FAUX ;
   'FINSI' ;
*   'MESSAGE' ('CHAINE' 'iaxi=' iaxi );
   'SI' ('EGA' idim 0) ;
      'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
   'FINSI' ;
   'SI' iaxi ;
      dprmt = '*' ('COORDONNEE' 1 _mt) ('*' PI 2.D0) ;
   'FINSI' ;
*
* Optimisation possible : construire la matrice par blocs
* qd valt et valq ne sont pas donnés
*
   numop = nincq ; numder = idim ; numvar = ninct ;
   numdat = ncof ; numcof = ncof ;
   A = ININLIN numop numvar numdat numcof numder ;
   'SI' lcof ;
      lvo = 'EGA' ('TYPE' valo) 'LISTREEL' ;
      'REPETER' iicof ncof ;
         icof = &iicof ;
         A . 'DAT' . icof . 'NOMDDL' = tcof . 'NOMPRI' . icof ;
         A . 'DAT' . icof . 'DISC'   = tcof . 'DISC' ;
         'SI' lvo ;
            A . 'DAT' . icof . 'VALEUR' = 'EXTRAIRE' valo icof ;
         'SINON' ;
            A . 'DAT' . icof . 'VALEUR' = valo ;
         'FINSI' ;
         A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
         A . 'COF' . icof . 'LDAT'   = 'LECT' icof ;
      'FIN' iicof ;
   'SINON' ;
      ll = 'LECT' ;
   'FINSI' ;
   lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
   'REPETER' iiincq nincq ;
      iincq = &iiincq ;
      iinct = 'MINIMUM' ('LECT' iincq ninct) ;
      A . 'VAR' . iinct . 'NOMDDL' = tnomt . 'NOMPRI' . iinct ;
      A . 'VAR' . iinct . 'DISC'   = tnomt . 'DISC'  ;
      'SI' lvalt ;
         'SI' lvt ;
            A . 'VAR' . iinct . 'VALEUR'   = 'EXTRAIRE' valt iinct ;
         'SINON' ;
            A . 'VAR' . iinct . 'VALEUR'   = valt ;
         'FINSI' ;
      'FINSI' ;
      'SI' lcof ;
         icof = 'MINIMUM' ('LECT' iincq ncof) ;
         A . iincq . iinct . 0  = 'LECT' icof ;
      'SINON' ;
         A . iincq . iinct . 0  = ll ;
      'FINSI' ;
   'FIN' iiincq ;
*
   'SI' iaxi ;
      numdat = 1 ;
      numcof = dim '+' 1 ;
   'SINON' ;
      numdat = 0 ;
      numcof = dim ;
   'FINSI' ;
   numvar = nincq ;
   B = ININLIN numop numvar numdat numcof numder ;
   icof = 0 ;
   'REPETER' iiidim idim ;
      icof = '+' icof 1 ;
      B . 'COF' . icof . 'COMPOR' = 'CHAINE' 'VNOR' &iiidim ;
      B . 'COF' . icof . 'LDAT'   = 'LECT' ;
   'FIN' iiidim ;
*
   'SI' iaxi ;
      icof = '+' icof 1 ;
      B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
      B . 'DAT' . 1 . 'DISC'   = discg ;
      B . 'DAT' . 1 . 'VALEUR' = dprmt ;
      B . 'COF' . icof . 'COMPOR' = 'IDEN' ;
      B . 'COF' . icof . 'LDAT'   = 'LECT' 1 ;
      ll = 'LECT' icof ;
   'SINON' ;
      ll = 'LECT' ;
   'FINSI' ;
   lvq = 'EGA' ('TYPE' valq) 'LISTREEL' ;
   'REPETER' iiincq nincq ;
      iincq = &iiincq ;
      B . 'VAR' . iincq . 'NOMDDL' = tnomq . 'NOMDUA' . iincq ;
      B . 'VAR' . iincq . 'DISC'   = tnomq . 'DISC'  ;
      'SI' lvalq ;
         'SI' lvq ;
            B . 'VAR' . iincq . 'VALEUR'   = 'EXTRAIRE' valq iincq ;
         'SINON' ;
            B . 'VAR' . iincq . 'VALEUR'   = valq ;
         'FINSI' ;
      'FINSI' ;
      B . iincq . iincq . 0  = ('LECT' iincq) 'ET' ll ;
   'FIN' iiincq ;
*
   mgnor = NLIN discg _mt A B methgau ;
*
'RESPRO' mgnor ;
'FINPROC' ;
*
* End of procedure file GNOR
*
*ENDPROCEDUR gnor
*BEGINPROCEDUR grespr
* GRESPR     PROCEDUR  GOUNAND   11/06/14    21:15:30     7005
************************************************************************
* NOM         : GRESPR
* DESCRIPTION : Résout un système par une méthode de projection
*               algébrique incrémentale
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 22/11/2007, version initiale
* HISTORIQUE : v1, 22/11/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' GRESPR ;
'ARGUMENT'  mat*'RIGIDITE' ;
*'ARGUMENT' ccli*'CHPOINT' ;
'ARGUMENT' smb*'CHPOINT' ;
'ARGUMENT' tdisc*'TABLE' ;
'ARGUMENT' rv*'TABLE' ;
*
* Nom des inconnues : vitesses et pressions
* Dans rvg . 'METHINV', on stocke l'éventuelle table d'inversion
* pour le laplacien de pression
*
debug = rv . 'debug' ;
rvg = rv ;
rvm = rv . 'METHINVQDM' ;
dblproj = rv . 'dblproj' ;
scaldia = rv . 'scaldia' ;
*
'SI' ('EXISTE' rvg 'noprec') ;
   noprec = rvg . 'noprec' ;
'SINON' ;
   noprec = faux ;
'FINSI' ;
'SI' debug ;
   'MESSAGE' ('CHAINE' 'noprec=' noprec) ;
'FINSI' ;
*
'SI' ('NON' noprec) ;
   'SI' ('NON' ('EXISTE' rvg 'preco')) ;
      rvg . 'preco' = 'TABLE' ;
   'FINSI' ;
   rvgp = rvg . 'preco' ;
* Précision relative utilisée pour le préconditionnement du calcul de la
* matrice de pression
*  precrel = 1.D-12 ;
*   precrel = 1.D-2 ;
   precrel = 1.D-1 ;
'SINON' ;
   'SI' ('EXISTE' rvg 'preco') ;
      rvgp = rvg . 'preco' ;
      'SI' ('EXISTE' rvgp 'mklapphi') ;
         mklapphi = rvgp . 'mklapphi' ;
         imkvito  = rvgp . 'imkvitd' ;
         'OUBLIER' rvgp 'mklapphi' ;
         'OUBLIER' rvgp 'imkvitd' ;
         'SI' debug ;
            'MESSAGE' 'Destruction mat pression' ;
         'FINSI' ;
*         'DETRUIT' mklapphi 'ELEMENTAIRE' ;
*         'DETRUIT' imkvito 'GEOMETRIE' ;
      'FINSI' ;
      'OUBLIER' rvg  'preco' ;
   'FINSI' ;
   rvg . 'preco' = 'TABLE' ;
   rvgp = rvg . 'preco' ;
'FINSI' ;
* On a mis chaine sinon plantage qd la pression s'appelle 'PRES' ;
nivit = @STBL (TDISC . 'UN' . 'NOMPRI') ;
nifor = @STBL (TDISC . 'XN' . 'NOMDUA') ;
nipre = @STBL (TDISC . 'PN' . 'NOMPRI') ;
niqre = @STBL (TDISC . 'PN' . 'NOMDUA') ;
*
'SI' ('NON' ('EXISTE' rv 'METHINVPRES')) ;
   rv . 'METHINVPRES' = rvm ;
'FINSI' ;
rvgm = rv . 'METHINVPRES' ;
*
* On commence par réduire le système en enlevant les conditions aux limites
*
matc smbc smb1 = 'KOPS' 'CONDENSE' mat smb ;
*
* Initialisation de la solution du système total
*
********************************************************************
*                                                                  *
* On résout la partie vitesse-pression par projection incrémentale *
*                                                                  *
********************************************************************
nivp = 'ET' nivit nipre ;
*
* Matrice : partie diagonale en vitesse et contrainte
*           on est obligé de reconstruire la transposée car KBBT
*           ne stocke que B et pas Bt
*
mkvit     = 'EXTRAIRE' matc nivit nifor ;
mkcontr   = 'EXTR' matc nivit niqre ;
mkcontrt  = 'EXTR' matc nipre nifor ;
*
* Construction du résidu * -1
*
residu = smbc ;
*
* Calcul de la matrice de pression et de la matrice masse diagonalisée
*
*  On lumpe la partie diagonale et on l'inverse
chpo1 = GCHPO _mtw 'UN' 1. ;
mkvitd  = '*' mkvit chpo1 ;
*'LISTE' mkvitd ;
*
* Choix 1 : si la matrice masse lumpée
*           a tous ses termes positifs, on en prend l'inverse
*           Avantage : simple, autorise le préconditionnement
*                      car ne change pas à chaque itération
*           Inconvénient : ne marche pas en axi quadratique sur l'axe
* Choix 2 : on prend la diagonale de la matrice totale
*           on la met à l'échelle pour avoir la même masse que mkvitd
*           Avantage : marche en axi, sur maillage déformés
*           Inconvénient : change à chaque itération quand le problème
*           est non linéaire (=> recalcul de la matrice de pression),
*           pas forcément positif.
* Choix 3 : valeur absolue du précédent + petite tolérance
*           Avantage : marche peut-être
*           Inconvénient : change à chaque itération.
*
*
*mmi = 'MINIMUM' mkvitd ;
*'SI' ('>' mmi 1.D-200) ;
*   imkvitd1 = 'INVERSE' mkvitd ;
*   'SI' debug ; 'MESSAGE' 'Choix 1 masse lumpée' ; 'FINSI' ;
*'SINON' ;
   mkdia = 'EXTR' mkvit 'DIAG' ;
*   'LISTE' mkdia ;
   mmi = 'MINIMUM' mkdia ;
   'SI' ('>' mmi 1.D-200) ;
      'SI' debug ; 'MESSAGE' 'Choix 2 masse lumpée' ; 'FINSI' ;
   'SINON' ;
      mkdia = 'ABS' mkdia ;
      mma = 'MAXIMUM' mkdia ;
      mkdia = '+' mkdia ('*' mma 1.D-8) ;
      'SI' debug ; 'MESSAGE' 'Choix 3 masse lumpée' ; 'FINSI' ;
   'FINSI' ;
   'SI' scaldia ;
      masvitd = 'MAXIMUM' ('RESULT' mkvitd) ;
      'SI' debug ;
         'MESSAGE' ('CHAINE' 'totmas=' masvitd) ;
      'FINSI' ;
      masdia = 'MAXIMUM' ('RESULT' mkdia) ;
      imkvitd = '*' ('INVERSE' mkdia) ('/' masdia masvitd) ;
   'SINON' ;
      imkvitd = 'INVERSE' mkdia ;
   'FINSI' ;
*'FINSI' ;
*
* Attention ! On préconditionne éventuellement le calcul de la matrice
* de pression
*
precok = FAUX ;
*'MESSAGE' ('CHAINE' 'precok = ' precok) ;
'SI' ('NON' noprec) ;
   exis1 = 'EXISTE' rvgp 'mklapphi' ;
   exis2 = 'EXISTE' rvgp 'imkvitd'   ;
   'SI' (exis1 'ET' exis2) ;
      imkvd = rvgp . 'imkvitd' ;
      mkl   = rvgp . 'mklapphi' ;
* On vérifie l'égalité des mkvitd à un facteur constant près
* et des ccli1
      ech  = 'MAXIMUM' imkvitd 'ABS' ;
      echp = 'MAXIMUM' imkvd 'ABS' ;
      alfa = '/' ech echp ;
*   'MESSAGE' ('CHAINE' 'alfa = ' alfa) ;
      dimkv = '-' imkvitd ('*' imkvd alfa) ;
      ndimkv = 'EXTRAIRE' dimkv 'COMP' ;
      dimkvs = '*' dimkv ('INVERSE' imkvitd) ndimkv ndimkv ndimkv ;
      mdimkvs = 'MAXIMUM' dimkvs 'ABS' ;
      egaimkv = 'EGA' mdimkvs 0. precrel ;
*
      'SI' debug ;
         'MESSAGE' ('CHAINE' 'egaimkv = ' egaimkv
                       '  a ' ('*' mdimkvs 100.) ' %' ) ;
      'FINSI' ;
      precok = egaimkv ;
   'FINSI' ;
'FINSI' ;
*  On surcharge par 0. l'inverse de la diagonale là où il y a des clims
*  de Dirichlet (pas besoin ici)
*
chdbg = 'CHAINE' 'Matrice de pression MP : ' ;
'SI' precok ;
   'SI' ('NON' ('EGA' alfa 1.D0 precrel)) ;
      ialfa = '/' 1. alfa ;
*      chdbg = 'CHAINE' chdbg 'MP(i) = MP(i-1) * ' ialfa ;
*      mklapphi = '*' (rvgp . 'mklapphi') ialfa ;
*      imkvitd  = '*' (rvgp . 'imkvitd')  ialfa ;
      chdbg = 'CHAINE' chdbg 'MP(i) = MP(i-1) * ' ialfa ;
      mklapphi = rvgp . 'mklapphi' ;
      imkvitd  = rvgp . 'imkvitd' ;
   'SINON' ;
      chdbg = 'CHAINE' chdbg 'MP(i) = MP(i-1)'  ;
      mklapphi = rvgp . 'mklapphi' ;
      imkvitd  = rvgp . 'imkvitd' ;
      ialfa = 1.D0 ;
   'FINSI' ;
'SINON' ;
* Dans le cas général où mkcontr et mkcontrt sont différents,
* il faudrait mettre :
*mklapphi = 'KOPS' 'CMCT' mkcontr imchd (kops mkcontrt 'TRANSPOS') ;
   chdbg = 'CHAINE' chdbg 'calcul...' ;
   'SI' ('EXISTE' rvgp 'mklapphi') ;
      mklapphi = rvgp . 'mklapphi' ;
      imkvito  = rvgp . 'imkvitd' ;
      'OUBLIER' rvgp 'mklapphi' ;
      'OUBLIER' rvgp 'imkvitd' ;
      'SI' debug ;
         'MESSAGE' 'Destruction de lancienne matpres' ;
      'FINSI' ;
*      'DETRUIT' mklapphi 'ELEMENTAIRE' ;
*      'DETRUIT' imkvito 'GEOMETRIE' ;
   'FINSI' ;
   mklapphi = 'KOPS' 'CMCT' mkcontr imkvitd mkcontr ;
   ialfa = 1.D0 ;
'FINSI' ;
'SI' debug ; 'MESSAGE' chdbg ; 'FINSI' ;
*'SI' ('NON' noprec) ;
   rvgp . 'mklapphi'  = mklapphi ;
   rvgp . 'imkvitd'   = imkvitd  ;
*'FINSI' ;
*
* Calcul éventuel d'une estimation de la pression (double projection)
*
'SI' dblproj ;
   desqdm = 'EXCO' nifor residu nifor 'NOID' ;
   idesq  = '*' desqdm imkvitd nifor nivit nivit ;
   didesq = '*' mkcontr idesq ;
   mtot  = mklapphi ;
   ftot  = didesq ;
*   'MESSAGE' 'Résolution Pression 1' ;
   dphi1 = MONRES mtot ftot rvgm ;
*   HCHPO dphi1 'dphi1' ;
   'SI' ('NEG' ialfa 1.D0 1.D-13) ;
      dphi1 = '/' dphi1 ialfa ;
   'FINSI' ;
   gdphi1 = '*' mkcontrt dphi1 ;
   residu = residu '-' gdphi1 ;
'FINSI' ;
*
* Calcul de l'incrément de vitesse dv* à partir de la qdm
* On ne splitte pas les équations par direction
*
dvitstar = 0. ;
mtot = mkvit ;
ftot = 'EXCO' nifor residu nifor 'NOID' ;
*   'MESSAGE' 'Résolution Vitesse' ;
ct   = MONRES  mtot ftot rvm ;
*HCHPO ct 'ct' ;
*TRVEC mtw nomvit ct 'ct' ;
dvitstar = '+' dvitstar ct ;
*
* Calcul de l'incrément de pression dphi pour vérifier la conservation
* de la masse
*
ftot1 = '*' mkcontr dvitstar ;
ftot2 = 'EXCO' niqre residu niqre 'NOID' ;
mtot  = mklapphi ;
ftot  = '-' ftot1 ftot2 ;
*   'MESSAGE' 'Résolution Pression 2' ;
dphi = MONRES mtot ftot rvgm ;
*HCHPO dphi 'dphi' ;
'SI' ('NEG' ialfa 1.D0 1.D-13) ;
   dphi = '/' dphi ialfa ;
'FINSI' ;
*
* Calcul de l'incrément final de vitesse dvit
*
gdphi = '*' mkcontrt dphi ;
gdphi = '*' gdphi imkvitd nifor nivit nivit ;
'SI' ('NEG' ialfa 1.D0 1.D-13) ;
   gdphi = '*' gdphi ialfa ;
'FINSI' ;
dvit  = dvitstar '-' gdphi ;
*
'SI' dblproj ;
   dphi = dphi '+' dphi1 ;
*   HCHPO dphi1 'dphi1' ;
'FINSI' ;
dsolc = dvit '+' dphi ;
*HCHPO dphi 'dphi' ;
*hchpo ('+' dphi dphi1) 'dphi+dphi1' ;
*
* Reconstruction de la solution du système total
*
dsol = 'KOPS' 'EVAPORE' dsolc mat smb smb1 ;
*
'RESPRO' dsol ;
*
* End of procedure file GRESPR
*
'FINPROC' ;
*ENDPROCEDUR grespr
*BEGINPROCEDUR grig2
************************************************************************
* NOM         : GRIG2
* DESCRIPTION :
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, ??/??/2007, version initiale
* HISTORIQUE : v1, ??/??/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' GRIG2 ;
'ARGUMENT'     _mt*'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
*
* Lectures
*
debug = FAUX ;
lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
                 'NCOF' 'FCOF' 'CCOF' 'LAPN' 'GMBT' ;
* Il faut initialiser valt et valq, sinon on peut capturer ceux de
* la procédure appelante
valt  = 'valt' ; valq = 'valq' ;
llapn = 0 ;
'REPETER' imotcle ;
   'ARGUMENT' motcle/'MOT' ;
   'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
   'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
      cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
   'FINSI' ;
   'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
   tst1 = 'EGA' motcle 'FPRI' ; tst2 = 'EGA' motcle 'FDUA' ;
   tst = tst1 'OU' tst2 ;
   'SI' tst ;
      'SI' tst1 ; tt = TDISC . nomt ; 'FINSI' ;
      'SI' tst2 ; tt = TDISC . nomq ; 'FINSI' ;
      isvec = ('>' ('DIME' (tt. 'NOMPRI')) 1) ;
      'SI' isvec ; 'ARGUMENT' valv*'LISTREEL' ; 'SINON' ;
         'ARGUMENT' valv*'FLOTTANT' ;
      'FINSI' ;
      'SI' tst1 ; valt = valv ; 'FINSI' ;
      'SI' tst2 ; valq = valv ; 'FINSI' ;
   'FINSI' ;
   'SI' ('EGA' motcle 'FCOF') ; 'ARGUMENT' valo*'FLOTTANT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'LAPN') ; llapn = 1 ; 'FINSI' ;
   'SI' ('EGA' motcle 'GMBT') ; llapn = 2 ; 'FINSI' ;
'FIN' imotcle ;
*
* Tests
*
discg = TDISC . 'GEOM' . 'DISC' ;
'SI' ('EXISTE' tdisc 'methgau') ;
   methgau = tdisc . 'methgau' . 'rigi' ;
'SINON' ;
   methgau = 'GAU7' ;
'FINSI' ;
tnomt = TDISC . nomt ;
lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
tnomq = TDISC . nomq ;
lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
* Scalaire ou vecteur
ninct = 'DIME' (tnomt . 'NOMPRI') ;
nincq = 'DIME' (tnomq . 'NOMPRI') ;
'SI' ('NEG' ninct nincq) ;
   cherr = 'CHAINE'
   'les primales et duales nont pas le meme nombre de composantes' ;
   'ERREUR' cherr ;
'FINSI' ;
'SI' ('NEG' ninct ('VALEUR' 'DIME')) ;
   cherr = 'CHAINE'
   'les inconnues doivent etre vectorielles' ;
   'ERREUR' cherr ;
'FINSI' ;
*
ninc = ninct ;
*
lcof = 'EXISTE' TDISC nomo ;
'SI' lcof ; ncof = 1 ; tcof = TDISC . nomo ;
'SINON' ; ncof = 0 ;
'FINSI' ;
*
'SI' debug ;
'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
'SINON' ; 'MESSAGE' 'pas de coef  detecte' ;
'FINSI' ;
'FINSI' ;
*
   vdim = 'VALEUR' 'DIME' ;
   vmod = 'VALEUR' 'MODE' ;
   idim = 0 ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
      idim = 2 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
      idim = 2 ;
      iaxi = VRAI ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
      idim = 3 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' vdim 1) ;
      idim = 1 ;
      iaxi = FAUX ;
   'FINSI' ;
*   'MESSAGE' ('CHAINE' 'iaxi=' iaxi );
   'SI' ('EGA' idim 0) ;
      'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
   'FINSI' ;
   'SI' iaxi ;
      rmt   = 'COORDONNEE' 1 _mt ;
      deupi = '*' PI 2.D0 ;
   'FINSI' ;
*
* Optimisation possible : construire la matrice par blocs
* qd valt et valq ne sont pas donnés
*
*
*Bug ?   numop  = ('**' ninc 2) '+' 1 ;
   numop  = '**' ninc 2 ;
   'SI' iaxi ; numop = '+' numop 1 ; 'FINSI' ;
   numder = idim ;
   numvar = ninc ;
   ncof = '+' ncof 1 ;
*delete   'SI' iaxi ; ncof = '+' ncof 1 ; 'FINSI' ;
   numdat = ncof ;
   numcof = ncof ;
*
   A = ININLIN numop numvar numdat numcof numder ;
*
   lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
   'REPETER' iiinc ninc ;
      iinc = &iiinc ;
      A . 'VAR' . iinc . 'NOMDDL' = tnomt . 'NOMPRI' . iinc ;
      A . 'VAR' . iinc . 'DISC'   = tnomt . 'DISC'  ;
      'SI' lvalt ;
         'SI' lvt ;
            A . 'VAR' . iinc . 'VALEUR'   = 'EXTRAIRE' valt iinc ;
         'SINON' ;
            A . 'VAR' . iinc . 'VALEUR'   = valt ;
         'FINSI' ;
      'FINSI' ;
   'FIN' iiinc ;
*
   icof = 0 ;
   ll = 'LECT' ;
   icof = '+' icof 1 ;
   A . 'DAT' . icof . 'NOMDDL' = 'MOTS' 'SCAL' ;
   A . 'DAT' . icof . 'DISC'   = 'CSTE' ;
   A . 'DAT' . icof . 'VALEUR' = 2.D0 ;
   A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
   A . 'COF' . icof . 'LDAT'   = 'LECT' icof ;
   ll  = 'LECT' ;
   ll2 = 'LECT' icof ;
   'SI' lcof ;
      icof = '+' icof 1 ;
      A . 'DAT' . icof . 'NOMDDL' = tcof . 'NOMPRI' . 1 ;
      A . 'DAT' . icof . 'DISC'   = tcof . 'DISC' ;
      A . 'DAT' . icof . 'VALEUR' = valo ;
      A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
      A . 'COF' . icof . 'LDAT'   = 'LECT' icof ;
      ll  = 'ET' ll  ('LECT' icof) ;
      ll2 = 'ET' ll2 ('LECT' icof) ;
   'FINSI' ;
*
   iop = 0 ;
   'REPETER' iidim idim ;
      'REPETER' jidim idim ;
         iop = '+' iop 1 ;
         'SI' ('EGA' &iidim &jidim) ;
            'SI' ('EGA' llapn 0) ;
               A . iop . &iidim . &jidim = ll2 ;
            'SINON' ;
               A . iop . &iidim . &jidim = ll ;
            'FINSI' ;
         'SINON' ;
            'SI' ('NEG' llapn 2) ;
               A . iop . &iidim . &jidim = ll ;
            'FINSI' ;
            'SI' ('NEG' llapn 1) ;
               A . iop . &jidim . &iidim = ll ;
            'FINSI' ;
         'FINSI' ;
      'FIN' jidim ;
   'FIN' iidim ;
   'SI' iaxi ;
      iop = '+' iop 1 ;
      'SI' ('EGA' llapn 0) ;
         A . iop . 1 . 0 = ll2 ;
      'SINON' ;
         A . iop . 1 . 0 = ll ;
      'FINSI' ;
   'FINSI' ;
*
   'SI' iaxi ;
      numdat = 2 ;
      numcof = 2 ;
   'SINON' ;
      numdat = 0 ;
      numcof = 0 ;
   'FINSI' ;
*
   B = ININLIN numop numvar numdat numcof numder ;
*
   icof = 0 ;
*
   'SI' iaxi ;
      icof = '+' icof 1 ;
      B . 'DAT' . icof . 'NOMDDL' = 'MOTS' 'SCAL' ;
      B . 'DAT' . icof . 'DISC'   = discg ;
      B . 'DAT' . icof . 'VALEUR' = rmt ;
      B . 'COF' . icof . 'COMPOR' = 'IDEN' ;
      B . 'COF' . icof . 'LDAT'   = 'LECT' icof ;
      icof = '+' icof 1 ;
      B . 'DAT' . icof . 'NOMDDL' = 'MOTS' 'SCAL' ;
      B . 'DAT' . icof . 'DISC'   = 'CSTE' ;
      B . 'DAT' . icof . 'VALEUR' = deupi ;
      B . 'COF' . icof . 'COMPOR' = 'IDEN' ;
      B . 'COF' . icof . 'LDAT'   = 'LECT' icof ;
      ll  = 'LECT' ('-' icof 1) icof ;
      mic = '*' ('-' icof 1) -1 ;
      llb = 'LECT' mic icof ;
   'FINSI' ;
*
   lvq = 'EGA' ('TYPE' valq) 'LISTREEL' ;
   'REPETER' iiinc ninc ;
      iinc = &iiinc ;
      B . 'VAR' . iinc . 'NOMDDL' = tnomq . 'NOMDUA' . iinc ;
      B . 'VAR' . iinc . 'DISC'   = tnomq . 'DISC'  ;
      'SI' lvalq ;
         'SI' lvq ;
            B . 'VAR' . iinc . 'VALEUR'   = 'EXTRAIRE' valq iinc ;
         'SINON' ;
            B . 'VAR' . iinc . 'VALEUR'   = valq ;
         'FINSI' ;
      'FINSI' ;
   'FIN' iiinc ;
*
   iop = 0 ;
   'REPETER' iidim idim ;
      'REPETER' jidim idim ;
         iop = '+' iop 1 ;
         B . iop . &iidim . &jidim = ll ;
      'FIN' jidim ;
   'FIN' iidim ;
   'SI' iaxi ;
      iop = '+' iop 1 ;
      B . iop . 1 . 0 = llb ;
*      B . iop . 1 . 0 = ll ;
   'FINSI' ;
*
   mgrig = NLINP discg _mt A B methgau ;
* Integration par parties
*   mgrig = '*' mgrig -1.D0 ;
*
'RESPRO' mgrig ;
'FINPROC' ;
*
* End of procedure file GRIG2
*
*ENDPROCEDUR grig2
*BEGINPROCEDUR gugrad2
************************************************************************
* NOM         : GUGRAD2
* DESCRIPTION : Une matrice de convection
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v2, 14/03/2006, mise à jour NLIN évolué
* VERSION    : v1, 13/05/2004, version initiale
* HISTORIQUE : v1, 13/05/2004, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' GUGRAD2 ;
'ARGUMENT'     _mt*'MAILLAGE' ;
'ARGUMENT'     _smt/'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
*
* Lectures
*
debug = FAUX ;
lmotcle = 'MOTS' 'NPRI' 'FPRI' 'CPRI' 'NDUA' 'FDUA' 'CDUA'
                 'NVIT' 'FVIT' 'CVIT' 'NCOF' 'FCOF' 'CCOF' ;
'REPETER' imotcle ;
   'ARGUMENT' motcle/'MOT' ;
   'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
   'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
      cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
   'FINSI' ;
   'SI' ('EGA' motcle 'NPRI') ; 'ARGUMENT' nomt*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NDUA') ; 'ARGUMENT' nomq*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NCOF') ; 'ARGUMENT' nomo*'MOT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'NVIT') ; 'ARGUMENT' nomu*'MOT' ; 'FINSI' ;
   tst1 = 'EGA' motcle 'FPRI' ; tst2 = 'EGA' motcle 'FDUA' ;
   tst = tst1 'OU' tst2 ;
   'SI' tst ;
      'SI' tst1 ; tt = TDISC . nomt ; 'FINSI' ;
      'SI' tst2 ; tt = TDISC . nomq ; 'FINSI' ;
      isvec = ('>' ('DIME' (tt. 'NOMINC')) 1) ;
      'SI' isvec ; 'ARGUMENT' valv*'LISTREEL' ; 'SINON' ;
         'ARGUMENT' valv*'FLOTTANT' ;
      'FINSI' ;
      'SI' tst1 ; valt = valv ; 'FINSI' ;
      'SI' tst2 ; valq = valv ; 'FINSI' ;
   'FINSI' ;
   'SI' ('EGA' motcle 'FCOF') ; 'ARGUMENT' valo*'FLOTTANT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'FVIT') ; 'ARGUMENT' valu*'LISTREEL' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CPRI') ; 'ARGUMENT' valt*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CDUA') ; 'ARGUMENT' valq*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CCOF') ; 'ARGUMENT' valo*'CHPOINT' ; 'FINSI' ;
   'SI' ('EGA' motcle 'CVIT') ; 'ARGUMENT' valu*'CHPOINT' ; 'FINSI' ;
'FIN' imotcle ;
*
* Tests
*
discg = TDISC . 'GEOM' . 'DISC' ;
'SI' ('EXISTE' tdisc 'methgau') ;
   methgau = tdisc . 'methgau' . 'amor' ;
'SINON' ;
   methgau = 'GAU7' ;
'FINSI' ;
tnomt = TDISC . nomt ;
lvalt = 'NEG' ('TYPE' valt) 'MOT' ;
tnomq = TDISC . nomq ;
lvalq = 'NEG' ('TYPE' valq) 'MOT' ;
tnomu = TDISC . nomu ;
* Scalaire ou vecteur
ninct = 'DIME' (tnomt . 'NOMPRI') ;
nincq = 'DIME' (tnomq . 'NOMPRI') ;
'SI' ('NEG' ninct nincq) ;
   cherr = 'CHAINE'
   'les primales et duales nont pas le meme nombre de composantes' ;
   'ERREUR' cherr ;
'FINSI' ;
ninc = ninct ;
*
lcof = 'EXISTE' TDISC nomo ;
'SI' lcof ; ncof = 1 ; tcof = TDISC . nomo ;
'SINON' ; ncof = 0 ;
'FINSI' ;
*
'SI' debug ;
'SI' lcof ; 'MESSAGE' 'Un coef a ete detecte' ;
'SINON' ; 'MESSAGE' 'pas de coef  detecte' ;
'FINSI' ;
'FINSI' ;
*
   vdim = 'VALEUR' 'DIME' ;
   vmod = 'VALEUR' 'MODE' ;
   idim = 0 ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'PLANDEFO')) ;
      idim = 2 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
      idim = 2 ;
      iaxi = VRAI ;
   'FINSI' ;
   'SI' ('ET' ('EGA' vdim 3) ('EGA' vmod 'TRID')) ;
      idim = 3 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' vdim 1) ;
      idim = 1 ;
      iaxi = FAUX ;
   'FINSI' ;
   'SI' ('EGA' idim 0) ;
      'ERREUR' ('CHAINE' 'vmod=' vmod ' et vdim=' vdim ' non prevu') ;
   'FINSI' ;
   'SI' iaxi ;
      dprmt = '*' ('COORDONNEE' 1 _mt) ('*' PI 2.D0) ;
   'FINSI' ;
* Scalaire ou vecteur
nincu = 'DIME' (tnomu . 'NOMPRI') ;
'SI' ('NEG' nincu idim) ;
   cherr = 'CHAINE'
   'la vitesse doit etre un vecteur' ;
   'ERREUR' cherr ;
'FINSI' ;
ncof = '+' ncof nincu ;
*
* Optimisation possible : construire la matrice par blocs
* qd valt et valq ne sont pas donnés
*
   numop = ninc ; numder = idim ; numvar = ninc ;
   numdat = ncof ; numcof = ncof ;
   A = ININLIN numop numvar numdat numcof numder ;
*
   lvt = 'EGA' ('TYPE' valt) 'LISTREEL' ;
   'REPETER' iiinc ninc ;
      iinc = &iiinc ;
      A . 'VAR' . iinc . 'NOMDDL' = tnomt . 'NOMPRI' . iinc ;
      A . 'VAR' . iinc . 'DISC'   = tnomt . 'DISC'  ;
      'SI' lvalt ;
         'SI' lvt ;
            A . 'VAR' . iinc . 'VALEUR'   = 'EXTRAIRE' valt iinc ;
         'SINON' ;
            A . 'VAR' . iinc . 'VALEUR'   = valt ;
         'FINSI' ;
      'FINSI' ;
   'FIN' iiinc ;
*
   lvu = 'EGA' ('TYPE' valu) 'LISTREEL' ;
   icof = 0 ;
   'REPETER' iiincu nincu ;
      iincu = &iiincu ;
      icof = icof '+' 1 ;
      A . 'DAT' . icof . 'NOMDDL' = tnomu . 'NOMPRI' . iincu ;
      A . 'DAT' . icof . 'DISC'   = tnomu . 'DISC' ;
      'SI' lvu ;
         A . 'DAT' . icof . 'VALEUR'   = 'EXTRAIRE' valu iincu ;
      'SINON' ;
         A . 'DAT' . icof . 'VALEUR'   = valu ;
      'FINSI' ;
      A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
      A . 'COF' . icof . 'LDAT'   = 'LECT' icof ;
   'FIN' iiincu ;
   'SI' lcof ;
      icof = icof '+' 1 ;
      A . 'DAT' . icof . 'NOMDDL' = tcof . 'NOMPRI' . 1 ;
      A . 'DAT' . icof . 'DISC'   = tcof . 'DISC' ;
      A . 'DAT' . icof . 'VALEUR' = valo ;
      A . 'COF' . icof . 'COMPOR' = 'IDEN' ;
      A . 'COF' . icof . 'LDAT'   = 'LECT' icof ;
      ll = 'LECT' icof ;
   'SINON' ;
      ll = 'LECT' ;
   'FINSI' ;
*
   'REPETER' iiinc ninc ;
      iinc = &iiinc ;
      'REPETER' iiidim idim ;
         iidim = &iiidim ;
         A . iinc . iinc . iidim = ('ET' ('LECT' iidim) ll) ;
      'FIN' iiidim ;
   'FIN' iiinc ;
*
   'SI' iaxi ;
      numdat = 1 ;
      numcof = 1 ;
   'SINON' ;
      numdat = 0 ;
      numcof = 0 ;
   'FINSI' ;
   B = ININLIN numop numvar numdat numcof numder ;
   'SI' iaxi ;
      B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
      B . 'DAT' . 1 . 'DISC'   = discg ;
      B . 'DAT' . 1 . 'VALEUR' = dprmt ;
      B . 'COF' . 1 . 'COMPOR' = 'IDEN' ;
      B . 'COF' . 1 . 'LDAT'   = 'LECT' 1 ;
      ll = 'LECT' 1 ;
   'SINON' ;
      ll = 'LECT' ;
   'FINSI' ;
   lvq = 'EGA' ('TYPE' valq) 'LISTREEL' ;
   'REPETER' iiinc ninc ;
      iinc = &iiinc ;
      B . 'VAR' . iinc . 'NOMDDL' = tnomq . 'NOMDUA' . iinc ;
      B . 'VAR' . iinc . 'DISC'   = tnomq . 'DISC'  ;
      'SI' lvalq ;
         'SI' lvq ;
            B . 'VAR' . iinc . 'VALEUR'   = 'EXTRAIRE' valq iinc ;
         'SINON' ;
            B . 'VAR' . iinc . 'VALEUR'   = valq ;
         'FINSI' ;
      'FINSI' ;
      B . iinc . iinc . 0  = ll ;
   'FIN' iiinc ;
*
   'SI' ('EXISTE' _smt) ;
       mgugrad2 = 'NLIN' discg _mt _smt A B methgau ;
   'SINON' ;
       mgugrad2 = NLINP discg _mt A B methgau ;
   'FINSI' ;
*
'RESPRO' mgugrad2 ;
'FINPROC' ;
*
* End of procedure file GUGRAD2
*
*ENDPROCEDUR gugrad2
*BEGINPROCEDUR gvol
************************************************************************
* NOM         : GVOL
* DESCRIPTION :
* Calcule le volume compris dans une surface fermée
* La normale doit être vers l'intérieur pour que le volume soit positif
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 22/04/2011, version initiale
* HISTORIQUE : v1, 22/04/2011, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
*
*
'DEBPROC' GVOL ;
'ARGUMENT' _surf*'MAILLAGE' ;
'ARGUMENT' tdisc*'TABLE' ;
'ARGUMENT' dbg/'LOGIQUE' ;
*
'SI' ('NON' ('EXISTE' dbg)) ;
   dbg = FAUX ;
'FINSI' ;
*
*  Vecteur position et calcul du volume
NOMVIT = @STBL (TDISC . 'XN' . 'NOMPRI') ;
DISCG  = TDISC . 'GEOM' . 'DISC' ;
vdim = 'VALEUR' 'DIME' ;
   'SI' ('ET' ('EGA' vdim 2) ('EGA' vmod 'AXIS')) ;
      fdim = 3 ;
   'SINON' ;
      fdim = vdim ;
   'FINSI' ;
   vposc = GETCOO _surf nomvit ;
*   'SI' iaxi ;
*      rs zs = 'COORDONNEE' _surf ;
*      nr = 'EXTRAIRE' nomvit 1 ;
*      nz = 'EXTRAIRE' nomvit 2 ;
*      vposc =
*   'FINSI' ;
*   fvol = GNOR _surf tdisc 'NPRI' ('CHAINE' discg 'V') 'CPRI' vpos
*                          'NDUA' 'CSTEV' ;
*   rfvol = 'RESULT' fvol ;
*   volx = 'MAXIMUM' ('EXCO' 'UX' rfvol) ;
*   voly = 'MAXIMUM' ('EXCO' 'UY' rfvol) ;
*   vol = '/' ('+' volx voly) vdim ;
   fvolc = GNOR _surf tdisc 'NPRI' discg
                          'NCOF' (chai discg 'V')
                          'CCOF' vposc
                          'NDUA' (chai discg 'V')
                          'FDUA' ('PROG' vdim * 1.) ;
   volc = '/' ('MAXIMUM' ('RESULT' fvolc))
              fdim ;
   vol = volc '*' -1. ;
'RESPRO' vol ;
*
* End of procedure file GVOL
*
'FINPROC' ;
*ENDPROCEDUR gvol
*BEGINPROCEDUR log10
************************************************************************
* NOM         : LOG10
* DESCRIPTION : Log_10
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 18/02/2003, version initiale
* HISTORIQUE : v1, 18/02/2003, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' LOG10 ;
'REPETER' bouc ;
   ok = FAUX ;
   'ARGUMENT' fl/'FLOTTANT' ;
   'ARGUMENT' lr/'LISTREEL' ;
   'ARGUMENT' cp/'CHPOINT ' ;
   'ARGUMENT' cm/'MCHAML  ' ;
   'SI' ('EXISTE' fl) ;
      ok = VRAI ;
      'RESPRO' ('/' ('LOG' fl) ('LOG' 10.D0)) ;
   'FINSI' ;
   'SI' ('EXISTE' lr) ;
      ok = VRAI ;
      'RESPRO' ('/' ('LOG' lr) ('LOG' 10.D0)) ;
   'FINSI' ;
   'SI' ('EXISTE' cp) ;
      ok = VRAI ;
      'RESPRO' ('/' ('LOG' cp) ('LOG' 10.D0)) ;
   'FINSI' ;
   'SI' ('EXISTE' cm) ;
      ok = VRAI ;
      'RESPRO' ('/' ('LOG' cm) ('LOG' 10.D0)) ;
   'FINSI' ;
   'SI' ('NON' ok) ;
      'QUITTER' bouc ;
   'FINSI' ;
'FIN' bouc ;
*
* End of procedure file LOG10
*
'FINPROC' ;
*ENDPROCEDUR log10
*BEGINPROCEDUR modenlin
************************************************************************
* NOM         : MODENLIN
* DESCRIPTION :
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, ??/??/2007, version initiale
* HISTORIQUE : v1, ??/??/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' MODENLIN ;
*
vdim = 'VALEUR' 'DIME' ;
vmod = 'VALEUR' 'MODE' ;
*
'SI' ('ET' ('NEG' vdim 2) ('NEG' vdim 3)) ;
   'ERREUR' ('CHAINE' 'vdim=' vdim) ;
'FINSI' ;
*
*  Noms par défaut
*
'SI' ('EGA' vdim 2) ;
   'SI' ('NEG' vmod 'AXIS') ;
      nomvec = 'MOTS' 'UX' 'UY' ;
      nomfor = 'MOTS' 'FX' 'FY' ;
   'SINON' ;
      nomvec = 'MOTS' 'UR' 'UZ' ;
      nomfor = 'MOTS' 'FR' 'FZ' ;
   'FINSI' ;
'FINSI' ;
'SI' ('EGA' vdim 3) ;
   nomvec = 'MOTS' 'UX' 'UY' 'UZ' ;
   nomfor = 'MOTS' 'FX' 'FY' 'FZ' ;
'FINSI' ;
*
nomsca = 'MOTS' 'SCAL' ;
*nomflu = 'MOTS' 'SCAL' ;
*
'ARGUMENT' tmode/'TABLE' ;
'SI' ('NON' ('EXISTE' tmode)) ;
   tmode = 'TABLE' ;
*
* Initialisation des inconnues par défaut
*
   lmdisc = 'MOTS' 'CSTE' 'LINE' 'QUAF' ;
   'REPETER' iidisc ('DIME' lmdisc) ;
      mdisc = 'EXTRAIRE' lmdisc &iidisc ;
      tmode . mdisc = 'TABLE' ;
      tmode . mdisc . 'DISC' = mdisc ;
      tmode . mdisc . 'NOMPRI' = 'TABLE' ;
      tmode . mdisc . 'NOMPRI' . 1 = nomsca ;
      tmode . mdisc . 'NOMDUA' = 'TABLE' ;
      tmode . mdisc . 'NOMDUA' . 1 = nomsca ;
      mdiscv = 'CHAINE' mdisc 'V' ;
      tmode . mdiscv = 'TABLE' ;
      tmode . mdiscv . 'DISC' = mdisc ;
      tmode . mdiscv . 'NOMPRI' = 'TABLE' ;
      'REPETER' idim vdim ;
         TMODE . mdiscv. 'NOMPRI' . &idim =
            'MOTS' ('EXTRAIRE' nomvec &idim) ;
      'FIN' idim ;
      tmode . mdiscv . 'NOMDUA' = 'TABLE' ;
      'REPETER' idim vdim ;
         TMODE . mdiscv. 'NOMDUA' . &idim =
            'MOTS' ('EXTRAIRE' nomfor &idim) ;
      'FIN' idim ;
   'FIN' iidisc ;
'FINSI' ;
*
* Lecture des mots clés et des inconnues
*
lmotcle = 'MOTS' 'GEOM' 'INCO' ;
ldiscdd = 'MOTS' 'LINM' 'CUBI' ;
ltypinc = 'MOTS' 'SCAL' 'VECT' ;
'REPETER' imotcle ;
   'ARGUMENT' motcle/'MOT' ;
   'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
   'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
      cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
   'FINSI' ;
   'SI' ('EGA' motcle 'GEOM') ;
      'ARGUMENT' discg*'MOT' ;
      TMODE . 'GEOM' = 'TABLE' ;
      TMODE . 'GEOM' . 'DISC' = discg ;
   'FINSI' ;
   'SI' ('EGA' motcle 'INCO') ;
      'ARGUMENT' nominc*'MOT' ;
      'SI' ('EGA' nominc 'GEOM') ;
         'ERREUR' 'GEOM nest pas un nom dinconnu acceptable' ;
      'FINSI' ;
      TMODE . nominc = 'TABLE' ;
      'ARGUMENT' disinc*'MOT' ;
      TMODE . nominc . 'DISC' = disinc ;
      ldd = 'EXISTE' ldiscdd disinc ;
      'ARGUMENT' typinc*'MOT' ;
      'SI' ('NON' ('EXISTE' ltypinc typinc)) ;
         cherr = 'CHAINE' 'Type=' typinc ' unknown.' ;
         'ERREUR' cherr ;
      'FINSI' ;
      'SI' ('EGA' typinc 'SCAL') ;
         nbinc = 1 ;
      'SINON' ;
         nbinc = vdim ;
      'FINSI' ;
      TMODE . nominc . 'NOMPRI' = 'TABLE' ;
      'REPETER' ibinc nbinc ;
         'SI' ldd ;
            'ARGUMENT' nomcom*'LISTMOTS' ;
         'SINON' ;
            'ARGUMENT' nomcom*'MOT' ;
            nomcom = 'MOTS' nomcom ;
         'FINSI' ;
         TMODE . nominc . 'NOMPRI' . &ibinc = nomcom ;
      'FIN' ibinc ;
      TMODE . nominc . 'NOMDUA' = 'TABLE' ;
      'REPETER' ibinc nbinc ;
         'SI' ldd ;
            'ARGUMENT' nomcom*'LISTMOTS' ;
         'SINON' ;
            'ARGUMENT' nomcom*'MOT' ;
            nomcom = 'MOTS' nomcom ;
         'FINSI' ;
         TMODE . nominc . 'NOMDUA' . &ibinc = nomcom ;
      'FIN' ibinc ;
   'FINSI' ;
'FIN' imotcle ;
*
'RESPRO' tmode ;
*
* End of procedure file MODENLIN
*
'FINPROC' ;
*ENDPROCEDUR modenlin
*BEGINPROCEDUR modulo
************************************************************************
* NOM         : MODULO
* DESCRIPTION : Calcule un entier modulo un autre...
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 15/10/2002, version initiale
* HISTORIQUE : v1, 15/10/2002, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' MODULO ;
'ARGUMENT' i*'ENTIER' j*'ENTIER' ;
'SI' ('EGA' j 0) ;
  'MESSAGE' 'Impossible de faire modulo 0' ;
  'ERREUR' 5 ;
'SINON' ;
  k=i '/' j ;
  mod=i '-' ( k '*'j ) ;
  'RESPRO' mod ;
'FINSI' ;
*
* End of procedure file MODULO
*
'FINPROC' ;
*ENDPROCEDUR modulo
*BEGINPROCEDUR myent
************************************************************************
* NOM         : MYENT
* DESCRIPTION :
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, ??/??/2007, version initiale
* HISTORIQUE : v1, ??/??/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' MYENT ;
'ARGUMENT' rr*'FLOTTANT' ;
'SI' ('>EG' rr 0.) ;
   ii = 'ENTIER' ('+' rr 0.5)  ;
'SINON' ;
   ii = 'ENTIER' ('-' rr 0.5)  ;
'FINSI' ;
'RESPRO' ii ;
*
* End of procedure file MYENT
*
'FINPROC' ;
*ENDPROCEDUR myent
*BEGINPROCEDUR projsysc
************************************************************************
* NOM         : PROJSYSC
* DESCRIPTION : Calcul matrice et second membre projetés suivant
*               un champ de directions données
*
*               Project a linear system with respect to a given
*               vector field
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 22/04/2011, version initiale
* HISTORIQUE : v1, 22/04/2011, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
*
*
'DEBPROC' PROJSYSC ;
'ARGUMENT' tdisc*'TABLE' ;
'ARGUMENT' vdirn*'CHPOINT' ;
'ARGUMENT' ktgra*'RIGIDITE' ;
'ARGUMENT' fpgrad*'CHPOINT' ;
'ARGUMENT' kvol/'CHPOINT' ;
lcnt = 'EXISTE' kvol ;
'SI' lcnt ;
   'ARGUMENT' dvol*'FLOTTANT' ;
'FINSI' ;
vdim = 'VALEUR' 'DIME' ;
NOMVIT = @STBL (TDISC . 'XN' . 'NOMPRI') ;
NOMFOR = @STBL (TDISC . 'XN' . 'NOMDUA') ;
*   fpgrad = 'PSCAL' fpgra vdirn nomfor nomvit ;
* Condensation de la matrice
   knord = 'MANUEL' 'RIGI' vdirn ;
   knord = 'CHANGER' 'INCO' knord nomvit nomvit
         nomfor nomvit  ;
*tc mise en commentaires du finsi cidessous
*   'FINSI' ;
   'SI' ('EGA' vdim 2) ;
      nomscal = 'MOTS' 'SCAL' 'SCAL' ;
   'SINON' ;
      nomscal = 'MOTS' 'SCAL' 'SCAL' 'SCAL' ;
   'FINSI' ;
   knor  = 'CHANGER' 'INCO' knord nomvit nomscal nomvit nomvit  ;
   knort = 'CHANGER' 'INCO' knord nomvit nomvit  nomvit nomscal ;
*   'MESSAGE' 'ktgra' ;
*   'LISTE' ('EXTRAIRE' ktgra 'COMP') ;
*   'LISTE' ('EXTRAIRE' ktgra 'COMP' 'DUAL') ;
*   'MESSAGE' 'knort' ;
*   'LISTE' ('EXTRAIRE' knort 'COMP') ;
*   'LISTE' ('EXTRAIRE' knort 'COMP' 'DUAL') ;
   ktg1  = 'KOPS' 'CMCT' ktgra knort ;
*   'MESSAGE' 'ktg1' ;
*   'LISTE' ('EXTRAIRE' ktg1 'COMP') ;
*   'LISTE' ('EXTRAIRE' ktg1 'COMP' 'DUAL') ;
   ktg2 = 'KOPS' 'TRANSPOS' ktg1 ;
*   'MESSAGE' 'ktg2' ;
*   'LISTE' ('EXTRAIRE' ktg2 'COMP') ;
*   'LISTE' ('EXTRAIRE' ktg2 'COMP' 'DUAL') ;
   ktg3 = 'KOPS' 'CMCT' knort ktg2 ;
*   'MESSAGE' 'ktg3' ;
*   'LISTE' ('EXTRAIRE' ktg3 'COMP') ;
*   'LISTE' ('EXTRAIRE' ktg3 'COMP' 'DUAL') ;
   ktot = ktg3 ;
   ftot = fpgrad ;
   'SI' lcnt ;
      ktvol = 'PSCAL' kvol vdirn nomvit nomvit ;
      ktv = rela ('NOMC' 'T' ktvol) ;
      smbvol = 'DEPIMPOSE' ktv dvol ;
      ktv1 = 'EXTRAIRE' ktv 'T' 'FLX' ;
      ktv2 = 'EXTRAIRE' ktv 'LX' 'Q' ;
      ktv = ktv1 'ET' ktv2 ;
      lpr1 = 'MOTS' 'T' 'LX' ; l2 = 'MOTS' 'SCAL' 'LXPD' ;
      ldu1 = 'MOTS' 'Q' 'FLX' ;
      ktv = 'CHANGER' 'INCO' ktv lpr1 l2 ldu1 l2 'MULT' ;
      smbvol = 'NOMC' ldu1 l2 smbvol ;
      ktot = ktot 'ET' ktv ;
      ftot = ftot '+' smbvol ;
   'FINSI' ;
'RESPRO' ktot ftot ;
*
* End of procedure file PROJSYSC
*
'FINPROC' ;
*ENDPROCEDUR projsysc
*BEGINPROCEDUR quafme
************************************************************************
* NOM         : QUAFME
* DESCRIPTION :
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 01/12/2004, version initiale
* HISTORIQUE : v1, 01/12/2004, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' QUAFME ;
'ARGUMENT' tmail/'TABLE' ;
'SI' ('NON' ('EXISTE' tmail)) ;
   tmail = 'TABLE' ;
'FINSI' ;
'REPETER' bcl ;
   'ARGUMENT' mquad/'MAILLAGE' ;
   'SI' ('EXISTE' mquad) ;
      'REPETER' bbcl 1 ;
         'ARGUMENT' nom*'MOT' ;
         tmail . nom = 'TABLE' ;
         dimm = DEADUTIL 'DIMM' mquad ;
         'SI' ('EGA' dimm 0) ;
            mquaf = mquad ;
            tmail . nom = mquad ;
            'QUITTER' bbcl ;
         'FINSI' ;
         typm = DEADUTIL 'TYPM' mquad ;
         'SI' ('EGA' typm 'QUAF') ;
            mquaf = mquad ;
            tmail . nom . 'QUAF' = mquaf ;
            'SI' ('EGA' dimm 2) ;
               tmail . nom . 'QUAI' = mquaf ;
               tmail . nom . 'QUAD' = mquaf ;
            'FINSI' ;
            'QUITTER' bbcl ;
         'FINSI' ;
         'SI' ('EGA' typm 'QUAI') ;
            mquaf = 'CHANGER' mquad 'QUAF' ;
            tmail . nom . 'QUAF' = mquaf ;
            tmail . nom . 'QUAI' = mquad ;
            tmail . nom . 'LINE' = 'CHANGER' 'LINEAIRE' mquad ;
            'QUITTER' bbcl ;
         'FINSI' ;
         'SI' ('EGA' typm 'LINE') ;
            mquai = 'CHANGER' 'QUADRATIQUE' mquad ;
            mquaf = 'CHANGER' mquai 'QUAF' ;
            tmail . nom . 'QUAF' = mquaf ;
            tmail . nom . 'QUAI' = mquai ;
            tmail . nom . 'LINE' = mquad ;
            'QUITTER' bbcl ;
         'FINSI' ;
      'FIN' bbcl ;
*      ielim = '+' ielim 1 ;
*      telim . ielim = mquaf ;
   'SINON' ;
      'QUITTER' bcl ;
   'FINSI' ;
'FIN' bcl ;
*
*
telim = 'TABLE' 'ESCLAVE' ;
ielim = 0 ;
tidx = 'INDEX' tmail ;
'REPETER' iidx ('DIME' tidx) ;
   idx = tidx . &iidx ;
   val = tmail . idx ;
   tval = 'TYPE' val ;
   'SI' ('EGA' tval 'MAILLAGE') ;
      ielim = '+' ielim 1 ;
      telim . ielim = val ;
   'FINSI' ;
   'SI' ('EGA' tval 'TABLE') ;
      'SI' ('EXISTE' val 'QUAF') ;
         ielim = '+' ielim 1 ;
         telim . ielim = val . 'QUAF' ;
      'FINSI' ;
   'FINSI' ;
'FIN' iidx ;
*
*
'ARGUMENT' tol*'FLOTTANT' ;
'ELIMINATION' ('ETG' telim) tol ;
*
'RESPRO' tmail ;
'FINPROC' ;
*
* End of procedure file QUAFME
*
*ENDPROCEDUR quafme
*BEGINPROCEDUR trvec
************************************************************************
* NOM         : TRVEC
* DESCRIPTION : Trace des champs de vecteurs.
*               Utile pour tracer des bilans de forces
*
*               Display vector fields.
*               Useful for visualization of force balance.
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 22/04/2011, version initiale
* HISTORIQUE : v1, 22/04/2011, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
*
*
'DEBPROC' TRVEC ;
*'ARGUMENT' tdisc*'TABLE' ;
*'ARGUMENT' motdom*'MOT' ;
'ARGUMENT' tdom*'MAILLAGE' ;
'ARGUMENT' nocomp*'LISTMOTS' ;
*
'ARGUMENT' tvec/'TABLE' ;
'SI' ('EXISTE' tvec) ;
   'ARGUMENT' ttit*'TABLE' ;
   i = 'DIME' tvec ;
   'SI' ('EXISTE' tvec 'SOUSTYPE') ; i = '-' i 1 ; 'FINSI' ;
   'SI' ('EXISTE' tvec 'CREATEUR') ; i = '-' i 1 ; 'FINSI' ;
'SINON' ;
   tvec = 'TABLE' ; ttit = 'TABLE' ;
   i = 0 ;
'FINSI' ;
*
'REPETER' livec ;
   'SI' ('EGA' i 0) ;
      'ARGUMENT' ccvec*'CHPOINT' ;
   'SINON' ;
      'ARGUMENT' ccvec/'CHPOINT' ;
   'FINSI' ;
   'SI' ('EXISTE' ccvec) ;
      'ARGUMENT' ttvec*'MOT' ;
   'SINON' ;
      'QUITTER' livec ;
   'FINSI' ;
   i = '+' i 1 ;
*   'MESSAGE' ('CHAINE' 'i=' i) ;
*   'LISTE' ccvec ;
*   'LISTE' tvec ;
   tvec . i = ccvec ;
   ttit . i = ttvec ;
'FIN' livec ;
'ARGUMENT' echv/'FLOTTANT' ;
'ARGUMENT' lnclk/'LOGIQUE' ;
'SI' ('NON' ('EXISTE' lnclk)) ;
   lnclk = faux ;
'FINSI' ;
lvec = faux ; lnor = faux ; lvno = faux ;
lmotcle = 'MOTS' 'VECT' 'NORM' 'VNOR' ;
nmotcle = 0 ;
'REPETER' imotcle ;
   'ARGUMENT' motcle/'MOT' ;
   'SI' ('NON' ('EXISTE' motcle)) ; 'QUITTER' imotcle ; 'FINSI' ;
   'SI' ('NON' ('EXISTE' lmotcle motcle)) ;
      cherr = 'CHAINE' 'Keyword ' motcle ' unknown.' ; 'ERREUR' cherr ;
   'FINSI' ;
   nmotcle = '+' nmotcle 1 ;
*   'MESSAGE' ('CHAINE' 'nmotcle=' nmotcle ' ' 'motcle=' motcle) ;
   'SI' ('EGA' motcle 'VECT') ; lvec = VRAI ; 'FINSI' ;
   'SI' ('EGA' motcle 'NORM') ; lnor = VRAI ; 'FINSI' ;
   'SI' ('EGA' motcle 'VNOR') ; lvno = VRAI ; 'FINSI' ;
*   'MESSAGE' ('CHAINE' '   lvec=' lvec) ;
*   'MESSAGE' ('CHAINE' '   lnor=' lnor) ;
*   'MESSAGE' ('CHAINE' '   lvno=' lvno) ;
'FIN' imotcle ;
'SI' ('EGA' nmotcle 0) ; lvec = vrai ; 'FINSI' ;
*'MESSAGE' ('CHAINE' 'lvec=' lvec) ;
*'MESSAGE' ('CHAINE' 'lnor=' lnor) ;
*'MESSAGE' ('CHAINE' 'lvno=' lvno) ;
lcoul = 'MOTS' 'JAUN' 'ROUG' 'BLAN' 'TURQ' 'VERT' 'OLIV'
 'AZUR' 'ORAN' 'VIOL' 'GRIS' 'OCEA' ;
*
*
'SI' ('OU' lnor lvno) ;
* Table des normes des champs
   tnvec = 'TABLE' ;
   dtvec = 'DIME' tvec ;
   'SI' ('EXISTE' tvec 'SOUSTYPE') ; dtvec = '-' dtvec 1 ; 'FINSI' ;
   'SI' ('EXISTE' tvec 'CREATEUR') ; dtvec = '-' dtvec 1 ; 'FINSI' ;
   'REPETER' ivec dtvec ;
      vec  = tvec . &ivec ;
*   nomdep = 'EXTRAIRE' vec 'COMP' ;
      nvec = '**' ('PSCAL' vec vec nocomp nocomp) 0.5 ;
      tnvec . &ivec = nvec ;
   'FIN' ivec ;
* Trouver la norme maxi locale
   nloc = tnvec . 1 ;
   'SI' ('>' dtvec 1) ;
      'REPETER' ivec ('-' dtvec 1) ;
         iivec = &ivec '+' 1 ;
         vec  = tnvec . iivec ;
         vec0  = '+' vec ('*' nloc 0.) ;
         nloc0 = '+' nloc ('*' vec 0.) ;
         maq = 'MASQUE' vec0 'SUPERIEUR' nloc0 ;
         imaq = ('-' maq 1.) '*' -1. ;
         nnloc = '+' ('*' maq vec) ('*' imaq nloc) ;
         nloc = nnloc ;
      'FIN' ivec ;
   'FINSI' ;
* Normer les champs
   nloc = '+' nloc 1.D-60 ;
   mm = 'MAXIMUM' nloc ;
*nloc = '+' nloc ('*' mm 1.D-4) ;
*lnloc = LOG10 nloc ;
   tvec2 = 'TABLE' ;
   'SI' lvno ;
      lnloc = nloc ;
      inloc = 'INVERSE' nloc ;
      'REPETER' ivec dtvec ;
         vec  = tvec . &ivec ;
         tvec2 . &ivec = '*' vec inloc ;
      'FIN' ivec ;
   'FINSI' ;
'SINON' ;
*
   lmax = 'PROG' ;
   'REPETER' ii i ;
      mm = 'MAXIMUM' (tvec . &ii) 'ABS' ;
      lmax = 'ET' lmax ('PROG' mm) ;
   'FIN' ii ;
   mm = '+' ('MAXIMUM' lmax) 1.D-60 ;
'FINSI' ;
*
'SI' ('NON' ('EXISTE' echv)) ;
   vmodo = 'EGA' ('VALEUR' 'MODE') 'AXIS' ;
   'SI' vmodo ; 'OPTI' 'MODE' 'PLAN' ; 'FINSI' ;
   mmas = gmass2 ('CHANGER' tdom 'QUAF') tdisc
      'NPRI' discv 'NDUA' discv ;
   mmasd = gmass2 ('CHANGER' tdom 'QUAF') tdisc
      'NPRI' 'CSTE' 'FPRI' 1.D0 'NDUA' discv ;
   mdiam = 'EXTR' mmas 'DIAG' ;
   mmdi = 'MAXIMUM' ('RESULT' mdiam) ;
   mmma = 'MAXIMUM' ('RESULT' mmasd) ;
   mdias = '*' mdiam ('/' mmma mmdi) ;
*   ctail = gmass2 ('CHANGER' tdom 'QUAF') tdisc
*      'NPRI' 'CSTE' 'FPRI' 1. 'NDUA' 'CSTE' 'FDUA' 1. ;
*   'LISTE' ctail ;
   'SI' vmodo ; 'OPTI' 'MODE' 'AXIS' ; 'FINSI' ;
   vdim = 'VALEUR' 'DIME' ;
   dimm = DEADUTIL 'DIMM' tdom ;
*   ctail = '**' ctail ('/' 1. ('-' vdim 1)) ;
   ctail = '**' mdias ('/' 1. dimm) ;
   tail = '**' ('*' ('MAXIMUM' ctail) ('MINIMUM' ctail)) 0.5 ;
*   tail =  'MINIMUM' ctail ;
*'LISTE' tail ;
*'LISTE' mm ;
   echv = '/' ('*' 1.4 tail) mm ;
*   echvn = ('*' 0.95 tail) ;
   echvn = 0.75 ;
   'SI' lvno ;
      'REPETER' ivec dtvec ;
         vec  = tvec2 . &ivec ;
         tvec2 . &ivec = '*' vec ctail ;
      'FIN' ivec ;
   'FINSI' ;
'FINSI' ;
*'MESSAGE' ('CHAINE' 'mm=' mm) ;
'SI' lvec ; tvec3 = 'TABLE' ; 'FINSI' ;
'REPETER' ii i ;
   cou = EXMOMOD lcoul &ii ;
   'SI' lvec ;
      tvec3 . &ii = 'VECT' (tvec . &ii) echv nocomp cou ;
   'FINSI' ;
   'SI' lvno ;
      tvec2 . &ii = 'VECT' (tvec2 . &ii) echvn nocomp cou ;
   'FINSI' ;
   'SI' ('EGA' &ii 1) ;
      tit = 'CHAINE' cou '=' (ttit . &ii) ;
   'SINON' ;
      tit = 'CHAINE' tit ' ' cou '=' (ttit . &ii) ;
   'FINSI' ;
'FIN' ii ;
idim = 'VALEUR' 'DIME' ;
vdim = DEADUTIL 'DIMM' tdom ;
ctdom = tdom ;
'SI' ('EGA' vdim 2) ;
   ctdom= 'CONTOUR' tdom ;
'FINSI' ;
'SI' ('EGA' vdim 3) ;
   ctdom= ARET tdom ;
'FINSI' ;
*
'SI' lvec ;
   titv  = 'CHAINE' 'Max=' (formar mm 1) ' ' tit ;
'FINSI' ;
'SI' lnor ;
   titn = tit ;
'FINSI' ;
'SI' lvno ;
   titvn = 'CHAINE' 'Norm' ' ' tit ;
'FINSI' ;
*
'SI' lnclk ;
   'SI' lvec ;
      'TRACER' (@stbl tvec3) tdom ctdom 'TITR' titv 'NCLK' ;
   'FINSI' ;
   'SI' lnor ;
      'TRACER' nloc tdom ctdom 'TITR' titn 'NCLK' ;
   'FINSI' ;
   'SI' lvno ;
      'TRACER' (@stbl tvec2) tdom ctdom 'TITR' titvn 'NCLK' ;
   'FINSI' ;
'SINON' ;
   'SI' lvec ;
      'TRACER' (@stbl tvec3) tdom ctdom 'TITR' titv ;
   'FINSI' ;
   'SI' lnor ;
      'TRACER' nloc tdom ctdom 'TITR' titn ;
   'FINSI' ;
   'SI' lvno ;
      'TRACER' (@stbl tvec2) tdom ctdom 'TITR' titvn ;
   'FINSI' ;
'FINSI' ;
*
* End of procedure file TRVEC
*
'FINPROC' ;
*ENDPROCEDUR trvec
*BEGINPROCEDUR tsurktan
************************************************************************
* NOM         : TSURKTAN
* DESCRIPTION : La matrice tangente pour la tension de surface
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 17/01/2007, version initiale
* HISTORIQUE : v1, 17/01/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' TSURKTAN ;
'ARGUMENT' _mt*'MAILLAGE' ;
'ARGUMENT' gdisc*'MOT' ;
'ARGUMENT' methgau*'MOT' ;
'ARGUMENT' dppri*'LISTMOTS' ;
'ARGUMENT' dpdua*'LISTMOTS' ;
*
dpdis = gdisc ;
*
idim = 'VALEUR' 'DIME' ;
vdim = DEADUTIL 'DIMM' _mt ;
laxi = DEADUTIL 'AXI?' ;
lsph = DEADUTIL 'SPH?' ;
*
loi = 'CHAINE' 'TSUJ' ;
loij = 'CHAINE' 'TSU' ;
*
'ARGUMENT' coef/'FLOTTANT' ;
   'SI' ('NON' ('EXISTE' coef)) ;
      'ARGUMENT' coef2/'CHPOINT ' ;
      'SI' ('NON' ('EXISTE' coef2)) ;
         'ERREUR' 'Il faut donner un coef FLOTTANT ou CHPOINT' ;
      'SINON' ;
         coef = coef2 ;
         'ARGUMENT' discc*'MOT     ' ;
      'FINSI' ;
   'SINON' ;
      discc = 'CSTE' ;
   'FINSI' ;
*
'ARGUMENT' jaco/'ENTIER' ;
'SI' ('NON' ('EXISTE' jaco)) ;
   jaco = 1 ;
   dir1 = VRAI ;
'FINSI' ;
'SI' ('OU' ('EGA' jaco 2) ('EGA' jaco 3)) ;
   'ARGUMENT' idir/'ENTIER' ;
   'SI' ('EXISTE' idir) ;
      ldir = 'LECT' idir ;
      'SI' ('EGA' idir 1) ;
         dir1 = VRAI ;
      'FINSI' ;
   'SINON' ;
      'ARGUMENT' ldir/'LISTENTI' ;
      'SI' ('NON' ('EXISTE' ldir)) ;
         ldir = 'LECT' 1 'PAS' 1 idim ;
         dir1 = VRAI ;
      'FINSI' ;
   'FINSI' ;
'FINSI' ;
'ARGUMENT' lterm/'LISTENTI' ;
llterm = 'EXISTE' lterm ;
'SI' llterm ;
   dlterm = 'DIME' lterm ;
'SINON' ;
   dlterm = 1 ;
'FINSI' ;
*
* Calcul du jacobien complet (jaco = 1)
*
'SI' ('EGA' jaco 1) ;
   numop  = '*' ('**' vdim 2) ('**' idim 2) ;
   numop = '*' numop dlterm ;
   'SI' ('OU' laxi lsph) ;
      numop  = '+' numop ('*' (vdim '*' idim) 2) ;
   'FINSI' ;
   'SI' lsph ;
      numop = '+' numop 1 ;
   'FINSI' ;
   numder = vdim ;
   numvar = idim ;
   numdat = 0 ;
   numcof = 0 ;
*
   A = ININLIN numop numvar numdat numcof numder ;
   numdat = 1 ;
   numcof = numop ;
   B = ININLIN numop numvar numdat numcof numder ;
   'REPETER' ivar numvar ;
      A . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dppri &ivar) ;
      A . 'VAR' . &ivar . 'DISC'   = dpdis ;
      B . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dpdua &ivar) ;
      B . 'VAR' . &ivar . 'DISC'   = dpdis ;
   'FIN' ivar ;
   iop = 0 ;
   'REPETER' h dlterm ;
      'REPETER' i idim ;
         'REPETER' j vdim ;
            'REPETER' k idim ;
               'REPETER' l vdim ;
                  iop = iop '+' 1 ;
                  A . iop . &i . &j  = 'LECT' ;
                  'SI' llterm ;
                     nl = 'EXTRAIRE' lterm &h ;
                     nomloi = 'CHAINE' loij nl &i &j &k &l ;
                  'SINON' ;
                     nomloi = 'CHAINE' loi &i &j &k &l ;
                  'FINSI' ;
*              'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
                  B . 'COF' . iop . 'COMPOR' = nomloi ;
                  B . 'COF' . iop . 'LDAT'   = 'LECT' 1 ;
                  B . iop . &k . &l  = 'LECT' iop ;
               'FIN' l ;
            'FIN' k ;
         'FIN' j ;
      'FIN' i ;
   'FIN' h ;
   'SI' ('OU' laxi lsph) ;
      'REPETER' i idim ;
         'REPETER' j vdim ;
            iop = iop '+' 1 ;
            A . iop . &i . &j  = 'LECT' ;
            nomloi = 'CHAINE' loi &i &j '10' ;
*              'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
            B . 'COF' . iop . 'COMPOR' = nomloi ;
            B . 'COF' . iop . 'LDAT'   = 'LECT' 1 ;
            B . iop . 1 . 0  = 'LECT' iop ;
         'FIN' j ;
      'FIN' i ;
      'REPETER' k idim ;
         'REPETER' l vdim ;
            iop = iop '+' 1 ;
            A . iop . 1 . 0  = 'LECT' ;
            nomloi = 'CHAINE' loi '10' &k &l ;
*              'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
            B . 'COF' . iop . 'COMPOR' = nomloi ;
            B . 'COF' . iop . 'LDAT'   = 'LECT' 1 ;
            B . iop . &k . &l  = 'LECT' iop ;
         'FIN' l ;
      'FIN' k ;
   'FINSI' ;
   'SI' lsph  ;
      iop = iop '+' 1 ;
      A . iop . 1 . 0  = 'LECT' ;
      nomloi = 'CHAINE' loi '1010' ;
*      'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
      B . 'COF' . iop . 'COMPOR' = nomloi ;
      B . 'COF' . iop . 'LDAT'   = 'LECT' 1 ;
      B . iop . 1 . 0  = 'LECT' iop ;
   'FINSI' ;
'FINSI' ;
'SI' ('EGA' jaco 2) ;
   nldir = 'DIME' ldir ;
   numop  = '*' nldir ('**' vdim 2) ;
   'SI' dir1 ;
      'SI' ('OU' laxi lsph) ;
         numop = '+' numop ('*' vdim 2) ;
      'FINSI' ;
      'SI' lsph ;
         numop = '+' numop 1 ;
      'FINSI' ;
   'FINSI' ;
*
   numder = vdim ;
   numvar = idim ;
   numdat = 0 ;
   numcof = 0 ;
*
   A = ININLIN numop numvar numdat numcof numder ;
   numdat = 1 ;
   numcof = numop ;
   B = ININLIN numop numvar numdat numcof numder ;
   'REPETER' ivar numvar ;
      A . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dppri &ivar) ;
      A . 'VAR' . &ivar . 'DISC'   = dpdis ;
      B . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dpdua &ivar) ;
      B . 'VAR' . &ivar . 'DISC'   = dpdis ;
   'FIN' ivar ;
   iop = 0 ;
   'REPETER' i nldir ;
      idir = 'EXTRAIRE' ldir &i ;
      'REPETER' j vdim ;
         'REPETER' l vdim ;
            iop = iop '+' 1 ;
            A . iop . idir . &j  = 'LECT' ;
            nomloi = 'CHAINE' loi idir &j idir &l ;
*               'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
            B . 'COF' . iop . 'COMPOR' = nomloi ;
            B . 'COF' . iop . 'LDAT'   = 'LECT' 1 ;
            B . iop . idir . &l  = 'LECT' iop ;
         'FIN' l ;
      'FIN' j ;
   'FIN' i ;
   'SI' dir1 ;
      'SI' ('OU' laxi lsph) ;
         'REPETER' j vdim ;
            iop = iop '+' 1 ;
            A . iop . 1 . &j  = 'LECT' ;
            nomloi = 'CHAINE' loi '1' &j '10' ;
*               'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
            B . 'COF' . iop . 'COMPOR' = nomloi ;
            B . 'COF' . iop . 'LDAT'   = 'LECT' 1 ;
            B . iop . 1 . 0  = 'LECT' iop ;
         'FIN' j ;
         'REPETER' l vdim ;
            iop = iop '+' 1 ;
            A . iop . 1 . 0 = 'LECT' ;
            nomloi = 'CHAINE' loi '101' &l ;
*               'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
            B . 'COF' . iop . 'COMPOR' = nomloi ;
            B . 'COF' . iop . 'LDAT'   = 'LECT' 1 ;
            B . iop . 1 . &l  = 'LECT' iop ;
         'FIN' l ;
      'FINSI' ;
      'SI' lsph ;
         iop = iop '+' 1 ;
         A . iop . 1 . 0  = 'LECT' ;
         nomloi = 'CHAINE' loi '1010' ;
*      'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
         B . 'COF' . iop . 'COMPOR' = nomloi ;
         B . 'COF' . iop . 'LDAT'   = 'LECT' 1 ;
         B . iop . 1 . 0  = 'LECT' iop ;
      'FINSI' ;
   'FINSI' ;
'FINSI' ;
*
'SI' ('EGA' jaco 3) ;
   nldir = 'DIME' ldir ;
*   numop  = '**' vdim 2 ;
   numop  = '*' nldir vdim ;
   'SI' ('ET' dir1 lsph) ;
      numop = '+' numop 1 ;
   'FINSI' ;
*
   numder = vdim ;
   numvar = idim ;
   numdat = 0 ;
   numcof = 0 ;
*
   A = ININLIN numop numvar numdat numcof numder ;
   numdat = 1 ;
   numcof = numop ;
   B = ININLIN numop numvar numdat numcof numder ;
   'REPETER' ivar numvar ;
      A . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dppri &ivar) ;
      A . 'VAR' . &ivar . 'DISC'   = dpdis ;
      B . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dpdua &ivar) ;
      B . 'VAR' . &ivar . 'DISC'   = dpdis ;
   'FIN' ivar ;
   iop = 0 ;
   'REPETER' i nldir ;
      idir = 'EXTRAIRE' ldir &i ;
      'REPETER' j vdim ;
         iop = iop '+' 1 ;
         A . iop . idir . &j  = 'LECT'  ;
         nomloi = 'CHAINE' loi idir &j idir &j ;
*               'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
         B . 'COF' . iop . 'COMPOR' = nomloi ;
         B . 'COF' . iop . 'LDAT'   = 'LECT' 1 ;
         B . iop . idir . &j  = 'LECT' iop ;
      'FIN' j ;
   'FIN' i ;
   'SI' ('ET' dir1 lsph) ;
      iop = iop '+' 1 ;
      A . iop . 1 . 0  = 'LECT' ;
      nomloi = 'CHAINE' loi '1010' ;
*      'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
      B . 'COF' . iop . 'COMPOR' = nomloi ;
      B . 'COF' . iop . 'LDAT'   = 'LECT' 1 ;
      B . iop . 1 . 0  = 'LECT' iop ;
   'FINSI' ;
'FINSI' ;
*
* Partie commune
*
   B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
   B . 'DAT' . 1 . 'DISC'   = discc ;
   B . 'DAT' . 1 . 'VALEUR' = coef ;
*
jac = 'NLIN' gdisc _mt A B 'EREF' methgau ;
*
'RESPRO' jac ;
*
* End of procedure file TSURKTAN
*
'FINPROC' ;
*ENDPROCEDUR tsurktan
*BEGINPROCEDUR tsurresi
************************************************************************
* NOM         : TSURRESI
* DESCRIPTION : Le résidu à annuler pour la tension de surface
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 17/01/2007, version initiale
* HISTORIQUE : v1, 17/01/2007, création
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
* 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' TSURRESI ;
'ARGUMENT' _mt*'MAILLAGE' ;
'ARGUMENT' gdisc*'MOT' ;
'ARGUMENT' methgau*'MOT' ;
'ARGUMENT' dpdua*'LISTMOTS' ;
*
dpdis = gdisc ;
*
idim  = 'VALEUR' 'DIME' ;
vdim  = DEADUTIL 'DIMM' _mt ;
laxi = DEADUTIL 'AXI?' ;
lsph = DEADUTIL 'SPH?' ;
*
loi = 'CHAINE' 'TSUR' ;
*
'ARGUMENT' coef/'FLOTTANT' ;
   'SI' ('NON' ('EXISTE' coef)) ;
      'ARGUMENT' coef2/'CHPOINT ' ;
      'SI' ('NON' ('EXISTE' coef2)) ;
         'ERREUR' 'Il faut donner un coef FLOTTANT ou CHPOINT' ;
      'SINON' ;
         coef = coef2 ;
         'ARGUMENT' discc*'MOT     ' ;
      'FINSI' ;
   'SINON' ;
      discc = 'CSTE' ;
   'FINSI' ;
*
dir1 = FAUX ;
'ARGUMENT' idir/'ENTIER' ;
'SI' ('EXISTE' idir) ;
   ldir = 'LECT' idir ;
   'SI' ('EGA' idir 1) ;
      dir1 = VRAI ;
   'FINSI' ;
'SINON' ;
   'ARGUMENT' ldir/'LISTENTI' ;
   'SI' ('NON' ('EXISTE' ldir)) ;
      ldir = 'LECT' 1 'PAS' 1 idim ;
      dir1 = VRAI ;
   'FINSI' ;
'FINSI' ;
*
* Calcul du résidu
*
nldir = 'DIME' ldir ;
*
numop  = '*' nldir vdim ;
term1 = ('OU' laxi lsph) 'ET' dir1 ;
'SI' term1 ;
   numop = '+' numop 1 ;
'FINSI' ;
numder = vdim ;
numvar = 1 ;
numdat = 0 ;
numcof = 0 ;
*
   A = ININLIN numop numvar numdat numcof numder ;
   A . 'VAR' . 1 . 'NOMDDL' = 'MOTS' 'DUMM' ;
   A . 'VAR' . 1 . 'DISC'   = 'CSTE' ;
   A . 'VAR' . 1 . 'VALEUR' = 1.D0 ;
*
   numvar = idim ;
   numdat = 1 ;
   numcof = numop ;
   B = ININLIN numop numvar numdat numcof numder ;
   'REPETER' ivar numvar ;
      B . 'VAR' . &ivar . 'NOMDDL' = 'MOTS' ('EXTRAIRE' dpdua &ivar) ;
      B . 'VAR' . &ivar . 'DISC'   = dpdis ;
   'FIN' ivar ;
*
   B . 'DAT' . 1 . 'NOMDDL' = 'MOTS' 'SCAL' ;
   B . 'DAT' . 1 . 'DISC'   = discc ;
   B . 'DAT' . 1 . 'VALEUR' = coef ;
*
   iop = 0 ;
   'REPETER' k nldir ;
      idir = 'EXTRAIRE' ldir &k ;
      'REPETER' l vdim ;
         iop = '+' iop 1 ;
         A . iop . 1 . 0  = 'LECT' ;
         nomloi = 'CHAINE' loi idir &l ;
*         'MESSAGE' ('CHAINE' 'Nomloi=' nomloi) ;
         B . 'COF' . iop . 'COMPOR' = nomloi ;
         B . 'COF' . iop . 'LDAT'   = 'LECT' 1 ;
         B . iop . idir . &l  = 'LECT' iop ;
      'FIN' l ;
   'FIN' k ;
*   'LISTE' A ; 'LISTE' iop ;
   'SI' term1 ;
      iop = '+' iop 1 ;
      A . iop . 1 . 0  = 'LECT' ;
      nomloi = 'CHAINE' loi '10' ;
      B . 'COF' . iop . 'COMPOR' = nomloi ;
      B . 'COF' . iop . 'LDAT'   = 'LECT' 1 ;
      B . iop . 1 . 0  = 'LECT' iop ;
   'FINSI' ;
*
res = 'NLIN' gdisc _mt A B 'EREF' methgau ;
*
'RESPRO' res ;
*
* End of procedure file TSURRESI
*
'FINPROC' ;
*ENDPROCEDUR tsurresi
*BEGINPROCEDUR gmail
************************************************************************
* NOM         : GMAIL
* DESCRIPTION : Extrait le maillage des points supports des degrés
*               de liberté d'une inconnue
*
*               Cette procédure est identique à GMAIL2 qui paraissait
*               mieux que l'ancienne GMAIL
*
*
*
* LANGAGE     : GIBIANE-CAST3M
* AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
*               mél : gounand@semt2.smts.cea.fr
**********************************************************************
* VERSION    : v1, 20/04/2012, version initiale
* HISTORIQUE :
* HISTORIQUE :
************************************************************************
*
*
'DEBPROC' GMAIL ;
'ARGUMENT'     _mt/'MAILLAGE' ;
'ARGU' tdisc/'TABLE' ;
xmt    = 'EXIS' _mt ;
xtdisc = 'EXIS' tdisc ;
*
'SI' ('NON' xmt) ;
   'SI' ('NON' xtdisc) ;
      'ARGUMENT' tdisc*'TABLE' ;
      xtdisc = vrai ;
   'FINS' ;
   'ARGUMENT' nmt*'MOT' ;
   'SI' ('NON' ('EXISTE' tdisc nmt)) ;
      cherr = 'CHAINE' 'Le domaine ' nmt ' nest pas defini dans'
      ' la table de discretisation' ;
      'ERREUR' cherr ;
   'FINSI' ;
   _mt = tdisc . nmt . 'QUAF' ;
'FINSI' ;
*
'ARGU' dua*'MOT' ;
'SI' ('NON' xtdisc) ;
   discdua = dua ;
*   nomdua  = dua ;
'SINO' ;
*   nomdua  = dua ;
   discdua = tdisc . dua . 'DISC' ;
'FINS' ;
*
*'ARGUMENT'    nomdua*'MOT     ' ;
*
* Raccourcis du 2019/01/25
* Est-ce vraiment astucieux ?...
*
'SI' ('EGA' discdua 'QUAF') ;
   'RESP' _mt ;
   'QUIT' GMAIL ;
'FINS' ;
'SI' ('EGA' discdua 'LINE') ;
   'RESP' ('CHAN' _mt 'LINE') ;
   'QUIT' GMAIL ;
'FINS' ;
'SI' ('EGA' discdua 'QUAI') ;
   'RESP' ('CHAN' _mt 'QUAD') ;
   'QUIT' GMAIL ;
'FINS' ;
*
mdim = DEADUTIL 'DIMM' _mt ;
*vdim = 'VALEUR' 'DIME' ;
'SI' ('EGA' discdua 'LINM') ;
   discdua = 'CSTE' ;
'FINSI' ;
*
numop = 1 ; numder = mdim ; numvar = 1 ; numdat = 0 ; numcof = 0 ;
A = ININLIN numop numvar numdat numcof numder ;
A . 'VAR' . 1 . 'VALEUR' = 0. ;
*
A . 1 . 1 . 0 = 0 ;
*
B = ININLIN numop numvar numdat numcof numder ;
B . 'VAR' . 1 . 'NOMDDL' = 'DUMM' ;
B . 'VAR' . 1 . 'DISC'   = discdua ;
*
B . 1 . 1 . 0 = 0 ;
*
gm = 'NLIN' 'LINE' _mt A B 'ERF1' 'GAU1' ;
dom = 'EXTRAIRE' gm 'MAIL' ;
*
'RESPRO' dom ;
'FINPROC' ;
*
* End of procedure file GMAIL
*
*ENDPROCEDUR gmail
**
************************************************************************
*
*
*   END OF PROCEDURES
*
*
************************************************************************
************************************************************************
*
*
*   MAIN : 1) MESH
*          2) COMPUTATIONAL LOOP
*          3) TESTs
*
************************************************************************
*
* Construction du "modèle" (maillage)
* et des paramètres de départ
*
'SI' complet ;
   raff = 2 ;
   nitermax = 25 ;
*   nitermax = 1 ;
   critnf = 5.D-5 ;
   critdf = 5.D-5 ;
*   critnf = 5.D-4 ;
*   critdf = 5.D-4 ;
'SINON' ;
   raff = 2 ;
   nitermax = 15 ;
   critnf = 5.D-3 ;
   critdf = 5.D-3 ;
'FINSI' ;
*
* Après ce 'SI' on met tous les paramètres constants
* et d'initialisation du calcul
*
lrest = 'EGA' ('TYPE' tdisc) 'TABLE' ;
'SI' ('NON' lrest) ;
*
* Paramètres de discrétisation spatiale
*
discg = 'QUAF' ; discv = 'QUAF' ; discp = 'LINM' ;
*
*
* Longueurs
*
*
rsou = 1. ; lav = 20. ; lap = 20. ; prof = 20. ;
diam = '*' rsou 2. ;
rzoom = 5. ; lrzoom = 'PROG' ('*' rzoom -1.) rzoom ;
*
*
* Création du maillage
*
tdisc = DEFMAIL rsou lav lap prof raff ;
tdisc . 'methgau' = 'TABLE' ;
tdisc . 'methgau' . 'mass' = 'GAM2' ;
tdisc . 'methgau' . 'amor'  = 'GAM2' ;
tdisc . 'methgau' . 'rigi'  = 'GAR2' ;
*tdisc . 'methgau' . 'mass' = 'GAU7' ;
*tdisc . 'methgau' . 'amor'  = 'GAU7' ;
*tdisc . 'methgau' . 'rigi'  = 'GAU7' ;
*
DISCG  = TDISC . 'GEOM' . 'DISC' ;
DISCU  = TDISC . 'XN' . 'DISC' ;
DISCV  = TDISC . 'UN' . 'DISC' ;
NOMDEP = @STBL (TDISC . 'XN' . 'NOMPRI') ;
NOMVIT = @STBL (TDISC . 'UN' . 'NOMPRI') ;
NOMPRE = @STBL (TDISC . 'PN' . 'NOMPRI') ;
NOMFOR = @STBL (TDISC . 'XN' . 'NOMDUA') ;
*
*
_mtw  = tdisc . 'mtw' . 'QUAF' ;
_cmtw  = tdisc . 'cmtw' . 'QUAF' ;
_hau = tdisc . 'hau' . 'QUAF' ;
_gau = tdisc . 'gau' . 'QUAF' ;
_bas = tdisc . 'bas' . 'QUAF' ;
_dro = tdisc . 'dro' . 'QUAF' ;
mtw = tdisc . 'mtw' . discg ;
cmtw = tdisc . 'cmtw' . discg ;
mtwi = 'DIFF' ('CHANGER' mtw 'POI1')
              ('CHANGER' cmtw 'POI1') ;
hau = tdisc . 'hau' . discg ;
ihau = 'INVERSE' hau ;
gau = tdisc . 'gau' . discg ;
bas = tdisc . 'bas' . discg ;
dro = tdisc . 'dro' . discg ;
haub = 'ET' ('POIN' hau 'INITIAL')
    ('POIN' hau 'FINAL') ;
'SI' ('EGA' vdim 3) ;
   _fro = tdisc . 'fro' .'QUAF' ;
    fro = tdisc . 'fro' . discg  ;
   _rea = tdisc . 'rea' .'QUAF' ;
    rea = tdisc . 'rea' . discg  ;
   _dhau = tdisc . 'dhau' .'QUAF' ;
    dhau = tdisc . 'dhau' . discg  ;
   _bhau = tdisc . 'bhau' .'QUAF' ;
    bhau = tdisc . 'bhau' . discg  ;
    amtw = 'ARETE' mtw ;
'FINSI' ;
*
* CB215821 : Ancienne syntaxe supprimee en 2015
*pA = 'POIN' hau  'PROC' ('POIN' ('PROG' vdim * 0.)) ;
*mpA = 'COULEUR' ('MANUEL' 'POI1' pA) 'ROUG' ;
* Maillage pour les blocages en vitesse
*'TRACER' 'CACH' mtw 'TITR' ('CHAINE' 'nbl=' ('NBEL' mtw))
*   'NCLK' ;
*
'SI' ('EGA' vdim 2) ;
   pA  = 'POIN' hau  'PROC' ('POIN' 0. 0. ) ;
   mpA = 'COULEUR' ('MANUEL' 'POI1' pA) 'ROUG' ;
 
   mblux = gau ;
   mblux = 'COULEUR' ('CHANGER' 'POI1' mblux) 'ROUG' ;
   mblun = 'COULEUR' ('CHANGER' 'POI1' hau) 'ROUG' ;
*   mbluy = gau 'ET' bas 'ET' dro ;
   mbluy = gau 'ET' bas ;
   mbluy = 'CHANGER' mbluy 'POI1' ;
   mpoi = 'INTERSECTION' mblun mbluy ;
   'SI' ('NEG' ('NBEL' mpoi) 0) ;
      mblun = 'DIFF' mblun mpoi ;
   'FINSI' ;
   mblun = 'COULEUR' mblun 'ROUG' ;
   mbluy = 'COULEUR' mbluy 'ROUG' ;
*mblpn = 'POIN' (GMAIL _mtr 'PN') 'INITIAL' ;
'FINSI' ;
'SI' ('EGA' vdim 3) ;
   pA  = 'POIN' hau  'PROC' ('POIN' 0. 0. 0.) ;
   mpA = 'COULEUR' ('MANUEL' 'POI1' pA) 'ROUG' ;
 
   mblux = gau ;
   mblux = 'COULEUR' ('CHANGER' 'POI1' mblux) 'ROUG' ;
   mbluy = gau 'ET' fro 'ET' rea ;
   mbluy = 'COULEUR' ('CHANGER' 'POI1' mbluy) 'ROUG' ;
   mblun = 'COULEUR' ('CHANGER' 'POI1' hau) 'ROUG' ;
*! On fait ceci pour permettre l'élimination des relations
*!   mbluz = gau 'ET' bas ;
   mbluz = gau 'ET' bas 'ET' bhau ;
   mbluz = 'CHANGER' mbluz 'POI1' ;
   mpoi = 'INTERSECTION' mblun mbluz ;
   'SI' ('NEG' ('NBEL' mpoi) 0) ;
      mblun = 'DIFF' mblun mpoi ;
   'FINSI' ;
   mblun = 'COULEUR' mblun 'ROUG' ;
   mbluz = 'COULEUR' mbluz 'ROUG' ;
*mblpn = 'POIN' (GMAIL _mtr 'PN') 'INITIAL' ;
'FINSI' ;
 
'SI' faux ;
   mdd =mtw ;
   'SI' ('EGA' vdim 3);  mdd = cmtw ; 'FINSI' ;
   'TRACER' ('ET' mdd mpA) 'TITR' 'Point dobservation pA' ;
   'TRACER' ('ET' mdd mblux) 'TITR' 'Blocage UX' ;
   'TRACER' ('ET' mdd mbluy) 'TITR' 'Blocage UY' ;
   'TRACER' ('ET' mdd mblun) 'TITR' 'Blocage UN' ;
'FINSI' ;
*   'OPTI' 'MODE' plan ;
met = calimet _mtw discg 'GAU1' ;
************************************************************************
*
* COMPUTATIONAL LOOP
*
************************************************************************
* Re : nombre de Reynolds (convection / viscosité)
* Eu : nombre d'Euler (pression imposée / convection)
* FrEu : nombre de Froude * nombre d'Euler (pression imposée / gravité)
* We : (pression imposée / tension de surface)
*
'DEBPROC' calcul ;
'ARGUMENT' Re*'FLOTTANT' ;
'ARGUMENT' Eu*'FLOTTANT' ;
'ARGUMENT' FrEu*'FLOTTANT' ;
'ARGUMENT' We*'FLOTTANT' ;
'ARGUMENT' direct*'LOGIQUE' ;
'ARGUMENT' miter*'LOGIQUE' ;
*
*** Test 1  Re=50. Eu=10000.  FrEu=1. We=10000.
***         La pression imposée est équilibrée par le poids
***         de fluide déplacé. L'écoulement n'a pas
***         d'influence sur la forme de la surface.
*** Test 2  Re=50. Eu=10000.  FrEu=1. We=1.
***         La pression imposée est équilibrée par le poids
***         de fluide déplacé et la tension de surface.
***         L'écoulement n'a pas d'influence sur la forme de la surface.
*** Test 3  Re=50. Eu=1.  FrEu=1. We=1.
***         La pression imposée est équilibrée par le poids
***         de fluide déplacé et la tension de surface.
***         L'écoulement influence la forme de la surface (aspiration).
*
*Re= 50. ; Eu=1. ; FrEu=1. ; We = 1. ;
*Re= 50. ; Eu=10000. ; FrEu=1. ; We = 10000. ;
*Re = 200. ; Eu=1. ; FrEu=1. ; We = 0.5 ;
*Re = 250. ; Eu=1. ; FrEu=1. ; We = 1. ;
*
*
* Numérique
*
* dblcond : Double condensation des relations
* maspvid : matrice de pression vide
* diagrig : Rajoute la diagonale de la rigidité dans la matrice tangente
* direct : vrai = solveur direct
*          faux = méthode de projection
* tbfor : faux = on prend en compte les forces de surface lorsqu'on
*                résout le problème sur le déplacement de la surface
*         vrai = on prend en compte les forces de surface lorsqu'on
*                résout l'écoulement
* tjaco : vrai = on teste que les matrices jacobiennes calculées
*                sont correctes en comparant avec une différence finie
* tsurb : faux = on supprime la composante 'UX' de la force de tension
*                de surface sur les bords gauches et droits
* lpdro : vrai = on met une force de pression sur la surface de sortie
*                chargée d'équilibrer l'incompressibilité du volume
* miter : vrai = méthodes itérative pour résoudre les systèmes linéaires
*         faux = méthode direct
* soldir : 1   = Solveur direct de KRES
*          0   = Solveur direct de RESOU
* solvimp : niveau d'impression pour les solveurs linéaires.
* ltime   : affichage des temps CPU pour les solveurs linéaires.
*
*direct  = vrai ; miter   = faux ;
dblcond = vrai ;
diagrig = ('NON' direct) ; adgr = 1. ;
solvimp = 0 ; ltime = faux ; soldir = 0 ;
maspvid = direct 'ET' miter ;
tbfor = faux ;
tjaco = FAUX ;
tsurb = faux ;
lpdro = vrai ;
idir = 1 ; jacoxf = 3 ; jacoxg = 0 ; jacoxv = 0 ;
           jacoxt = 'LECT' 1 2 3 4;
*
* Paramètres du schéma non linéaire
*
*omegn = 0.5 ; nitn = 20 ; critn = 5.D-3 ;
*omegd = 0.5 ; nitd = 20 ; critd = 1.D-1 ; nbck = 16 ;
*omegn = 0.5 ; nitn = 5 ; critn = 1.D-3 ;
*omegd = 0.75 ; nitd = 20 ; critd = 1.D-3 ; nbck = 1 ;
omegn = 0.5 ; nitn = 2 ; critn = 1.D-4 ;
omegd = 0.75 ; nitd = 20 ; critd = 1.D-4 ; nbck = 1 ;
*omegn = 0.75 ; nitn = 2 ; critn = 1.D-4 ;
*omegd = 0.85 ; nitd = 20 ; critd = 1.D-4 ; nbck = 1 ;
*
* Paramètres de résolution des problèmes linéaires (cf. notice KRES)
*  * 1 Le problème d'écoulement en méthode couplée
'SI' direct ;
   tqp = 'TABLE' 'METHINV' ;
   tqp . 'IMPR' = solvimp ; tqp . 'LTIME' =  ltime ;
   'SI' miter ;
** Le solveur de RESO ne fonctionne pas pour celui-là
      tqp . 'TYPINV' = 3 ;
   'SINON' ;
      tqp . 'TYPINV' = 1 ;
   'FINSI' ;
*  * 1bis Le problème d'écoulement en méthode de projection
'SINON' ;
   tqproj = 'TABLE' ;
   tqproj . 'debug' = faux ;
   tqproj . 'dblproj' = vrai ;
   tqproj . 'scaldia' = vrai ;
   tqproj . 'METHINVQDM'  = 'TABLE' 'METHINV' ;
   tqdm  = tqproj . 'METHINVQDM' ;
   tqdm . 'IMPR' = solvimp ;   tqdm . 'LTIME' =  ltime ;
   tqproj . 'METHINVPRES' = 'TABLE' 'METHINV' ;
   tpres = tqproj . 'METHINVPRES' ;
   tpres . 'IMPR' = solvimp ;   tpres . 'LTIME' =  ltime ;
   'SI' miter ;
      tqdm  . 'TYPINV' = 3 ; tqdm . 'NITMAX' =500 ;
      tpres . 'TYPINV' = 3 ; tpres . 'NITMAX' =500 ;
   'SINON' ;
      tqdm  . 'TYPINV' = soldir ;
      tpres . 'TYPINV' = soldir ;
   'FINSI' ;
'FINSI' ;
*  * 2 Le problème de déplacement de la surface
tds = 'TABLE' 'METHINV' ;
tds . 'IMPR' = solvimp ; tds . 'LTIME' = ltime ;
tds . 'LDEPE' = FAUX ;
'SI' ('NON' miter) ;
** Le solveur de RESO ne fonctionne pas forcément pour celui-là
   tds . 'TYPINV' = 1 ;
'SINON' ;
   tds . 'TYPINV' = 3 ; tds . 'IMPR' = solvimp ;
   tds . 'NITMAX' =500 ;
'FINSI' ;
*  * 3 Le problème de déplacement du maillage (DEDU ADAP)
tdm = 'TABLE' 'METHINV' ;
tdm . 'IMPR' = solvimp ; tdm . 'LTIME' = ltime ;
'SI' ('NON' miter) ;
   tdm . 'TYPINV' = soldir ;
'SINON' ;
   tdm . 'TYPINV' = 2 ;   tdm . 'NITMAX' =500 ;
*   tdm . 'TYPINV' = 7 ;   tdm . 'NITMAX' =500 ;
'FINSI' ;
*
ckonv = '/' 1. Eu ;
*ckonv = ckonv '*' 0. ;
crig  = '/' 1. ('*' Re Eu) ;
cgra  = '/' 1. FrEu ;
ctsur = '/' 1. We. ;
*
nomdir = 'CHAINE' '/test4/gounand/kong/' ;
nomfic = 'CHAINE' 'defila2' vdim 'd'
                        'Re' (MYENT Re)
                        'Eu' (MYENT Eu)
                        'FE' (MYENT FrEu)
                        'We' (MYENT We)
                        'r' raff ;
*************************************************************
*         Début de l'algorithme de Résolution               *
vol = GVOL _cmtw tdisc faux ;
'SI' debug ;
'MESSAGE' ('CHAINE' 'Volume  initial = ' (formar vol 2)) ;
'FINSI' ;
pini = 0. ;
delp = 0. ; delp2 = 0. ;
lza  = 'PROG' ;
*
*sol   = vide 'CHPO' ;
*sol   = 'EXCO' 'UX' (GCHPO _mtw 'UN' 1.) 'UX' ;
mailun = GMAIL _mtw tdisc 'UN' ;
sol = 'MANU' 'CHPO' mailun 2 'UX' 1. 'UY' 0. ;
itcou = 0 ;
*
* Ce FINSI correspond à 'SI' ('NON' lrest) ;
*tc mise en commentaire du finsi ci dessous
*
*'FINSI' ;
*DEFSUMM ;
*
* Ici, on calcule les termes qui ne dépendent pas de la déformation
* du maillage et de la solution
*
*    1 Eventuelle matrice de pression vide
'SI' maspvid ;
    mpvid = GMASS2 _mtw tdisc 'NPRI' 'PN' 'NDUA' 'PN'
      'NCOF' 'CSTE' 'FCOF' 0. ;
'FINSI' ;
*    2 Matrices de blocage en x et y (voire z et P)
cblox = '+' ('MANUEL' 'CHPO' mblux 1 'UX' 0.)
            ('MANUEL' 'CHPO' gau 1 'UX' 1.) ;
mblox = 'BLOQUE' 'UX' mblux ;
smblox = 'DEPIMPOSE' mblox cblox ;
*
mbloy = 'BLOQUE' 'UY' mbluy ;
'SI' ('EGA' vdim 3) ;
   mbloz = 'BLOQUE' 'UZ' mbluz ;
'FINSI' ;
*  mp = GMAIL _mtw 'PN' ;
* mblop = 'BLOQUE' 'T' ('POIN' mp 'INITIAL') ;
*
* Boucle d'itérations globales
*
'REPETER' it nitermax ;
*      'MESSAGE' ('CHAINE' 'itcou = ' ('+' itcou 1)) ;
*
*  1      Résolution Navier-Stokes Incompressible (surface fixée)
*
   'SI' lpdro ; delp = '+' delp delp2 ; 'FINSI' ;
*
*  1.1  Calcul des termes LINEAIRES de Navier-Stokes
*
*  1.1.1 Pression imposée en surface
*
   'SI' ('EGA' vdim 2) ;
      rhau = 'COORDONNEE' 1 hau ;
   'SINON' ;
*       rhau = 'COORDONNEE' 1 hau ;
      xhau = 'COORDONNEE' 1 hau ;
      yhau = 'COORDONNEE' 2 hau ;
      rhau = '**' ('+' ('**' xhau 2) ('**' yhau 2)) 0.5 ;
   'FINSI' ;
   pfor = 'EXP' ('*' ('**' ('/' rhau rsou) 2) -1.) ;
   fpfor = GFORC _hau tdisc pfor ;
   ktforx = GKFORC _hau tdisc pfor jacoxf ;
*   ktforx0 = GKFORC _hau tdisc pfor 0 ;
   'SI' ('EGA' jacoxf 3) ;
      ktforx = '*' ktforx 0. ;
   'FINSI' ;
*   TRVEC hau fpfor 'FPFOR' ;
*
*  1.1.2 Potentiel de gravité imposé en surface
*
   fpgra = GGRAVI _hau tdisc cgra 0. ;
   ktgrax = GKGRAVI _hau tdisc jacoxg ('*' cgra -1.) 0. ;
*    ktgrax0 = GKGRAVI _hau tdisc 0 ('*' cgra -1.) 0. ;
*
*  1.1.3 Forces de tension de surface
*
   mgtens = tdisc . 'methgau' . 'rigi' ;
   ftsuri = TSURRESI _hau discg mgtens nomfor ('*' -1. ctsur) ;
   'SI' ('NON' tsurb) ;
      ftsurx = 'EXCO' 'FX' ftsuri 'FX' ;
      ftsura = 'REDU' ftsurx haub ;
      ftsur = '-' ftsuri ftsura ;
   'SINON' ;
      ftsur = ftsuri ;
   'FINSI' ;
   ktsurx = TSURKTAN _hau discg mgtens nomvit nomfor ('*' -1. ctsur)
      jacoxt ;
*   ktsurx0 = TSURKTAN _hau discg mgtens nomvit nomvit
*      ('*' -1. ctsur) ('LECT' 1 2 3 4) ;
*
*  1.1.4 Eventuel terme de pression sur le bord droit
*
*   'MESSAGE' ('CHAINE' 'delp=' delp) ;
*   'MESSAGE' ('CHAINE' 'delp2=' delp2) ;
   'SI' lpdro ;
      lpn = 'PROG' vdim * 0. ;
      'REMPLACER' lpn 1 ('*' delp -1.) ;
      ftpdr = GMASS2 _dro tdisc 'NPRI' 'UN'
         'FPRI' lpn
         'NDUA' 'UN' ;
   'FINSI' ;
*
*  1.1.5 Rigidité
*
   gri = GRIG2 _mtw tdisc 'NPRI' 'UN' 'NDUA' 'UN' ;
   gri = '*' gri crig ;
*
*  1.1.6 Couplage vitesse-pression (conservation de la masse)
*
   gdiv = GDIV2 _mtw tdisc 'NPRI' 'UN' 'NDUA' 'PN' 'GBBT' ;
   gdiv = '*' gdiv -1. ;
*
*  1.1.7 Eventuelle diagonale de la rigidité
*
   'SI' diagrig ;
      gridia = 'EXTRAIRE' gri 'DIAG' ;
      grid   = 'MANUEL' 'RIGI' gridia ;
      grid = grid '*' adgr ;
*      nel = 'NBEL' _mtw ;
*      dx  = '/' 1. ('**' nel ('/' 1. vdim)) ;
*      dt = '/' 0.01 dx ;
*      lidt = 'PROG' vdim * ('/' 1. dt) ;
*      'MESSAGE' ('CHAINE' 'Nel=' nel ' dx=' (formar dx 2)
*                                     ' dt=' (formar dt 2)) ;
*      gma = GMASS2 _mtw tdisc 'NPRI' 'UN' 'NDUA' 'UN'
*                              'NCOF' 'CSTE' 'FCOF' ('/' 1. dt) ;
   'FINSI' ;
*
*  1.2 Boucle d'itérations NON linéaires pour Navier-Stokes
*
   'REPETER' itn nitn ;
      sol = 'EXCO' (nomvit 'ET' nompre) sol 'NOID' ;
*
*  1.2.1 Matrice de convection
*
      gugr = gugrad2 _mtw tdisc 'NPRI' 'UN' 'NDUA' 'UN'
                               'NVIT' 'UN' 'CVIT' sol ;
      gugr = '*' gugr ckonv ;
*
*
*  1.2.2 Blocage de la surface dans la direction normale
*
      vnor = DEFDD tdisc 0 ;
*      TRVEC hau nomdep vnor 'Vnor' ;
      mblon = 'BLOQUE' 'DEPL' 'DIRE' vnor mblun ;
      mblot = mblox 'ET' mbloy 'ET' mblon ;
      'SI' ('EGA' vdim 3) ;
         mblot = mblot 'ET' mbloz ;
      'FINSI' ;
*
* 1.2.3 Matrices tangente et déséquilibre
*
      ktot = gugr 'ET' gri 'ET' gdiv 'ET' mblot ;
      ktg  = ('*' (gugr 'ET' gri) ('/' 1. omegn)) 'ET' gdiv
         'ET' mblot ;
      'SI' diagrig ;
         ktg  = ktg 'ET' ('*' grid ('/' 1. omegn)) ;
      'FINSI' ;
      'SI' maspvid ;
         ktg = ktg 'ET' mpvid ;
      'FINSI' ;
      ftot = smblox ;
      'SI' tbfor ;
         ftot = ftot '+' fpgra '+' ftsur '+' fpfor ;
      'FINSI' ;
      'SI' lpdro ;
         ftot = '+' ftot ftpdr ;
      'FINSI' ;
      ftot = '-' ftot ('*' ktot sol) ;
*
* 1.2.4 Résolution (si dblcond=vrai, on effectue une deuxième fois
*       l'élimination des contraintes, la première étant faite
*       dans KRES ou GRESPR. Ceci permet d'éliminer effectivement
*       toutes les relations de blocage suivant la normale)
*
      'OPTI' impi 0 ;
      'SI' dblcond ;
         ktgc ftotc ftot1 = 'KOPS' 'CONDENSE' ktg ftot ;
      'SINON' ;
         ktgc ftotc = ktg ftot ;
      'FINSI' ;
      'SI' direct ;
         dsolc = MONRES ktgc ftotc tqp ;
      'SINON' ;
         tqproj . 'noprec' =  'EGA' &itn 1 ;
         dsolc = GRESPR ktgc ftotc tdisc tqproj ;
      'FINSI' ;
*
      'SI' dblcond ;
         dsol = 'KOPS' 'EVAPORE' dsolc ktg ftot ftot1 ;
      'SINON' ;
         dsol = dsolc ;
      'FINSI' ;
*      HCHPO dsol 'dsol' ;
      'OPTI' impi 0 ;
      sol = '+' sol dsol ;
      vit = 'EXCO' nomdep sol ;
      pre = 'EXCO' nompre sol ;
      mdv  = 'MAXIMUM' dsol 'ABS' 'AVEC' nomvit ;
      mdvr = '/' mdv ('MAXIMUM' vit 'ABS') ;
      'SI' debug ;
      'MESSAGE' ('CHAINE' '  ' 'itn=' &itn
         ' mdv=' (formar mdv 2) ' mdvr=' (formar mdvr 2)) ;
*         HCHPO sol 'sol' ;
      'FINSI' ;
*
      'SI' ('EGA' vdim 2) ; rvit = vit ;
      'SINON' ; rvit = 'REDU' vit cmtw ; 'FINSI' ;
      'SI' (graph 'ET' interact) ;
         TRVEC mtw nomvit rvit 'Vit' vrai ;
      'FINSI' ;
      'SI' ('<' mdvr critn) ;
         'QUITTER' itn ;
      'FINSI' ;
   'FIN' itn ;
*
* 1.3 Post-traitement après calcul de la vitesse
*
* 1.3.1 Variables primales
*
   'SI' graphqdm ;
      DEFVIT ;
      DEFVSURF ;
*      DEFVSURF rzoom ;
   'FINSI' ;
*
* 1.3.2 Bilan des forces
*
   fpre  = 'EXCO' nomfor ('*' ('*' gdiv -1.) sol) ;
   fblot = 'REACTION' mblot sol ;
   frig = '*' ('*' gri -1.) sol ;
   fugr = '*' ('*' gugr -1.) sol ;
   tqf = 'TABLE' 'ESCLAVE' ; tqt = 'TABLE' 'ESCLAVE' ; i=1 ;
   tqf . i = fpre  ; tqt . i = 'pr' ; i = i '+' 1 ;
   tqf . i = frig  ; tqt . i = 'ri' ; i = i '+' 1 ;
   tqf . i = fugr  ; tqt . i = 'co' ; i = i '+' 1 ;
   tqf . i = fblot ; tqt . i = 'bl' ; i = i '+' 1 ;
   'SI' lpdro ;
      tqf . i = ftpdr ; tqt . i = 'ps' ; i = i '+' 1 ;
   'FINSI' ;
   'SI' tbfor ;
      tqf . i = fpgra ; tqt . i = 'gr' ; i = i '+' 1 ;
      tqf . i = ftsur ; tqt . i = 'ts' ; i = i '+' 1 ;
      tqf . i = fpfor ; tqt . i = 'pi' ; i = i '+' 1 ;
   'FINSI' ;
   desfq = 'ETG' tqf ;
   'SI' ('ET' graphqdm ('EGA' vdim 2)) ;
      tqf . i = '*' desfq -1. ; tqt . i = 'ds' ; i = i '+' 1 ;
      TRVEC mtw nomfor tqf tqt 'VECT' 'VNOR' ;
      tqfr = 'ASSI' 'TOUS' 'REDU' tqf mtwi ;
      TRVEC mtw nomfor tqfr tqt 'VECT' 'VNOR' ;
*
      xhau = 'COORDONNEE' 1 hau ;
      mhau = gmass2 _hau tdisc 'NPRI' discv 'NDUA' discv ;
      mvnor = '*' vnor -1. ;
      vtan = '*' ('PVEC' vnor nomdep nomdep) -1. ;
      DESSVEC ihau xhau 'x' 'Fnor' nomfor mvnor tqf tqt mhau ;
*      DESSVEC ihau xhau 'x' 'Fnor' nomfor mvnor tqf tqt mhau lrzoom ;
      DESSVEC ihau xhau 'x' 'Ftan' nomfor vtan tqf tqt mhau ;
*      DEFQFORV tqf tqt ;
*      DEFQFORS rzoom ;
   'FINSI' ;
*   'OPTION' 'DONN' 5 ;
*
*  2      Résolution Déplacement de la surface dans une direction
*         donnée (vdir) à partir du déséquilibre calculé
*         à l'étape 1 (normal à la surface).
*
*  2.1  Déséquilibre normal à la surface
*
   fblon  = 'REDU' fblot hau ;
   mfblon = '*' fblon -1. ;
*
*  2.2  Suivant tbfor, les forces normales "dans le continu" qui ne
*       travaillent donc pas quand la surface est fixée peuvent
*       être prises en compte, soit à l'étape 1, soit à l'étape 2
*       Il s'agit : de la tension de surface, de la pression imposée
*                   et du potentiel lié à la gravité
*       Il semble que la prise en compte à l'étape 2 soit plus stable
*       par contre, il y a léger déséquilibre tangentiel à la fin
*       du calcul FAUX ? tbfor marche maintenant
*
   vdir = '*' (DEFDD tdisc idir) -1. ;
*   TRVEC hau nomdep vdir 'Vdir' ;
   mfblond = 'PSCAL' mfblon vdir nomfor nomvit ;
   mfblond = 'CHANGER' 'ATTRIBUT' mfblond 'NATURE' 'DISCRET' ;
   'SI' tbfor ;
      desnd = mfblond ;
   'SINON' ;
      fpgrad = 'PSCAL' fpgra vdir nomfor nomvit ;
      fpgrad = 'CHANGER' 'ATTRIBUT' fpgrad 'NATURE' 'DISCRET' ;
      ftsurd = 'PSCAL' ftsur vdir nomfor nomvit ;
      ftsurd = 'CHANGER' 'ATTRIBUT' ftsurd 'NATURE' 'DISCRET' ;
      fpford = 'PSCAL' fpfor vdir nomfor nomvit ;
      fpford = 'CHANGER' 'ATTRIBUT' fpford 'NATURE' 'DISCRET' ;
      desnd  = mfblond '+' fpgrad '+' ftsurd '+' fpford ;
   'FINSI' ;
*   hchpo desn 'desn' ;
   ktotx = ktgrax 'ET' ktforx 'ET' ktsurx ;
   ktotx = '*' ktotx -1. ;
*
* 2.3 Contrainte sur le volume
*
   volc = GVOL _cmtw tdisc ;
   dvol = ('-' vol volc) '*' -1. ;
   kvol = GKVOL _cmtw tdisc jacoxv ;
*
   ok = FAUX ;
   omd = omegd ;
*
* 2.5 Résolution avec backtracking
*
   'REPETER' bck nbck ;
* 2.5.1 Réduction du système sur l'inconnue déplacement dans une direction
      ktot ftot = PROJSYSC tdisc vdir ('*' ktotx ('/' 1. omd))
         desnd kvol dvol ;
      sold = MONRES ktot ftot tds ;
      dep  = '*' ('EXCO' 'SCAL' sold) vdir ;
      mulag   = 'EXCO' 'LXPD' sold 'LXPD' ;
      delp2    = 'MAXIMUM' mulag ;
*
      mdep = 'MAXIMUM' dep 'ABS' ;
      mdepr = '/' mdep prof ;
      odep = dep ;
      'SI' debug ;
         chmes = 'CHAINE' ' dvol=' (formar dvol 1) ;
         chmes2 = 'CHAINE' '  mdep=' (formar mdep 1)
                       ' mdepr=' (formar mdepr 1)
                       ' delp2=' (formar delp2 1) ;
         'SI' ('>' &bck 1) ;
            chmes = 'CHAINE' '! bck=' &bck ' ' chmes ;
         'FINSI' ;
         'MESSAGE' chmes ; 'MESSAGE' chmes2 ;
      'FINSI' ;
      'SI' ('<' mdep critd) ;
         'QUITTER' bck ;
      'SINON' ;
         omd = '/' omd 2. ;
      'FINSI' ;
   'FIN' bck ;
*
* 2.6 Post-traitement après calcul du déplacement
*
*  Forces
   fpvolnd = '*' ktot ('*' mulag -1.) ;
*
   tdf = 'TABLE' 'ESCLAVE' ; tdt = 'TABLE' 'ESCLAVE' ; i=1 ;
   tdf . i = mfblond ; tdt . i = 'Fint'    ; i = i '+' 1 ;
   tdf . i = fpvolnd  ; tdt . i = 'Incomp.' ; i = i '+' 1 ;
   'SI' ('NON' tbfor) ;
      tdf . i = fpgrad ; tdt . i = 'RhoG'  ; i = i '+' 1 ;
      tdf . i = ftsurd ; tdt . i = 'Tsurf' ; i = i '+' 1 ;
      tdf . i = fpford ; tdt . i = 'SouP' ; i = i '+' 1 ;
   'FINSI' ;
   desfd = 'ETG' tdf ;
   'SI' ('ET' graphsur ('EGA' vdim 2)) ;
      tdf . i = desfd ; tdt . i = 'Des' ; i = i '+' 1 ;
*
      xhau = 'COORDONNEE' 1 hau ;
      mhau = gmass2 _hau tdisc 'NPRI' discv 'NDUA' discv ;
      DESSVEC ihau xhau 'x' 'Fdir' tdf tdt mhau ;
   'FINSI' ;
*
* 3       Déplacement du maillage
*
   'SI' ('EGA' vdim 2) ;
      muy = bas 'ET' hau ;
*   bux = 'BLOQUE' 'UX' mtw ;
      buy = 'BLOQUE' 'UY' muy ;
      btot = buy ;
      cblo =  ('MANUEL' 'CHPO' cmtw 2 'UX' 0. 'UY' 0.) '+'
         odep ;
      ftot = 'DEPIMPOSE' buy cblo ;
      dxv = 'DEDU' adap mtw btot ftot 'METR' met 'CSTE'
         'NITM' 1 'IDIR' vdim 'TINV' tdm ;
   'SINON' ;
      muz = bas 'ET' hau ;
*   bux = 'BLOQUE' 'UX' mtw ;
*   buy = 'BLOQUE' 'UY' mtw ;
      buz = 'BLOQUE' 'UZ' muz ;
*   btot = bux 'ET' buy 'ET' buz ;
      btot = buz ;
      cblo =  ('MANUEL' 'CHPO' muz 1 'UZ' 0.) '+'
         odep ;
      ftot = 'DEPIMPOSE' buz cblo ;
      dxv = 'DEDU' adap mtw btot ftot 'METR' met 'CSTE'
         'NITM' 1 'IDIR' vdim 'TINV' tdm ;
   'FINSI' ;
   'FORME' dxv ;
*
* 3.2 Post-traitement après déplacement du maillage
*
   za = 'COORDONNEE' vdim pA ;
   zcmt = 'COORDONNEE' vdim hau ;
   miz = 'MINIMUM' zcmt ; maz = 'MAXIMUM' zcmt ;
   miv = 'MINIMUM' vit ; mav = 'MAXIMUM' vit ;
   'MESSAGE' ('CHAINE' 'It=' ('+' itcou 1)
      ' zA=' (formar zA 2) ' miz=' (formar miz 2)
         ' maz=' (formar maz 2)  ' miv=' (formar miv 2)
         ' mav=' (formar mav 2) ) ;
   lza = 'ET' lza ('PROG' zA) ;
*
   'SI' graphsur ;
      'SI' ('EGA' vdim 2) ;
         lhau = ihau ;
      'SINON' ;
         lhau = dhau ;
      'FINSI' ;
      xhau = 'COORDONNEE' 1 lhau ;
      zhau = 'COORDONNEE' vdim lhau ;
      DESSVEC lhau xhau 'x' 'z(x)' zhau 'z' ;
      DESSVEC lhau xhau 'x' 'z(x)' zhau 'z' lrzoom ;
   'FINSI' ;
*
* 4       Critère de convergence
*
  'SI' ('ET' ('<' mdepr critdf) ('<' mdvr critnf)) ;
      itcou = '+' itcou 1 ;
      'QUITTER' it ;
   'FINSI' ;
   itcou = '+' itcou 1 ;
'FIN' it ;
*
* 5       Post traitement en fin de calcul
*
* 5.1 Altitude de A en fonction des itérations
*
'SI' graph ;
   dza = 'DIME' lza ;
   'SI' (> dza 1) ;
      lit = 'PROG' 1. 'PAS' 1 dza ;
      tabev = 'TABLE' ;   tabt  = 'TABLE' ;
      tabev . 1 = 'EVOL' 'MANU' lit lza ;
      tabt  . 1 = 'CHAINE' 'Z pA' ;
      tix = 'iter' ; tiy = 'z' ; tit = 'CHAINE' tiy '(' tix ')' ;
      dessevol (@STBL tabev) tabt tit tix tiy ;
   'FINSI' ;
*
* 5.2 Quantités primales
*
   DEFVIT ;
   DEFVSURF ;
*
* 5.3 Bilan global des forces :
*     Ce serait mieux de les recalculer ?
*     Mais bon, si on a convergé, c'est ok.
*
   fblor = fblot '+' mfblon ;
   tgf = 'TABLE' 'ESCLAVE' ; tgt = 'TABLE' 'ESCLAVE' ; i=1 ;
   tgf . i = fpre  ; tgt . i = 'pr' ; i = i '+' 1 ;
   tgf . i = frig  ; tgt . i = 'ri' ; i = i '+' 1 ;
   tgf . i = fugr  ; tgt . i = 'co' ; i = i '+' 1 ;
   tgf . i = fblor ; tgt . i = 'bl' ; i = i '+' 1 ;
   'SI' lpdro ;
      tgf . i = ftpdr ; tgt . i = 'ps' ; i = i '+' 1 ;
   'FINSI' ;
   tgf . i = fpgra ; tgt . i = 'gr' ; i = i '+' 1 ;
   tgf . i = ftsur ; tgt . i = 'ts' ; i = i '+' 1 ;
   tgf . i = fpfor ; tgt . i = 'pi' ; i = i '+' 1 ;
*
   desfg = 'ETG' tgf ;
   'SI' ('ET' graph ('EGA' vdim 2)) ;
      tgf . i = '*' desfg -1. ; tgt . i = 'ds' ; i = i '+' 1 ;
      TRVEC mtw nomfor tgf tgt 'VECT' 'VNOR' ;
      tgfr = 'ASSI' 'TOUS' 'REDU' tgf mtwi ;
      TRVEC mtw nomfor tgfr tgt 'VECT' 'VNOR' ;
*
      xhau = 'COORDONNEE' 1 hau ;
      mhau = gmass2 _hau tdisc 'NPRI' discv 'NDUA' discv ;
      mvnor = '*' vnor -1. ;
      vtan = '*' ('PVEC' vnor nomdep nomdep) -1. ;
      DESSVEC ihau xhau 'x' 'Fnor' nomfor mvnor tgf tgt mhau ;
      DESSVEC ihau xhau 'x' 'Ftan' nomfor vtan tgf tgt mhau ;
   'FINSI' ;
'FINSI' ;
*
* 6       Sauvegarde
*
'SI' lsauv ;
   'OUBL' KTG ; 'OUBL' KTOT ; 'OUBLIER' MHAU ;
   'OUBLIER' btot ;'OUBLIER' buz ; 'OUBLIER' mbloy ; 'OUBLIER' mblox ;
   'OUBL'  BUY  ; 'OUBL' BUX ; 'OUBL'     KTOTX ; 'OUBL'   KTPREX0 ;
   'OUBLIER' ktotr ; 'OUBLIER' ktgc ;
   'OUBL' KTPREX ;  'OUBL' KTUGRX ; 'OUBL' KTRIGX ;
   'OUBL'  MBLOT ; 'OUBL'   MBLON ; 'OUBL'   MBLOZ  ;
   'OUBL'  GUGR  ; 'OUBL'   GRI     ; 'OUBL'   GRID     ;
   'OUBL' KTSURX ; 'OUBL'  KTGRAX   ; 'OUBLIER' ktforx ;
   'OUBL'  GDIV    ; 'OUBLIER' mpvid ;
   'SI' ('NON' direct) ; 'OUBLIER' tqproj 'preco' ; 'FINSI' ;
    list *MATRIK ;
    list *rigidite ;
   fic = 'CHAINE' nomdir nomfic 'it' itcou '.sauv' ;
   'MESSAGE' ('CHAINE' 'Saving ' fic '...') ;
   'OPTI' 'SAUV' fic ;
   'SAUV' ;
   'TEMPS' 'IMPR' ;
'FINSI' ;
*'MESSAGE' ('CHAINE' 'zA=' zA) ;
*'MESSAGE' ('CHAINE' 'miz=' miz) ;
*'MESSAGE' ('CHAINE' 'maz=' maz) ;
*'MESSAGE' ('CHAINE' 'miv=' miv) ;
*'MESSAGE' ('CHAINE' 'mav=' mav) ;
*'SI' interact ;
*   'OPTION' 'DONN' 5 'ECHO' 1 ;
*'SINON' ;
*   'FIN' ;
*'FINSI' ;
'RESPRO' zA miz maz ;
'FINPROC' ;
************************************************************************
*
* END OF COMPUTATIONAL LOOP
*
************************************************************************
*
************************************************************************
*
* TEST PART
*
************************************************************************
*
lpass = VRAI ;
'SAUTER' 2 'LIGNE' ; 'OPTI' 'ECHO' 1 ;
***
*** Test 1  Re=50. Eu=10000.  FrEu=1. We=10000.
***
***         La pression imposée est équilibrée par le poids
***         de fluide déplacé. L'écoulement n'a pas
***         d'influence sur la forme de la surface.
***
***
'OPTI' 'ECHO' 0 ;
Re = 50. ;
Eu= 10000. ;
FrEu =  1. ;
We = 10000.  ;
direct = faux ; miter = faux ;
za miz maz = calcul Re Eu FrEu We direct miter ;
dz = '-' maz miz ;
* Reference values
zaref = -0.956 ; dzref = 1. ;
* Tests
AFFVAR zaref 'zaref' za 'za' ;
AFFVAR dzref 'dzref' dz 'dz' ;
errv = 2.D-2 ;
err1 = errrel za zaref ; tst1 = '<' err1 errv ;
err2 = errrel dz dzref ; tst2 = '<' err2 errv ;
tst = tst1 'ET' tst2 ;
'MESSAGE' ('CHAINE' 'Test 1 :')   ;
'MESSAGE' ('CHAINE' '         err1 = ' err1)  ;
'MESSAGE' ('CHAINE' '         err2 = ' err2)  ;
'SI' tst ;
   'MESSAGE' 'Test 1 OK' ;
'SINON' ;
   'MESSAGE' '!!! Test 1 not passed  ' ;
'FINSI' ;
lpass = lpass 'ET' tst ;
'SAUTER' 2 'LIGNE' ; 'OPTI' 'ECHO' 1 ;
*'OPTION' 'DONN' 5 ;
***
*** Test 2  Re=50. Eu=10000.  FrEu=1. We=1.
***
***         La pression imposée est équilibrée par le poids
***         de fluide déplacé et la tension de surface.
***         L'écoulement n'a pas d'influence sur la forme de la surface.
***
'OPTI' 'ECHO' 0 ;
Re = 50. ;
Eu= 10000. ;
FrEu =  1. ;
We = 1.  ;
direct = faux ; miter = vrai ;
za miz maz = calcul Re Eu FrEu We direct miter ;
dz = '-' maz miz ;
* Reference values
zaref = -0.505 ; dzref = 0.55 ;
* Tests
AFFVAR zaref  'zaref' za 'za' ;
AFFVAR dzref 'dzref' dz 'dz' ;
errv = 2.D-2 ;
err1 = errrel za zaref ; tst1 = '<' err1 errv ;
err2 = errrel dz dzref ; tst2 = '<' err2 errv ;
tst = tst1 'ET' tst2 ;
'MESSAGE' ('CHAINE' 'Test 2 :')   ;
'MESSAGE' ('CHAINE' '         err1 = ' err1)  ;
'MESSAGE' ('CHAINE' '         err2 = ' err2)  ;
'SI' tst ;
   'MESSAGE' 'Test 2 OK' ;
'SINON' ;
   'MESSAGE' '!!! Test 2 not passed  ' ;
'FINSI' ;
lpass = lpass 'ET' tst ;
'SAUTER' 2 'LIGNE' ; 'OPTI' 'ECHO' 1 ;
***
*** Test 3  Re=50. Eu=1.  FrEu=1. We=1.
***
***         La pression imposée est équilibrée par le poids
***         de fluide déplacé et la tension de surface.
***         L'écoulement influence la forme de la surface (aspiration).
***
***
'OPTI' 'ECHO' 0 ;
Re = 50. ;
Eu= 1. ;
FrEu =  1. ;
We = 1.  ;
direct = vrai ; miter = faux ;
za miz maz = calcul Re Eu FrEu We direct miter ;
dz = '-' maz miz ;
* Reference values
zaref = -1.07 ; dzref = 1.21 ;
* Tests
AFFVAR zaref  'zaref' za 'za' ;
AFFVAR dzref  'dzref' dz 'dz' ;
errv = 2.D-2 ;
err1 = errrel za zaref ; tst1 = '<' err1 errv ;
err2 = errrel dz dzref ; tst2 = '<' err2 errv ;
tst = tst1 'ET' tst2 ;
'MESSAGE' ('CHAINE' 'Test 3 :')   ;
'MESSAGE' ('CHAINE' '         err1 = ' err1)  ;
'MESSAGE' ('CHAINE' '         err2 = ' err2)  ;
'SI' tst ;
   'MESSAGE' 'Test 3 OK' ;
'SINON' ;
   'MESSAGE' '!!! Test 3 not passed  ' ;
'FINSI' ;
lpass = lpass 'ET' tst ;
 
 
 
'SAUTER' 2 'LIGNE' ;
'SI' lpass ;
   'MESSAGE' 'Tout sest bien passe' ;
'SINON' ;
   'MESSAGE' 'Il y a eu des erreurs' ;
'FINSI' ;
'SAUTER' 2 'LIGNE' ;
'SI' ('NON' lpass) ;
   'ERREUR' 5 ;
'FINSI' ;
************************************************************************
*
* END OF TEST PART
*
************************************************************************
'SI' interact ;
   'OPTION' 'DONN' 5 'ECHO' 1 ;
'FINSI' ;
*
'FIN' ;