Download lapn_impl_centre.dgibi

Back to the list

* fichier :  lapn_impl_centre.dgibi
***********************************************************
***********************************************************
**** APPROCHE VF "Cell-Centred Formulation" pour la    ****
**** solution des                                      **** 
**** Equations d'Euler pour un gaz parfait             ****
**** OPERATEURS PRIM, PENT, LAPN                       ****
**** Implicit: calcul du jacobien du residu            ****
****                                                   ****
**** Cas gaz monoespece, "calorically perfect"         ****
****                                                   ****
**** Methodes: DIAMANT                                 ****
****                                                   ****
**** A. BECCANTINI DRN/DMT/SEMT/LTMF    AOUT  2000     ****
**** S. GOUNAND    DEN/DM2S/SFME/LTMF   AOUT  2001     ****
***********************************************************
***********************************************************

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

*
*** 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) ;

 DOM5 = DOM5 'COULEUR' 'ROUG ' ;
 DOM6 = DOM6 'COULEUR' 'VERT' ;
 DOM7 = DOM7 'COULEUR' 'JAUN' ;
 
 DOMTOT = DOM1 'ET' DOM2 'ET' DOM3 'ET' DOM4 'ET' DOM5 'ET'
          DOM6 'ET' DOM7 'ET' DOM8 'ET' DOM9 'ELIMINATION' 0.0001 ;
 
*
**** On calcule le JACOBIAN dans le centre de DOM5;
*    qui doit dependre de la valeur en DOM5, DOM6, mais il ne doit
*    pas dependre de DOM7
*

 $DOMTOT = 'MODELISER' DOMTOT 'EULER';

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

 TDOMTOT = 'DOMA'  $DOMTOT 'VF';
 
 TDOM6 = 'DOMA'  $DOM6 'VF';
 TDOM7 = 'DOMA'  $DOM7 'VF';
 TDOM5 = 'DOMA'  $DOM5 'VF';

 MDOM6 = TDOM6 . 'QUAF' ;
 MDOM7 = TDOM7 . 'QUAF' ;
 MDOM5 = TDOM5 . 'QUAF' ;          
 
**** old stuff  $DOMTOT = 'DOMA' DOMTOT ;

 MDOMTOT = TDOMTOT . 'QUAF' ;

 'ELIMINATION'  (MDOMTOT ET MDOM6) 0.0001 ;
 'ELIMINATION'  (MDOMTOT ET MDOM7) 0.0001 ;
 'ELIMINATION'  (MDOMTOT ET MDOM5) 0.0001 ;

 RN0 = 'BRUI' 'BLAN' 'UNIF' 1.11 0.5 ('DOMA' $DOMTOT 'CENTRE') ;
 RN0 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' RN0 ;
 PN0 = 'BRUI' 'BLAN' 'UNIF' 1234.1 800 ('DOMA' $DOMTOT 'CENTRE') ;
 GAMMAN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 1.4 ;

 CSONN = (GAMMAN '*' PN0) '/' RN0 ;
 CSONN = 'KOPS' CSONN '**' 0.5 ;
 UXN0  = 0.2 *  CSONN ;
 UYN0  = 0.3  * CSONN ;
 
 GN0 = 'KCHT' $DOMTOT 'VECT' 'CENTRE' (
   ('NOMC' (RN0 '*' UXN0) 'UX' 'NATU' 'DISCRET') 'ET'
   ('NOMC' (RN0 '*' UYN0) 'UY' 'NATU' 'DISCRET')) ;
 
*
*** retgd
*

 ECIN = 0.5D0 '*' RN0 '*' ((UXN0  '*' UXN0) '+' (UYN0  '*' UYN0)) ;
 RETN0  = 'KCHT' $DOMTOT 'SCAL' 'CENTRE'
        ((PN0  '/' (GAMMAN  '-' 1.0)) '+' ECIN) ;

 CV=PI ;
 TN0 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE'
        ((PN0  '/' (GAMMAN  '-' 1.0)) '/' (RN0 '*' CV)) ;
         
 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 valeurs en $DOM5 'CENTRE'
