* fichier :  domall.dgibi
*
**** Compatibility check
*
 GRAPH = FAUX ;
* CHECKERR to see the error messages. 
 CHECKERR = FAUX ;

*
**** Nothing from here
*
 'SI' FAUX ;
* 
**** Till here
*
 'FINSI' ;
 
*****************************
**** SEG3 in 2D *************
*****************************

 'OPTION' 'DIME' 2 'ELEM' 'TRI3' ;

 P1 = 1.11 0.0 ;
 P2 = 1.0 1.13 ;

 MAIL = P1 'DROIT' 5 P2 ;
*
*
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM 'VF' ;
 QMAIL = TABM . 'QUAF' ;


 'SI' GRAPH ;
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
    'TITRE' 'CENTRE et MAILLAGE';
 'FINSI' ;

 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'SEG3' ;
    'ERREUR' 5 ;
 'FINSI' ; 

* Volume

  VOL0 = 'PSCAL' (P2 'MOIN' P1) (P2 'MOIN' P1) ;
  VOL0 = VOL0 '**' 0.5 ;
  AA = 'MANUEL' 'CHPO' (TABM .'CENTRE') 1 'SCAL' 1.0 'NATU'
       'DISCRETE' ;
  VOL1 = 'XTY' (TABM . 'XXVOLUM')  AA ('MOTS' 'SCAL')
         ('MOTS' 'SCAL') ;

  ERRO = (VOL0 '-' VOL1) 'ABS' ;
  'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'SEG3' ;
    'ERREUR' 5 ;
 'FINSI' ;
    
  
*****************************
**** SEG3 in 3D *************
*****************************

 'SI' CHECKERR ;
    'OPTION' 'DIME' 3 'ELEM' 'TRI3' ;

    P1 = 0.0 0.0 0.0 ;
    P2 = 1.0 0.0 0.0 ;

    MAIL = P1 'DROIT' 1 P2 ;
*
*
* QMAIL1 = QMAIL 'PLUS' (0.0 0.0) ;
    MODM  = 'MODELISER' MAIL 'EULER' ;
    TABM  = 'DOMA' MODM 'VF' ;
 'FINSI' ;
 
*****************************
**** TRI3 in 3D *************
*****************************

 'OPTION' 'DIME' 3 'ELEM' 'TET4' ;
 
 P1 = 0.0 0.0 0.0 ;
 P2 = 1.0 0.0 0.0 ;
 P3 = 1.11 1.13 0.0 ;
 P4 = 0.0 1.0 0.0 ;

 P1 P2 P3 P4 = P1 P2 P3 P4 'TOURNER' 15.0
   (0.0 0.0 2.0) (1.1 1.1 2.0) ;
  
 MAIL = 'DALL' (P1 'DROIT' 1 P2) (P2 'DROIT' 1 P3) (P3 'DROIT' 1 P4)
               (P4 'DROIT' 1 P1) 'PLAN'  ;
*
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM 'VF' ;
 QMAIL = TABM . 'QUAF' ;
*
 'SI' GRAPH ;
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
    'TITRE' 'CENTRE et MAILLAGE';
 'FINSI' ;
*
 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'TRI3' ;
    'ERREUR' 5 ;
 'FINSI' ; 

* Volume 
 
 NN = 'NBEL' MAIL ;
 'REPETER' BL1 NN ;
    AA = 'MESURE' ('ELEM' &BL1 MAIL) ;
    BB = 'EXTRAIRE' (TABM . 'XXVOLUM') 'SCAL'
         ('POIN' &BL1 (TABM . 'CENTRE')) ;
    'LISTE' BB ;
    ERRO = (AA '-' BB) 'ABS' ;
    'SI' (ERRO > 1.0D-6) ;
      'MESSAGE' 'TRI3' ;
      'ERREUR' 5 ;
    'FINSI' ;
 'FIN' BL1 ;

*****************************
**** QUA4 in 3D *************
*****************************

 'OPTION' 'DIME' 3 'ELEM' 'CUB8' ;

 P1 = 0.0 0.0 0.0 ;
 P2 = 1.0 0.0 0.0 ;
 P3 = 1.11 1.13 0.0 ;
 P4 = 0.0 1.0 0.0 ;

 P1 P2 P3 P4 = P1 P2 P3 P4 'TOURNER' 15.0
   (0.0 0.0 2.0) (1.1 1.1 2.0) ;
 
 MAIL = 'DALL' (P1 'DROIT' 2 P2) (P2 'DROIT' 2 P3) (P3 'DROIT' 2 P4)
               (P4 'DROIT' 2 P1) 'PLAN'  ;
*
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM 'VF' ;
 QMAIL = TABM . 'QUAF' ;

 'SI' GRAPH ; 
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
       'TITRE' 'CENTRE et MAILLAGE';
 'FINSI' ;
*
 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'QUA4' ;
    'ERREUR' 5 ;
 'FINSI' ; 

* Volume 
 
 NN = 'NBEL' MAIL ;
 'REPETER' BL1 NN ;
*
*   NB If the elts are not in the same plane
*      'MESURE' gives not the same results
*
    AA = 'MESURE' ('ELEM' &BL1 MAIL) ;
    BB = 'EXTRAIRE' (TABM . 'XXVOLUM') 'SCAL'
         ('POIN' &BL1 (TABM . 'CENTRE')) ;
    ERRO = (AA '-' BB) 'ABS' ;
    'SI' (ERRO > 1.0D-6) ;
      'MESSAGE' 'QUA4' ;
      'ERREUR' 21 ;
    'FINSI' ;
 'FIN' BL1 ;

 
*****************************
**** TRI3 and QUA4  in 3D ***
*****************************


 'OPTION' 'DIME' 3 'ELEM' 'TET4' ;
 
 P1 = 0.0 0.0 0.0 ;
 P2 = 1.0 0.0 0.0 ;
 P3 = 1.11 1.13 0.0 ;
 P4 = 0.0 1.0 0.0 ;
 P5 = 0.0 2.0 0.0 ;
 P6 = 1.0 2.0 0.0 ;
 
 MAIL1 = 'DALL' (P1 'DROIT' 2 P2) (P2 'DROIT' 2 P3) (P3 'DROIT' 2 P4)
               (P4 'DROIT' 2 P1) 'PLAN'  ;

 'OPTION' 'ELEM' 'CUB8' ;

 MAIL2 = 'SURF' ((P4 'DROIT' 3 P5) 'ET' (P5 'DROIT' 2 P6) 'ET'
                 (P6 'DROIT' 2 P3) 'ET' (P3 'DROIT' 2 P4)) 'PLAN'  ;

 MAIL = MAIL1 'ET' MAIL2 ;
 'ELIMINATION' MAIL 0.0001 ;
 
*
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM 'VF' ;
 QMAIL = TABM . 'QUAF' ;
 
* 
 'SI' GRAPH ; 
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
       'TITRE' 'CENTRE et MAILLAGE';
 'FINSI' ;
*
 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 5 ;
 'FINSI' ; 

* Volume 

 MAIL = (TABM . 'MAILLAGE') ;
 TYPEL = MAIL 'ELEM' 'TYPE' ;
 NN0 = 'DIME' TYPEL ;
 NNE=0 ;
 'REPETER' BL0 NN0 ;
    MAILA = MAIL 'ELEM' ('EXTRAIRE' TYPEL &BL0) ;
    NN = 'NBEL' MAILA ;
    'REPETER' BL1 NN ;
       NNE = NNE '+' 1 ;
*
*   NB If the elts are not on the same plane
*      'MESURE' gives not the same results
*
       AA = 'MESURE' ('ELEM' &BL1 MAILA) ;
       BB = 'EXTRAIRE' (TABM . 'XXVOLUM') 'SCAL'
            ('POIN' NNE (TABM . 'CENTRE')) ;
       ERRO = (AA '-' BB) 'ABS' ;
      'SI' (ERRO > 1.0D-6) ;
         'MESSAGE' 'TRI3 et QUA4' ;
         'ERREUR' 21 ;
      'FINSI' ;
    'FIN' BL1 ;
 'FIN' BL0 ;


*****************************
**** TRI3 in 2D *************
*****************************

 'OPTION' 'DIME' 2 'ELEM' 'TRI3' ;
 
 P1 = 0.0 0.0  ;
 P2 = 1.0 0.0  ;
 P3 = 1.11 1.13  ;
 P4 = 0.0 1.0  ;
 
 MAIL = 'DALL' (P1 'DROIT' 1 P2) (P2 'DROIT' 1 P3) (P3 'DROIT' 1 P4)
               (P4 'DROIT' 1 P1) 'PLAN'  ;
*
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM 'VF' ;
 QMAIL = TABM . 'QUAF' ;
 
*
 'SI' GRAPH ;
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
    'TITRE' 'CENTRE et MAILLAGE';
    VEC1 = 'VECTEUR' (TABM . 'XXNORMAF') 0.1 'UX' 'UY' 'ROUG' ;
    'TRACER' MAIL VEC1 'TITRE' 'Normals' ;
 'FINSI' ;
*
 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'TRI3' ;
    'ERREUR' 5 ;
 'FINSI' ; 

*****************************
**** QUA4 bizarre in 2D *****
*****************************

 'SI' CHECKERR ;

* Center of gravity outside the domain
 
    'OPTION' 'DIME' 2 'ELEM' 'QUA4' ;

    P1 = 0.0 0.0  ;
    P2 = 1.0 0.0  ;
    P3 = 1.0 1.0  ;
    P4 = 0.75 0.1 ;
 
    MAIL = 'DALL' (P1 'DROIT' 1 P2) (P2 'DROIT' 1 P3) (P3 'DROIT' 1 P4)
               (P4 'DROIT' 1 P1) 'PLAN'  ;
*
    MODM  = 'MODELISER' MAIL 'EULER' ;
    TABM  = 'DOMA' MODM 'VF' ;
 
 'FINSI' ;

* 
* Element with singular Jacobian allowed since
* the center of gravity is inside of the elements
* Clockwise orientation of the noeuds.
*

 'OPTION' 'DIME' 2 'ELEM' 'QUA4' ;

 P1 = 1.0 1.0  ;
 P2 = 1.0 0.0  ;
 P3 = 0.0 0.0  ;
 P4 = 0.5 0.35  ;
 
 MAIL = 'MANUEL' 'QUA4' P1 P2 P3 P4 ;
*
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM 'VF' ;
 QMAIL = TABM . 'QUAF' ;

 'SI' GRAPH ; 
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et
             ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
             'TITRE' 'CENTRE et MAILLAGE';
    VEC1 = 'VECTEUR' (TABM . 'XXNORMAF') 0.1 'UX' 'UY' 'ROUG' ;
    'TRACER' MAIL VEC1 'TITRE' 'Normals' ;
 'FINSI' ;
*
 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'QUA4' ;
    'ERREUR' 5 ;
 'FINSI' ; 

* Position of the center of gravity

 VOL1 = 'MESURE' ('MANUEL' 'TRI3' P1 P2 P4) ;
 VOL2 = 'MESURE' ('MANUEL' 'TRI3' P2 P3 P4) ;
 PCEN1 = (P1 'PLUS' P2 'PLUS' P4) '/' 3. ;
 PCEN2 = (P2 'PLUS' P3 'PLUS' P4) '/' 3. ;
 PCEN = (VOL1 * PCEN1) 'PLUS' (VOL2 * PCEN2) ;
 PCEN = PCEN '/' (VOL1 '+' VOL2) ;
 'ELIMINATION' ((TABM . CENTRE) 'ET' PCEN) 0.001 ;
 AA = 'DIFF' (TABM . CENTRE) ('MANUEL' 'POI1' PCEN) ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'QUA4' ;
    'ERREUR' 5 ;
 'FINSI' ;

