* fichier :  konv_impl2ord2.dgibi
***********************************************************
***********************************************************
**** APPROCHE VF "Cell-Centred Formulation" pour la    ****
**** solution des                                      **** 
**** Equations d'Euler pour un gaz parfait             ****
**** OPERATEURS PRIM, PRET, KONV                       ****
**** Implicit: calcul du jacobien du residu            ****
****           (2nd order consistent)                  ****
****                                                   ****
**** Cas gaz monoespece, "calorically perfect"         ****
****                                                   ****
**** Methodes: AUSMPLUS                                ****
****                                                   ****
**** A. BECCANTINI DRN/DMT/SEMT/LTMF    SEPT  2000     ****
***********************************************************
***********************************************************

 'OPTION'  'DIME' 2 ;
 'OPTION'  'ELEM' QUA4 ;
 'OPTION'  'ECHO'  0 ;
 'OPTION'  'TRAC' 'X' ;

 'MESSAGE' 'A eliminer' ;
 'FIN' ;

*
*** GRAPH
*

 GRAPH = FAUX ;
* GRAPH = VRAI ;

 ERRTOL = 1.0D-3 ;

***************************
***** DOMAINE SPATIAL  ****
***************************


 A0 = 0.0D0 0.0D0;
 A1 = 1.0D0 0.0D0;
 A2 = 2.0D0 0.0D0;
 A3 = 3.0D0 0.0D0;

 A0A1 = A0 'DROIT' 1 A1;
 A1A2 = A1 'DROIT' 1 A2;
 A2A3 = A2 'DROIT'  1 A3;


 DOM1 = 'TRANSLATION' A0A1 1 (0.0 1.0) ;
 DOM2 = DOM1 'PLUS' (0.0 1.0) ;
 DOM3 = DOM2 'PLUS' (0.0 1.0) ;
 DOM4 = 'TRANSLATION' A1A2 1 (0.0 1.0) ;
 DOM5 = DOM4 'PLUS' (0.0 1.0) ;
 DOM6 = DOM5 'PLUS' (0.0 1.0) ;
 DOM7 = 'TRANSLATION' A2A3 1 (0.0 1.0) ; 
 DOM8 = DOM7 'PLUS' (0.0 1.0) ;
 DOM9 = DOM8 'PLUS' (0.0 1.0) ;
 DOM1 = 'CHANGER' DOM1 'TRI3' ;

 
 DOMTOT = DOM1 'ET' DOM2 'ET' DOM3 'ET' DOM4 'ET' DOM5 'ET'
          DOM6 'ET' DOM7 'ET' DOM8 'ET' DOM9 'ELIMINATION' 0.0001 ;
 DOMCON = 'CONTOUR' DOMTOT ;

 CHPBRU = ('NOMC' 'UX' ('BRUI' 'BLAN' 'UNIF' 0.0  0.4 DOMTOT)
            'NATU' 'DISCRET') 'ET'
            ('NOMC' 'UY' ('BRUI' 'BLAN' 'UNIF' 0.0  0.2 DOMTOT)
            'NATU' 'DISCRET') ;
            
 'FORME' CHPBRU ;

 'SI' GRAPH ;  
   'TRACER' DOMTOT 'TITRE' 'Maillage' ; 
 'FINSI' ;

 $DOMTOT = 'MODELISER' DOMTOT 'EULER';

 $DOMCON = 'MODELISER'  DOMCON 'EULER';
 $DOM5 = 'MODELISER'  DOM5 'EULER';
 $DOM6 = 'MODELISER'  DOM6 'EULER';
 $DOM7 = 'MODELISER'  DOM7 'EULER';
 $DOM9 = 'MODELISER'  DOM9 'EULER';

 TDOMTOT = 'DOMA'  $DOMTOT 'VF';
 
 TDOMCON = 'DOMA'  $DOMCON 'VF';
 TDOM5 = 'DOMA'  $DOM5 'VF';
 TDOM6 = 'DOMA'  $DOM6 'VF';
 TDOM7 = 'DOMA'  $DOM7 'VF';
 TDOM9 = 'DOMA'  $DOM9 'VF';

 MDOMCON = TDOMCON . 'QUAF' ;
 MDOM5 = TDOM5 . 'QUAF' ;
 MDOM6 = TDOM6 . 'QUAF' ;
 MDOM7 = TDOM7 . 'QUAF' ;
 MDOM9 = TDOM9 . 'QUAF' ;          
 
 MDOMTOT = TDOMTOT . 'QUAF' ;

 'ELIMINATION'  (MDOMTOT ET MDOMCON) 0.0001 ;
 'ELIMINATION'  (MDOMTOT ET MDOM5) 0.0001 ;
 'ELIMINATION'  (MDOMTOT ET MDOM6) 0.0001 ;
 'ELIMINATION'  (MDOMTOT ET MDOM7) 0.0001 ;
 'ELIMINATION'  (MDOMTOT ET MDOM9) 0.0001 ;
 
 VNLIM = 'MANUEL' 'CHPO' ('DOMA' $DOMCON 'CENTRE') 2 'UX' 1.0
   'UY' 0.0 ;
   