*

 PCON = ('DOMA' $DOM5 'CENTRE') 'POIN' 1 ;
 ro0 = 'EXTRAIRE' RN0 PCON 'SCAL' ;
 cson0 = 'EXTRAIRE' CSONN PCON 'SCAL' ;
 ret0 = 'EXTRAIRE' RETN0 PCON 'SCAL' ;
 

***********************
**** Les CHPOINTs  ****
***********************

 RN  =  'COPIER' RN0 ;

 GN  =  'COPIER' GN0  ;

 RETN = 'COPIER' RETN0 ;

 LISTINCO = 'MOTS' 'RN' 'RUXN' 'RUYN' 'RETN' ;

 V0 PRES   = 'PRIM' 'PERFMONO'  RN GN RETN GAMMAN ;

 CV     = PI ;
 MU     = '*' PI PI ;
 KAPPA  = '*' ('*' PI PI) PI ;
 T0 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE'
        ((PRES  '/' (GAMMAN  '-' 1.0)) '/' (RN '*' CV)) ;
 GRADT0 MCHAMT = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' T0 ;
 GRADV0 MCHAMV = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' V0  ;
*
 JACO CHPRES DT = 'LAPN' 'VF' 'PROPCOST' 'RESI' 'IMPL'
    $DOMTOT  MU KAPPA CV RN V0 T0 GRADV0  GRADT0 MCHAMV MCHAMT
    LISTINCO  ;

***********************
***** Le residu *******
***********************

 DEBRN0   = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 1) ;
 DEBGNX0  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 2) ;
 DEBGNY0  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 3) ;
 DEBRETN0 = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 4) ;


*****************************************************
*****************************************************
******* TEST1 ***************************************
*****************************************************
*****************************************************
*
* On compare le jacobien et la variation du residu
* en $DOM5 'CENTRE' par rapport à une variation
* infinitésimal en $DOM5 'CENTRE'
*
 
 UNRO = ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 
       4 'RN' 0.0 'RUXN' 0.0 'RUYN' 0.0 'RETN' 0.0
        'NATURE' 'DISCRET') 'ET' 
      ('MANUEL' 'CHPO' ('DOMA' $DOM5 'CENTRE') 4 'RN' 1.0
      'RUXN' 0.0 'RUYN' 0.0 'RETN' 0.0
        'NATURE' 'DISCRET') ;
 UNGX = ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 
       4 'RN' 0.0 'RUXN' 0.0 'RUYN' 0.0 'RETN' 0.0
        'NATURE' 'DISCRET') 'ET'
       ('MANUEL' 'CHPO' ('DOMA' $DOM5 'CENTRE') 4 'RN' 0.0
      'RUXN' 1.0 'RUYN' 0.0 'RETN' 0.0 
        'NATURE' 'DISCRET') ;
 UNGY = ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 
       4 'RN' 0.0 'RUXN' 0.0 'RUYN' 0.0 'RETN' 0.0
        'NATURE' 'DISCRET') 'ET'
          ('MANUEL' 'CHPO' ('DOMA' $DOM5 'CENTRE') 4 'RN' 0.0
      'RUXN' 0.0 'RUYN' 1.0 'RETN' 0.0
        'NATURE' 'DISCRET')  ;
 UNRET = ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 
       4 'RN' 0.0 'RUXN' 0.0 'RUYN' 0.0 'RETN' 0.0
        'NATURE' 'DISCRET') 'ET'
        ('MANUEL' 'CHPO' ('DOMA' $DOM5 'CENTRE') 4 'RN' 0.0
      'RUXN' 0.0 'RUYN' 0.0 'RETN' 1.0
        'NATURE' 'DISCRET') ;

   
 DDRHO = 'KOPS' JACO 'MULT' UNRO ;