* Volume

 VOL = 'EXTRAIRE' (TABM . 'XXVOLUM') PCEN 'SCAL' ;
 ERRO = (VOL '-' (VOL1 '+' VOL2)) 'ABS' ;

 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

* Normals

 XP1 = 'COORDONNEE' 1 P1 ;
 YP1 = 'COORDONNEE' 2 P1 ;
 XP2 = 'COORDONNEE' 1 P2 ;
 YP2 = 'COORDONNEE' 2 P2 ;
 PCEN3 = (P1 'PLUS' P2) '/' 2 ;
 'ELIMINATION' ((TABM . 'FACE') 'ET' PCEN3) 0.001 ;
 
 VAL = (((XP1 '-' XP2) '**' 2) '+' ((YP2 '-' YP1) '**' 2)) '**' 0.5 ;
 VALNOR1 = 'MANUEL' 'CHPO' ('MANUEL' 'POI1' PCEN3) 2 'UX'
           ((YP1 '-' YP2) '/' VAL) 'UY' ((XP2 '-' XP1) '/' VAL) ;
 VALNOR2 = 'REDU' (TABM . 'XXNORMAF') ('MANUEL' 'POI1' PCEN3) ;
 ERRO = 'MAXIMUM' (VALNOR1 '-' VALNOR2) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;

*
* We check that the sum of (normals times surfaces) is zero
*

 VALNOR = TABM . 'XXNORMAF' ;
 VALSUR = TABM . 'XXSURFAC' ;
 
 AA = VALNOR '*' VALSUR ;
 BB = 'MANUEL' 'CHPO' (TABM . 'FACE') 1 'SCAL' 1.0 ;
 ERRO = 'XTY' AA BB ('MOTS' 'UX' 'UY') ('MOTS' 'SCAL' 'SCAL') ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
 
*
* We check that the normals modulus are 1
*

 AA = 'PSCAL' VALNOR VALNOR ('MOTS' 'UX' 'UY') ('MOTS' 'UX' 'UY') ;
 'SI' ((('MAXIMUM' AA) 'NEG' 1.0 1.0D-12) 'OU'
       (('MINIMUM' AA) 'NEG' 1.0 1.0D-12)) ;
    'ERREUR' ('CHAIN' 'Problem QUA4') ;
 'FINSI' ;

*
* Matrot
*

 AA = 'EXCO' (TABM . 'MATROT') ('MOTS' 'MX' 'MY')
      ('MOTS' 'UX' 'UY') ;
 ERRO ='MAXIMUM' (AA '-' (TABM . 'XXNORMAF')) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 BB = 'EXCO' (TABM . 'MATROT') ('MOTS' 'RX' 'RY')
      ('MOTS' 'UX' 'UY') ;
 ERRO = 'MAXIMUM' ('PSCAL' AA BB ('MOTS' 'UX' 'UY') ('MOTS' 'UX' 'UY'))
        'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
* Vector product
 COMPZ = 'PSCAL' AA BB  ('MOTS' 'UX') ('MOTS' 'UY') ;
 COMPZ = COMPZ '-' ('PSCAL' AA BB  ('MOTS' 'UY') ('MOTS' 'UX')) ;
 ERRO = 'MAXIMUM' (COMPZ '-' 1.) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;


*****************************
**** QUA4 in 2D *************
*****************************

 'OPTION' 'DIME' 2 'ELEM' 'QUA4' ;

 P1 = 0.0 0.0  ;
 P2 = 1.0 0.0  ;
 P3 = 1.11 1.13  ;
 P4 = 0.0 1.0 ;
 
 MAIL = 'DALL' (P1 'DROIT' 2 P2) (P2 'DROIT' 2 P3) (P3 'DROIT' 2 P4)
               (P4 'DROIT' 2 P1) 'PLAN'  ;
*
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM 'VF' ; QMAIL = TABM . 'QUAF' ;
 

 'SI' GRAPH ; 
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
       'TITRE' 'CENTRE et MAILLAGE';
    VEC1 = 'VECTEUR' (TABM . 'XXNORMAF') 0.1 'UX' 'UY' 'ROUG' ;
    'TRACER' MAIL VEC1 'TITRE' 'Normals' ;
 'FINSI' ;
*
 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'QUA4' ;
    'ERREUR' 5 ;
 'FINSI' ; 

*****************************
**** TRI3 and QUA4  in 2D ***
*****************************


 'OPTION' 'DIME' 2 'ELEM' 'TRI3' ;
 
 P1 = 0.0 0.0  ;
 P2 = 1.0 0.0  ;
 P3 = 1.11 1.13  ;
 P4 = 0.0 1.0  ;
 P5 = 0.0 2.0  ;
 P6 = 1.0 2.0  ;
 
 MAIL1 = 'DALL' (P1 'DROIT' 2 P2) (P2 'DROIT' 2 P3) (P3 'DROIT' 2 P4)
               (P4 'DROIT' 2 P1) 'PLAN'  ;
 'OPTION' 'ELEM' 'QUA4' ;
 MAIL2 = 'SURF' ((P4 'DROIT' 3 P5) 'ET' (P5 'DROIT' 2 P6) 'ET'
                 (P6 'DROIT' 2 P3) 'ET' (P3 'DROIT' 2 P4)) 'PLAN'  ;

 MAIL = MAIL1 'ET' MAIL2 ;
 'ELIMINATION' MAIL 0.0001 ;

 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM 'VF  ' ;
 QMAIL = TABM . 'QUAF' ;
* 
 'SI' GRAPH ; 
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
       'TITRE' 'CENTRE et MAILLAGE';
    VEC1 = 'VECTEUR' (TABM . 'XXNORMAF') 0.1 'UX' 'UY' 'ROUG' ;
    'TRACER' MAIL VEC1 'TITRE' 'Normals' ;
 'FINSI' ;
*
 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 5 ;
 'FINSI' ; 

* Volume 

 MAIL = (TABM . 'MAILLAGE') ;
 TYPEL = MAIL 'ELEM' 'TYPE' ;
 NN0 = 'DIME' TYPEL ;
 NNE=0 ;
 'REPETER' BL0 NN0 ;
    MAILA = MAIL 'ELEM' ('EXTRAIRE' TYPEL &BL0) ;
    NN = 'NBEL' MAILA ;
    'REPETER' BL1 NN ;
       NNE = NNE '+' 1 ;
       AA = 'MESURE' ('ELEM' &BL1 MAILA) ;
       BB = 'EXTRAIRE' (TABM . 'XXVOLUM') 'SCAL'
            ('POIN' NNE (TABM . 'CENTRE')) ;
       ERRO = (AA '-' BB) 'ABS' ;
      'SI' (ERRO > 1.0D-6) ;
         'MESSAGE' 'TRI3 et QUA4' ;
         'ERREUR' 21 ;
      'FINSI' ;
    'FIN' BL1 ;
 'FIN' BL0 ;

* Sommet

 MAILS = 'CHANGER' MAIL 'POI1' ;
 AA = 'DIFF' MAILS (TABM . 'SOMMET') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

* Comparaison complete 'NAVIER_STOKES' - 'EULER'

 MOD2 = 'MODELISER' QMAIL 'NAVIER_STOKES' 'LINE' ;
 TAB2 = 'DOMA' MOD2 'TABLE' ;
 TAB2 . 'PRECONDI' = 1 ;

 AA = 'DIFF' (TAB2 . 'QUAF') (TABM . 'QUAF') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
 AA = 'DIFF' (TAB2 . 'MAILLAGE') (TABM . 'MAILLAGE') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
 AA = 'DIFF' (TAB2 . 'CENTRE') (TABM . 'CENTRE') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
 AA = 'DIFF' (TAB2 . 'SOMMET') (TABM . 'SOMMET') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
 AA = 'DIFF' (TAB2 . 'FACEL') (TABM . 'FACEL') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 AA = 'DIFF' (TAB2 . 'FACEP') (TABM . 'FACEP') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 AA = 'DIFF' (TAB2 . 'ELTFA') (TABM . 'ELTFA') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

* XXSURFAC, XXNORMAF, XXDIEMIN  (comparaison 'NAVIER_STOKES' - 'EULER')

 CHPO1 = 'DOMA' TAB2 'XXSURFAC' ;
 CHPO2 = TABM . 'XXSURFAC' ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 CHPO1 = 'DOMA' TAB2 'XXVOLUM' ;
 CHPO2 = TABM . 'XXVOLUM' ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 CHPO1 = 'ABS' (TABM . 'XXNORMAF') ;
 CHPO2 = 'ABS' ('DOMA' TAB2  'XXNORMAF') ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 CHPO1 = (TABM . 'XXDIEMIN') ;
 CHPO2 = ('DOMA' TAB2  'XXDIEMIN') ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

*
* We check that the sum of (normals times surfaces) is zero
* on the border
*

 AA = 'CONTOUR' MAIL ;
 CC = 'MODELISER' AA 'EULER' ;
 TT = 'DOMA' CC 'VF' ;
 BB = TT. 'QUAF' ;

 'ELIMINATION' BB (TABM . 'QUAF') 0.0001 ;
 DD = TT . 'CENTRE' ;

 VALNOR = 'REDU' (TABM . 'XXNORMAF') DD ;
 VALSUR = 'REDU' (TABM . 'XXSURFAC') DD ;
 
 AA = VALNOR '*' VALSUR ;
 BB = 'MANUEL' 'CHPO' DD 1 'SCAL' 1.0 ;
 ERRO = 'XTY' AA BB ('MOTS' 'UX' 'UY') ('MOTS' 'SCAL' 'SCAL') ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' 'TRI3 et QUA4' ;
 'FINSI' ;

*
* We check that each element satisfies the property  
* (normals times surfaces) = zero
*
* ELTFA and MAIL have the same structure
*
 LISTYP =  (TABM . 'ELTFA') 'ELEM' 'TYPE' ;
 LISTYP1 = (TABM . 'QUAF')  'ELEM' 'TYPE' ;