* VNLIM PIPI = 'KOPS' MATRIK ;

***************************************************
***** Densité, pression, vitesse, gamma ***********
***************************************************

 RN0 = 'BRUI' 'BLAN' 'UNIF' 1.11 0.5 ('DOMA' $DOMTOT 'CENTRE') ;
 PN0 = 'BRUI' 'BLAN' 'UNIF' 1234.1 800 ('DOMA' $DOMTOT 'CENTRE') ;
 GAMMAN = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'SCAL' 1.4
   'NATU' 'DISCRET' ;
 CSONN = (GAMMAN '*' PN0) '/' RN0 ;
 CSONN = 'KOPS' CSONN '**' 0.5 ;
 UXN0  = 0.2 *  CSONN ;
 UYN0  = 0.3  * CSONN ;
 
 GN0 = ('NOMC' (RN0 '*' UXN0) 'UX' 'NATU' 'DISCRET') 'ET'
       ('NOMC' (RN0 '*' UYN0) 'UY' 'NATU' 'DISCRET') ;

 ECIN = 0.5D0 '*' RN0 '*' ((UXN0  '*' UXN0) '+' (UYN0  '*' UYN0)) ;
 RETN0  = 'NOMC' 'SCAL' ((PN0  '/' (GAMMAN  '-' 1.0)) '+' ECIN)
   'NATU' 'DISCRET' ;

 VIT PRES = PRIM 'PERFMONO' RN0 GN0 RETN0 GAMMAN ;

 ERR1 = 'MAXIMUM' (PRES '-' PN0) 'ABS' ;
 ERR2 = 'MAXIMUM' (GN0 '-' (RN0 '*' VIT ('MOTS' 'SCAL' 'SCAL')
         ('MOTS' 'UX' 'UY') ('MOTS' 'UX' 'UY'))) 'ABS' ;

 'SI' (('MAXIMUM' ('PROG' ERR1 ERR2 ) 'ABS' ) > 1.0D-6) ;
      'MESSAGE' 'Problem 0' ;
      'ERREUR' 5 ;
 'FINSI' ;

 'SI' GRAPH;
      'TRACER' (('DOMA' $DOMTOT 'MAILLAGE') 'ET' 
      ('DOMA' $DOMTOT 'CENTRE')) 'TITRE' 'Domaine et centre' ;
 'FINSI' ;

*
**** Les variables conservative
*
*    RN0 (densité)
*    GN0 (quantité de mouvement)
*    RETN0 (énergie totale par unité de volume)
*
*    sont definiés
*

****************************************************
****************************************************
******* Calcul du jacobien et du residu   **********
****************************************************
****************************************************
*
* JACO est le jacobien
*
* DEBRN0 le residu concernant la densité
* DEBGNX0 le residu concernant la quantité de mouvement (axe x)
* DEBGNY0 le residu concernant la quantité de mouvement (axe y)
* DEBRETN0 le residu concernant l'enrgie totale par unité de volume
*

* Noms des variables

 NOMDEN = 'RN  ' ;
 NOMMOX = 'RUXN' ;
 NOMMOY = 'RUYN' ;
 NOMRET = 'RETN' ;

* Metode
 
 METO = 'AUSMPLUS' ;
 
  
 LISTINCO = 'MOTS'  NOMDEN  NOMMOX  NOMMOY  NOMRET ;
 NOMRN = 'MOTS' 'SCAL' ;
 NOMVN = 'MOTS' 'UX' 'UY'  ;

 VITESSE PRES   = 'PRIM' 'PERFMONO'  RN0 GN0 RETN0 GAMMAN ;


 GRN IARN CHAMRN = 'PENT' $DOMTOT 'CENTRE' 'EULESCAL' 'NOLIMITE'
    NOMRN RN0 ;
 GPN IAPN CHAMPN = 'PENT' $DOMTOT 'CENTRE' 'EULESCAL' 'NOLIMITE'
    NOMRN PRES ;
 GVN IAVN CHAMVN = 'PENT' $DOMTOT 'CENTRE' 'EULEVECT' 'NOLIMITE'
    NOMVN VITESSE 'CLIM' VNLIM ;


 IARN = IARN * 0.5 ;
 IAVN = ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE')
         2 'P1' 0.70 'P2' 0.6) ;
 IAPN = (IAPN *  0.8) ;  
  
 ROF VITF PF GAMF   =  'PRET' 'PERFMONO'  2 1
                       $DOMTOT
                        RN0   GRN IARN
                        VITESSE GVN IAVN
                        PRES GPN  IAPN
                        GAMMAN ;

  CHPRES0 DT =  'KONV' 'VF' 'PERFMONO' 'RESI' METO
          $DOMTOT ROF VITF PF GAMF LISTINCO ;

  JACO = 'KONV' 'VF' 'PERFMONO' 'JACOCON2' METO 
          $DOMTOT LISTINCO ROF VITF PF GAMF 
          RN0 VITESSE PRES GAMMAN
          GRN GVN GPN
          IARN IAVN IAPN           
          CHAMRN CHAMVN CHAMPN VNLIM ;

 DEBRN0   = 'EXCO' NOMDEN CHPRES0 ;
 DEBGNX0  = 'EXCO' NOMMOX CHPRES0 ;
 DEBGNY0  = 'EXCO' NOMMOY CHPRES0 ;
 DEBRETN0 = 'EXCO' NOMRET CHPRES0 ;