*   
*
***** DDRHO contient: dRES_RN '/' dRN ('RN') ;
*                     dRES_GXN '/' dRN ('RUXN') ; 
*                     dRES_GYN '/' dRN ('RUYN') ;
*                     dRES_RETN '/' dRN ('RETN') ;
*    DDGX contient ...
*
*
   
 DDGX = 'KOPS' JACO 'MULT' UNGX ;
 DDGY = 'KOPS' JACO 'MULT' UNGY ;    
 DDRET = 'KOPS' JACO 'MULT' UNRET ;

 DRR  = 'EXTRAIRE' DDRHO PCON  'RN' ;
 DGXR = 'EXTRAIRE' DDRHO PCON 'RUXN' ;
 DGYR = 'EXTRAIRE' DDRHO PCON 'RUYN' ;
 DRETR = 'EXTRAIRE' DDRHO PCON 'RETN' ;

 DRGX   = 'EXTRAIRE' DDGX PCON  'RN' ;
 DGXGX  = 'EXTRAIRE' DDGX PCON 'RUXN' ;
 DGYGX  = 'EXTRAIRE' DDGX PCON 'RUYN' ;
 DRETGX = 'EXTRAIRE' DDGX PCON 'RETN' ;

 DRGY   = 'EXTRAIRE' DDGY PCON  'RN' ;
 DGXGY = 'EXTRAIRE' DDGY PCON 'RUXN' ;
 DGYGY = 'EXTRAIRE' DDGY PCON 'RUYN' ;
 DRETGY = 'EXTRAIRE' DDGY PCON 'RETN' ;

 DRRET   = 'EXTRAIRE' DDRET PCON  'RN' ;
 DGXRET = 'EXTRAIRE' DDRET PCON 'RUXN' ;
 DGYRET = 'EXTRAIRE' DDRET PCON 'RUYN' ;
 DRETRET = 'EXTRAIRE' DDRET PCON 'RETN' ;


*********************************************************************
*********************************************************************
***** On calcule les residues pour ro1 = ro0 * (1.'+' DELTA)  *******
*********************************************************************
*********************************************************************

 DELTA = 1.0D-6 ;
 ro1 = ro0 '*' (1 '+' DELTA) ;
 RNDOM5 = 'KCHT' $DOM5 'SCAL' 'CENTRE' ro1 ;
    
 RN  = 'KCHT' $DOMTOT   'SCAL' 'CENTRE' RN0 RNDOM5 ;

 GN  = 'COPIER' GN0 ;

 RETN = 'COPIER' RETN0 ;

 V0 PRES   = 'PRIM' 'PERFMONO'  RN GN RETN GAMMAN ;
 CV     = PI ;
 MU     = '*' PI PI ;
 KAPPA  = '*' ('*' PI PI) PI ;
 T0 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE'
        ((PRES  '/' (GAMMAN  '-' 1.0)) '/' (RN '*' CV)) ;
 GRADT0 MCHAMT = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' T0 ;
 GRADV0 MCHAMV = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' V0  ;
*
*
*
 LISTINCO = 'MOTS' 'RN' 'RUXN' 'RUYN' 'RETN' ;
 IJACEL CHPRES DT = 'LAPN' 'VF' 'PROPCOST' 'RESI' 'IMPL'
    $DOMTOT  MU KAPPA CV RN V0 T0 GRADV0  GRADT0 MCHAMV MCHAMT
    LISTINCO  ;

***********************
***** Le residu *******
***********************

 DEBRN1   = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 1) ;
 DEBGNX1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 2) ;
 DEBGNY1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 3) ;
 DEBRETN1 = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 4) ;

*
**** On calcule le jacobien numeriquement
*
   
 DRRN = (DEBRN1 '-' DEBRN0) '/' (DELTA '*' ro0) ;
 DGXRN = (DEBGNX1 '-' DEBGNX0) '/' (DELTA '*' ro0) ;
 DGYRN = (DEBGNY1 '-' DEBGNY0) '/' (DELTA '*' ro0) ;
 DRETRN = (DEBRETN1 '-' DEBRETN0) '/' (DELTA '*' ro0) ;

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

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

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

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

   
***************************************************************************
***************************************************************************
***** On calcule les residues pour gnx1 = gnx0 '+' (DELTA ro0 cson0)  *****
***************************************************************************
***************************************************************************


 DELTA = 1.0D-6 ;
 gnx0 = 'EXTRAIRE' GN0 PCON 'UX' ;
 gny0 = 'EXTRAIRE' GN0 PCON 'UY' ;
 gnx1 = gnx0 '+' (ro0 '*' cson0 '*'  DELTA) ; 
 GNDOM5 = 'KCHT' $DOM5 'VECT' 'CENTRE' (gnx1 gny0) ;

    
 RN  = 'COPIER' RN0 ;

 GN  = 'KCHT' $DOMTOT 'VECT' 'CENTRE' GN0 GNDOM5 ;

 RETN = 'COPIER' RETN0 ;

 V0 PRES   = 'PRIM' 'PERFMONO'  RN GN RETN GAMMAN ;
 CV     = PI ;
 MU     = '*' PI PI ;
 KAPPA  = '*' ('*' PI PI) PI ;
 T0 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE'
        ((PRES  '/' (GAMMAN  '-' 1.0)) '/' (RN '*' CV)) ;
 GRADT0 MCHAMT = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' T0 ;
 GRADV0 MCHAMV = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' V0  ;