*
 NELEM = 0 ;
 'REPETER' BL1 ('DIME' LISTYP) ;
    MAIL1 = 'ELEM' (TABM . 'ELTFA') ('EXTRAIRE' LISTYP &BL1) ;
    MAIL2 = 'ELEM' (TABM . 'QUAF')  ('EXTRAIRE' LISTYP1 &BL1) ;
    NM1 = 'NBEL' MAIL1 ;
    NM2 = 'NBEL' MAIL2 ;
    'SI' ('NEG' NM1 NM2) ;
       'ERREUR' 'QUAF and ELTFA non compatible' ;
       'ERREUR' 21 ;
    'FINSI' ;
    'REPETER' BL2 ('NBEL' MAIL1) ;
       NELEM = NELEM '+' 1 ;
       ELE1 = 'ELEM' MAIL1 &BL2 ;
       ELE2 = 'CHANGER' ELE1 'POI1' ;
       NM2 = 'NBNO' ELE2 ;
       ELE3 = 'ELEM' MAIL2 &BL2 ;
       ELE4 = 'CHANGER' ELE3 'POI1' ;
       NM4 = 'NBNO' ELE4 ;
       ELE5 = 'DIFF' ELE4 ELE2 ;
       NM5 = 'NBNO' ELE5 ;
       'SI' ('NEG' (NM2 '+' NM5) NM4) ;
          'ERREUR' 'QUAF and ELTFA non compatible' ;
          'ERREUR' 21 ;
       'FINSI' ;
       PCEN = 'POIN' (TABM . 'CENTRE') NELEM ;
       DSURF = 0.0 0.0 ;
       'REPETER' BL3 ('NBEL' ELE2) ;
          PFAC = 'POIN' &BL3 ELE2 ;
          SURFF = 'EXTRAIRE' (TABM . 'XXSURFAC') 'SCAL' PFAC ;
          PNORM = ('EXTRAIRE' (TABM . 'XXNORMAF') 'UX' PFAC)
                  ('EXTRAIRE' (TABM . 'XXNORMAF') 'UY' PFAC);
          PCF = PFAC 'MOIN' PCEN ;
          DORIE = 'SIGNE' ('PSCAL' PNORM PCF) ;
          DSURF = DSURF 'PLUS' (DORIE * PNORM * SURFF) ;
       'FIN' BL3 ;
       AA = 'MAXIMUM' ('ABS' ('PROG' ('COORDONNEE' 1 DSURF)
             ('COORDONNEE' 2 DSURF) ('COORDONNEE' 3 DSURF))) ;
       'SI' (AA > 1.0D-12) ;
          'MESSAGE' 'GENE' ;
          'ERREUR' 21 ;
       'FINSI' ;
    'FIN' BL2 ;
 'FIN' BL1 ;

* 
******************************
**** TET4 ********************
******************************
*
 'OPTION' 'DIME' 3 'ELEM' 'TET4' ;
* P1, P2, P3 clockwise 
 P1 = 1.0 0.0 0.0 ;
 P2 = 0.0 0.0 0.0 ;
 P3 = 1.11 1.13 0.0 ;
 P4 = -0.11 1.73 1.11 ;
 P1 P2 P3 P4 = P1 P2 P3 P4  'TOURNER' 15.0
   (0.0 0.0 2.0) (1.1 1.1 2.0) ;
 PCEN = P1 'PLUS' P2 'PLUS' P3 'PLUS' P4 ;
 PCEN = PCEN '/' 4. ;
 MAIL = 'MANUEL' 'TET4' P1 P2 P3 P4 ;
*
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM  'VF' ; QMAIL = TABM . 'QUAF' ;
 VALNOR = TABM . 'XXNORMAF' ;
 VALSUR = TABM . 'XXSURFAC' ;
* 
 'SI' GRAPH ;
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
    'TITRE' 'CENTRE et MAILLAGE';
    VEC1 = 'VECTEUR' (TABM . 'XXNORMAF') 0.1 'UX' 'UY' 'UZ' 'ROUG' ;
    'TRACER' MAIL VEC1 'TITRE' 'Normals' ;
 'FINSI' ;