*****************************************************
*****************************************************
******** PROCEDURES *********************************
*****************************************************
*****************************************************

*
* Derivé partielle du residu en un point par rapport
* aux variable en un autre point
*

 'DEBPROC' JACOVA JACO*'MATRIK' $MODE*'MMODEL' LISTINCO*'LISTMOTS'
    PPRIM*'POINT' PDUAL*'POINT' MOTPRI*'MOT' MOTDUA*'MOT';
    
*   PPRIM = point ou est localisé la variable primale
*   PDUAL = point ou est localisé la variable duale
*   MOTPRI = nom de la composante concernante la variable primale
*   MOTDUA = nom de la composante concernante la variable duale

    ELT1 = 'MANUEL' 'POI1' PPRIM ;
    NDIM = 'DIME' LISTINCO ;
    CHPUN = 'MANUEL' 'CHPO' ELT1 1 MOTPRI 1.0 'NATURE' 'DISCRET' ;
    'REPETER' BL1 NDIM ;
       MOTCEL = 'EXTRAIRE' LISTINCO &BL1 ;
       CHPUN = CHPUN 'ET' ('MANUEL' 'CHPO' ('DOMA' $MODE 'CENTRE')
          1 MOTCEL 0.0 'NATURE' 'DISCRET') ;
    'FIN' BL1 ;
    D_DMOT = 'KOPS' JACO 'MULT' CHPUN ;
    SCAL = 'EXTRAIRE' D_DMOT PDUAL MOTDUA ;

 'FINPROC' SCAL  ;


 
 'DEBPROC' JACNUM $MODE*'MMODEL' NOMMET*'MOT' RN*'CHPOINT' GN*'CHPOINT'
    RETN*'CHPOINT' GAMN*'CHPOINT' VNLIM*'CHPOINT'
    IARN*'CHPOINT' IAVN*'CHPOINT' IAPN*'CHPOINT'    
    LISTINCO*'LISTMOTS' PPRIM*'POINT' PDUAL*'POINT'
    MOTPRI*'MOT' MOTDUA*'MOT' EPSILON*'FLOTTANT';

*   PPRIM = point ou est localisé la variable primale
*   PDUAL = point ou est localisé la variable duale
*   MOTPRI = nom de la composante concernante la variable primale
*   MOTDUA = nom de la composante concernante la variable duale

*   Le valeur dans l'état non-perturbé en PDUAL ;
    
    RNCEL   = 'COPIER'  RN ;
    GNCEL   = 'COPIER'  GN ;
    RETNCEL = 'COPIER'  RETN ;

    VITESSE PRES   = 'PRIM' 'PERFMONO'  RNCEL GNCEL RETNCEL GAMN ;


    GRN IARN0 CHAMRN = 'PENT' $MODE 'CENTRE' 'EULESCAL' 'NOLIMITE'
      NOMRN  RNCEL ;
    GPN IAPN0 CHAMPN = 'PENT' $MODE 'CENTRE' 'EULESCAL' 'NOLIMITE'
      NOMRN  PRES ;
    GVN IAVN0 CHAMVN = 'PENT' $MODE 'CENTRE' 'EULEVECT' 'NOLIMITE'
      NOMVN  VITESSE 'CLIM' VNLIM ;

* IARN0, IAPN0, IAVN0 ne sont pas utilisés
        
    ROF VITF PF GAMF   =  'PRET' 'PERFMONO'  2 1
                           $MODE
                           RNCEL   GRN IARN
                           VITESSE GVN IAVN
                           PRES GPN  IAPN
                           GAMN ;
    
    CHPRES0 DT = 'KONV' 'VF' 'PERFMONO' 'RESI' NOMMET
          $MODE ROF VITF PF GAMF LISTINCO ;

    VAL0 = 'EXTRAIRE' CHPRES0 PDUAL MOTDUA ;
    