*
*
*
 LISTINCO = 'MOTS' 'RN' 'RUXN' 'RUYN' 'RETN' ;
 IJACEL CHPRES DT = 'LAPN' 'VF' 'PROPCOST' 'RESI' 'IMPL'
    $DOMTOT  MU KAPPA CV RN V0 T0 GRADV0  GRADT0 MCHAMV MCHAMT
    LISTINCO  ;

***********************
***** Le residu *******
***********************

 DEBRN1   = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 1) ;
 DEBGNX1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 2) ;
 DEBGNY1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 3) ;
 DEBRETN1 = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 4) ;

*
**** On calcule le jacobien numeriquement
*
   
 DRGXN = (DEBRN1 '-' DEBRN0) '/' (DELTA '*' ro0 '*' cson0) ;
 DGXGXN = (DEBGNX1 '-' DEBGNX0) '/' (DELTA '*' ro0 '*' cson0) ;
 DGYGXN = (DEBGNY1 '-' DEBGNY0) '/'  (DELTA '*' ro0 '*' cson0) ;
 DRETGXN = (DEBRETN1 '-' DEBRETN0) '/'  (DELTA '*' ro0 '*' cson0) ;

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

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

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

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


***************************************************************************
***************************************************************************
***** On calcule les residues pour gny1 = gny0 '+' (DELTA ro0 cson0)  ******
***************************************************************************
***************************************************************************


 DELTA = 1.0D-6 ;
    
 gnx0 = 'EXTRAIRE' GN0 PCON 'UX' ;
 gny0 = 'EXTRAIRE' GN0 PCON 'UY' ;
 gny1 = gny0 '+' (ro0 '*' cson0 '*'  DELTA) ; 
 GNDOM5 = 'KCHT' $DOM5 'VECT' 'CENTRE' (gnx0 gny1) ;

    
 RN  = 'COPIER' RN0 ;

 GN  = 'KCHT' $DOMTOT 'VECT' 'CENTRE' GN0 GNDOM5 ;

 RETN = 'COPIER' RETN0 ;

 V0 PRES   = 'PRIM' 'PERFMONO'  RN GN RETN GAMMAN ;
 CV     = PI ;
 MU     = '*' PI PI ;
 KAPPA  = '*' ('*' PI PI) PI ;
 T0 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE'
        ((PRES  '/' (GAMMAN  '-' 1.0)) '/' (RN '*' CV)) ;
 GRADT0 MCHAMT = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' T0 ;
 GRADV0 MCHAMV = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' V0  ;
*
*
*
 LISTINCO = 'MOTS' 'RN' 'RUXN' 'RUYN' 'RETN' ;
 IJACEL CHPRES DT = 'LAPN' 'VF' 'PROPCOST' 'RESI' 'IMPL'
    $DOMTOT  MU KAPPA CV RN V0 T0 GRADV0  GRADT0 MCHAMV MCHAMT
    LISTINCO  ;

***********************
***** Le residu *******
***********************

 DEBRN1   = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 1) ;
 DEBGNX1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 2) ;
 DEBGNY1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 3) ;
 DEBRETN1 = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 4) ;

*
**** On calcule le jacobien numeriquement
*
   
 DRGYN = (DEBRN1 '-' DEBRN0) '/' (DELTA '*' ro0 '*' cson0) ;
 DGXGYN = (DEBGNX1 '-' DEBGNX0) '/' (DELTA '*' ro0 '*' cson0) ;
 DGYGYN = (DEBGNY1 '-' DEBGNY0) '/'  (DELTA '*' ro0 '*' cson0) ;
 DRETGYN = (DEBRETN1 '-' DEBRETN0) '/'  (DELTA '*' ro0 '*' cson0) ;

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

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

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


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


