Télécharger la plaquette Cast3M


fr en

Télécharger konv_implbm.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier :  konv_implbm.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. ***********************************************************
  5. ***********************************************************
  6. **** APPROCHE VF "Cell-Centred Formulation" pour la    ****
  7. **** solution des                                      **** 
  8. **** Equations d'Euler pour un gaz parfait             ****
  9. **** OPERATEURS PRIM, PRET, KONV                       ****
  10. **** Implicit: calcul du jacobien du residu            ****
  11. ****                                                   ****
  12. **** Cas gaz monoespece, "calorically perfect"         ****
  13. ****                                                   ****
  14. **** Methodes: AUSMLM                                  ****
  15. ****                                                   ****
  16. **** S.V. KUDRIAKOV  DM2S/SFME/LTMF    JUNE 2002       ****
  17. ***********************************************************
  18. ***********************************************************
  19.  
  20.  'OPTION'  'DIME' 2 ;
  21.  'OPTION'  'ELEM' QUA4 ;
  22.  'OPTION'  'ECHO'  0 ;
  23.  'OPTION'  'TRAC' 'X' ;
  24.  
  25. *
  26. *** GRAPH
  27. *
  28.  
  29.  GRAPH = FAUX ;
  30. * GRAPH = VRAI ;
  31.  
  32.  ERRTOL = 1.0D-4 ;
  33.  
  34. ***************************
  35. ***** DOMAINE SPATIAL  ****
  36. ***************************
  37.  
  38.  
  39.  A0 = 0.0D0 0.0D0;
  40.  A1 = 1.0D0 0.0D0;
  41.  A2 = 2.0D0 0.0D0;
  42.  A3 = 3.0D0 0.0D0;
  43.  
  44.  A0A1 = A0 'DROIT' 1 A1;
  45.  A1A2 = A1 'DROIT' 1 A2;
  46.  A2A3 = A2 'DROIT'  1 A3;
  47.  
  48.  
  49.  DOM1 = 'TRANSLATION' A0A1 1 (0.0 1.0) ;
  50.  DOM2 = DOM1 'PLUS' (0.0 1.0) ;
  51.  DOM3 = DOM2 'PLUS' (0.0 1.0) ;
  52.  DOM4 = 'TRANSLATION' A1A2 1 (0.0 1.0) ;
  53.  DOM5 = DOM4 'PLUS' (0.0 1.0) ;
  54.  DOM6 = DOM5 'PLUS' (0.0 1.0) ;
  55.  DOM7 = 'TRANSLATION' A2A3 1 (0.0 1.0) ; 
  56.  DOM8 = DOM7 'PLUS' (0.0 1.0) ;
  57.  DOM9 = DOM8 'PLUS' (0.0 1.0) ;
  58.  
  59.  
  60.  DOMTOT = DOM1 'ET' DOM2 'ET' DOM3 'ET' DOM4 'ET' DOM5 'ET'
  61.           DOM6 'ET' DOM7 'ET' DOM8 'ET' DOM9 'ELIMINATION' 0.0001 ;
  62.  
  63.  $DOMTOT = 'MODELISER' DOMTOT 'EULER';
  64.  
  65.  $DOM6 = 'MODELISER'  DOM6 'EULER';
  66.  $DOM7 = 'MODELISER'  DOM7 'EULER';
  67.  $DOM9 = 'MODELISER'  DOM9 'EULER';
  68.  
  69.  TDOMTOT = 'DOMA'  $DOMTOT 'VF';
  70.  
  71.  TDOM6 = 'DOMA'  $DOM6 'VF';
  72.  TDOM7 = 'DOMA'  $DOM7 'VF';
  73.  TDOM9 = 'DOMA'  $DOM9 'VF';
  74.  
  75.  MDOM6 = TDOM6 . 'QUAF' ;
  76.  MDOM7 = TDOM7 . 'QUAF' ;
  77.  MDOM9 = TDOM9 . 'QUAF' ;          
  78.  
  79. **** old stuff  $DOMTOT = 'DOMA' DOMTOT ;
  80.  
  81.  MDOMTOT = TDOMTOT . 'QUAF' ;
  82.  
  83.  'ELIMINATION'  (MDOMTOT ET MDOM6) 0.0001 ;
  84.  'ELIMINATION'  (MDOMTOT ET MDOM7) 0.0001 ;
  85.  'ELIMINATION'  (MDOMTOT ET MDOM9) 0.0001 ;
  86.  
  87. ***************************************************
  88. ***** Densité, pression, vitesse, gamma ***********
  89. ***************************************************
  90.  
  91.  RN0 = 'BRUI' 'BLAN' 'UNIF' 1.11 0.5 ('DOMA '$DOMTOT 'CENTRE') ;
  92.  PN0 = 'BRUI' 'BLAN' 'UNIF' 1234.1 800 ('DOMA '$DOMTOT 'CENTRE') ;
  93.  GAMMAN = 'MANUEL' 'CHPO' ('DOMA '$DOMTOT 'CENTRE') 1 'SCAL' 1.4
  94.    'NATU' 'DISCRET' ;
  95.  CSONN = (GAMMAN '*' PN0) '/' RN0 ;
  96.  CSONN = 'KOPS' CSONN '**' 0.5 ;
  97.  UXN0  = 0.05 *  CSONN ;
  98.  UYN0  = 0.01  * CSONN ;
  99. ***************************************************
  100.  UINF = 'MANUEL' 'CHPO' ('DOMA '$DOMTOT 'CENTRE') 1 'SCAL' 0.01
  101.    'NATU' 'DISCRET' ;
  102.  UPRI = 'MANUEL' 'CHPO' ('DOMA '$DOMTOT 'CENTRE') 1 'SCAL' 0.0001
  103.    'NATU' 'DISCRET' ;   
  104. *************************************************** 
  105.  GN0 = ('NOMC' (RN0 '*' UXN0) 'UX' 'NATU' 'DISCRET') 'ET'
  106.        ('NOMC' (RN0 '*' UYN0) 'UY' 'NATU' 'DISCRET') ;
  107.  
  108.  ECIN = 0.5D0 '*' RN0 '*' ((UXN0  '*' UXN0) '+' (UYN0  '*' UYN0)) ;
  109.  RETN0  = 'NOMC' 'SCAL' ((PN0  '/' (GAMMAN  '-' 1.0)) '+' ECIN)
  110.    'NATU' 'DISCRET' ;
  111.  
  112.  VIT PRES = PRIM 'PERFMONO' RN0 GN0 RETN0 GAMMAN ;
  113.  
  114.  ERR1 = 'MAXIMUM' (PRES '-' PN0) 'ABS' ;
  115.  ERR2 = 'MAXIMUM' (GN0 '-' (RN0 '*' VIT ('MOTS' 'SCAL' 'SCAL')
  116.          ('MOTS' 'UX' 'UY') ('MOTS' 'UX' 'UY'))) 'ABS' ;
  117.  
  118.  'SI' (('MAXIMUM' ('PROG' ERR1 ERR2 ) 'ABS' ) > 1.0D-6) ;
  119.       'MESSAGE' 'Problem 0' ;
  120.       'ERREUR' 5 ;
  121.  'FINSI' ;
  122.  
  123.  
  124.  'SI' GRAPH;
  125.       'TRACER' (('DOMA' $DOMTOT 'MAILLAGE' ) 'ET' 
  126.            ('DOMA' $DOMTOT 'CENTRE')) 'TITRE' 'Domaine et centre' ;
  127.  'FINSI' ;
  128.  
  129. *
  130. **** Les variables conservative
  131. *
  132. *    RN0 (densité)
  133. *    GN0 (quantité de mouvement)
  134. *    RETN0 (énergie totale par unité de volume)
  135. *
  136. *    sont definiés
  137. *
  138.  
  139. ****************************************************
  140. * Creating a mesh WITHOUT mesh 
  141. **************************************************** 
  142.   P1 = 0.5 0.5 ;
  143.   ELP1 = 'MANU' 'POI1' P1 ;
  144.   MELLIM = 'DIFF' ELP1 ELP1 ;
  145. ****************************************************
  146. ****************************************************
  147. ******* Calcul du jacobien et du residu   **********
  148. ****************************************************
  149. ****************************************************
  150. *
  151. * JACO est le jacobien
  152. *
  153. * DEBRN0 le residu concernant la densité
  154. * DEBGNX0 le residu concernant la quantité de mouvement (axe x)
  155. * DEBGNY0 le residu concernant la quantité de mouvement (axe y)
  156. * DEBRETN0 le residu concernant l'enrgie totale par unité de volume
  157. *
  158.  
  159. * Noms des variables
  160.  
  161.  NOMDEN = 'RN  ' ;
  162.  NOMMOX = 'RUXN' ;
  163.  NOMMOY = 'RUYN' ;
  164.  NOMRET = 'RETN' ;
  165.  
  166. * Metode
  167.  
  168.  METO = 'AUSMPLM' ;
  169.  
  170.  
  171.  LISTINCO = 'MOTS'  NOMDEN  NOMMOX  NOMMOY  NOMRET ;
  172.  
  173.  VITESSE PRES   = 'PRIM' 'PERFMONO'  RN0 GN0 RETN0 GAMMAN ;
  174.  
  175.  ROF VITF PF GAMF   =  'PRET' 'PERFMONO'  1 1
  176.                            $DOMTOT RN0 VITESSE PRES GAMMAN ;
  177.  
  178.  CHPRES0 DT = 'KONV' 'VF' 'PERFMONO' 'RESI' METO
  179.           $DOMTOT ROF VITF PF GAMF LISTINCO
  180.           MELLIM UINF UPRI ;
  181.  
  182.  JACO = 'KONV' 'VF' 'PERFMONO' 'JACOCONS'
  183.          $DOMTOT LISTINCO MELLIM 
  184.          METO  RN0 VITESSE PRES GAMMAN
  185.           UINF UPRI ;
  186.  
  187.  
  188.  DEBRN0   = 'EXCO' NOMDEN CHPRES0 ;
  189.  DEBGNX0  = 'EXCO' NOMMOX CHPRES0 ;
  190.  DEBGNY0  = 'EXCO' NOMMOY CHPRES0 ;
  191.  DEBRETN0 = 'EXCO' NOMRET CHPRES0 ;
  192.  
  193. *****************************************************
  194. *****************************************************
  195. ******** PROCEDURES *********************************
  196. *****************************************************
  197. *****************************************************
  198.  
  199. *
  200. * Derivé partielle du residu en un point par rapport
  201. * aux variable en un autre point
  202. *
  203.  
  204.  'DEBPROC' JACOVA JACO*'MATRIK' $MODE*'MMODEL' LISTINCO*'LISTMOTS'
  205.     PPRIM*'POINT' PDUAL*'POINT' MOTPRI*'MOT' MOTDUA*'MOT';
  206.  
  207. *   PPRIM = point ou est localisé la variable primale
  208. *   PDUAL = point ou est localisé la variable duale
  209. *   MOTPRI = nom de la composante concernante la variable primale
  210. *   MOTDUA = nom de la composante concernante la variable duale
  211.  
  212.     ELT1 = 'MANUEL' 'POI1' PPRIM ;
  213.     NDIM = 'DIME' LISTINCO ;
  214.     CHPUN = 'MANUEL' 'CHPO' ELT1 1 MOTPRI 1.0 'NATURE' 'DISCRET' ;
  215.     'REPETER' BL1 NDIM ;
  216.        MOTCEL = 'EXTRAIRE' LISTINCO &BL1 ;
  217.        CHPUN = CHPUN 'ET' ('MANUEL' 'CHPO' ('DOMA' $MODE 'CENTRE')
  218.           1 MOTCEL 0.0 'NATURE' 'DISCRET') ;
  219.     'FIN' BL1 ;
  220.     D_DMOT = 'KOPS' JACO 'MULT' CHPUN ;
  221.     SCAL = 'EXTRAIRE' D_DMOT PDUAL MOTDUA ;
  222.  
  223.  'FINPROC' SCAL  ;
  224.  
  225.  
  226.  
  227.  'DEBPROC' JACNUM $MODE*'MMODEL' NOMMET*'MOT' RN*'CHPOINT' GN*'CHPOINT'
  228.     RETN*'CHPOINT' GAMN*'CHPOINT'
  229.     LISTINCO*'LISTMOTS' MELLIM*'MAILLAGE' UINF*'CHPOINT' UPRI*'CHPOINT'
  230.     PPRIM*'POINT' PDUAL*'POINT'
  231.     MOTPRI*'MOT' MOTDUA*'MOT' EPSILON*'FLOTTANT';
  232.  
  233. *   PPRIM = point ou est localisé la variable primale
  234. *   PDUAL = point ou est localisé la variable duale
  235. *   MOTPRI = nom de la composante concernante la variable primale
  236. *   MOTDUA = nom de la composante concernante la variable duale
  237.  
  238. *   Le valeur dans l'état non-perturbé en PDUAL ;
  239.  
  240.     RNCEL   = 'COPIER'  RN ;
  241.     GNCEL   = 'COPIER'  GN ;
  242.     RETNCEL = 'COPIER'  RETN ;
  243.  
  244.     VITESSE PRES   = 'PRIM' 'PERFMONO'  RNCEL GNCEL RETNCEL GAMN ;
  245.  
  246.     ROF VITF PF GAMF   =  'PRET' 'PERFMONO'  1 1
  247.                            $MODE RNCEL VITESSE PRES GAMN ;
  248.  
  249.     CHPRES0 DT = 'KONV' 'VF' 'PERFMONO' 'RESI' NOMMET
  250.           $MODE ROF VITF PF GAMF LISTINCO MELLIM UINF UPRI ;
  251.  
  252. *    'LISTE' DT ;      
  253.  
  254.     VAL0 = 'EXTRAIRE' CHPRES0 PDUAL MOTDUA ;
  255.  
  256. *   EPSILON = perturbation
  257.  
  258. *   Adimensionalisation
  259.  
  260.     dens0 = 'EXTRAIRE' RN PPRIM 'SCAL' ;
  261.  
  262.     VN PN = 'PRIM' 'PERFMONO' RN GN RETN GAMN ;
  263.     CN2 = GAMN '*' (PN '/' RN) ;
  264.     cson0 = ('EXTRAIRE' CN2 PPRIM 'SCAL') '**' 0.5 ;
  265.  
  266.     ret0 = ('EXTRAIRE' RETN PPRIM 'SCAL') '**' 0.5 ;
  267.  
  268. *   On etabli la variable à perturber
  269.  
  270.     NDIM = 'DIME' LISTINCO ;
  271.     'REPETER' BL1 NDIM ;
  272.        MOTCEL = 'EXTRAIRE' LISTINCO &BL1 ;
  273.        'SI' ('EGA' MOTCEL MOTPRI) ;
  274.           ICEL = &BL1 ;
  275.           'QUITTER' BL1 ;
  276.        'FINSI' ;
  277.     'FIN' BL1 ;
  278.  
  279.     'SI' (ICEL > NDIM) ;
  280.        'MESSAGE' 'Procedure JACNUM' ;
  281.        'MESSAGE' 'MOTPRI = ??? ';
  282.        'ERREUR' 21 ;
  283.     'FINSI' ;
  284.  
  285.     ELT1 = 'MANUEL' 'POI1' PPRIM ;
  286.  
  287. *   ICEL = 1 -> On perturbe la densité
  288.  
  289.     'SI' ('EGA' ICEL 1) ;
  290.        DELTATOT = (EPSILON * dens0) ;
  291.        RNCEL = ('MANUEL' 'CHPO' ELT1 1 'SCAL' DELTATOT
  292.           'NATURE' 'DISCRET') 'ET' RN ;
  293.        GNCEL = 'COPIER' GN ;
  294.        RETNCEL = 'COPIER' RETN ;
  295.     'FINSI' ;
  296.  
  297. *   ICEL = 2 -> On perturbe la q.d.m. long  l'ax x
  298.  
  299.     'SI' ('EGA' ICEL 2) ;
  300.        DELTATOT = (EPSILON * dens0 * cson0) ;
  301.        GNCEL = ('MANUEL' 'CHPO' ELT1 1 'UX' DELTATOT
  302.           'NATURE' 'DISCRET') 'ET' GN ;
  303.        RNCEL = 'COPIER' RN ;
  304.        RETNCEL = 'COPIER' RETN ;
  305.     'FINSI' ;
  306.  
  307. *   ICEL = 3 -> On perturbe la q.d.m. long  l'ax y
  308.  
  309.     'SI' ('EGA' ICEL 3) ;
  310.        DELTATOT = (EPSILON * dens0 * cson0) ;
  311.        GNCEL = ('MANUEL' 'CHPO' ELT1 1 'UY' DELTATOT
  312.           'NATURE' 'DISCRET') 'ET' GN ;
  313.        RNCEL = 'COPIER' RN ;
  314.        RETNCEL = 'COPIER' RETN ;
  315.     'FINSI' ;
  316.  
  317. *   ICEL = 4 -> On perturbe l'énergie totale
  318.  
  319.     'SI' ('EGA' ICEL 4) ;
  320.        DELTATOT = (EPSILON * ret0) ;
  321.        RETNCEL = ('MANUEL' 'CHPO' ELT1 1 'SCAL' DELTATOT
  322.           'NATURE' 'DISCRET') 'ET' RETN ;
  323.        RNCEL = 'COPIER' RN ;
  324.        GNCEL = 'COPIER' GN ;
  325.     'FINSI' ;
  326.  
  327.     VITESSE PRES   = 'PRIM' 'PERFMONO'  RNCEL GNCEL RETNCEL GAMN ;
  328.  
  329.     ROF VITF PF GAMF   =  'PRET' 'PERFMONO'  1 1
  330.                            $MODE RNCEL VITESSE PRES GAMN ;
  331.  
  332.     CHPRES1 DT = 'KONV' 'VF' 'PERFMONO' 'RESI' NOMMET
  333.           $MODE ROF VITF PF GAMF LISTINCO MELLIM UINF UPRI;
  334.  
  335.     VAL1 = 'EXTRAIRE' CHPRES1 PDUAL MOTDUA ;
  336.  
  337.  'FINPROC' ((VAL1 '-' VAL0) '/' DELTATOT)  ;
  338.  
  339. *****************************************************
  340. *****************************************************
  341. ******** FIN PROCEDURES *****************************
  342. *****************************************************
  343. *****************************************************
  344. *****************************************************
  345.  
  346.  
  347. *****************************************************
  348. ******* TEST1 ***************************************
  349. *****************************************************
  350. *
  351. * On compare le jacobien et la variation du residu
  352. * en $DOM9 'CENTRE' par rapport à une variation
  353. * infinitésimal en $DOM9 'CENTRE'
  354. *
  355.  
  356.  PCEN9 = 'POIN' 1 ('DOMA' $DOM9 'CENTRE') ;
  357.  
  358. * Les grandeurs pour adimesionner les erreurs
  359.  
  360.  ro0 = 'EXTRAIRE' RN0 PCEN9 'SCAL' ;
  361.  
  362.  VN0 PN0 = 'PRIM' 'PERFMONO' RN0 GN0 RETN0 GAMMAN ;
  363.  CN20 = GAMMAN '*' (PN0 '/' RN0) ;
  364.  cson0 = ('EXTRAIRE' CN20 PCEN9 'SCAL') '**' 0.5 ;
  365.  
  366.  ret0 = ('EXTRAIRE' RETN0 PCEN9 'SCAL') '**' 0.5 ;
  367.  
  368. *
  369. * Le jacobien exact. 
  370. * DRR  = d(RES_ro)/dro (variable primale en PCEN9, variable duale en PCEN9) ;
  371. * DGXR = d(RES_gx)/dro (variable primale en PCEN9, variable duale en PCEN9) ;
  372. * ...
  373.  
  374.  DRR   = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMDEN NOMDEN ;
  375.  DGXR  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMDEN NOMMOX ;
  376.  DGYR  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMDEN NOMMOY ;
  377.  DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMDEN NOMRET ;
  378.  
  379.  DRGX   = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOX NOMDEN ;
  380.  DGXGX  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOX NOMMOX ;
  381.  DGYGX  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOX NOMMOY ;
  382.  DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOX NOMRET ;
  383.  
  384.  DRGY   = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOY NOMDEN ;
  385.  DGXGY  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOY NOMMOX ;
  386.  DGYGY  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOY NOMMOY ;
  387.  DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOY NOMRET ;
  388.  
  389.  DRRET   = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMRET NOMDEN ;
  390.  DGXRET  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMRET NOMMOX ;
  391.  DGYRET  = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMRET NOMMOY ;
  392.  DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMRET NOMRET ;
  393.  
  394. * Le jacobien numerique
  395.  
  396.  
  397.  DELTA = 1.0D-6 ;
  398.  
  399.  DRRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  400.         UPRI PCEN9 PCEN9 NOMDEN NOMDEN DELTA ;
  401.  DGXRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  402.         UPRI PCEN9 PCEN9 NOMDEN NOMMOX DELTA ;
  403.  DGYRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  404.         UPRI PCEN9 PCEN9 NOMDEN NOMMOY DELTA ;
  405.  DRETRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  406.         UPRI PCEN9 PCEN9 NOMDEN NOMRET DELTA ;
  407.  
  408.  DRGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  409.         UPRI PCEN9 PCEN9 NOMMOX NOMDEN DELTA ;
  410.  DGXGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  411.         UPRI PCEN9 PCEN9 NOMMOX NOMMOX DELTA ;
  412.  DGYGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  413.         UPRI PCEN9 PCEN9 NOMMOX NOMMOY DELTA ;
  414.  DRETGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  415.         UPRI PCEN9 PCEN9 NOMMOX NOMRET DELTA ;
  416.  
  417.  DRGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  418.         UPRI PCEN9 PCEN9 NOMMOY NOMDEN DELTA ;
  419.  DGXGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  420.         UPRI PCEN9 PCEN9 NOMMOY NOMMOX DELTA ;
  421.  DGYGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  422.         UPRI PCEN9 PCEN9 NOMMOY NOMMOY DELTA ;
  423.  DRETGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  424.         UPRI PCEN9 PCEN9 NOMMOY NOMRET DELTA ;
  425.  
  426.  DRRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  427.         UPRI PCEN9 PCEN9 NOMRET NOMDEN DELTA ;
  428.  DGXRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  429.         UPRI PCEN9 PCEN9 NOMRET NOMMOX DELTA ;
  430.  DGYRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  431.         UPRI PCEN9 PCEN9 NOMRET NOMMOY DELTA ;
  432.  DRETRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM
  433.  UINF UPRI PCEN9 PCEN9 NOMRET NOMRET DELTA ;
  434.  
  435. * Test des comparaisons jacobien exact-jacobien numerique
  436. **************************************************************
  437.  'SI' (('ABS'(DRR)) < DELTA) ;
  438.     ERR1 = ('ABS' (DRRN '-' DRR)) ;
  439.  'SINON' ;  
  440.     ERR1 = ('ABS' (DRRN '-' DRR)) '/' ('ABS' (DRR)) ;
  441.  'FINSI' ; 
  442. ***********************
  443.  'SI' (ERR1 > ERRTOL) ;
  444.       'MESSAGE' 'Problem 1';
  445.       'ERREUR' 5 ;
  446.  'FINSI' ;
  447. *************************************************************
  448.  'SI' (('ABS'(DGXR)) < DELTA) ;
  449.     ERR1 = ('ABS' (DGXRN '-' DGXR)) ;
  450.  'SINON' ;  
  451.     ERR1 = ('ABS' (DGXRN '-' DGXR)) '/' ('ABS' (DGXR)) ;
  452.  'FINSI' ; 
  453. ***********************
  454.  'SI' (ERR1 > ERRTOL) ;
  455.       'MESSAGE' 'Problem 2';
  456.       'ERREUR' 5 ;
  457.  'FINSI' ;
  458. **************************************************************
  459.  'SI' (('ABS'(DGYR)) < DELTA) ;
  460.     ERR1 = ('ABS' (DGYRN '-' DGYR)) ;
  461.  'SINON' ;  
  462.     ERR1 = ('ABS' (DGYRN '-' DGYR)) '/' ('ABS' (DGYR)) ;
  463.  'FINSI' ; 
  464. ************************
  465.  'SI' (ERR1 > ERRTOL) ;
  466.       'MESSAGE' 'Problem 3';
  467.        'ERREUR' 5 ;
  468.  'FINSI' ;
  469. **************************************************************
  470.  'SI' (('ABS'(DRETR)) < DELTA) ;
  471.     ERR1 = ('ABS' (DRETRN '-' DRETR)) ;
  472.  'SINON' ;  
  473.     ERR1 = ('ABS' (DRETRN '-' DRETR)) '/' ('ABS' (DRETR)) ;
  474.  'FINSI' ;
  475. **************************
  476.  'SI' (ERR1 > ERRTOL) ;
  477.       'MESSAGE' 'Problem 4';
  478.       'ERREUR' 5 ;
  479.  'FINSI' ;
  480. **************************************************************
  481.  'SI' (('ABS'(DRGX)) < DELTA) ;
  482.     ERR1 = ('ABS' (DRGXN '-' DRGX)) ;
  483.  'SINON' ;   
  484.     ERR1 = ('ABS' (DRGXN '-' DRGX)) '/' ('ABS' (DRGX))  ;
  485.  'FINSI' ;
  486. ***************************
  487.  'SI' (ERR1 > ERRTOL) ;
  488.       'MESSAGE' 'Problem 5';
  489.       'ERREUR' 5 ;
  490.  'FINSI' ;
  491.  
  492. *************************************************************
  493.  'SI' (('ABS'(DGXGX)) < DELTA) ;
  494.    ERR1 = ('ABS' (DGXGXN '-' DGXGX)) ;
  495.  'SINON' ;
  496.    ERR1 = ('ABS' (DGXGXN '-' DGXGX)) '/' ('ABS' (DGXGX)) ;
  497.  'FINSI' ;
  498. *****************************   
  499.  'SI' (ERR1 > ERRTOL) ;
  500.         'MESSAGE' 'Problem 6';
  501.         'ERREUR' 5 ;
  502.  'FINSI' ;
  503. *************************************************************
  504.  'SI' (('ABS'(DGYGX)) < DELTA) ;
  505.     ERR1 = ('ABS' (DGYGXN '-' DGYGX)) ;
  506.  'SINON' ;
  507.     ERR1 = ('ABS' (DGYGXN '-' DGYGX)) '/' ('ABS' (DGYGX)) ;
  508.  'FINSI' ;
  509. ******************************
  510.  'SI' (ERR1 > ERRTOL) ;
  511.         'MESSAGE' 'Problem 7';
  512.         'ERREUR' 5 ;
  513.  'FINSI' ;
  514. *************************************************************
  515.  'SI' (('ABS'(DRETGX)) < DELTA) ;
  516.     ERR1 = ('ABS' (DRETGXN '-' DRETGX)) ;
  517.  'SINON' ;
  518.     ERR1 = ('ABS' (DRETGXN '-' DRETGX)) '/' ('ABS' (DRETGX)) ;
  519.  'FINSI' ;
  520. *******************************    
  521.  'SI' (ERR1 > ERRTOL) ;
  522.         'MESSAGE' 'Problem 8';
  523.         'ERREUR' 5 ;
  524.  'FINSI' ;
  525. **************************************************************
  526.  'SI' (('ABS'(DRGY)) < DELTA) ;
  527.     ERR1 = ('ABS' (DRGYN '-' DRGY)) ;
  528.  'SINON' ;
  529.     ERR1 = ('ABS' (DRGYN '-' DRGY)) '/' ('ABS' (DRGY)) ;
  530.  'FINSI' ;
  531. *******************************
  532.  'SI' (ERR1 > ERRTOL) ;
  533.         'MESSAGE' 'Problem 9';
  534.         'ERREUR' 5 ;
  535.  'FINSI' ;
  536. **************************************************************
  537.  'SI' (('ABS'(DGXGY)) < DELTA) ;
  538.     ERR1 = ('ABS' (DGXGYN '-' DGXGY)) ;
  539.  'SINON' ;
  540.     ERR1 = ('ABS' (DGXGYN '-' DGXGY))  '/' ('ABS' (DGXGY)) ;
  541.  'FINSI' ;
  542. ********************************
  543.  'SI' (ERR1 > ERRTOL) ;
  544.         'MESSAGE' 'Problem 10';
  545.         'ERREUR' 5 ;
  546.  'FINSI' ;
  547. **************************************************************
  548.  'SI' (('ABS'(DGYGY)) < DELTA) ;
  549.     ERR1 = ('ABS' (DGYGYN '-' DGYGY)) ;
  550.  'SINON' ;
  551.     ERR1 = ('ABS' (DGYGYN '-' DGYGY)) '/' ('ABS' (DGYGY)) ;
  552.  'FINSI' ;
  553. ***************
  554.  'SI' (ERR1 > ERRTOL) ;
  555.         'MESSAGE' 'Problem 11';
  556.         'ERREUR' 5 ;
  557.  'FINSI' ;
  558. *************************************************************
  559.  'SI' (('ABS'(DRETGY)) < DELTA) ;
  560.     ERR1 = ('ABS' (DRETGYN '-' DRETGY)) ;
  561.  'SINON' ;
  562.     ERR1 = ('ABS' (DRETGYN '-' DRETGY)) '/' ('ABS' (DRETGY)) ;
  563.  'FINSI' ;
  564. ****************     
  565.  'SI' (ERR1 > ERRTOL) ;
  566.         'MESSAGE' 'Problem 12';
  567.         'ERREUR' 5 ;
  568.  'FINSI' ;
  569. **************************************************************
  570.  'SI' (('ABS'(DRRET)) < DELTA) ;
  571.     ERR1 = ('ABS' (DRRETN '-' DRRET)) ;
  572.  'SINON' ;
  573.     ERR1 = ('ABS' (DRRETN '-' DRRET)) '/' ('ABS' (DRRET)) ;
  574.  'FINSI' ;
  575. *****************
  576.  'SI' (ERR1 > ERRTOL) ;
  577.         'MESSAGE' 'Problem 13';
  578.         'ERREUR' 5 ;
  579.  'FINSI' ;
  580. **************************************************************
  581.  'SI' (('ABS'(DGXRET)) < DELTA) ;
  582.     ERR1 = ('ABS' (DGXRETN '-' DGXRET)) ;
  583.  'SINON' ;
  584.     ERR1 = ('ABS' (DGXRETN '-' DGXRET)) '/' ('ABS' (DGXRET)) ;
  585.  'FINSI' ;
  586. *****************
  587.  'SI' (ERR1 > ERRTOL) ;
  588.         'MESSAGE' 'Problem 14';
  589.         'ERREUR' 5 ;
  590.  'FINSI' ;
  591. *************************************************************
  592.  'SI' (('ABS'(DGYRET)) < DELTA) ;
  593.     ERR1 = ('ABS' (DGYRETN '-' DGYRET)) ; 
  594.  'SINON' ;
  595.     ERR1 = ('ABS' (DGYRETN '-' DGYRET)) '/' ('ABS' (DGYRET)) ;
  596.  'FINSI' ;
  597. *****************    
  598.  'SI' (ERR1 > ERRTOL) ;
  599.         'MESSAGE' 'Problem 15';
  600.         'ERREUR' 5 ;
  601.  'FINSI' ;
  602. **************************************************************
  603.  'SI' (('ABS'(DRETRET)) < DELTA) ;
  604.     ERR1 = ('ABS' (DRETRETN '-' DRETRET)) ;
  605.  'SINON' ;
  606.     ERR1 = ('ABS' (DRETRETN '-' DRETRET)) '/' ('ABS' (DRETRET)) ;
  607.  'FINSI' ;
  608. *****************
  609.  'SI' (ERR1 > ERRTOL) ;
  610.         'MESSAGE' 'Problem 16';
  611.         'ERREUR' 5 ;
  612.  'FINSI' ;
  613.  
  614. *****************************************************
  615. *****************************************************
  616. ******* TEST2 ***************************************
  617. *****************************************************
  618. *****************************************************
  619. *
  620. * On compare le jacobien et la variation du residu
  621. * en $DOM9 'CENTRE' par rapport à une variation
  622. * infinitésimal en $DOM6 'CENTRE'
  623. *
  624.  
  625.  PCEN6 = ('DOMA' $DOM6 'CENTRE') 'POIN' 1 ;
  626.  
  627.  ro0 = 'EXTRAIRE' RN0 PCEN6 'SCAL' ;
  628.  cson0 = ('EXTRAIRE' CN20 PCEN6 'SCAL') '**' 0.5 ;
  629.  ret0 = 'EXTRAIRE' RETN0 PCEN6 'SCAL' ;
  630.  
  631. *
  632. * Le jacobien exact. 
  633. *
  634.  
  635.  DRR   = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMDEN NOMDEN ;
  636.  DGXR  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMDEN NOMMOX ;
  637.  DGYR  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMDEN NOMMOY ;
  638.  DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMDEN NOMRET ;
  639.  
  640.  DRGX   = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOX NOMDEN ;
  641.  DGXGX  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOX NOMMOX ;
  642.  DGYGX  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOX NOMMOY ;
  643.  DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOX NOMRET ;
  644.  
  645.  DRGY   = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOY NOMDEN ;
  646.  DGXGY  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOY NOMMOX ;
  647.  DGYGY  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOY NOMMOY ;
  648.  DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOY NOMRET ;
  649.  
  650.  DRRET   = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMRET NOMDEN ;
  651.  DGXRET  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMRET NOMMOX ;
  652.  DGYRET  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMRET NOMMOY ;
  653.  DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMRET NOMRET ;
  654.  
  655. * Le jacobien numerique
  656.  
  657.  
  658.  DELTA = 1.0D-6 ;
  659.  
  660.  DRRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  661.    UPRI PCEN6 PCEN9  NOMDEN NOMDEN DELTA ;
  662.  DGXRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  663.    UPRI PCEN6 PCEN9 NOMDEN NOMMOX DELTA ;
  664.  DGYRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  665.    UPRI PCEN6 PCEN9  NOMDEN NOMMOY DELTA ;
  666.  DRETRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  667.    UPRI PCEN6 PCEN9  NOMDEN NOMRET DELTA ;
  668.  
  669.  DRGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  670.    UPRI PCEN6 PCEN9  NOMMOX NOMDEN DELTA ;
  671.  DGXGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  672.    UPRI PCEN6 PCEN9  NOMMOX NOMMOX DELTA ;
  673.  DGYGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  674.    UPRI PCEN6 PCEN9  NOMMOX NOMMOY DELTA ;
  675.  DRETGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  676.    UPRI PCEN6 PCEN9  NOMMOX NOMRET DELTA ;
  677.  
  678.  DRGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  679.    UPRI PCEN6 PCEN9  NOMMOY NOMDEN DELTA ;
  680.  DGXGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  681.    UPRI PCEN6 PCEN9  NOMMOY NOMMOX DELTA ;
  682.  DGYGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  683.    UPRI PCEN6 PCEN9  NOMMOY NOMMOY DELTA ;
  684.  DRETGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  685.    UPRI PCEN6 PCEN9  NOMMOY NOMRET DELTA ;
  686.  
  687.  DRRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  688.    UPRI PCEN6 PCEN9  NOMRET NOMDEN DELTA ;
  689.  DGXRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  690.    UPRI PCEN6 PCEN9  NOMRET NOMMOX DELTA ;
  691.  DGYRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM UINF
  692.    UPRI PCEN6 PCEN9  NOMRET NOMMOY DELTA ;
  693.  DRETRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO MELLIM
  694.    UINF UPRI PCEN6   PCEN9  NOMRET NOMRET DELTA ;
  695.  
  696. * Test des comparaisons jacobien exact-jacobien numerique
  697. **************************************************************
  698.  'SI' (('ABS'(DRR)) < DELTA) ;
  699.     ERR1 = ('ABS' (DRRN '-' DRR)) ;
  700.  'SINON' ;  
  701.     ERR1 = ('ABS' (DRRN '-' DRR)) '/' ('ABS' (DRR)) ;
  702.  'FINSI' ; 
  703. ***********************
  704.  'SI' (ERR1 > ERRTOL) ;
  705.       'MESSAGE' 'Problem 1';
  706.       'ERREUR' 5 ;
  707.  'FINSI' ;
  708. *************************************************************
  709.  'SI' (('ABS'(DGXR)) < DELTA) ;
  710.     ERR1 = ('ABS' (DGXRN '-' DGXR)) ;
  711.  'SINON' ;  
  712.     ERR1 = ('ABS' (DGXRN '-' DGXR)) '/' ('ABS' (DGXR)) ;
  713.  'FINSI' ; 
  714. ***********************
  715.  'SI' (ERR1 > ERRTOL) ;
  716.       'MESSAGE' 'Problem 2';
  717.       'ERREUR' 5 ;
  718.  'FINSI' ;
  719. **************************************************************
  720.  'SI' (('ABS'(DGYR)) < DELTA) ;
  721.     ERR1 = ('ABS' (DGYRN '-' DGYR)) ;
  722.  'SINON' ;  
  723.     ERR1 = ('ABS' (DGYRN '-' DGYR)) '/' ('ABS' (DGYR)) ;
  724.  'FINSI' ; 
  725. ************************
  726.  'SI' (ERR1 > ERRTOL) ;
  727.       'MESSAGE' 'Problem 3';
  728.        'ERREUR' 5 ;
  729.  'FINSI' ;
  730. **************************************************************
  731.  'SI' (('ABS'(DRETR)) < DELTA) ;
  732.     ERR1 = ('ABS' (DRETRN '-' DRETR)) ;
  733.  'SINON' ;  
  734.     ERR1 = ('ABS' (DRETRN '-' DRETR)) '/' ('ABS' (DRETR)) ;
  735.  'FINSI' ;
  736. **************************
  737.  'SI' (ERR1 > ERRTOL) ;
  738.       'MESSAGE' 'Problem 4';
  739.       'ERREUR' 5 ;
  740.  'FINSI' ;
  741. **************************************************************
  742.  'SI' (('ABS'(DRGX)) < DELTA) ;
  743.     ERR1 = ('ABS' (DRGXN '-' DRGX)) ;
  744.  'SINON' ;   
  745.     ERR1 = ('ABS' (DRGXN '-' DRGX)) '/' ('ABS' (DRGX))  ;
  746.  'FINSI' ;
  747. ***************************
  748.  'SI' (ERR1 > ERRTOL) ;
  749.       'MESSAGE' 'Problem 5';
  750.       'ERREUR' 5 ;
  751.  'FINSI' ;
  752. *************************************************************
  753.  'SI' (('ABS'(DGXGX)) < DELTA) ;
  754.    ERR1 = ('ABS' (DGXGXN '-' DGXGX)) ;
  755.  'SINON' ;
  756.    ERR1 = ('ABS' (DGXGXN '-' DGXGX)) '/' ('ABS' (DGXGX)) ;
  757.  'FINSI' ;  
  758. *****************************   
  759.  'SI' (ERR1 > ERRTOL) ;
  760.         'MESSAGE' 'Problem 6';
  761.         'ERREUR' 5 ;
  762.  'FINSI' ;
  763. *************************************************************
  764.  'SI' (('ABS'(DGYGX)) < DELTA) ;
  765.     ERR1 = ('ABS' (DGYGXN '-' DGYGX)) ;
  766.  'SINON' ;
  767.     ERR1 = ('ABS' (DGYGXN '-' DGYGX)) '/' ('ABS' (DGYGX)) ;
  768.  'FINSI' ;
  769. ******************************
  770.  'SI' (ERR1 > ERRTOL) ;
  771.         'MESSAGE' 'Problem 7';
  772.         'ERREUR' 5 ;
  773.  'FINSI' ;
  774. *************************************************************
  775.  'SI' (('ABS'(DRETGX)) < DELTA) ;
  776.     ERR1 = ('ABS' (DRETGXN '-' DRETGX)) ;
  777.  'SINON' ;
  778.     ERR1 = ('ABS' (DRETGXN '-' DRETGX)) '/' ('ABS' (DRETGX)) ;
  779.  'FINSI' ;
  780. *******************************    
  781.  'SI' (ERR1 > ERRTOL) ;
  782.         'MESSAGE' 'Problem 8';
  783.         'ERREUR' 5 ;
  784.  'FINSI' ;
  785. **************************************************************
  786.  'SI' (('ABS'(DRGY)) < DELTA) ;
  787.     ERR1 = ('ABS' (DRGYN '-' DRGY)) ;
  788.  'SINON' ;
  789.     ERR1 = ('ABS' (DRGYN '-' DRGY)) '/' ('ABS' (DRGY)) ;
  790.  'FINSI' ;
  791. *******************************
  792.  'SI' (ERR1 > ERRTOL) ;
  793.         'MESSAGE' 'Problem 9';
  794.         'ERREUR' 5 ;
  795.  'FINSI' ;
  796. **************************************************************
  797.  'SI' (('ABS'(DGXGY)) < DELTA) ;
  798.     ERR1 = ('ABS' (DGXGYN '-' DGXGY)) ;
  799.  'SINON' ;
  800.     ERR1 = ('ABS' (DGXGYN '-' DGXGY))  '/' ('ABS' (DGXGY)) ;
  801.  'FINSI' ;
  802. ********************************
  803.  'SI' (ERR1 > ERRTOL) ;
  804.         'MESSAGE' 'Problem 10';
  805.         'ERREUR' 5 ;
  806.  'FINSI' ;
  807. **************************************************************
  808.  'SI' (('ABS'(DGYGY)) < DELTA) ;
  809.     ERR1 = ('ABS' (DGYGYN '-' DGYGY)) ;
  810.  'SINON' ;
  811.     ERR1 = ('ABS' (DGYGYN '-' DGYGY)) '/' ('ABS' (DGYGY)) ;
  812.  'FINSI' ;
  813. ***************
  814.  'SI' (ERR1 > ERRTOL) ;
  815.         'MESSAGE' 'Problem 11';
  816.         'ERREUR' 5 ;
  817.  'FINSI' ;
  818. *************************************************************
  819.  'SI' (('ABS'(DRETGY)) < DELTA) ;
  820.     ERR1 = ('ABS' (DRETGYN '-' DRETGY)) ;
  821.  'SINON' ;
  822.     ERR1 = ('ABS' (DRETGYN '-' DRETGY)) '/' ('ABS' (DRETGY)) ;
  823.  'FINSI' ;
  824. ****************     
  825.  'SI' (ERR1 > ERRTOL) ;
  826.         'MESSAGE' 'Problem 12';
  827.         'ERREUR' 5 ;
  828.  'FINSI' ;
  829. **************************************************************
  830.  'SI' (('ABS'(DRRET)) < DELTA) ;
  831.     ERR1 = ('ABS' (DRRETN '-' DRRET)) ;
  832.  'SINON' ;
  833.     ERR1 = ('ABS' (DRRETN '-' DRRET)) '/' ('ABS' (DRRET)) ;
  834.  'FINSI' ;
  835. *****************
  836.  'SI' (ERR1 > ERRTOL) ;
  837.         'MESSAGE' 'Problem 13';
  838.         'ERREUR' 5 ;
  839.  'FINSI' ;
  840. **************************************************************
  841.  'SI' (('ABS'(DGXRET)) < DELTA) ;
  842.     ERR1 = ('ABS' (DGXRETN '-' DGXRET)) ;
  843.  'SINON' ;
  844.     ERR1 = ('ABS' (DGXRETN '-' DGXRET)) '/' ('ABS' (DGXRET)) ;
  845.  'FINSI' ;
  846. *****************
  847.  'SI' (ERR1 > ERRTOL) ;
  848.         'MESSAGE' 'Problem 14';
  849.         'ERREUR' 5 ;
  850.  'FINSI' ;
  851. *************************************************************
  852.  'SI' (('ABS'(DGYRET)) < DELTA) ;
  853.     ERR1 = ('ABS' (DGYRETN '-' DGYRET)) ; 
  854.  'SINON' ;
  855.     ERR1 = ('ABS' (DGYRETN '-' DGYRET)) '/' ('ABS' (DGYRET)) ;
  856.  'FINSI' ;
  857. *****************    
  858.  'SI' (ERR1 > ERRTOL) ;
  859.         'MESSAGE' 'Problem 15';
  860.         'ERREUR' 5 ;
  861.  'FINSI' ;
  862. **************************************************************
  863.  'SI' (('ABS'(DRETRET)) < DELTA) ;
  864.     ERR1 = ('ABS' (DRETRETN '-' DRETRET)) ;
  865.  'SINON' ;
  866.     ERR1 = ('ABS' (DRETRETN '-' DRETRET)) '/' ('ABS' (DRETRET)) ;
  867.  'FINSI' ;
  868. *****************
  869.  'SI' (ERR1 > ERRTOL) ;
  870.         'MESSAGE' 'Problem 16';
  871.         'ERREUR' 5 ;
  872.  'FINSI' ;
  873.  
  874.  
  875.  
  876. *****************************************************
  877. *****************************************************
  878. ******* TEST3 ***************************************
  879. *****************************************************
  880. *****************************************************
  881. *
  882. * On observe que la variation du residu
  883. * en $DOM9 'CENTRE' par rapport à une variation
  884. * infinitésimal en $DOM7 'CENTRE' doit etre nulle.
  885. * 
  886.  
  887.  PCEN7 = ('DOMA' $DOM7 'CENTRE') 'POIN' 1 ;
  888.  
  889.  ro0 = 'EXTRAIRE' RN0 PCEN7 'SCAL' ;
  890.  cson0 = ('EXTRAIRE' CN20 PCEN7 'SCAL') '**' 0.5 ;
  891.  ret0 = 'EXTRAIRE' RETN0 PCEN7 'SCAL' ;
  892.  
  893. *
  894. * Le jacobien exact. 
  895. *
  896.  
  897.  DRR   = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMDEN NOMDEN ;
  898.  DGXR  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMDEN NOMMOX ;
  899.  DGYR  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMDEN NOMMOY ;
  900.  DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMDEN NOMRET ;
  901.  
  902.  DRGX   = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOX NOMDEN ;
  903.  DGXGX  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOX NOMMOX ;
  904.  DGYGX  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOX NOMMOY ;
  905.  DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOX NOMRET ;
  906.  
  907.  DRGY   = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOY NOMDEN ;
  908.  DGXGY  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOY NOMMOX ;
  909.  DGYGY  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOY NOMMOY ;
  910.  DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOY NOMRET ;
  911.  
  912.  DRRET   = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMRET NOMDEN ;
  913.  DGXRET  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMRET NOMMOX ;
  914.  DGYRET  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMRET NOMMOY ;
  915.  DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMRET NOMRET ;
  916.  
  917.  
  918.  'SI' (('MAXIMUM'
  919.     ('PROG'  DRR  DGXR  DGYR  DRETR
  920.              DRGX DGXGX DGYGX DRETGX
  921.              DRGY DGXGY DGYGY DRETGY
  922.              DRRET DGXRET DGYRET DRETRET ) 'ABS' )
  923.       > ERRTOL) ;
  924.       'MESSAGE' 'Problem final' ;
  925.       'ERREUR' 5 ;
  926.  'FINSI' ;    
  927.  
  928. *****************************************************
  929. *****************************************************
  930. ******* TEST4 ***************************************
  931. *****************************************************
  932. *****************************************************
  933. * AUSMPLUS = AUSMPBM in the following cases:
  934. *            1) VN0 > CN
  935. *            2) UINF1 > CN
  936. *            3) UINF2 > CN
  937. *
  938. * Let us start from RN0, VN0, PN0, GAMMAN
  939. *
  940.  
  941.  ROF VITF PF GAMF   =  'PRET' 'PERFMONO'  1 1
  942.                            $DOMTOT RN0 VN0 PN0 GAMMAN ;
  943.  
  944.  ECIN = 0.5D0 * RN0 * ('PSCAL' VN0 VN0 ('MOTS' 'UX' 'UY')
  945.         ('MOTS' 'UX' 'UY')) ;
  946.  RETN0  = 'NOMC' 'SCAL' ((PN0  '/' (GAMMAN  '-' 1.0)) '+' ECIN)
  947.    'NATU' 'DISCRET' ;
  948.  GN0 = VN0 '*' RN0 ('MOTS' 'UX' 'UY') ('MOTS' 'SCAL' 'SCAL')
  949.     ('MOTS' 'UX' 'UY');
  950.  CN0 = (GAMMAN '*' (PN0 '/' RN0)) '**' 0.5 ;
  951.  
  952.  UINF2 = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'SCAL' 0.0 ;
  953.  UINF1 = CN0 * 1.17 ;
  954.  
  955.  CHPRES DT = 'KONV' 'VF' 'PERFMONO' 'RESI' 'AUSMPLM'
  956.           $DOMTOT ROF VITF PF GAMF LISTINCO
  957.           MELLIM UINF1 UINF2 ;
  958.  
  959.  JACO = 'KONV' 'VF' 'PERFMONO' 'JACOCONS'
  960.          $DOMTOT LISTINCO MELLIM 
  961.          'AUSMPLM'  RN0 VN0 PN0 GAMMAN
  962.          UINF2 UINF1 ;
  963.  
  964.  CHPRES1 DT = 'KONV' 'VF' 'PERFMONO' 'RESI' 'AUSMPLUS'
  965.           $DOMTOT ROF VITF PF GAMF LISTINCO
  966.           MELLIM ;
  967.  
  968.  JACO1 = 'KONV' 'VF' 'PERFMONO' 'JACOCONS'
  969.          $DOMTOT LISTINCO MELLIM 
  970.          'AUSMPLUS' RN0 VN0 PN0 GAMMAN ;
  971.  
  972.  ERRO = 'MAXIMUM' (CHPRES1 '-' CHPRES) 'ABS' ;
  973.  
  974.  'SI' (ERRO > ERRTOL) ;
  975.       'MESSAGE' 'Problem final 2' ;
  976.       'ERREUR' 5 ;
  977.  'FINSI' ;
  978.  
  979.  UNTEST = ('NOMC' RN0 NOMDEN 'NATU' 'DISCRET') 'ET'
  980.           ('NOMC' GN0 ('MOTS' 'UX' 'UY') ('MOTS' NOMMOX NOMMOY)
  981.            'NATU' 'DISCRET') 'ET'
  982.           ('NOMC' RETN0 NOMRET 'NATU' 'DISCRET') ;
  983.  
  984.  ERRO = 'MAXIMUM' (('KOPS' 'MULT' JACO UNTEST) '-'
  985.                    ('KOPS' 'MULT' JACO1 UNTEST)) 'ABS' ;
  986.  
  987.  SI (ERRO > ERRTOL) ;           
  988.       'MESSAGE' 'Problem final 2' ;
  989.       'ERREUR' 5 ;
  990.  'FINSI' ;
  991.  
  992. * Redefinition of VN0  
  993.  
  994.  CN0 = (GAMMAN '*' (PN0 '/' RN0)) '**' 0.5 ;
  995.  VN0 = ('NOMC' (CN0 * 1.11) 'UX' 'NATU' 'DISCRET') 'ET'
  996.        ('NOMC' (CN0 * 2.11) 'UY' 'NATU' 'DISCRET') ;
  997.  
  998.  ROF VITF PF GAMF   =  'PRET' 'PERFMONO'  1 1
  999.                            $DOMTOT RN0 VN0 PN0 GAMMAN ;
  1000.  
  1001.  UINF2 = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'SCAL' 0.0 ;
  1002.  UINF1 = UINF2 ;
  1003.  
  1004.  CHPRES DT = 'KONV' 'VF' 'PERFMONO' 'RESI' 'AUSMPLM'
  1005.           $DOMTOT ROF VITF PF GAMF LISTINCO
  1006.           MELLIM UINF1 UINF2 ;
  1007.  
  1008.  JACO = 'KONV' 'VF' 'PERFMONO' 'JACOCONS'
  1009.          $DOMTOT LISTINCO MELLIM 
  1010.          'AUSMPLM'  RN0 VN0 PN0 GAMMAN
  1011.          UINF2 UINF1 ;
  1012.  
  1013.  CHPRES1 DT = 'KONV' 'VF' 'PERFMONO' 'RESI' 'AUSMPLUS'
  1014.           $DOMTOT ROF VITF PF GAMF LISTINCO
  1015.           MELLIM ;
  1016.  
  1017.  JACO1 = 'KONV' 'VF' 'PERFMONO' 'JACOCONS'
  1018.          $DOMTOT LISTINCO MELLIM 
  1019.          'AUSMPLUS' RN0 VN0 PN0 GAMMAN ;
  1020.  
  1021.  ERRO = 'MAXIMUM' (CHPRES1 '-' CHPRES) 'ABS' ;
  1022.  
  1023.  'SI' (ERRO > ERRTOL) ;
  1024.       'MESSAGE' 'Problem final 2' ;
  1025.       'ERREUR' 5 ;
  1026.  'FINSI' ;
  1027.  
  1028.  UNTEST = ('NOMC' RN0 NOMDEN 'NATU' 'DISCRET') 'ET'
  1029.           ('NOMC' GN0 ('MOTS' 'UX' 'UY') ('MOTS' NOMMOX NOMMOY)
  1030.            'NATU' 'DISCRET') 'ET'
  1031.           ('NOMC' RETN0 NOMRET 'NATU' 'DISCRET') ;
  1032.  
  1033.  ERRO = 'MAXIMUM' (('KOPS' 'MULT' JACO UNTEST) '-'
  1034.                    ('KOPS' 'MULT' JACO1 UNTEST)) 'ABS' ;
  1035.  
  1036.  SI (ERRO > ERRTOL) ;           
  1037.       'MESSAGE' 'Problem final 3' ;
  1038.       'ERREUR' 5 ;
  1039.  'FINSI' ;
  1040.  
  1041. *****************************************************
  1042. *****************************************************
  1043. ******* TEST5 ***************************************
  1044. *****************************************************
  1045. *****************************************************
  1046. * Primitive variables
  1047. * Redefinition of VN0  
  1048. *
  1049.  UINF2 = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'SCAL' 0.0 ;
  1050.  UINF1 = UINF2 ;
  1051.  
  1052.  RN0 = 'BRUI' 'BLAN' 'UNIF' 1.11 0.5 ('DOMA' $DOMTOT 'CENTRE') ;
  1053.  PN0 = 'BRUI' 'BLAN' 'UNIF' 1234.1 800 ('DOMA' $DOMTOT 'CENTRE') ;
  1054.  GAMMAN = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'SCAL' 1.4
  1055.    'NATU' 'DISCRET' ;
  1056.  CN0 = (GAMMAN '*' PN0) '/' RN0 ;
  1057.  CN0 = 'KOPS' CSONN '**' 0.5 ;
  1058.  VN0 = ('NOMC' (CN0 * 0.11) 'UX' 'NATU' 'DISCRET') 'ET'
  1059.        ('NOMC' (CN0 * 2.11) 'UY' 'NATU' 'DISCRET') ;
  1060.  
  1061.  LISTCONS = 'MOTS' 'RN' 'RUX' 'RUY' 'RETN' ;
  1062.  
  1063.  JACO  = 'KONV' 'VF' 'PERFMONO' 'JACOCONS'
  1064.           $DOMTOT LISTCONS METO RN0 VN0 PN0 GAMMAN
  1065.           UINF1 UINF2;
  1066.  
  1067.  LISTPRIM = 'MOTS' 'RN' 'UX' 'UY' 'PN' ;
  1068.  
  1069.  DUDV = 'PRIM' 'CONSPRIM' ('DOMA' $DOMTOT 'CENTRE')
  1070.          LISTCONS LISTPRIM RN0 VN0 PN0 GAMMAN ;
  1071.  
  1072.  JACOP = 'KONV' 'VF' 'PERFMONO' 'JACOPRIM'
  1073.          $DOMTOT LISTCONS LISTPRIM METO RN0 VN0 PN0 GAMMAN
  1074.          UINF1 UINF2; 
  1075.  
  1076.  UNP = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 4  'RN' 0.11
  1077.        'UX' 0.15 'UY' 0.19 'PN' 0.21 ;
  1078.  
  1079. *  JACOP * UNP = JACO * DUDV * UNP, \forall UNP       
  1080.  
  1081.  ERRO = 'MAXIMUM' (
  1082.      ('KOPS' 'MULT'  JACOP  UNP) '-' 
  1083.      ('KOPS' 'MULT'  JACO  ('KOPS' 'MULT' DUDV  UNP))
  1084.      ) ;
  1085.  
  1086.  'SI' (ERRO > ERRTOL) ;
  1087.       'MESSAGE' 'Problem final' ;
  1088.       'ERREUR' 5 ;
  1089.  'FINSI' ;
  1090.  
  1091. *
  1092.  
  1093.  UINF2 = CN0 ;
  1094.  UINF1 = UINF2 ;
  1095.  
  1096.  JACO  = 'KONV' 'VF' 'PERFMONO' 'JACOCONS'
  1097.           $DOMTOT LISTCONS METO RN0 VN0 PN0 GAMMAN
  1098.           UINF1 UINF2;
  1099.  
  1100.  JACOP = 'KONV' 'VF' 'PERFMONO' 'JACOPRIM'
  1101.          $DOMTOT LISTCONS LISTPRIM METO RN0 VN0 PN0 GAMMAN
  1102.          UINF1 UINF2; 
  1103.  
  1104. *  JACOP * UNP = JACO * DUDV * UNP, \forall UNP       
  1105.  
  1106.  ERRO = 'MAXIMUM' (
  1107.      ('KOPS' 'MULT'  JACOP  UNP) '-' 
  1108.      ('KOPS' 'MULT'  JACO  ('KOPS' 'MULT' DUDV  UNP))
  1109.      ) ;
  1110.  
  1111.  'SI' (ERRO > ERRTOL) ;
  1112.       'MESSAGE' 'Problem final' ;
  1113.       'ERREUR' 5 ;
  1114.  'FINSI' ;
  1115.  
  1116.  
  1117.  'FIN' ;
  1118.  
  1119.  
  1120.  
  1121.  
  1122.  
  1123.  
  1124.  
  1125.  
  1126.  

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