*   EPSILON = perturbation

*   Adimensionalisation

    dens0 = 'EXTRAIRE' RN PPRIM 'SCAL' ;
    
    VN PN = 'PRIM' 'PERFMONO' RN GN RETN GAMN ;
    CN2 = GAMN '*' (PN '/' RN) ;
    cson0 = ('EXTRAIRE' CN2 PPRIM 'SCAL') '**' 0.5 ;

    ret0 = ('EXTRAIRE' RETN PPRIM 'SCAL') '**' 0.5 ;

*   On etabli la variable à perturber

 
    NDIM = 'DIME' LISTINCO ;
    ICEL = NDIM '+' 1 ;
    'REPETER' BL1 NDIM ;
       MOTCEL = 'EXTRAIRE' LISTINCO &BL1 ;
       'SI' ('EGA' MOTCEL MOTPRI) ;
          ICEL = &BL1 ;
          'QUITTER' BL1 ;
       'FINSI' ;
    'FIN' BL1 ;
    
    'SI' (ICEL > NDIM) ;
       'MESSAGE' 'Procedure JACNUM' ;
       'MESSAGE' 'MOTPRI = ??? ';
       'ERREUR' 21 ;
    'FINSI' ;
    
    ELT1 = 'MANUEL' 'POI1' PPRIM ;
    
*   ICEL = 1 -> On perturbe la densité

    'SI' ('EGA' ICEL 1) ;
       DELTATOT = (EPSILON * dens0) ;
       RNCEL = ('MANUEL' 'CHPO' ELT1 1 'SCAL' DELTATOT
          'NATURE' 'DISCRET') 'ET' RN ;
       GNCEL = 'COPIER' GN ;
       RETNCEL = 'COPIER' RETN ;
    'FINSI' ;

*   ICEL = 2 -> On perturbe la q.d.m. long  l'ax x

    'SI' ('EGA' ICEL 2) ;
       DELTATOT = (EPSILON * dens0 * cson0) ;
       GNCEL = ('MANUEL' 'CHPO' ELT1 1 'UX' DELTATOT
          'NATURE' 'DISCRET') 'ET' GN ;
       RNCEL = 'COPIER' RN ;
       RETNCEL = 'COPIER' RETN ;
    'FINSI' ;

*   ICEL = 3 -> On perturbe la q.d.m. long  l'ax y

    'SI' ('EGA' ICEL 3) ;
       DELTATOT = (EPSILON * dens0 * cson0) ;
       GNCEL = ('MANUEL' 'CHPO' ELT1 1 'UY' DELTATOT
          'NATURE' 'DISCRET') 'ET' GN ;
       RNCEL = 'COPIER' RN ;
       RETNCEL = 'COPIER' RETN ;
    'FINSI' ;

*   ICEL = 4 -> On perturbe l'énergie totale

    'SI' ('EGA' ICEL 4) ;
       DELTATOT = (EPSILON * ret0) ;
       RETNCEL = ('MANUEL' 'CHPO' ELT1 1 'SCAL' DELTATOT
          'NATURE' 'DISCRET') 'ET' RETN ;
       RNCEL = 'COPIER' RN ;
       GNCEL = 'COPIER' GN ;
    'FINSI' ;

    VITESSE PRES   = 'PRIM' 'PERFMONO'  RNCEL GNCEL RETNCEL GAMN ;
    

    GRN IARN0 CHAMRN = 'PENT' $MODE 'CENTRE' 'EULESCAL' 'NOLIMITE'
      NOMRN  RNCEL ;
    GPN IAPN0 CHAMPN = 'PENT' $MODE 'CENTRE' 'EULESCAL' 'NOLIMITE'
      NOMRN   PRES ;
    GVN IAVN0 CHAMVN = 'PENT' $MODE 'CENTRE' 'EULEVECT' 'NOLIMITE'
      NOMVN  VITESSE 'CLIM' VNLIM ;

* IARN0, IAPN0, IAVN0 ne sont pas utilisés

    ROF VITF PF GAMF   =  'PRET' 'PERFMONO'  2 1
                           $MODE
                           RNCEL   GRN IARN
                           VITESSE GVN IAVN
                           PRES GPN  IAPN
                           GAMN ;
    
    CHPRES1 DT = 'KONV' 'VF' 'PERFMONO' 'RESI' NOMMET
          $MODE ROF VITF PF GAMF LISTINCO ;

    VAL1 = 'EXTRAIRE' CHPRES1 PDUAL MOTDUA ;

 'FINPROC' ((VAL1 '-' VAL0) '/' DELTATOT)  ;

*****************************************************
*****************************************************
******** FIN PROCEDURES *****************************
*****************************************************
*****************************************************
*****************************************************