***************************************************************************
***************************************************************************
***** On calcule les residues pour ret1 = ret0 '*' (1. '+' DELTA)    ******
***************************************************************************
***************************************************************************


 DELTA = 1.0D-6 ;
 ret1 = ret0 '*' (1. '+' DELTA) ;
 RETNDOM5 = 'KCHT' $DOM5 'SCAL' 'CENTRE' ret1 ;

 RN  = 'COPIER' RN0 ;

 GN  = 'COPIER' GN0 ;

 RETN = 'KCHT' $DOMTOT   'SCAL' 'CENTRE' RETN0 RETNDOM5 ;

 V0 PRES   = 'PRIM' 'PERFMONO'  RN GN RETN GAMMAN ;
 CV     = PI ;
 MU     = '*' PI PI ;
 KAPPA  = '*' ('*' PI PI) PI ;
 T0 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE'
        ((PRES  '/' (GAMMAN  '-' 1.0)) '/' (RN '*' CV)) ;
 GRADT0 MCHAMT = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' T0 ;
 GRADV0 MCHAMV = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' V0  ;
*
*
*
 LISTINCO = 'MOTS' 'RN' 'RUXN' 'RUYN' 'RETN' ;
 IJACEL CHPRES DT = 'LAPN' 'VF' 'PROPCOST' 'RESI' 'IMPL'
    $DOMTOT  MU KAPPA CV RN V0 T0 GRADV0  GRADT0 MCHAMV MCHAMT
    LISTINCO  ;

***********************
***** Le residu *******
***********************

 DEBRN1   = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 1) ;
 DEBGNX1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 2) ;
 DEBGNY1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 3) ;
 DEBRETN1 = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 4) ;

*
**** On calcule le jacobien numeriquement
*
   
 DRRETN = (DEBRN1 '-' DEBRN0) '/' (DELTA '*' ret0) ;
 DGXRETN = (DEBGNX1 '-' DEBGNX0) '/' (DELTA '*' ret0) ;
 DGYRETN = (DEBGNY1 '-' DEBGNY0) '/'  (DELTA '*' ret0) ;
 DRETRETN = (DEBRETN1 '-' DEBRETN0) '/'  (DELTA '*' ret0) ;

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

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

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


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

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

 PCELL = ('DOMA' $DOM6 'CENTRE') 'POIN' 1 ;
 ro0 = 'EXTRAIRE' RN0 PCELL 'SCAL' ;
 cson0 = 'EXTRAIRE' CSONN PCELL 'SCAL' ;
 ret0 = 'EXTRAIRE' RETN0 PCELL 'SCAL' ;

 
 UNRO = ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 
       4 'RN' 0.0 'RUXN' 0.0 'RUYN' 0.0 'RETN' 0.0
        'NATURE' 'DISCRET') 'ET' 
      ('MANUEL' 'CHPO' ('DOMA' $DOM6 'CENTRE') 4 'RN' 1.0
      'RUXN' 0.0 'RUYN' 0.0 'RETN' 0.0
        'NATURE' 'DISCRET') ;
 UNGX = ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 
       4 'RN' 0.0 'RUXN' 0.0 'RUYN' 0.0 'RETN' 0.0
        'NATURE' 'DISCRET') 'ET'
       ('MANUEL' 'CHPO' ('DOMA' $DOM6 'CENTRE') 4 'RN' 0.0
      'RUXN' 1.0 'RUYN' 0.0 'RETN' 0.0 
        'NATURE' 'DISCRET') ;
 UNGY = ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 
       4 'RN' 0.0 'RUXN' 0.0 'RUYN' 0.0 'RETN' 0.0
        'NATURE' 'DISCRET') 'ET'
          ('MANUEL' 'CHPO' ('DOMA' $DOM6 'CENTRE') 4 'RN' 0.0
      'RUXN' 0.0 'RUYN' 1.0 'RETN' 0.0
        'NATURE' 'DISCRET')  ;
 UNRET = ('MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 
       4 'RN' 0.0 'RUXN' 0.0 'RUYN' 0.0 'RETN' 0.0
        'NATURE' 'DISCRET') 'ET'
        ('MANUEL' 'CHPO' ('DOMA' $DOM6 'CENTRE') 4 'RN' 0.0
      'RUXN' 0.0 'RUYN' 0.0 'RETN' 1.0
        'NATURE' 'DISCRET') ;

   
 DDRHO = 'KOPS' JACO 'MULT' UNRO ;