*
 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'TET4' ;
    'ERREUR' 5 ;
 'FINSI' ; 
 'ELIMINATION' ((TABM. 'CENTRE') 'ET' PCEN) 0.0001 ;
 AA = 'DIFF' ('MANUEL'  'POI1' PCEN) (TABM. 'CENTRE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'TET4' ;
    'ERREUR' 5 ;
 'FINSI' ;
 
* Volume 
 
 NN = 'NBEL' MAIL ;
 'REPETER' BL1 NN ;
    AA = 'MESURE' ('ELEM' &BL1 MAIL) ; 
    BB = 'EXTRAIRE' (TABM . 'XXVOLUM') 'SCAL' 
         ('POIN' &BL1 (TABM . 'CENTRE')) ; 
    ERRO = (AA '-' BB) 'ABS' ;
    'SI' (ERRO > 1.0D-6) ;
      'MESSAGE' 'TET4' ;
      'ERREUR' 21 ;
    'FINSI' ;
 'FIN' BL1 ;

* Sommet

 MAILS = 'CHANGER' MAIL 'POI1' ;
 AA = 'DIFF' MAILS (TABM . 'SOMMET') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TET4' ;
    'ERREUR' 21 ;
 'FINSI' ;

* Comparaison complete 'NAVIER_STOKES' - 'EULER'

 MOD2 = 'MODELISER' QMAIL 'NAVIER_STOKES' 'LINE' ;
 TAB2 = 'DOMA' MOD2 'TABLE' ;
 TAB2 . 'PRECONDI' = 1 ;

 AA = 'DIFF' (TAB2 . 'QUAF') (TABM . 'QUAF') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
 AA = 'DIFF' (TAB2 . 'MAILLAGE') (TABM . 'MAILLAGE') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
 AA = 'DIFF' (TAB2 . 'CENTRE') (TABM . 'CENTRE') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
 AA = 'DIFF' (TAB2 . 'SOMMET') (TABM . 'SOMMET') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
 AA = 'DIFF' (TAB2 . 'FACEL') (TABM . 'FACEL') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 AA = 'DIFF' (TAB2 . 'FACEP') (TABM . 'FACEP') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 AA = 'DIFF' (TAB2 . 'ELTFA') (TABM . 'ELTFA') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

* XXSURFAC, XXNORMAF, XXDIEMIN  (comparaison 'NAVIER_STOKES' - 'EULER')

 CHPO1 = 'DOMA' TAB2 'XXSURFAC' ;
 CHPO2 = TABM . 'XXSURFAC' ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 CHPO1 = 'DOMA' TAB2 'XXVOLUM' ;
 CHPO2 = TABM . 'XXVOLUM' ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 CHPO1 = 'ABS' (TABM . 'XXNORMAF') ;
 CHPO2 = 'ABS' ('DOMA' TAB2  'XXNORMAF') ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 CHPO1 = (TABM . 'XXDIEMIN') ;
 CHPO2 = ('DOMA' TAB2  'XXDIEMIN') ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

*
* We check the interfaces
*

 PCEN1 = (P1 'PLUS' P2 'PLUS' P3) '/' 3. ;
 PCEN2 = (P2 'PLUS' P3 'PLUS' P4) '/' 3. ;
 PCEN3 = (P3 'PLUS' P4 'PLUS' P1) '/' 3. ;
 PCEN4 = (P4 'PLUS' P1 'PLUS' P2) '/' 3. ;
*
 A1 = 'PVECT' (P2 'MOIN' P1) (P3 'MOIN' P1) ;
 A2 = 'PVECT' (P3 'MOIN' P2) (P4 'MOIN' P2) ;
 A3 = 'PVECT' (P4 'MOIN' P3) (P1 'MOIN' P3) ;
 A4 = 'PVECT' (P1 'MOIN' P4) (P2 'MOIN' P4) ;
 VALSUR1 = 'PROG' ('PSCAL' A1 A1) ('PSCAL' A2 A2)
                  ('PSCAL' A3 A3) ('PSCAL' A4 A4) ;
 VALSUR1 = VALSUR1 '**' 0.5 ;
 VALSUR1 = VALSUR1 '*' 0.5 ;
 VALSURM = 'MANUEL' 'QUA4' PCEN1 PCEN2 PCEN3 PCEN4 ;
 'ELIMINATION' VALSURM 0.001 (TABM . 'FACE') ;
 VALSUR2 = 'MANUEL' 'CHPO' VALSURM 1 'SCAL' VALSUR1 ;   
 ERRO = 'MAXI' (VALSUR '-' VALSUR2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
* We check that the normals modulus are 1
*
 AA = 'PSCAL' VALNOR VALNOR ('MOTS' 'UX' 'UY' 'UZ')
      ('MOTS' 'UX' 'UY' 'UZ') ;
 'SI' ((('MAXIMUM' AA) 'NEG' 1.0 1.0D-12) 'OU'
       (('MINIMUM' AA) 'NEG' 1.0 1.0D-12)) ;
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
* We check the normal in one face (P1, P2, P4)
*
 VECNOR1 = 'PVEC' (P4 'MOIN' P1) (P2 'MOIN' P1) ;
 VECNOR1 = VECNOR1 '/' (('PSCAL' VECNOR1 VECNOR1) '**' 0.5) ;
 VALNOR1 = 'MANUEL' 'CHPO' ('MANUEL' 'POI1' PCEN4) 3
           'UX' ('COORDONNEE' 1 VECNOR1)
           'UY' ('COORDONNEE' 2 VECNOR1)
           'UZ' ('COORDONNEE' 3 VECNOR1)
           ;
 ERRO = 'MAXI' (VALNOR1 '-' ('REDU' VALNOR ('MANUEL' 'POI1' PCEN4)))
        'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
* We check that the sum of (normals times surfaces) is zero
*
 AA = VALNOR '*' VALSUR ;
 BB = 'MANUEL' 'CHPO' (TABM . 'FACE') 1 'SCAL' 1.0 ;
 ERRO = 'XTY' AA BB ('MOTS' 'UX' 'UY' 'UZ')
        ('MOTS' 'SCAL' 'SCAL' 'SCAL') ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;

*
* Matrot
*
 AA = 'EXCO' (TABM . 'MATROT') ('MOTS' 'MX' 'MY' 'MZ')
      ('MOTS' 'UX' 'UY' 'UZ') ;
 ERRO ='MAXIMUM' (AA '-' (TABM . 'XXNORMAF')) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 BB = 'EXCO' (TABM . 'MATROT') ('MOTS' 'RX' 'RY' 'RZ')
      ('MOTS' 'UX' 'UY' 'UZ') ;
 ERRO = 'MAXIMUM' ('PSCAL' AA BB ('MOTS' 'UX' 'UY' 'UZ')
                                 ('MOTS' 'UX' 'UY' 'UZ'))
        'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 COMPZ = 'PSCAL' AA BB  ('MOTS' 'UX') ('MOTS' 'UY') ;
 COMPZ = COMPZ '-' ('PSCAL' AA BB  ('MOTS' 'UY') ('MOTS' 'UX')) ;
 COMPZ1 = 'EXCO' (TABM . 'MATROT') ('MOTS' 'UZ')
          ('MOTS' 'SCAL') ;
 ERRO = 'MAXIMUM' (COMPZ '-' COMPZ1) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 COMPX = 'PSCAL' AA BB  ('MOTS' 'UY') ('MOTS' 'UZ') ;
 COMPX = COMPX '-' ('PSCAL' AA BB  ('MOTS' 'UZ') ('MOTS' 'UY')) ;
 COMPX1 = 'EXCO' (TABM . 'MATROT') ('MOTS' 'UX')
          ('MOTS' 'SCAL') ;
 ERRO = 'MAXIMUM' (COMPX '-' COMPX1) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 COMPY = 'PSCAL' AA BB  ('MOTS' 'UZ') ('MOTS' 'UX') ;
 COMPY = COMPY '-' ('PSCAL' AA BB  ('MOTS' 'UX') ('MOTS' 'UZ')) ;
 COMPY1 = 'EXCO' (TABM . 'MATROT') ('MOTS' 'UY')
          ('MOTS' 'SCAL') ;
 ERRO = 'MAXIMUM' (COMPY '-' COMPY1) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;

******************************
**** PYR5 ********************
******************************
 

*
* General case
*
 'OPTION' 'DIME' 3 'ELEM' 'PYR5' ;
 P1 = 0.0 0.0 0.11 ;
 P2 = 1.1 0.0 0.13 ;
 P3 = 1.0 1.3 0.15 ;
 P4 = 0.0 1.0 0.0;
 P5 = 0.5 0.5 1.0 ;
 
 P1 P2 P3 P4 P5 = P1 P2 P3 P4 P5 'TOURNER' 15.0
   (0.0 0.0 4.0) (1.1 1.1 4.0) ;
 MAIL = 'MANUEL' 'PYR5' P1 P2 P3 P4 P5 ;
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM  'VF' ; QMAIL = TABM . 'QUAF' ;
 MAIL1 = 'MANUEL' 'QUA4' P1 P2 P3 P4 ;
 MODM1  = 'MODELISER' MAIL1 'EULER' ;
 TABM1  = 'DOMA' MODM1  'VF' ; QMAIL1 = TABM1 . 'QUAF' ;

 'ELIMINATION' QMAIL QMAIL1 0.001 ;

 PCENF = 'POIN' 1 (TABM1 . 'CENTRE') ;

* 
 'SI' GRAPH ;
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
    'TITRE' 'CENTRE et MAILLAGE';
    VEC1 = 'VECTEUR' (TABM . 'XXNORMAF') 0.1 'UX' 'UY' 'UZ' 'ROUG' ;
    'TRACER' MAIL VEC1 'TITRE' 'Normals' ;
 'FINSI' ;
*
 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'PYR5' ;
    'ERREUR' 5 ;
 'FINSI' ; 

* Volume as the volume of 4 tetrahedra

 VOL1 = 'MESURE' ('MANUEL' 'TET4' P1 P2 PCENF P5) ;
 VOL2 = 'MESURE' ('MANUEL' 'TET4' P2 P3 PCENF P5) ;
 VOL3 = 'MESURE' ('MANUEL' 'TET4' P3 P4 PCENF P5) ;
 VOL4 = 'MESURE' ('MANUEL' 'TET4' P4 P1 PCENF P5) ;
 VOL = VOL1 '+' VOL2 '+' VOL3 '+' VOL4 ;

 VOL0 = 'MAXIMUM' (TABM . 'XXVOLUM') ;

 ERRO = (VOL0 '-' VOL) 'ABS' ;
 
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'PYR5' ;
    'ERREUR' 21 ;
 'FINSI' ;

* Centre
 PCEN1 =  P1 'PLUS' P2 'PLUS' PCENF 'PLUS' P5 ;
 PCEN1 = PCEN1 '/' 4. ;
 PCEN2 =  P2 'PLUS' P3 'PLUS' PCENF 'PLUS' P5 ;
 PCEN2 = PCEN2 '/' 4. ;
 PCEN3 =  P3 'PLUS' P4 'PLUS' PCENF 'PLUS' P5 ;
 PCEN3 = PCEN3 '/' 4. ;
 PCEN4 =  P4 'PLUS' P1 'PLUS' PCENF 'PLUS' P5 ;
 PCEN4 = PCEN4 '/' 4. ;

 PCEN = (VOL1 * PCEN1) 'PLUS' (VOL2 * PCEN2) 'PLUS'
        (VOL3 * PCEN3) 'PLUS' (VOL4 * PCEN4) ;
 PCEN = PCEN '/' VOL ;
 'ELIMINATION'  (PCEN 'ET' (TABM . 'CENTRE')) 0.0001 ;

 AA = 'DIFF' ('MANUEL' 'POI1' PCEN) (TABM . 'CENTRE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'PYR5' ;
    'ERREUR' 5 ;
 'FINSI' ;

* Normals and surface

 PNOR1 = 0.5 '*' ('PVEC' (P1 'MOIN' PCENF) (P4 'MOIN' PCENF)) ;
 PNOR2 = 0.5 '*' ('PVEC' (P2 'MOIN' PCENF) (P1 'MOIN' PCENF)) ;
 PNOR3 = 0.5 '*' ('PVEC' (P3 'MOIN' PCENF) (P2 'MOIN' PCENF)) ;
 PNOR4 = 0.5 '*' ('PVEC' (P4 'MOIN' PCENF) (P3 'MOIN' PCENF)) ;

 PNOR = PNOR1 'PLUS' PNOR2 'PLUS' PNOR3 'PLUS' PNOR4 ;
 XSURF = (('PSCAL' PNOR PNOR) '**' 0.5) ;
 PNOR = PNOR '/' XSURF ;
 
 AA = ('REDU' (tabm . 'XXNORMAF') ('MANUEL' 'POI1' PCENF)) ;
 BB = 'MANUEL' 'CHPO' ('MANU' 'POI1' PCENF) 3
      'UX' ('COORDONNEE' 1 PNOR)
      'UY' ('COORDONNEE' 2 PNOR)
      'UZ' ('COORDONNEE' 3 PNOR);
 ERRO = 'MAXIMUM' (AA '-' BB) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'PYR5' ;
    'ERREUR' 21 ;
 'FINSI' ;
      
 AA = ('REDU' (tabm . 'XXSURFAC') ('MANUEL' 'POI1' PCENF)) ;
 BB = 'MANUEL' 'CHPO' ('MANU' 'POI1' PCENF) 1 'SCAL' XSURF ;
 ERRO = 'MAXIMUM' (AA '-' BB) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'PYR5' ;
    'ERREUR' 21 ;
 'FINSI' ;
      
*
* Particular case (base bizarre)
* Clockwise orientation.

 'OPTION' 'DIME' 3 'ELEM' 'PYR5' ;
 P1 = 0.0 0.0 0.12  ;
 P2 = 1.0 0.0 -0.1 ;
 P3 = 1.0 1.0 0.31 ;
 P4 = 0.5 0.4 -0.21 ;
 P5 = 0.5 0.5 2.0 ;

 P1 P2 P3 P4 P5 = P1 P2 P3 P4 P5 'TOURNER' 15.0
   (0.0 0.0 4.0) (1.1 1.1 4.0) ;
 MAIL = 'MANUEL' 'PYR5' P1 P2 P3 P4 P5 ;
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM  'VF' ; QMAIL = TABM . 'QUAF' ;
 MAIL1 = 'MANUEL' 'QUA4' P1 P2 P3 P4 ;
 MODM1  = 'MODELISER' MAIL1 'EULER' ;
 TABM1  = 'DOMA' MODM1  'VF' ; QMAIL1 = TABM1 . 'QUAF' ;

 'ELIMINATION' QMAIL QMAIL1 0.001 ;

 PCENF = 'POIN' 1 (TABM1 . 'CENTRE') ;

* 
 'SI' GRAPH ;
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
    'TITRE' 'CENTRE et MAILLAGE';
    VEC1 = 'VECTEUR' (TABM . 'XXNORMAF') 0.1 'UX' 'UY' 'UZ' 'ROUG' ;
    'TRACER' MAIL VEC1 'TITRE' 'Normals' ;
 'FINSI' ;
*
 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'PYR5' ;
    'ERREUR' 5 ;
 'FINSI' ; 

* Volume as the volume of 4 tetrahedra

 VOL1 = 'MESURE' ('MANUEL' 'TET4' P1 P2 PCENF P5) ;
 VOL2 = 'MESURE' ('MANUEL' 'TET4' P2 P3 PCENF P5) ;
 VOL3 = 'MESURE' ('MANUEL' 'TET4' P3 P4 PCENF P5) ;
 VOL4 = 'MESURE' ('MANUEL' 'TET4' P4 P1 PCENF P5) ;
 VOL = VOL1 '+' VOL2 '+' VOL3 '+' VOL4 ;

 VOL0 = 'MAXIMUM' (TABM . 'XXVOLUM') ;

 ERRO = (VOL0 '-' VOL) 'ABS' ;
 
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'PYR5' ;
    'ERREUR' 21 ;
 'FINSI' ;

* Centre
 PCEN1 =  P1 'PLUS' P2 'PLUS' PCENF 'PLUS' P5 ;
 PCEN1 = PCEN1 '/' 4. ;
 PCEN2 =  P2 'PLUS' P3 'PLUS' PCENF 'PLUS' P5 ;
 PCEN2 = PCEN2 '/' 4. ;
 PCEN3 =  P3 'PLUS' P4 'PLUS' PCENF 'PLUS' P5 ;
 PCEN3 = PCEN3 '/' 4. ;
 PCEN4 =  P4 'PLUS' P1 'PLUS' PCENF 'PLUS' P5 ;
 PCEN4 = PCEN4 '/' 4. ;

 PCEN = (VOL1 * PCEN1) 'PLUS' (VOL2 * PCEN2) 'PLUS'
        (VOL3 * PCEN3) 'PLUS' (VOL4 * PCEN4) ;
 PCEN = PCEN '/' VOL ;
 'ELIMINATION'  (PCEN 'ET' (TABM . 'CENTRE')) 0.0001 ;

 AA = 'DIFF' ('MANUEL' 'POI1' PCEN) (TABM . 'CENTRE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'PYR5' ;
    'ERREUR' 5 ;
 'FINSI' ;

* Normals and surface

 PNOR1 = 0.5 '*' ('PVEC' (P1 'MOIN' PCENF) (P4 'MOIN' PCENF)) ;
 PNOR2 = 0.5 '*' ('PVEC' (P2 'MOIN' PCENF) (P1 'MOIN' PCENF)) ;
 PNOR3 = 0.5 '*' ('PVEC' (P3 'MOIN' PCENF) (P2 'MOIN' PCENF)) ;
 PNOR4 = 0.5 '*' ('PVEC' (P4 'MOIN' PCENF) (P3 'MOIN' PCENF)) ;

 PNOR = PNOR1 'PLUS' PNOR2 'PLUS' PNOR3 'PLUS' PNOR4 ;
 XSURF = (('PSCAL' PNOR PNOR) '**' 0.5) ;
 PNOR = PNOR '/' XSURF ;
 
 AA = ('REDU' (tabm . 'XXNORMAF') ('MANUEL' 'POI1' PCENF)) ;
 BB = 'MANUEL' 'CHPO' ('MANU' 'POI1' PCENF) 3
      'UX' ('COORDONNEE' 1 PNOR)
      'UY' ('COORDONNEE' 2 PNOR)
      'UZ' ('COORDONNEE' 3 PNOR);
 ERRO = 'MAXIMUM' (AA '-' BB) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'PYR5' ;
    'ERREUR' 21 ;
 'FINSI' ;
      
 AA = ('REDU' (tabm . 'XXSURFAC') ('MANUEL' 'POI1' PCENF)) ;
 BB = 'MANUEL' 'CHPO' ('MANU' 'POI1' PCENF) 1 'SCAL' XSURF ;
 ERRO = 'MAXIMUM' (AA '-' BB) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'PYR5' ;
    'ERREUR' 21 ;
 'FINSI' ;

*
* We check that the normals modulus are 1
*
 VALNOR = TABM . 'XXNORMAF' ;
 VALSUR = TABM . 'XXSURFAC' ;

 AA = 'PSCAL' VALNOR VALNOR ('MOTS' 'UX' 'UY' 'UZ')
      ('MOTS' 'UX' 'UY' 'UZ') ;
 'SI' ((('MAXIMUM' AA) 'NEG' 1.0 1.0D-12) 'OU'
       (('MINIMUM' AA) 'NEG' 1.0 1.0D-12)) ;
    'ERREUR' 'PYR5' ;
 'FINSI' ;
 
*
* We check that the sum of (normals times surfaces) is zero
*
 AA = VALNOR '*' VALSUR ;
 BB = 'MANUEL' 'CHPO' (TABM . 'FACE') 1 'SCAL' 1.0 ;
 ERRO = 'XTY' AA BB ('MOTS' 'UX' 'UY' 'UZ')
        ('MOTS' 'SCAL' 'SCAL' 'SCAL') ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' 'PYR5' ;
 'FINSI' ;
*
* We check the triangular faces.
*

 PCEN1 = (P1 'PLUS' P2 'PLUS' P5) '/' 3. ;
 PCEN2 = (P2 'PLUS' P3 'PLUS' P5) '/' 3. ;
 PCEN3 = (P3 'PLUS' P4 'PLUS' P5) '/' 3. ;
 PCEN4 = (P4 'PLUS' P1 'PLUS' P5) '/' 3. ;
*
 A1 = 'PVECT' (P1 'MOIN' P5) (P2 'MOIN' P5) ;
 A2 = 'PVECT' (P2 'MOIN' P5) (P3 'MOIN' P5) ;
 A3 = 'PVECT' (P3 'MOIN' P5) (P4 'MOIN' P5) ;
 A4 = 'PVECT' (P4 'MOIN' P5) (P1 'MOIN' P5) ;
 VALSUR1 = 'PROG' ('PSCAL' A1 A1) ('PSCAL' A2 A2)
                  ('PSCAL' A3 A3) ('PSCAL' A4 A4) ;
 VALSUR1 = VALSUR1 '**' 0.5 ;
 VALSUR1 = VALSUR1 '*' 0.5 ;
 VALSURM = 'MANUEL' 'QUA4' PCEN1 PCEN2 PCEN3 PCEN4 ;
 'ELIMINATION' VALSURM 0.001 (TABM . 'FACE') ;
 VALSUR2 = 'MANUEL' 'CHPO' VALSURM 1 'SCAL' VALSUR1 ;   
 ERRO = 'MAXI' (('REDU' VALSUR VALSURM) '-' VALSUR2) 'ABS' ;
 'SI' (ERRO > 1.0D-7) ;
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;

*
* Matrot
*
 AA = 'EXCO' (TABM . 'MATROT') ('MOTS' 'MX' 'MY' 'MZ')
      ('MOTS' 'UX' 'UY' 'UZ') ;
 ERRO ='MAXIMUM' (AA '-' (TABM . 'XXNORMAF')) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 BB = 'EXCO' (TABM . 'MATROT') ('MOTS' 'RX' 'RY' 'RZ')
      ('MOTS' 'UX' 'UY' 'UZ') ;
 ERRO = 'MAXIMUM' ('PSCAL' AA BB ('MOTS' 'UX' 'UY' 'UZ')
                                 ('MOTS' 'UX' 'UY' 'UZ'))
        'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 COMPZ = 'PSCAL' AA BB  ('MOTS' 'UX') ('MOTS' 'UY') ;
 COMPZ = COMPZ '-' ('PSCAL' AA BB  ('MOTS' 'UY') ('MOTS' 'UX')) ;
 COMPZ1 = 'EXCO' (TABM . 'MATROT') ('MOTS' 'UZ')
          ('MOTS' 'SCAL') ;
 ERRO = 'MAXIMUM' (COMPZ '-' COMPZ1) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 COMPX = 'PSCAL' AA BB  ('MOTS' 'UY') ('MOTS' 'UZ') ;
 COMPX = COMPX '-' ('PSCAL' AA BB  ('MOTS' 'UZ') ('MOTS' 'UY')) ;
 COMPX1 = 'EXCO' (TABM . 'MATROT') ('MOTS' 'UX')
          ('MOTS' 'SCAL') ;
 ERRO = 'MAXIMUM' (COMPX '-' COMPX1) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 COMPY = 'PSCAL' AA BB  ('MOTS' 'UZ') ('MOTS' 'UX') ;
 COMPY = COMPY '-' ('PSCAL' AA BB  ('MOTS' 'UX') ('MOTS' 'UZ')) ;
 COMPY1 = 'EXCO' (TABM . 'MATROT') ('MOTS' 'UY')
          ('MOTS' 'SCAL') ;
 ERRO = 'MAXIMUM' (COMPY '-' COMPY1) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;

*
**** Plusieur pyr5 et CUB8 ;
*

 P1 = 1.0 0.0 0.0 ;
 P2 = 0.0 0.0 0.0 ;
 P3 = 0.0 1.0 0.0 ;
 P4 = 1.0 1.0 0.0 ;
 P5 = 1.0 0.0 1.0 ;
 P6 = 0.0 0.0 1.0 ;
 P7 = 0.0 1.0 1.0 ;
 P8 = 1.0 1.0 1.0 ;
 P9 = 0.5 0.5 0.5 ;
 P10= 1.0 0.0 2.0 ;
 P11= 0.0 0.0 2.0 ;
 P12= 0.0 1.0 2.0 ;
 P13= 1.0 1.0 2.0 ;
 

*

 MAIL1 = 'MANUEL' 'PYR5' P1 P2 P6 P5 P9 ;
 MAIL2 = 'MANUEL' 'PYR5' P1 P4 P8 P5 P9 ;
 MAIL3 = 'MANUEL' 'PYR5' P8 P7 P3 P4 P9 ;
 MAIL4 = 'MANUEL' 'PYR5' P7 P6 P2 P3 P9 ;
 MAIL5 = 'MANUEL' 'PYR5' P1 P2 P3 P4 P9 ;
 MAIL6 = 'MANUEL' 'PYR5' P5 P6 P7 P8 P9 ;
 MAIL7 = 'MANUEL' 'CUB8' P5 P6 P7 P8 P10 P11 P12 P13 ;
 
 
 MAILT1 = MAIL1 'ET' MAIL2 'ET' MAIL3 'ET' MAIL4 'ET' MAIL5
          'ET' MAIL6 'ET' MAIL7 ;
 MAILT2 = MAILT1 'PLUS' (1.0 0.0 0.0) ;
 MAILT3 = MAILT2 'PLUS' (1.0 0.0 0.0) ;
 MAILT1 = MAILT1 'ET' MAILT2 'ET' MAILT3 ;
 'ELIMINATION' MAILT1 0.001 ;
 MAILT2 = MAILT1 'PLUS' (0.0 1.0 0.0) ;
 MAILT3 = MAILT2 'PLUS' (0.0 1.0 0.0) ;
 MAILT1 = MAILT1 'ET' MAILT2 'ET' MAILT3 ;
 'ELIMINATION' MAILT1 0.001 ;
 MAILT2 = MAILT1 'PLUS' (0.0 0.0 2.0) ;
 MAILT3 = MAILT2 'PLUS' (0.0 0.0 2.0) ;
 MAIL   = MAILT1 'ET' MAILT2 'ET' MAILT3 ;
 'ELIMINATION' MAIL   0.001 ;
  
 ENVMAI = 'ENVELOPPE' MAIL ;
 INTMAIL = 'DIFF' ('CHANGER' 'POI1' MAIL)
           ('CHANGER' 'POI1' ENVMAI) ;
           
 CHPO1 = 'BRUI' 'BLAN' 'UNIF' 0.0 0.2 INTMAIL ;
 CHPO2 = 'BRUI' 'BLAN' 'UNIF' 0.0 0.2 INTMAIL ;
 CHPO3 = 'BRUI' 'BLAN' 'UNIF' 0.0 0.2 INTMAIL ;
 CHPO = ('NOMC' CHPO1 'UX' 'NATU' 'DISCRETE') 'ET'
        ('NOMC' CHPO2 'UY' 'NATU' 'DISCRETE') 'ET'
        ('NOMC' CHPO3 'UZ' 'NATU' 'DISCRETE') ;
        
 'FORME' CHPO ;

 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM  'VF' ; QMAIL = TABM . 'QUAF' ;

 VOLU0 = 'MESURE' MAIL ;
 VOLU1 = 'XTY' (TABM . 'XXVOLUM') ('MANUEL' 'CHPO' (TABM . 'CENTRE')
    1 'SCAL' 1.0) ('MOTS' 'SCAL') ('MOTS' 'SCAL') ;

 ERRO = (VOLU0 '-' VOLU1) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem PYR5') ;
 'FINSI' ;
   
*
******************************
**** PRI6 ********************
******************************
*
 
*
* Clockwise direction
* 

 'OPTION' 'DIME' 3 'ELEM' 'PRI6' ;
 P1 = 0.1 0.0 0.1 ;
 P3 = 1.0 1.1 -0.13 ;
 P5 = 1.11 0.11 0.13 ;
 P10 = 0.11 0.2 3.1 ;
 P12 = 1.11 1. 2.7 ;
 P14 = 1.1 0.11 3.23 ;
* 
 P1 P3 P5 P10 P12 P14 = P1 P3 P5 P10 P12 P14
    'TOURNER' 15.0 (0.0 0.0 5.0) (1.1 1.1 6.0) ;
*
 MAIL = 'MANUEL' 'PRI6' P1 P3 P5 P10 P12 P14 ;
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM  'VF' ; QMAIL = TABM . 'QUAF' ;
 QMAIL1 = 'CHANGER' QMAIL 'POI1' ;
*
 'SI' GRAPH ;
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
        'TITRE' 'CENTRE et MAILLAGE';
    VEC1 = 'VECTEUR' (TABM . 'XXNORMAF') 0.1 'UX' 'UY' 'UZ' 'ROUG' ;
    'TRACER' MAIL VEC1 'TITRE' 'Normals' ;
 'FINSI' ;
* 
 P16 = 'POIN' 16 QMAIL1 ;
 P17 = 'POIN' 17 QMAIL1 ;
 P18 = 'POIN' 18 QMAIL1 ;
 P19 = 'POIN' 19 QMAIL1 ;
 P20 = 'POIN' 20 QMAIL1 ;
*
 P21 = 'POIN' 21 QMAIL1 ; 
*
* Volume as the volume of 14 tetrahedra
*
 VOL1 = 'MESURE' ('MANUEL' 'TET4' P1  P3  P5  P21) ;
 VOL2 = 'MESURE' ('MANUEL' 'TET4' P10 P12 P14 P21) ;
* 
 VOL3 = 'MESURE' ('MANUEL' 'TET4' P1 P3 P16 P21) ;
 VOL4 = 'MESURE' ('MANUEL' 'TET4' P3 P12 P16 P21) ;
 VOL5 = 'MESURE' ('MANUEL' 'TET4' P12 P10 P16 P21) ;
 VOL6 = 'MESURE' ('MANUEL' 'TET4' P10 P1  P16 P21) ;
* 
 VOL7  = 'MESURE' ('MANUEL' 'TET4' P3 P5  P17 P21) ;
 VOL8  = 'MESURE' ('MANUEL' 'TET4' P5 P14 P17 P21) ;
 VOL9  = 'MESURE' ('MANUEL' 'TET4' P14 P12 P17 P21) ;
 VOL10 = 'MESURE' ('MANUEL' 'TET4' P12 P3  P17 P21) ;
*
 VOL11  = 'MESURE' ('MANUEL' 'TET4' P5 P14  P18 P21) ;
 VOL12  = 'MESURE' ('MANUEL' 'TET4' P14 P10 P18 P21) ;
 VOL13  = 'MESURE' ('MANUEL' 'TET4' P10 P1  P18 P21) ;
 VOL14  = 'MESURE' ('MANUEL' 'TET4' P1  P5  P18 P21) ;
*
 VOL = VOL1 '+' VOL2 '+' VOL3 '+' VOL4 '+' VOL5 '+' VOL6
   '+' VOL7 '+' VOL8 '+' VOL9 '+' VOL10 '+' VOL11 '+' VOL12
   '+' VOL13 '+' VOL14 ;
*   
 VOL0 = 'MAXIMUM' (TABM . 'XXVOLUM') ;
*
 ERRO = (VOL0 '-' VOL) 'ABS' ;
* 
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'PRI6' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
* Sommet

 MAILS = 'CHANGER' MAIL 'POI1' ;
 AA = 'DIFF' MAILS (TABM . 'SOMMET') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'PRI6' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
* FACEL, FACEP, ELTFA (comparaison 'NAVIER_STOKES' - 'EULER')

 MOD2 = 'MODELISER' QMAIL 'NAVIER_STOKES' 'LINE' ;
 TAB2 = 'DOMA' MOD2 'TABLE' ; 

 AA = 'DIFF' (TAB2 . 'FACEL') (TABM . 'FACEL') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'PRI6' ;
    'ERREUR' 21 ;
 'FINSI' ;

 AA = 'DIFF' (TAB2 . 'FACEP') (TABM . 'FACEP') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'PRI6' ;
    'ERREUR' 21 ;
 'FINSI' ;
  
 MAIL1 = (TAB2 . 'ELTFA') ;
 MAIL2 = (TABM . 'ELTFA') ;
 AA = 'DIFF' MAIL1 MAIL2 ;
 'SI' (('NBNO' AA) > 0 ) ;
      'MESSAGE' 'PRI6' ;
      'ERREUR' 21 ;
 'FINSI' ;

* XXDIEMIN  (comparaison 'NAVIER_STOKES' - 'EULER')

 CHPO1 = 'ABS' (TABM . 'XXDIEMIN') ;
 CHPO2 = 'ABS' ('DOMA' TAB2  'XXDIEMIN') ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

*
* Here we deplace the center of the gravity and we put
* it outside.
* But 'DOMA' rechanges its position.
*
* 'SI' CHECKERR ;
*    'MESSAGE' 'PRI6' ;
*    'MESSAGE' 'On deplace le centre de l''elt dehors' ;
*    'MESSAGE' ;
*    AA = 'CHANGER' 'POI1' QMAIL ;
*    DP = ('POIN' 19 AA) 'MOIN' ('POIN' 21 AA) ;
*    BB = 'MANUEL' 'CHPO' ('MANUEL' 'POI1' ('POIN' 21 AA)) 3
*         'UX' (2 '*' ('COORDONNEE' 1 DP)) 
*         'UY' (2 '*' ('COORDONNEE' 2 DP)) 
*         'UZ' (2 '*' ('COORDONNEE' 3 DP)) ;
*    'FORME' BB ;
*    'SI' GRAPH ;
*       'TRACER' QMAIL  'TITRE' 'Centre dehors' ;
*    'FINSI' ;
*    MODM = 'MODELISER' QMAIL 'EULER' ;
*    TABM = 'DOMA' MODM 'VF' ;
*    'OPTION' DONN 5 ;
* 'FINSI' ;
*

* Normals and surface

* We consider the QUA5 P1 P3 P12 P10 P16
* The point are clockwise oriented

 PNOR1 = 0.5 '*' ('PVEC' (P1 'MOIN' P16) (P3 'MOIN' P16)) ;
 PNOR2 = 0.5 '*' ('PVEC' (P3 'MOIN' P16) (P12 'MOIN' P16)) ;
 PNOR3 = 0.5 '*' ('PVEC' (P12 'MOIN' P16) (P10 'MOIN' P16)) ;
 PNOR4 = 0.5 '*' ('PVEC' (P10 'MOIN' P16) (P1 'MOIN' P16)) ;

 PNOR = PNOR1 'PLUS' PNOR2 'PLUS' PNOR3 'PLUS' PNOR4 ;
 XSURF = (('PSCAL' PNOR PNOR) '**' 0.5) ;
 PNOR = -1 '*' (PNOR '/' XSURF) ;
 
 AA = ('REDU' (tabm . 'XXNORMAF') ('MANUEL' 'POI1' P16)) ;
 BB = 'MANUEL' 'CHPO' ('MANU' 'POI1' P16) 3
      'UX' ('COORDONNEE' 1 PNOR)
      'UY' ('COORDONNEE' 2 PNOR)
      'UZ' ('COORDONNEE' 3 PNOR);
 ERRO = 'MAXIMUM' (AA '-' BB) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'PRI6' ;
    'ERREUR' 21 ;
 'FINSI' ;
      
 AA = ('REDU' (tabm . 'XXSURFAC') ('MANUEL' 'POI1' P16)) ;
 BB = 'MANUEL' 'CHPO' ('MANU' 'POI1' P16) 1 'SCAL' XSURF ;
 ERRO = 'MAXIMUM' (AA '-' BB) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'PRI6' ;
    'ERREUR' 21 ;
 'FINSI' ;

 
*
* We check that the normals modulus are 1
*
 VALNOR = TABM . 'XXNORMAF' ;
 VALSUR = TABM . 'XXSURFAC' ;
 
 AA = 'PSCAL' VALNOR VALNOR ('MOTS' 'UX' 'UY' 'UZ')
      ('MOTS' 'UX' 'UY' 'UZ') ;
 'SI' ((('MAXIMUM' AA) 'NEG' 1.0 1.0D-12) 'OU'
       (('MINIMUM' AA) 'NEG' 1.0 1.0D-12)) ;
    'ERREUR' 'PRI6' ;
 'FINSI' ;
 
*
* We check that the sum of (normals times surfaces) is zero
*

 AA = VALNOR '*' VALSUR ;
 BB = 'MANUEL' 'CHPO' (TABM . 'FACE') 1 'SCAL' 1.0 ;
 ERRO = 'XTY' AA BB ('MOTS' 'UX' 'UY' 'UZ')
        ('MOTS' 'SCAL' 'SCAL' 'SCAL') ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' 'PRI6';
 'FINSI' ;
 
*
* We check one triangular face.
*

 PCEN1 = (P1 'PLUS' P3 'PLUS' P5) '/' 3. ;
 'ELIMINATION' ((TABM . 'FACE') 'ET' PCEN1) 0.0001 ;
 
 PNOR = 0.5 '*' ('PVEC' (P3 'MOIN' P1) (P5 'MOIN' P1)) ;
 XSURF = (('PSCAL' PNOR PNOR) '**' 0.5) ;
 PNOR = (PNOR '/' XSURF) ;
 
 AA = ('REDU' (tabm . 'XXNORMAF') ('MANUEL' 'POI1' PCEN1)) ;
 BB = 'MANUEL' 'CHPO' ('MANU' 'POI1' PCEN1) 3
      'UX' ('COORDONNEE' 1 PNOR)
      'UY' ('COORDONNEE' 2 PNOR)
      'UZ' ('COORDONNEE' 3 PNOR);
 ERRO = 'MAXIMUM' (AA '-' BB) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'PRI6' ;
    'ERREUR' 21 ;
 'FINSI' ;
      
 AA = ('REDU' (tabm . 'XXSURFAC') ('MANUEL' 'POI1' PCEN1)) ;
 BB = 'MANUEL' 'CHPO' ('MANU' 'POI1' PCEN1) 1 'SCAL' XSURF ;
 ERRO = 'MAXIMUM' (AA '-' BB) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'PRI6' ;
    'ERREUR' 21 ;
 'FINSI' ;

  
******************************
**** CUB8 ********************
******************************
*
* Particular case
*
 'OPTION' 'DIME' 3 'ELEM' 'CUB8' ;

 P1 = 0.0 0.0 0.0 ;
 P2 = 1.0 0.0 0.0 ;
 P3 = 1.0 1.0 0.0 ;
 P4 = 0.0 1.0 0.0 ;
 
 MAILB = 'DALLER' (P1 'DROIT' 4 P2) (P2 'DROIT' 4 P3)
                 (P3 'DROIT' 4 P4) (P4 'DROIT' 4 P1)
        'PLAN' ;
 MAIL  = 'VOLUME' MAILB 'TRANSLATION' 4 (0.0 0.0 4.0) ;
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM 'VF' ;
 QMAIL = TABM . 'QUAF' ;

* Comparaison complete 'NAVIER_STOKES' - 'EULER'

 MOD2 = 'MODELISER' QMAIL 'NAVIER_STOKES' 'LINE' ;
 TAB2 = 'DOMA' MOD2 'TABLE' ;
 TAB2 . 'PRECONDI' = 1 ;

 AA = 'DIFF' (TAB2 . 'QUAF') (TABM . 'QUAF') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
 AA = 'DIFF' (TAB2 . 'MAILLAGE') (TABM . 'MAILLAGE') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
 AA = 'DIFF' (TAB2 . 'CENTRE') (TABM . 'CENTRE') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
 AA = 'DIFF' (TAB2 . 'SOMMET') (TABM . 'SOMMET') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
 AA = 'DIFF' (TAB2 . 'FACEL') (TABM . 'FACEL') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 AA = 'DIFF' (TAB2 . 'FACEP') (TABM . 'FACEP') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 AA = 'DIFF' (TAB2 . 'ELTFA') (TABM . 'ELTFA') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

* XXSURFAC, XXNORMAF, XXDIEMIN  (comparaison 'NAVIER_STOKES' - 'EULER')

 CHPO1 = 'DOMA' TAB2 'XXSURFAC' ;
 CHPO2 = TABM . 'XXSURFAC' ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 CHPO1 = 'DOMA' TAB2 'XXVOLUM' ;
 CHPO2 = TABM . 'XXVOLUM' ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 CHPO1 = 'ABS' (TABM . 'XXNORMAF') ;
 CHPO2 = 'ABS' ('DOMA' TAB2  'XXNORMAF') ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

* Normales dans la meme direction que FACEL

 NNO = 'NBEL' (TABM . 'FACEL') ;
 'REPETER' BL1 NNO ;
    AA = 'ELEM' &BL1 (TABM . 'FACEL') ;
    BB = 'CHANGER' AA 'POI1' ;
    P1 = 'POIN' 1 BB ;
    P2 = 'POIN' 2 BB ;
    VALX = ('COORDONNEE' 1 P2) '-' ('COORDONNEE' 1 P1) ;
    VALY = ('COORDONNEE' 2 P2) '-' ('COORDONNEE' 2 P1) ;
    VALZ = ('COORDONNEE' 3 P2) '-' ('COORDONNEE' 3 P1) ;
    VALX1 = 'EXTRAIRE' (TABM . 'XXNORMAF') 'UX' P2 ;
    VALY1 = 'EXTRAIRE' (TABM . 'XXNORMAF') 'UY' P2 ;
    VALZ1 = 'EXTRAIRE' (TABM . 'XXNORMAF') 'UZ' P2 ;
    VAL = (VALX '*' VALX1) '+' (VALY '*' VALY1)
       '+' (VALZ '*' VALZ1) ;
    'SI' (VAL < 0) ;
       'MESSAGE' 'CUB8' ;
       'ERREUR' 21 ;
    'FINSI' ;
 'FIN' BL1 ;

 CHPO1 = (TABM . 'XXDIEMIN') ;
 CHPO2 = ('DOMA' TAB2  'XXDIEMIN') ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;


*
**** MATROT dans le cas NAVIER_STOKES
*

  

*
* General case
* Clockwise orientation
*
 'OPTION' 'DIME' 3 'ELEM' 'CUB8' ;
 P1  = 0.10 0.20 0.11 ;
 P3  = 0.01 1.03 -0.11 ;
 P5  = 1.03 1.04 0.05 ;
 P7  = 1.01 0.02 0.03 ;
 P13 = 0.05 0.04 2.11 ;
 P15 = 0.08 1.01 1.91 ;
 P17 = 1.07 1.02 2.15 ;
 P19 = 1.06 0.03 2.13 ;
*
* P1 P3 P5 P7 P13 P15 P17 P19 =  P1 P3 P5 P7 P13 P15 P17 P19
*   'TOURNER' 15.0  (0.0 0.0 5.0) (1.1 1.1 6.0) ;   
*
 MAIL = 'MANUEL' 'CUB8' P1 P3 P5 P7 P13 P15 P17 P19 ;
 MODM  = 'MODELISER' MAIL 'EULER' ; QMAIL = TABM . 'QUAF' ;
 TABM  = 'DOMA' MODM 'VF' ;
 QMAIL = TABM . 'QUAF' ;
*
 QMAIL1 = 'CHANGER' QMAIL 'POI1' ;
 P21 = 'POIN' QMAIL1 21 ;
 P22 = 'POIN' QMAIL1 22 ;
 P23 = 'POIN' QMAIL1 23 ;
 P24 = 'POIN' QMAIL1 24 ;
 P25 = 'POIN' QMAIL1 25 ;
 P26 = 'POIN' QMAIL1 26 ;
 P27 = 'POIN' QMAIL1 27 ;
*
 'SI' GRAPH ;
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
    'TITRE' 'CENTRE et MAILLAGE';
 'FINSI' ;
*
 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'CUB8' ;
    'ERREUR' 5 ;
 'FINSI' ;
* 
*
* Volume
*
 VOL = 0.0D0 ;
 QC = P27 ;
 'REPETER' BL1 6 ;
    'SI' ('EGA' &BL1 1) ;
       Q1 = P1 ;  Q2 = P3 ;  Q3 = P5 ;   Q4 = P7 ;   QF = P25 ;
    'FINSI' ;
    'SI' ('EGA' &BL1 2) ;
       Q1 = P13 ;  Q2 = P15 ;  Q3 = P17 ;   Q4 = P19 ;   QF = P26 ;
    'FINSI' ;
    'SI' ('EGA' &BL1 3) ;
       Q1 = P1 ;  Q2 = P3 ;  Q3 = P15 ;   Q4 = P13 ;   QF = P21 ;
    'FINSI' ;
    'SI' ('EGA' &BL1 4) ;
       Q1 = P3 ;  Q2 = P5 ;  Q3 = P17 ;   Q4 = P15 ;   QF = P22 ;
    'FINSI' ;
    'SI' ('EGA' &BL1 5) ;
       Q1 = P7 ;  Q2 = P5 ;  Q3 = P17 ;   Q4 = P19 ;   QF = P23 ;
    'FINSI' ;
    'SI' ('EGA' &BL1 6) ;
       Q1 = P1 ;  Q2 = P13 ;  Q3 = P19 ;   Q4 = P7 ;   QF = P24 ;
    'FINSI' ;
    VOL1 = 'MESURE' ('MANUEL' 'TET4' Q1 Q2 QF QC) ;
    VOL2 = 'MESURE' ('MANUEL' 'TET4' Q2 Q3 QF QC) ;
    VOL3 = 'MESURE' ('MANUEL' 'TET4' Q3 Q4 QF QC) ;
    VOL4 = 'MESURE' ('MANUEL' 'TET4' Q4 Q1 QF QC) ;
    VOL = VOL '+' VOL1 '+' VOL2 '+' VOL3 '+' VOL4 ;
 'FIN' BL1 ;
*
 VOL0 = 'MAXIMUM' (TABM . 'XXVOLUM') ;
*
 ERRO = (VOL0 '-' VOL) 'ABS' ;
* 
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'CUB8' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
* Sommet

 MAILS = 'CHANGER' MAIL 'POI1' ;
 AA = 'DIFF' MAILS (TABM . 'SOMMET') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'CUB8' ;
    'ERREUR' 21 ;
 'FINSI' ;
 
* FACEL, FACEP, ELTFA (comparaison 'NAVIER_STOKES' - 'EULER')

 MOD2 = 'MODELISER' QMAIL 'NAVIER_STOKES' 'LINE' ;
 TAB2 = 'DOMA' MOD2 'TABLE' ; 

 AA = 'DIFF' (TAB2 . 'FACEL') (TABM . 'FACEL') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'CUB8' ;
    'ERREUR' 21 ;
 'FINSI' ;

 AA = 'DIFF' (TAB2 . 'FACEP') (TABM . 'FACEP') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'CUB8' ;
    'ERREUR' 21 ;
 'FINSI' ;
  
*MAIL1 = 'CHANGER' (TAB2 . 'ELTFA') 'POI1' ;
*MAIL2 = 'CHANGER' (TABM . 'ELTFA') 'POI1' ;
 MAIL1 = (TAB2 . 'ELTFA') ;
 MAIL2 = (TABM . 'ELTFA') ;
 AA = 'DIFF' MAIL1 MAIL2 ;
 'SI' (('NBNO' AA) > 0 ) ;
     'MESSAGE' 'CUB8' ;
      'ERREUR' 21 ;
 'FINSI' ;

* XXDIEMIN  (comparaison 'NAVIER_STOKES' - 'EULER')

 CHPO1 = 'ABS' (TABM . 'XXDIEMIN') ;
 CHPO2 = 'ABS' ('DOMA' TAB2  'XXDIEMIN') ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

* Normals and surface
*
* We consider the QUA5 P1 P3 P5 P7 P25
* The point are clockwise oriented
*

 PNOR1 = 0.5 '*' ('PVEC' (P1 'MOIN' P25) (P3 'MOIN' P25)) ;
 PNOR2 = 0.5 '*' ('PVEC' (P3 'MOIN' P25) (P5 'MOIN' P25)) ;
 PNOR3 = 0.5 '*' ('PVEC' (P5 'MOIN' P25) (P7 'MOIN' P25)) ;
 PNOR4 = 0.5 '*' ('PVEC' (P7 'MOIN' P25) (P1 'MOIN' P25)) ;

 PNOR = PNOR1 'PLUS' PNOR2 'PLUS' PNOR3 'PLUS' PNOR4 ;
 XSURF = (('PSCAL' PNOR PNOR) '**' 0.5) ;
 PNOR = (PNOR '/' XSURF) ;
 
 AA = ('REDU' (tabm . 'XXNORMAF') ('MANUEL' 'POI1' P25)) ;
 BB = 'MANUEL' 'CHPO' ('MANU' 'POI1' P25) 3
      'UX' ('COORDONNEE' 1 PNOR)
      'UY' ('COORDONNEE' 2 PNOR)
      'UZ' ('COORDONNEE' 3 PNOR);
 ERRO = 'MAXIMUM' (AA '-' BB) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'CUB8' ;
    'ERREUR' 21 ;
 'FINSI' ;
      
 AA = ('REDU' (tabm . 'XXSURFAC') ('MANUEL' 'POI1' P25)) ;
 BB = 'MANUEL' 'CHPO' ('MANU' 'POI1' P25) 1 'SCAL' XSURF ;
 ERRO = 'MAXIMUM' (AA '-' BB) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'CUB8' ;
    'ERREUR' 21 ;
 'FINSI' ;

*
* We check that the normals modulus are 1
*
 VALNOR = TABM . 'XXNORMAF' ;
 VALSUR = TABM . 'XXSURFAC' ;
 
 AA = 'PSCAL' VALNOR VALNOR ('MOTS' 'UX' 'UY' 'UZ')
      ('MOTS' 'UX' 'UY' 'UZ') ;
 'SI' ((('MAXIMUM' AA) 'NEG' 1.0 1.0D-12) 'OU'
       (('MINIMUM' AA) 'NEG' 1.0 1.0D-12)) ;
    'ERREUR' 'CUB8' ;
 'FINSI' ;
 
*
* We check that the sum of (normals times surfaces) is zero
*

 AA = VALNOR '*' VALSUR ;
 BB = 'MANUEL' 'CHPO' (TABM . 'FACE') 1 'SCAL' 1.0 ;
 ERRO = 'XTY' AA BB ('MOTS' 'UX' 'UY' 'UZ')
        ('MOTS' 'SCAL' 'SCAL' 'SCAL') ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' 'CUB8';
 'FINSI' ;

 
****************************************
***** GENERAL **************************
****************************************


 'OPTION' 'DIME' 3 'ELEM' 'CUB8' ;

 P1 = 0.0 0.0 0.0 ;
 P2 = 1.0 0.0 0.0 ;
 P3 = 1.0 1.0 0.0 ;
 P4 = 0.0 1.0 0.0 ;
 P5 = 0.0 0.0 1.0 ;
 P6 = 1.0 0.0 1.0 ;
 P7 = 1.0 1.0 1.0 ;
 P8 = 0.0 1.0 1.0 ;

 P1P2 = P1 'DROIT' 4 P2 ;
 P2P3 = P2 'DROIT' 5 P3 ;
 P3P4 = P3 'DROIT' 4 P4 ;
 P4P1 = P4 'DROIT' 5 P1 ;
 P5P6 = P5 'DROIT' 4 P6 ;
 P6P7 = P6 'DROIT' 5 P7 ;
 P7P8 = P7 'DROIT' 4 P8 ;
 P8P5 = P8 'DROIT' 5 P5 ;
 P1P5 = P1 'DROIT' 4 P5 ;
 P2P6 = P2 'DROIT' 5 P6 ;
 P3P7 = P3 'DROIT' 4 P7 ;
 P4P8 = P4 'DROIT' 5 P8 ;

 SURFB = 'SURFACE' (P1P2 'ET' P2P3 'ET' P3P4 'ET' P4P1) 'PLAN' ;
 SURFH = 'SURFACE' (P5P6 'ET' P6P7 'ET' P7P8 'ET' P8P5) 'PLAN' ;
 SURFL1 = 'SURFACE' (P1P5 'ET' ('INVERSE' P8P5) 'ET'
   ('INVERSE' P4P8) 'ET' P4P1) 'PLAN' ;
 SURFL2 = 'SURFACE' (P1P2 'ET' P2P6 'ET' ('INVERSE' P5P6) 'ET'
   ('INVERSE' P1P5)) 'PLAN' ;
 SURFL3 = 'SURFACE' (P2P6 'ET' P6P7 'ET' ('INVERSE' P3P7) 'ET'
   ('INVERSE' P2P3)) 'PLAN' ;
 SURFL4 = 'SURFACE' (('INVERSE' P3P4) 'ET' P3P7 'ET' P7P8 'ET'
 ('INVERSE' P4P8)) 'PLAN' ;

 SURFT = SURFB 'ET' SURFH 'ET' SURFL1 'ET' SURFL2 'ET' SURFL3
       'ET' SURFL4 ;

 MAIL = 'VOLUME' SURFT ;
* 
 MODM  = 'MODELISER' MAIL 'EULER' ;
 TABM  = 'DOMA' MODM 'VF' ; QMAIL = TABM . 'QUAF' ;
* 
 'SI' GRAPH ;
    'TRACER' (MAIL 'ET' QMAIL) 'TITRE' 'QUAF' ;
    'TRACER' ((tabm . 'MAILLAGE') et ((tabm . 'CENTRE') 'COUL' 'ROUGE'))
    'TITRE' 'CENTRE et MAILLAGE';
 'FINSI' ;
*
 AA = 'DIFF' MAIL (TABM . 'MAILLAGE') ;
 'SI' ('NEG' ('NBEL' AA) 0) ;
    'MESSAGE' 'CUB8' ;
    'ERREUR' 5 ;
 'FINSI' ; 
*
*** Border Surfaces
*
 MODMS  = 'MODELISER' SURFT 'EULER' ;
 TABMS  = 'DOMA' MODMS 'VF' ;
 QSURFT = TABMS . 'QUAF' ;
 'ELIMINATION' (TABMS . 'QUAF') (TABM . 'QUAF') 0.0001 ;

 VALNOR = 'REDU' (TABM . 'XXSURFAC') (TABMS . 'CENTRE') ;
 VALSUR = 'REDU' (TABM . 'XXNORMAF') (TABMS . 'CENTRE') ;
 AA = VALNOR '*' VALSUR ;
 BB = 'MANUEL' 'CHPO' (TABM . 'FACE') 1 'SCAL' 1.0 ;
 ERRO = 'XTY' AA BB ('MOTS' 'UX' 'UY' 'UZ')
        ('MOTS' 'SCAL' 'SCAL' 'SCAL') ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' 'GENE';
 'FINSI' ;
*
 LISTYP =  (TABMS . 'MAILLAGE') 'ELEM' 'TYPE' ;
*
  NELEM = 0 ;
 'REPETER' BL1 ('DIME' LISTYP) ;
    MAIL1 = 'ELEM' (TABMS . 'MAILLAGE') ('EXTRAIRE' LISTYP &BL1) ;
    'REPETER' BL2 ('NBEL' MAIL1) ;
       ELE1 = 'ELEM' MAIL1 &BL2 ;
       AA  = 'MESURE' ELE1 ;
       NELEM = NELEM '+' 1;
       P1 = 'POIN' (TABMS . 'CENTRE') NELEM ;
       BB  = 'EXTRAIRE' (TABMS . 'XXVOLUM') 'SCAL' P1 ;
       ERRO = (AA '-' BB) 'ABS' ;
       'SI' (ERRO > 1.0D-6) ;
          'MESSAGE' 'GENE' ;
          'ERREUR' 21 ;
       'FINSI' ;
    'FIN' BL2 ;
 'FIN' BL1 ;
 
* Sommet
   
 MAILS = 'CHANGER' MAIL 'POI1' ;
 AA = 'DIFF' MAILS (TABM . 'SOMMET') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'GENE' ;
    'ERREUR' 21 ;
 'FINSI' ;

  
* FACEL, FACEP, ELTFA (comparaison 'NAVIER_STOKES' - 'EULER')

 MOD2 = 'MODELISER' QMAIL 'NAVIER_STOKES' 'LINE' ;
 TAB2 = 'DOMA' MOD2 'TABLE' ; 

 AA = 'DIFF' (TAB2 . 'FACE') (TABM . 'FACE') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'GENE' ;
    'ERREUR' 21 ;
 'FINSI' ;

 AA = 'DIFF' (TAB2 . 'FACEL') (TABM . 'FACEL') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;

 AA = 'DIFF' (TAB2 . 'FACEP') (TABM . 'FACEP') ;
 'SI' (('NBNO' AA) > 0 ) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
  
 MAIL1 = (TAB2 . 'ELTFA')  ;
 MAIL2 = (TABM . 'ELTFA')  ;
 AA = 'DIFF' MAIL1 MAIL2 ;
 'SI' (('NBNO' AA) > 0 ) ;
      'MESSAGE' 'TRI3 et QUA4' ;
      'ERREUR' 21 ;
 'FINSI' ;

* XXDIEMIN  (comparaison 'NAVIER_STOKES' - 'EULER')

 CHPO1 = 'ABS' (TABM . 'XXDIEMIN') ;
 CHPO2 = 'ABS' ('DOMA' TAB2  'XXDIEMIN') ;
 ERRO = 'MAXIMUM' (CHPO1 '-' CHPO2) 'ABS' ;
 'SI' (ERRO > 1.0D-6) ;
    'MESSAGE' 'TRI3 et QUA4' ;
    'ERREUR' 21 ;
 'FINSI' ;
*
* We check that the normals modulus are 1
*
 VALNOR = TABM . 'XXNORMAF' ;
 VALSUR = TABM . 'XXSURFAC' ;
 
 AA = 'PSCAL' VALNOR VALNOR ('MOTS' 'UX' 'UY' 'UZ')
      ('MOTS' 'UX' 'UY' 'UZ') ;
 'SI' ((('MAXIMUM' AA) 'NEG' 1.0 1.0D-12) 'OU'
       (('MINIMUM' AA) 'NEG' 1.0 1.0D-12)) ;
    'ERREUR' 'CUB8' ;
 'FINSI' ;

*
* We check that each element satisfies the property  
* (normals times surfaces) = zero
*
* ELTFA and MAIL have the same structure
*
 LISTYP =  (TABM . 'ELTFA')     'ELEM' 'TYPE' ;
 LISTYP1 = (TABM . 'QUAF')      'ELEM' 'TYPE' ;
 LISTYP0 = (TABM . 'MAILLAGE')  'ELEM' 'TYPE' ;
*
 NELEM = 0 ;
* Volume as sum of single volumes 
 VOL = 0.0 ;
 'REPETER' BL1 ('DIME' LISTYP) ;
    MAIL0 = 'ELEM' (TABM . 'MAILLAGE') ('EXTRAIRE' LISTYP0 &BL1) ;
    MAIL1 = 'ELEM' (TABM . 'ELTFA') ('EXTRAIRE' LISTYP &BL1) ;
    MAIL2 = 'ELEM' (TABM . 'QUAF')  ('EXTRAIRE' LISTYP1 &BL1) ;
    NM1 = 'NBEL' MAIL1 ;
    NM2 = 'NBEL' MAIL2 ;
    'SI' ('NEG' NM1 NM2) ;
       'ERREUR' 21 ;
    'FINSI' ;
    'REPETER' BL2 ('NBEL' MAIL1) ;
       NELEM = NELEM '+' 1 ;
       ELEM1 = 'ELEM' MAIL0 &BL2 ;
       ELE1 = 'ELEM' MAIL1 &BL2 ;
       ELE2 = 'CHANGER' ELE1 'POI1' ;
       NM2 = 'NBNO' ELE2 ;
       ELE3 = 'ELEM' MAIL2 &BL2 ;
       ELE4 = 'CHANGER' ELE3 'POI1' ;
       NM4 = 'NBNO' ELE4 ;
       ELE5 = 'DIFF' ELE4 ELE2 ;
       NM5 = 'NBNO' ELE5 ;
       'SI' ('NEG' (NM2 '+' NM5) NM4) ;
          'ERREUR' 'QUAF and ELTFA non compatible' ;
          'ERREUR' 21 ;
       'FINSI' ;
       MOD2 = 'MODELISER' ELEM1 'EULER' ;
       TAB2 = 'DOMA' MOD2 'VF' ;
       ELE7 = TAB2 . 'QUAF' ;
       'ELIMINATION' ELE3 ELE7 0.0001 ;
       AA = 'DIFF' ELE7 ELE3 ;
       'SI' ('NEG' ('NBEL' AA) 0) ;
         'MESSAGE' 'GENE' ;
         'ERREUR' 'Anormale' ;
       'FINSI' ;
       VOL = VOL '+' ('MAXIMUM' (TAB2 . 'XXVOLUM')) ;
       PCEN = 'POIN' (TABM . 'CENTRE') NELEM ;
       DSURF = 0.0 0.0 0.0 ;
       'REPETER' BL3 ('NBEL' ELE2) ;
          PFAC = 'POIN' &BL3 ELE2 ;
          SURFF = 'EXTRAIRE' (TABM . 'XXSURFAC') 'SCAL' PFAC ;
          PNORM = ('EXTRAIRE' (TABM . 'XXNORMAF') 'UX' PFAC)
                  ('EXTRAIRE' (TABM . 'XXNORMAF') 'UY' PFAC)
                  ('EXTRAIRE' (TABM . 'XXNORMAF') 'UZ' PFAC) ;
          PCF = PFAC 'MOIN' PCEN ;
          DORIE = 'SIGNE' ('PSCAL' PNORM PCF) ;
          DSURF = DSURF 'PLUS' (DORIE * PNORM * SURFF) ;
       'FIN' BL3 ;
       AA = 'MAXIMUM' ('ABS' ('PROG' ('COORDONNEE' 1 DSURF)
             ('COORDONNEE' 2 DSURF) ('COORDONNEE' 3 DSURF))) ;
       'SI' (AA > 1.0D-12) ;
          'MESSAGE' 'GENE' ;
          'ERREUR' 21 ;
       'FINSI' ;
    'FIN' BL2 ;
 'FIN' BL1 ;

*
* Volume
*
 AA    = 'MESURE' MAIL ;
 CHPUN = 'MANUEL' 'CHPO' (TABM . 'CENTRE') 1 'SCAL' 1.0 ;
 BB = 'XTY' (TABM . 'XXVOLUM') CHPUN ('MOTS' 'SCAL') ('MOTS' 'SCAL') ;
 ERRO = (AA '-' BB) 'ABS' ;
 'SI' (ERRO > 5.0D-6) ;
    'MESSAGE' 'GENE' ;
    'ERREUR' 21 ;
 'FINSI' ;
 ERRO = VOL '-' BB ;
  'SI' (ERRO > 5.0D-6) ;
    'MESSAGE' 'GENE' ;
    'ERREUR' 21 ;
 'FINSI' ;

*
* Matrot
*
 AA = 'EXCO' (TABM . 'MATROT') ('MOTS' 'MX' 'MY' 'MZ')
      ('MOTS' 'UX' 'UY' 'UZ') ;
 ERRO ='MAXIMUM' (AA '-' (TABM . 'XXNORMAF')) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 BB = 'EXCO' (TABM . 'MATROT') ('MOTS' 'RX' 'RY' 'RZ')
      ('MOTS' 'UX' 'UY' 'UZ') ;
 ERRO = 'MAXIMUM' ('PSCAL' AA BB ('MOTS' 'UX' 'UY' 'UZ')
                                 ('MOTS' 'UX' 'UY' 'UZ'))
        'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 COMPZ = 'PSCAL' AA BB  ('MOTS' 'UX') ('MOTS' 'UY') ;
 COMPZ = COMPZ '-' ('PSCAL' AA BB  ('MOTS' 'UY') ('MOTS' 'UX')) ;
 COMPZ1 = 'EXCO' (TABM . 'MATROT') ('MOTS' 'UZ')
          ('MOTS' 'SCAL') ;
 ERRO = 'MAXIMUM' (COMPZ '-' COMPZ1) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 COMPX = 'PSCAL' AA BB  ('MOTS' 'UY') ('MOTS' 'UZ') ;
 COMPX = COMPX '-' ('PSCAL' AA BB  ('MOTS' 'UZ') ('MOTS' 'UY')) ;
 COMPX1 = 'EXCO' (TABM . 'MATROT') ('MOTS' 'UX')
          ('MOTS' 'SCAL') ;
 ERRO = 'MAXIMUM' (COMPX '-' COMPX1) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
*
 COMPY = 'PSCAL' AA BB  ('MOTS' 'UZ') ('MOTS' 'UX') ;
 COMPY = COMPY '-' ('PSCAL' AA BB  ('MOTS' 'UX') ('MOTS' 'UZ')) ;
 COMPY1 = 'EXCO' (TABM . 'MATROT') ('MOTS' 'UY')
          ('MOTS' 'SCAL') ;
 ERRO = 'MAXIMUM' (COMPY '-' COMPY1) 'ABS' ;
 'SI' (ERRO > 1.0D-12) ; 
    'ERREUR' ('CHAIN' 'Problem ' ('VALEUR' 'ELEM')) ;
 'FINSI' ;
 
 'FIN' ;