*****************************************************
******* TEST0 ***************************************
*****************************************************
*
* On compare le jacobien et la variation du residu
* en $DOM5 'CENTRE' par rapport à une variation
* infinitésimal en $DOM5 'CENTRE'
*

 PCEN5 = 'POIN' 1 ('DOMA' $DOM5 'CENTRE') ;

* Les grandeurs pour adimesionner les erreurs

 ro0 = 'EXTRAIRE' RN0 PCEN5 'SCAL' ;
    
 VN0 PN0 = 'PRIM' 'PERFMONO' RN0 GN0 RETN0 GAMMAN ;
 CN20 = GAMMAN '*' (PN0 '/' RN0) ;
 cson0 = ('EXTRAIRE' CN20 PCEN5 'SCAL') '**' 0.5 ;

 ret0 = ('EXTRAIRE' RETN0 PCEN5 'SCAL') '**' 0.5 ;

*
* Le jacobien exact. 
* DRR  = d(RES_ro)/dro (variable primale en PCEN5, variable duale en PCEN5) ;
* DGXR = d(RES_gx)/dro (variable primale en PCEN5, variable duale en PCEN5) ;
* ...

 DRR   = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMDEN NOMDEN ;
 DGXR  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMDEN NOMMOX ;
 DGYR  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMDEN NOMMOY ;
 DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMDEN NOMRET ;

 DRGX   = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOX NOMDEN ;
 DGXGX  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOX NOMMOX ;
 DGYGX  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOX NOMMOY ;
 DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOX NOMRET ;

 DRGY   = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOY NOMDEN ;
 DGXGY  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOY NOMMOX ;
 DGYGY  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOY NOMMOY ;
 DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOY NOMRET ;

 DRRET   = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMRET NOMDEN ;
 DGXRET  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMRET NOMMOX ;
 DGYRET  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMRET NOMMOY ;
 DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMRET NOMRET ;

* Le jacobien numerique


 DELTA = 1.0D-4 ;
 
 DRRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM 
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMDEN NOMDEN DELTA ;
 DGXRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMDEN NOMMOX DELTA ;
 DGYRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMDEN NOMMOY DELTA ;
 DRETRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMDEN NOMRET DELTA ;
   
 DRGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMMOX NOMDEN DELTA ;
 DGXGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMMOX NOMMOX DELTA ;
 DGYGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMMOX NOMMOY DELTA ;
 DRETGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMMOX NOMRET DELTA ;

 DRGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMMOY NOMDEN DELTA ;
 DGXGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMMOY NOMMOX DELTA ;
 DGYGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMMOY NOMMOY DELTA ;
 DRETGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMMOY NOMRET DELTA ;

 DRRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMRET NOMDEN DELTA ;
 DGXRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMRET NOMMOX DELTA ;
 DGYRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMRET NOMMOY DELTA ;
 DRETRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN5 PCEN5
   NOMRET NOMRET DELTA ;