*   
*
***** DDRHO contient: dRES_RN '/' dRN ('RN') ;
*                     dRES_GXN '/' dRN ('RUXN') ; 
*                     dRES_GYN '/' dRN ('RUYN') ;
*                     dRES_RETN '/' dRN ('RETN') ;
*    DDGX contient ...
*
*
   
 DDGX = 'KOPS' JACO 'MULT' UNGX ;
 DDGY = 'KOPS' JACO 'MULT' UNGY ;    
 DDRET = 'KOPS' JACO 'MULT' UNRET ;

 DRR  = 'EXTRAIRE' DDRHO PCON  'RN' ;
 DGXR = 'EXTRAIRE' DDRHO PCON 'RUXN' ;
 DGYR = 'EXTRAIRE' DDRHO PCON 'RUYN' ;
 DRETR = 'EXTRAIRE' DDRHO PCON 'RETN' ;

 DRGX   = 'EXTRAIRE' DDGX PCON  'RN' ;
 DGXGX  = 'EXTRAIRE' DDGX PCON 'RUXN' ;
 DGYGX  = 'EXTRAIRE' DDGX PCON 'RUYN' ;
 DRETGX = 'EXTRAIRE' DDGX PCON 'RETN' ;

 DRGY   = 'EXTRAIRE' DDGY PCON  'RN' ;
 DGXGY = 'EXTRAIRE' DDGY PCON 'RUXN' ;
 DGYGY = 'EXTRAIRE' DDGY PCON 'RUYN' ;
 DRETGY = 'EXTRAIRE' DDGY PCON 'RETN' ;

 DRRET   = 'EXTRAIRE' DDRET PCON  'RN' ;
 DGXRET = 'EXTRAIRE' DDRET PCON 'RUXN' ;
 DGYRET = 'EXTRAIRE' DDRET PCON 'RUYN' ;
 DRETRET = 'EXTRAIRE' DDRET PCON 'RETN' ;



*********************************************************************
*********************************************************************
***** On calcule les residues pour ro1 = ro0 * (1.'+' DELTA)  *******
*********************************************************************
*********************************************************************

 DELTA = 1.0D-6 ;
 ro1 = ro0 '*' (1 '+' DELTA) ;
 RNDOM6 = 'KCHT' $DOM6 'SCAL' 'CENTRE' ro1 ;
    
 RN  = 'KCHT' $DOMTOT   'SCAL' 'CENTRE' RN0 RNDOM6 ;

 GN  = 'COPIER' GN0 ;

 RETN = 'COPIER' RETN0 ;

 V0 PRES   = 'PRIM' 'PERFMONO'  RN GN RETN GAMMAN ;
 CV     = PI ;
 MU     = '*' PI PI ;
 KAPPA  = '*' ('*' PI PI) PI ;
 T0 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE'
        ((PRES  '/' (GAMMAN  '-' 1.0)) '/' (RN '*' CV)) ;
 GRADT0 MCHAMT = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' T0 ;
 GRADV0 MCHAMV = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' V0  ;
*
*
*
 LISTINCO = 'MOTS' 'RN' 'RUXN' 'RUYN' 'RETN' ;
 IJACEL CHPRES DT = 'LAPN' 'VF' 'PROPCOST' 'RESI' 'IMPL'
    $DOMTOT  MU KAPPA CV RN V0 T0 GRADV0  GRADT0 MCHAMV MCHAMT
    LISTINCO  ;

***********************
***** Le residu *******
***********************

 DEBRN1   = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 1) ;
 DEBGNX1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 2) ;
 DEBGNY1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 3) ;
 DEBRETN1 = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 4) ;