* Test des comparaisons jacobien exact-jacobien numerique

 
 ERR1 = ('ABS' (DRRN '-' DRR)) '*' (ro0 '/' (ro0 '*' cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 1';
      'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXRN '-' DGXR)) '*'
          (ro0 '/' (ro0 '*' cson0 '*' cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 2';
      'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYRN '-' DGYR)) '*'
          (ro0 '/' (ro0 '*' cson0 '*' cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 3';
       'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRETRN '-' DRETR)) '/'
        ('MAXIMUM' ('PROG'  ((ret0 '*' cson0) '/' ro0)
         DRETR) 'ABS') ; 
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 4';
      'ERREUR' 5 ;
 'FINSI' ;

   
 ERR1 = ('ABS' (DRGXN '-' DRGX))  ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 5';
      'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXGXN '-' DGXGX)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 6';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYGXN '-' DGYGX)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 7';
        'ERREUR' 5 ;
 'FINSI' ;


 ERR1 = ('ABS' (DRETGXN '-' DRETGX)) '/'
        ('MAXIMUM' ('PROG' (ret0 '/' ro0) DRETGX)
        'ABS' ) ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 8';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRGYN '-' DRGY))  ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 9';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXGYN '-' DGXGY))  '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 10';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYGYN '-' DGYGY)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 11';
        'ERREUR' 5 ;
 'FINSI' ;


 ERR1 = ('ABS' (DRETGYN '-' DRETGY))  '/'
        ('MAXIMUM' ('PROG' (ret0 / ro0)
         DRETGYN) 'ABS' ) ;
         
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 12';
        'ERREUR' 5 ;
 'FINSI' ;

 
 ERR1 = ('ABS' (DRRETN '-' DRRET)) '*' (ret0 '/' (ro0 '*' cson0))   ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 13';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXRETN '-' DGXRET)) '*'
          (ret0 '/' (ro0 '*' cson0 * cson0)) ; 
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 14';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYRETN '-' DGYRET)) '*' 
             (ret0 '/' (ro0 '*' cson0 * cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 15';
        'ERREUR' 5 ;
 'FINSI' ;


 ERR1 = ('ABS' (DRETRETN '-' DRETRET)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 16';
        'ERREUR' 5 ;
 'FINSI' ;


*****************************************************
******* TEST1 ***************************************
*****************************************************
*
* On compare le jacobien et la variation du residu
* en $DOM9 'CENTRE' par rapport à une variation
* infinitésimal en $DOM9 'CENTRE'
*

 PCEN9 = 'POIN' 1 ('DOMA' $DOM9 'CENTRE') ;

* Les grandeurs pour adimesionner les erreurs

 ro0 = 'EXTRAIRE' RN0 PCEN9 'SCAL' ;
    
 VN0 PN0 = 'PRIM' 'PERFMONO' RN0 GN0 RETN0 GAMMAN ;
 CN20 = GAMMAN '*' (PN0 '/' RN0) ;
 cson0 = ('EXTRAIRE' CN20 PCEN9 'SCAL') '**' 0.5 ;

 ret0 = ('EXTRAIRE' RETN0 PCEN9 'SCAL') '**' 0.5 ;

*
* Le jacobien exact. 
* DRR  = d(RES_ro)/dro (variable primale en PCEN9, variable duale en PCEN9) ;
* DGXR = d(RES_gx)/dro (variable primale en PCEN9, variable duale en PCEN9) ;
* ...

 DRR   = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMDEN NOMDEN ;
 DGXR  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMDEN NOMMOX ;
 DGYR  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMDEN NOMMOY ;
 DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMDEN NOMRET ;

 DRGX   = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOX NOMDEN ;
 DGXGX  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOX NOMMOX ;
 DGYGX  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOX NOMMOY ;
 DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOX NOMRET ;

 DRGY   = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOY NOMDEN ;
 DGXGY  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOY NOMMOX ;
 DGYGY  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOY NOMMOY ;
 DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOY NOMRET ;

 DRRET   = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMRET NOMDEN ;
 DGXRET  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMRET NOMMOX ;
 DGYRET  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMRET NOMMOY ;
 DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMRET NOMRET ;

* Le jacobien numerique


 DELTA = 1.0D-4 ;
 
 DRRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMDEN NOMDEN DELTA ;
 DGXRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMDEN NOMMOX DELTA ;
 DGYRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMDEN NOMMOY DELTA ;
 DRETRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMDEN NOMRET DELTA ;
   
 DRGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMMOX NOMDEN DELTA ;
 DGXGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMMOX NOMMOX DELTA ;
 DGYGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMMOX NOMMOY DELTA ;
 DRETGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMMOX NOMRET DELTA ;

 DRGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMMOY NOMDEN DELTA ;
 DGXGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMMOY NOMMOX DELTA ;
 DGYGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMMOY NOMMOY DELTA ;
 DRETGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMMOY NOMRET DELTA ;

 DRRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMRET NOMDEN DELTA ;
 DGXRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMRET NOMMOX DELTA ;
 DGYRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMRET NOMMOY DELTA ;
 DRETRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN9 PCEN9
   NOMRET NOMRET DELTA ;


* Test des comparaisons jacobien exact-jacobien numerique

 ERR1 = ('ABS' (DRRN '-' DRR)) '*' (ro0 '/' (ro0 '*' cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 1';
      'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXRN '-' DGXR)) '*'
          (ro0 '/' (ro0 '*' cson0 '*' cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 2';
      'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYRN '-' DGYR)) '*'
          (ro0 '/' (ro0 '*' cson0 '*' cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 3';
       'ERREUR' 5 ;
 'FINSI' ;

 
 ERR1 = ('ABS' (DRETRN '-' DRETR)) '/'
        ('MAXIMUM' ('PROG' ((ret0 '*' cson0) '/' ro0 )
          DRETRN )'ABS' ) ; 
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 4';
      'ERREUR' 5 ;
 'FINSI' ;

   
 ERR1 = ('ABS' (DRGXN '-' DRGX))  ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 5';
      'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXGXN '-' DGXGX)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 6';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYGXN '-' DGYGX)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 7';
        'ERREUR' 5 ;
 'FINSI' ;


 ERR1 = ('ABS' (DRETGXN '-' DRETGX)) '*'
        (ro0  / ret0) ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 8';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRGYN '-' DRGY))  ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 9';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXGYN '-' DGXGY))  '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 10';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYGYN '-' DGYGY)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 11';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRETGYN '-' DRETGY))  '/'
        ('MAXIMUM' ('PROG' (ret0 '/' ro0) DRETGY) 'ABS' ) ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 12';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRRETN '-' DRRET)) '*' (ret0 '/' (ro0 '*' cson0))   ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 13';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXRETN '-' DGXRET)) '*'
          (ret0 '/' (ro0 '*' cson0 * cson0))   ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 14';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYRETN '-' DGYRET)) '*' 
             (ret0 '/' (ro0 '*' cson0 * cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 15';
        'ERREUR' 5 ;
 'FINSI' ;


 ERR1 = ('ABS' (DRETRETN '-' DRETRET)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 16';
        'ERREUR' 5 ;
 'FINSI' ;
 
*****************************************************
*****************************************************
******* TEST2 ***************************************
*****************************************************
*****************************************************
*
* On compare le jacobien et la variation du residu
* en $DOM9 'CENTRE' par rapport à une variation
* infinitésimal en $DOM6 'CENTRE'
*

 PCEN6 = ('DOMA' $DOM6 'CENTRE') 'POIN' 1 ;
 
 ro0 = 'EXTRAIRE' RN0 PCEN6 'SCAL' ;
 cson0 = ('EXTRAIRE' CN20 PCEN6 'SCAL') '**' 0.5 ;
 ret0 = 'EXTRAIRE' RETN0 PCEN6 'SCAL' ;

*
* Le jacobien exact. 
*

 DRR   = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMDEN NOMDEN ;
 DGXR  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMDEN NOMMOX ;
 DGYR  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMDEN NOMMOY ;
 DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMDEN NOMRET ;

 DRGX   = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOX NOMDEN ;
 DGXGX  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOX NOMMOX ;
 DGYGX  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOX NOMMOY ;
 DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOX NOMRET ;

 DRGY   = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOY NOMDEN ;
 DGXGY  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOY NOMMOX ;
 DGYGY  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOY NOMMOY ;
 DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOY NOMRET ;

 DRRET   = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMRET NOMDEN ;
 DGXRET  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMRET NOMMOX ;
 DGYRET  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMRET NOMMOY ;
 DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMRET NOMRET ;

* Le jacobien numerique


 DELTA = 1.0D-4 ;
 
 DRRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMDEN NOMDEN DELTA ;
 DGXRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMDEN NOMMOX DELTA ;
 DGYRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMDEN NOMMOY DELTA ;
 DRETRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMDEN NOMRET DELTA ;
   
 DRGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMMOX NOMDEN DELTA ;
 DGXGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMMOX NOMMOX DELTA ;
 DGYGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMMOX NOMMOY DELTA ;
 DRETGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMMOX NOMRET DELTA ;

 DRGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMMOY NOMDEN DELTA ;
 DGXGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMMOY NOMMOX DELTA ;
 DGYGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMMOY NOMMOY DELTA ;
 DRETGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMMOY NOMRET DELTA ;

 DRRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMRET NOMDEN DELTA ;
 DGXRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMRET NOMMOX DELTA ;
 DGYRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6 PCEN9 NOMRET NOMMOY DELTA ;
 DRETRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN6  PCEN9
   NOMRET NOMRET DELTA ;

* Test des comparaisons jacobien exact-jacobien numerique

 ERR1 = ('ABS' (DRRN '-' DRR)) '*' (ro0 '/' (ro0 '*' cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 1';
      'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXRN '-' DGXR)) '*'
          (ro0 '/' (ro0 '*' cson0 '*' cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 2';
      'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYRN '-' DGYR)) '*'
          (ro0 '/' (ro0 '*' cson0 '*' cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 3';
       'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRETRN '-' DRETR)) '*'
        (ro0 '/' (ret0 '*' cson0)) ; 
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 4';
      'ERREUR' 5 ;
 'FINSI' ;

   
 ERR1 = ('ABS' (DRGXN '-' DRGX))  ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 5';
      'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXGXN '-' DGXGX)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 6';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYGXN '-' DGYGX)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 7';
        'ERREUR' 5 ;
 'FINSI' ;


 ERR1 = ('ABS' (DRETGXN '-' DRETGX)) '*'
        (ro0  / ret0) ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 8';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRGYN '-' DRGY))  ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 9';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXGYN '-' DGXGY))  '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 10';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYGYN '-' DGYGY)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 11';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRETGYN '-' DRETGY))  '/'
        ('MAXIMUM' ('PROG' (ret0 '/' ro0) DRETGY) 'ABS' ) ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 12';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRRETN '-' DRRET)) '*' (ret0 '/' (ro0 '*' cson0))   ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 13';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXRETN '-' DGXRET)) '*'
          (ret0 '/' (ro0 '*' cson0 * cson0))   ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 14';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYRETN '-' DGYRET)) '*' 
             (ret0 '/' (ro0 '*' cson0 * cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 15';
        'ERREUR' 5 ;
 'FINSI' ;


 ERR1 = ('ABS' (DRETRETN '-' DRETRET)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 16';
        'ERREUR' 5 ;
 'FINSI' ;
 
  
*****************************************************
*****************************************************
******* TEST3 ***************************************
*****************************************************
*****************************************************
*
* On compare le jacobien et la variation du residu
* en $DOM9 'CENTRE' par rapport à une variation
* infinitésimal en $DOM7 'CENTRE'
*

 PCEN7 = ('DOMA' $DOM7 'CENTRE') 'POIN' 1 ;
 
 ro0 = 'EXTRAIRE' RN0 PCEN7 'SCAL' ;
 cson0 = ('EXTRAIRE' CN20 PCEN7 'SCAL') '**' 0.5 ;
 ret0 = 'EXTRAIRE' RETN0 PCEN7 'SCAL' ;

*
* Le jacobien exact. 
*

 DRR   = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMDEN NOMDEN ;
 DGXR  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMDEN NOMMOX ;
 DGYR  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMDEN NOMMOY ;
 DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMDEN NOMRET ;

 DRGX   = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOX NOMDEN ;
 DGXGX  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOX NOMMOX ;
 DGYGX  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOX NOMMOY ;
 DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOX NOMRET ;

 DRGY   = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOY NOMDEN ;
 DGXGY  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOY NOMMOX ;
 DGYGY  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOY NOMMOY ;
 DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOY NOMRET ;

 DRRET   = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMRET NOMDEN ;
 DGXRET  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMRET NOMMOX ;
 DGYRET  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMRET NOMMOY ;
 DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMRET NOMRET ;

* Le jacobien numerique


 DELTA = 1.0D-4 ;
 
 DRRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMDEN NOMDEN DELTA ;
 DGXRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMDEN NOMMOX DELTA ;
 DGYRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMDEN NOMMOY DELTA ;
 DRETRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMDEN NOMRET DELTA ;
   
 DRGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMMOX NOMDEN DELTA ;
 DGXGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMMOX NOMMOX DELTA ;
 DGYGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMMOX NOMMOY DELTA ;
 DRETGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMMOX NOMRET DELTA ;

 DRGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMMOY NOMDEN DELTA ;
 DGXGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMMOY NOMMOX DELTA ;
 DGYGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMMOY NOMMOY DELTA ;
 DRETGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMMOY NOMRET DELTA ;

 DRRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMRET NOMDEN DELTA ;
 DGXRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMRET NOMMOX DELTA ;
 DGYRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7 PCEN9 NOMRET NOMMOY DELTA ;
 DRETRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN VNLIM   
   IARN IAVN IAPN LISTINCO
   PCEN7  PCEN9
   NOMRET NOMRET DELTA ;

* Test des comparaisons jacobien exact-jacobien numerique

 
 ERR1 = ('ABS' (DRRN '-' DRR)) '*' (ro0 '/' (ro0 '*' cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 1';
      'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXRN '-' DGXR)) '*'
          (ro0 '/' (ro0 '*' cson0 '*' cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 2';
      'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYRN '-' DGYR)) '*'
          (ro0 '/' (ro0 '*' cson0 '*' cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 3';
       'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRETRN '-' DRETR)) '*'
        (ro0 '/' (ret0 '*' cson0)) ; 
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 4';
      'ERREUR' 5 ;
 'FINSI' ;

   
 ERR1 = ('ABS' (DRGXN '-' DRGX))  ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 5';
      'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXGXN '-' DGXGX)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 6';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYGXN '-' DGYGX)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 7';
        'ERREUR' 5 ;
 'FINSI' ;


 ERR1 = ('ABS' (DRETGXN '-' DRETGX)) '/'
        ('MAXIMUM' ('PROG' (ret0 '/' ro0) DRETGX)
        'ABS' ) ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 8';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRGYN '-' DRGY))  ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 9';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXGYN '-' DGXGY))  '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 10';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYGYN '-' DGYGY)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 11';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRETGYN '-' DRETGY))  '*'
        (ro0 / ret0)  ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 12';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DRRETN '-' DRRET)) '*' (ret0 '/' (ro0 '*' cson0))   ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 13';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGXRETN '-' DGXRET)) '*'
          (ret0 '/' (ro0 '*' cson0 * cson0))   ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 14';
        'ERREUR' 5 ;
 'FINSI' ;

 ERR1 = ('ABS' (DGYRETN '-' DGYRET)) '*' 
             (ret0 '/' (ro0 '*' cson0 * cson0)) ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 15';
        'ERREUR' 5 ;
 'FINSI' ;


 ERR1 = ('ABS' (DRETRETN '-' DRETRET)) '/' cson0 ;
 'SI' (ERR1 > ERRTOL) ;
        'MESSAGE' 'Problem 16';
        'ERREUR' 5 ;
 'FINSI' ;
 
  
 'FIN' ;