*
**** On calcule le jacobien numeriquement
*
   
 DRRN = (DEBRN1 '-' DEBRN0) '/' (DELTA '*' ro0) ;
 DGXRN = (DEBGNX1 '-' DEBGNX0) '/' (DELTA '*' ro0) ;
 DGYRN = (DEBGNY1 '-' DEBGNY0) '/' (DELTA '*' ro0) ;
 DRETRN = (DEBRETN1 '-' DEBRETN0) '/' (DELTA '*' ro0) ;

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

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

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

 ERR1 = ('ABS' (DRETRN '-' DRETR)) '*'
        (ro0 '/' (ret0 '*' cson0)) ; 
 'MESSAGE' ('CHAINE' 'DRETRN=' DRETRN ' DRETR=' DRETR) ;
 'SI' (ERR1 > ERRTOL) ;
      'MESSAGE' 'Problem 4';
      'ERREUR' 5 ;
 'FINSI' ;
  
***************************************************************************
***************************************************************************
***** On calcule les residues pour gnx1 = gnx0 '+' (DELTA ro0 cson0)  ******
***************************************************************************
***************************************************************************

 DELTA = 1.0D-6 ;
 gnx0 = 'EXTRAIRE' GN0 PCELL 'UX' ;
 gny0 = 'EXTRAIRE' GN0 PCELL 'UY' ;
 gnx1 = gnx0 '+' (ro0 '*' cson0 '*'  DELTA) ; 
 GNDOM6 = 'KCHT' $DOM6 'VECT' 'CENTRE' (gnx1 gny0) ;

    
 RN  = 'COPIER' RN0 ;

 GN  = 'KCHT' $DOMTOT 'VECT' 'CENTRE' GN0 GNDOM6 ;

 RETN = 'COPIER' RETN0 ;

 V0 PRES   = 'PRIM' 'PERFMONO'  RN GN RETN GAMMAN ;
 CV     = PI ;
 MU     = '*' PI PI ;
 KAPPA  = '*' ('*' PI PI) PI ;
 T0 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE'
        ((PRES  '/' (GAMMAN  '-' 1.0)) '/' (RN '*' CV)) ;
 GRADT0 MCHAMT = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' T0 ;
 GRADV0 MCHAMV = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' V0  ;
*
*
*
 LISTINCO = 'MOTS' 'RN' 'RUXN' 'RUYN' 'RETN' ;
 IJACEL CHPRES DT = 'LAPN' 'VF' 'PROPCOST' 'RESI' 'IMPL'
    $DOMTOT  MU KAPPA CV RN V0 T0 GRADV0  GRADT0 MCHAMV MCHAMT
    LISTINCO  ;

***********************
***** Le residu *******
***********************

 DEBRN1   = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 1) ;
 DEBGNX1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 2) ;
 DEBGNY1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 3) ;
 DEBRETN1 = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 4) ;

*
**** On calcule le jacobien numeriquement
*
   
 DRGXN = (DEBRN1 '-' DEBRN0) '/' (DELTA '*' ro0 '*' cson0) ;
 DGXGXN = (DEBGNX1 '-' DEBGNX0) '/' (DELTA '*' ro0 '*' cson0) ;
 DGYGXN = (DEBGNY1 '-' DEBGNY0) '/'  (DELTA '*' ro0 '*' cson0) ;
 DRETGXN = (DEBRETN1 '-' DEBRETN0) '/'  (DELTA '*' ro0 '*' cson0) ;

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

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

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


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


***************************************************************************
***************************************************************************
***** On calcule les residues pour gny1 = gny0 '+' (DELTA ro0 cson0)  ******
***************************************************************************
***************************************************************************


 DELTA = 1.0D-6 ;
    
 gnx0 = 'EXTRAIRE' GN0 PCELL 'UX' ;
 gny0 = 'EXTRAIRE' GN0 PCELL 'UY' ;
 gny1 = gny0 '+' (ro0 '*' cson0 '*'  DELTA) ; 
 GNDOM6 = 'KCHT' $DOM6 'VECT' 'CENTRE' (gnx0 gny1) ;

    
 RN  = 'COPIER' RN0 ;

 GN  = 'KCHT' $DOMTOT 'VECT' 'CENTRE' GN0 GNDOM6 ;

 RETN = 'COPIER' RETN0 ;

 V0 PRES   = 'PRIM' 'PERFMONO'  RN GN RETN GAMMAN ;
 CV     = PI ;
 MU     = '*' PI PI ;
 KAPPA  = '*' ('*' PI PI) PI ;
 T0 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE'
        ((PRES  '/' (GAMMAN  '-' 1.0)) '/' (RN '*' CV)) ;
 GRADT0 MCHAMT = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' T0 ;
 GRADV0 MCHAMV = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' V0  ;
*
*
*
 LISTINCO = 'MOTS' 'RN' 'RUXN' 'RUYN' 'RETN' ;
 IJACEL CHPRES DT = 'LAPN' 'VF' 'PROPCOST' 'RESI' 'IMPL'
    $DOMTOT  MU KAPPA CV RN V0 T0 GRADV0  GRADT0 MCHAMV MCHAMT
    LISTINCO  ;

***********************
***** Le residu *******
***********************

 DEBRN1   = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 1) ;
 DEBGNX1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 2) ;
 DEBGNY1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 3) ;
 DEBRETN1 = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 4) ;


*
**** On calcule le jacobien numeriquement
*
   
 DRGYN = (DEBRN1 '-' DEBRN0) '/' (DELTA '*' ro0 '*' cson0) ;
 DGXGYN = (DEBGNX1 '-' DEBGNX0) '/' (DELTA '*' ro0 '*' cson0) ;
 DGYGYN = (DEBGNY1 '-' DEBGNY0) '/'  (DELTA '*' ro0 '*' cson0) ;
 DRETGYN = (DEBRETN1 '-' DEBRETN0) '/'  (DELTA '*' ro0 '*' cson0) ;

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

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

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

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


***************************************************************************
***************************************************************************
***** On calcule les residues pour ret1 = ret0 '*' (1. '+' DELTA)    ******
***************************************************************************
***************************************************************************


 DELTA = 1.0D-6 ;
 ret1 = ret0 '*' (1. '+' DELTA) ;
 RETNDOM6 = 'KCHT' $DOM6 'SCAL' 'CENTRE' ret1 ;

 RN  = 'COPIER' RN0 ;

 GN  = 'COPIER' GN0 ;

 RETN = 'KCHT' $DOMTOT   'SCAL' 'CENTRE' RETN0 RETNDOM6 ;

 V0 PRES   = 'PRIM' 'PERFMONO'  RN GN RETN GAMMAN ;
 CV     = PI ;
 MU     = '*' PI PI ;
 KAPPA  = '*' ('*' PI PI) PI ;
 T0 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE'
        ((PRES  '/' (GAMMAN  '-' 1.0)) '/' (RN '*' CV)) ;
 GRADT0 MCHAMT = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' T0 ;
 GRADV0 MCHAMV = 'PENT'  $DOMTOT  'FACE' 'DIAMANT' V0  ;
*
*
*
 LISTINCO = 'MOTS' 'RN' 'RUXN' 'RUYN' 'RETN' ;
 IJACEL CHPRES DT = 'LAPN' 'VF' 'PROPCOST' 'RESI' 'IMPL'
    $DOMTOT  MU KAPPA CV RN V0 T0 GRADV0  GRADT0 MCHAMV MCHAMT
    LISTINCO  ;

***********************
***** Le residu *******
***********************

 DEBRN1   = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 1) ;
 DEBGNX1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 2) ;
 DEBGNY1  = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 3) ;
 DEBRETN1 = 'EXTRAIRE' CHPRES PCON ('EXTRAIRE' LISTINCO 4) ;

*
**** On calcule le jacobien numeriquement
*
   
 DRRETN = (DEBRN1 '-' DEBRN0) '/' (DELTA '*' ret0) ;
 DGXRETN = (DEBGNX1 '-' DEBGNX0) '/' (DELTA '*' ret0) ;
 DGYRETN = (DEBGNY1 '-' DEBGNY0) '/'  (DELTA '*' ret0) ;
 DRETRETN = (DEBRETN1 '-' DEBRETN0) '/'  (DELTA '*' ret0) ;

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


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

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


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

 
 'FIN' ;

 

 

 

 

 

© Cast3M 2003 - All rights reserved.
Disclaimer