Télécharger la plaquette Cast3M


fr en

Télécharger konv_impl_centre2.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier :  konv_impl_centre2.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: AUSM+                                     ****
  15. ****                                                     ****
  16. **** A. BECCANTINI, S. KUDRIAKOV SFME/LTMF, JUILLET 2001 ****
  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-3 ;
  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.  
  64.  $DOMTOT = 'MODELISER' DOMTOT 'EULER';
  65.  
  66.  $DOM6 = 'MODELISER'  DOM6 'EULER';
  67.  $DOM7 = 'MODELISER'  DOM7 'EULER';
  68.  $DOM5 = 'MODELISER'  DOM5 'EULER';
  69.  
  70.  TDOMTOT = 'DOMA'  $DOMTOT 'VF';
  71.  
  72.  TDOM6 = 'DOMA'  $DOM6 'VF';
  73.  TDOM7 = 'DOMA'  $DOM7 'VF';
  74.  TDOM5 = 'DOMA'  $DOM5 'VF';
  75.  
  76.  MDOM6 = TDOM6 . 'QUAF' ;
  77.  MDOM7 = TDOM7 . 'QUAF' ;
  78.  MDOM5 = TDOM5 . 'QUAF' ;          
  79.  
  80. **** old stuff  $DOMTOT = 'DOMA' DOMTOT ;
  81.  
  82.  MDOMTOT = TDOMTOT . 'QUAF' ;
  83.  
  84.  'ELIMINATION'  (MDOMTOT ET MDOM6) 0.0001 ;
  85.  'ELIMINATION'  (MDOMTOT ET MDOM7) 0.0001 ;
  86.  'ELIMINATION'  (MDOMTOT ET MDOM5) 0.0001 ;
  87.  
  88. ***************************************************
  89. ***** Densité, pression, vitesse, gamma ***********
  90. ***************************************************
  91.  
  92.  RN0 = 'BRUI' 'BLAN' 'UNIF' 1.11 0.5 ('DOMA' $DOMTOT 'CENTRE') ;
  93.  PN0 = 'BRUI' 'BLAN' 'UNIF' 1234.1 800 ('DOMA' $DOMTOT 'CENTRE') ;
  94.  GAMMAN = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'SCAL' 1.4
  95.    'NATU' 'DISCRET' ;
  96.  CSONN = (GAMMAN '*' PN0) '/' RN0 ;
  97.  CSONN = 'KOPS' CSONN '**' 0.5 ;
  98.  UXN0  = 0.2 *  CSONN ;
  99.  UYN0  = 0.3  * CSONN ;
  100.  
  101.  GN0 = ('NOMC' (RN0 '*' UXN0) 'UX' 'NATU' 'DISCRET') 'ET'
  102.        ('NOMC' (RN0 '*' UYN0) 'UY' 'NATU' 'DISCRET') ;
  103.  
  104.  ECIN = 0.5D0 '*' RN0 '*' ((UXN0  '*' UXN0) '+' (UYN0  '*' UYN0)) ;
  105.  RETN0  = 'NOMC' 'SCAL' ((PN0  '/' (GAMMAN  '-' 1.0)) '+' ECIN)
  106.    'NATU' 'DISCRET' ;
  107.  
  108.  VIT PRES = PRIM 'PERFMONO' RN0 GN0 RETN0 GAMMAN ;
  109.  
  110.  ERR1 = 'MAXIMUM' (PRES '-' PN0) 'ABS' ;
  111.  ERR2 = 'MAXIMUM' (GN0 '-' (RN0 '*' VIT ('MOTS' 'SCAL' 'SCAL')
  112.          ('MOTS' 'UX' 'UY') ('MOTS' 'UX' 'UY'))) 'ABS' ;
  113.  
  114.  'SI' (('MAXIMUM' ('PROG' ERR1 ERR2 ) 'ABS' ) > 1.0D-6) ;
  115.       'MESSAGE' 'Problem 0' ;
  116.       'ERREUR' 5 ;
  117.  'FINSI' ;
  118.  
  119.  
  120.  'SI' GRAPH;
  121.       'TRACER' (('DOMA' $DOMTOT 'MAILLAGE' ) 'ET' 
  122.           ('DOMA' $DOMTOT 'CENTRE')) 'TITRE' 'Domaine et centre' ;
  123.  'FINSI' ;
  124.  
  125. *
  126. **** Les variables conservative
  127. *
  128. *    RN0 (densité)
  129. *    GN0 (quantité de mouvement)
  130. *    RETN0 (énergie totale par unité de volume)
  131. *
  132. *    sont definiés
  133. *
  134.  
  135. ****************************************************
  136. ****************************************************
  137. ******* Calcul du jacobien et du residu   **********
  138. ****************************************************
  139. ****************************************************
  140. *
  141. * JACO est le jacobien
  142. *
  143. * DEBRN0 le residu concernant la densité
  144. * DEBGNX0 le residu concernant la quantité de mouvement (axe x)
  145. * DEBGNY0 le residu concernant la quantité de mouvement (axe y)
  146. * DEBRETN0 le residu concernant l'enrgie totale par unité de volume
  147. *
  148.  
  149. * Noms des variables
  150.  
  151.  NOMDEN = 'RN  ' ;
  152.  NOMMOX = 'RUXN' ;
  153.  NOMMOY = 'RUYN' ;
  154.  NOMRET = 'RETN' ;
  155.  
  156. * Metode
  157.  
  158.  METO = 'AUSMPLUS' ;
  159.  
  160.  
  161.  LISTINCO = 'MOTS'  NOMDEN  NOMMOX  NOMMOY  NOMRET ;
  162.  
  163.  VITESSE PRES   = 'PRIM' 'PERFMONO'  RN0 GN0 RETN0 GAMMAN ;
  164.  
  165.  ROF VITF PF GAMF   =  'PRET' 'PERFMONO'  1 1
  166.                            $DOMTOT RN0 VITESSE PRES GAMMAN ;
  167.  
  168.  CHPRES0 DT = 'KONV' 'VF' 'PERFMONO' 'RESI' METO
  169.           $DOMTOT ROF VITF PF GAMF LISTINCO ;
  170.  
  171.  JACO = 'KONV' 'VF' 'PERFMONO' 'JACOCONS' 
  172.      $DOMTOT LISTINCO
  173.     METO RN0 VITESSE PRES GAMMAN ;
  174.  
  175.  DEBRN0   = 'EXCO' NOMDEN CHPRES0 ;
  176.  DEBGNX0  = 'EXCO' NOMMOX CHPRES0 ;
  177.  DEBGNY0  = 'EXCO' NOMMOY CHPRES0 ;
  178.  DEBRETN0 = 'EXCO' NOMRET CHPRES0 ;
  179.  
  180. *****************************************************
  181. *****************************************************
  182. ******** PROCEDURES *********************************
  183. *****************************************************
  184. *****************************************************
  185.  
  186. *
  187. * Derivé partielle du residu en un point par rapport
  188. * aux variable en un autre point
  189. *
  190.  
  191.  'DEBPROC' JACOVA JACO*'MATRIK' $MODE*'MMODEL' LISTINCO*'LISTMOTS'
  192.     PPRIM*'POINT' PDUAL*'POINT' MOTPRI*'MOT' MOTDUA*'MOT';
  193.  
  194. *   PPRIM = point ou est localisé la variable primale
  195. *   PDUAL = point ou est localisé la variable duale
  196. *   MOTPRI = nom de la composante concernante la variable primale
  197. *   MOTDUA = nom de la composante concernante la variable duale
  198.  
  199.     ELT1 = 'MANUEL' 'POI1' PPRIM ;
  200.     NDIM = 'DIME' LISTINCO ;
  201.     CHPUN = 'MANUEL' 'CHPO' ELT1 1 MOTPRI 1.0 'NATURE' 'DISCRET' ;
  202.     'REPETER' BL1 NDIM ;
  203.        MOTCEL = 'EXTRAIRE' LISTINCO &BL1 ;
  204.        CHPUN = CHPUN 'ET' ('MANUEL' 'CHPO' ('DOMA' $MODE 'CENTRE')
  205.           1 MOTCEL 0.0 'NATURE' 'DISCRET') ;
  206.     'FIN' BL1 ;
  207.     D_DMOT = 'KOPS' JACO 'MULT' CHPUN ;
  208.     SCAL = 'EXTRAIRE' D_DMOT PDUAL MOTDUA ;
  209.  
  210.  'FINPROC' SCAL  ;
  211.  
  212.  
  213.  
  214.  'DEBPROC' JACNUM $MODE*'MMODEL' NOMMET*'MOT' RN*'CHPOINT' GN*'CHPOINT'
  215.     RETN*'CHPOINT' GAMN*'CHPOINT'
  216.     LISTINCO*'LISTMOTS' PPRIM*'POINT' PDUAL*'POINT'
  217.     MOTPRI*'MOT' MOTDUA*'MOT' EPSILON*'FLOTTANT';
  218.  
  219. *   PPRIM = point ou est localisé la variable primale
  220. *   PDUAL = point ou est localisé la variable duale
  221. *   MOTPRI = nom de la composante concernante la variable primale
  222. *   MOTDUA = nom de la composante concernante la variable duale
  223.  
  224. *   Le valeur dans l'état non-perturbé en PDUAL ;
  225.  
  226.     RNCEL   = 'COPIER'  RN ;
  227.     GNCEL   = 'COPIER'  GN ;
  228.     RETNCEL = 'COPIER'  RETN ;
  229.  
  230.     VITESSE PRES   = 'PRIM' 'PERFMONO'  RNCEL GNCEL RETNCEL GAMN ;
  231.  
  232.     ROF VITF PF GAMF   =  'PRET' 'PERFMONO'  1 1
  233.                            $MODE RNCEL VITESSE PRES GAMN ;
  234.  
  235.     CHPRES0 DT = 'KONV' 'VF' 'PERFMONO' 'RESI' NOMMET
  236.           $MODE ROF VITF PF GAMF LISTINCO ;
  237.  
  238.     VAL0 = 'EXTRAIRE' CHPRES0 PDUAL MOTDUA ;
  239.  
  240. *   EPSILON = perturbation
  241.  
  242. *   Adimensionalisation
  243.  
  244.     dens0 = 'EXTRAIRE' RN PPRIM 'SCAL' ;
  245.  
  246.     VN PN = 'PRIM' 'PERFMONO' RN GN RETN GAMN ;
  247.     CN2 = GAMN '*' (PN '/' RN) ;
  248.     cson0 = ('EXTRAIRE' CN2 PPRIM 'SCAL') '**' 0.5 ;
  249.  
  250.     ret0 = ('EXTRAIRE' RETN PPRIM 'SCAL') '**' 0.5 ;
  251.  
  252. *   On etabli la variable à perturber
  253.  
  254.     NDIM = 'DIME' LISTINCO ;
  255.     'REPETER' BL1 NDIM ;
  256.        MOTCEL = 'EXTRAIRE' LISTINCO &BL1 ;
  257.        'SI' ('EGA' MOTCEL MOTPRI) ;
  258.           ICEL = &BL1 ;
  259.           'QUITTER' BL1 ;
  260.        'FINSI' ;
  261.     'FIN' BL1 ;
  262.  
  263.     'SI' (ICEL > NDIM) ;
  264.        'MESSAGE' 'Procedure JACNUM' ;
  265.        'MESSAGE' 'MOTPRI = ??? ';
  266.        'ERREUR' 21 ;
  267.     'FINSI' ;
  268.  
  269.     ELT1 = 'MANUEL' 'POI1' PPRIM ;
  270.  
  271. *   ICEL = 1 -> On perturbe la densité
  272.  
  273.     'SI' ('EGA' ICEL 1) ;
  274.        DELTATOT = (EPSILON * dens0) ;
  275.        RNCEL = ('MANUEL' 'CHPO' ELT1 1 'SCAL' DELTATOT
  276.           'NATURE' 'DISCRET') 'ET' RN ;
  277.        GNCEL = 'COPIER' GN ;
  278.        RETNCEL = 'COPIER' RETN ;
  279.     'FINSI' ;
  280.  
  281. *   ICEL = 2 -> On perturbe la q.d.m. long  l'ax x
  282.  
  283.     'SI' ('EGA' ICEL 2) ;
  284.        DELTATOT = (EPSILON * dens0 * cson0) ;
  285.        GNCEL = ('MANUEL' 'CHPO' ELT1 1 'UX' DELTATOT
  286.           'NATURE' 'DISCRET') 'ET' GN ;
  287.        RNCEL = 'COPIER' RN ;
  288.        RETNCEL = 'COPIER' RETN ;
  289.     'FINSI' ;
  290.  
  291. *   ICEL = 3 -> On perturbe la q.d.m. long  l'ax y
  292.  
  293.     'SI' ('EGA' ICEL 3) ;
  294.        DELTATOT = (EPSILON * dens0 * cson0) ;
  295.        GNCEL = ('MANUEL' 'CHPO' ELT1 1 'UY' DELTATOT
  296.           'NATURE' 'DISCRET') 'ET' GN ;
  297.        RNCEL = 'COPIER' RN ;
  298.        RETNCEL = 'COPIER' RETN ;
  299.     'FINSI' ;
  300.  
  301. *   ICEL = 4 -> On perturbe l'énergie totale
  302.  
  303.     'SI' ('EGA' ICEL 4) ;
  304.        DELTATOT = (EPSILON * ret0) ;
  305.        RETNCEL = ('MANUEL' 'CHPO' ELT1 1 'SCAL' DELTATOT
  306.           'NATURE' 'DISCRET') 'ET' RETN ;
  307.        RNCEL = 'COPIER' RN ;
  308.        GNCEL = 'COPIER' GN ;
  309.     'FINSI' ;
  310.  
  311.     VITESSE PRES   = 'PRIM' 'PERFMONO'  RNCEL GNCEL RETNCEL GAMN ;
  312.  
  313.     ROF VITF PF GAMF   =  'PRET' 'PERFMONO'  1 1
  314.                            $MODE RNCEL VITESSE PRES GAMN ;
  315.  
  316.     CHPRES1 DT = 'KONV' 'VF' 'PERFMONO' 'RESI' NOMMET
  317.           $MODE ROF VITF PF GAMF LISTINCO ;
  318.  
  319.     VAL1 = 'EXTRAIRE' CHPRES1 PDUAL MOTDUA ;
  320.  
  321.  'FINPROC' ((VAL1 '-' VAL0) '/' DELTATOT)  ;
  322.  
  323. *****************************************************
  324. *****************************************************
  325. ******** FIN PROCEDURES *****************************
  326. *****************************************************
  327. *****************************************************
  328. *****************************************************
  329.  
  330.  
  331. *****************************************************
  332. ******* TEST1 ***************************************
  333. *****************************************************
  334. *
  335. * On compare le jacobien et la variation du residu
  336. * en $DOM5 . 'CENTRE' par rapport à une variation
  337. * infinitésimal en $DOM5 . 'CENTRE'
  338.  
  339.  
  340.  PCEN5 = 'POIN' 1 ('DOMA' $DOM5 'CENTRE') ;
  341.  
  342. * Les grandeurs pour adimesionner les erreurs
  343.  
  344.  ro0 = 'EXTRAIRE' RN0 PCEN5 'SCAL' ;
  345.  
  346.  VN0 PN0 = 'PRIM' 'PERFMONO' RN0 GN0 RETN0 GAMMAN ;
  347.  CN20 = GAMMAN '*' (PN0 '/' RN0) ;
  348.  cson0 = ('EXTRAIRE' CN20 PCEN5 'SCAL') '**' 0.5 ;
  349.  
  350.  ret0 = ('EXTRAIRE' RETN0 PCEN5 'SCAL') '**' 0.5 ;
  351.  
  352. *
  353. * Le jacobien exact. 
  354. * DRR  = d(RES_ro)/dro (variable primale en PCEN5, variable duale en PCEN5) ;
  355. * DGXR = d(RES_gx)/dro (variable primale en PCEN5, variable duale en PCEN5) ;
  356. * ...
  357.  
  358.  DRR   = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMDEN NOMDEN ;
  359.  DGXR  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMDEN NOMMOX ;
  360.  DGYR  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMDEN NOMMOY ;
  361.  DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMDEN NOMRET ;
  362.  
  363.  DRGX   = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOX NOMDEN ;
  364.  DGXGX  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOX NOMMOX ;
  365.  DGYGX  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOX NOMMOY ;
  366.  DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOX NOMRET ;
  367.  
  368.  DRGY   = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOY NOMDEN ;
  369.  DGXGY  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOY NOMMOX ;
  370.  DGYGY  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOY NOMMOY ;
  371.  DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMMOY NOMRET ;
  372.  
  373.  DRRET   = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMRET NOMDEN ;
  374.  DGXRET  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMRET NOMMOX ;
  375.  DGYRET  = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMRET NOMMOY ;
  376.  DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN5 PCEN5 NOMRET NOMRET ;
  377.  
  378. * Le jacobien numerique
  379.  
  380.  
  381.  DELTA = 1.0D-5 ;
  382.  
  383.  DRRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  384.    NOMDEN NOMDEN DELTA ;
  385.  DGXRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  386.    NOMDEN NOMMOX DELTA ;
  387.  DGYRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  388.    NOMDEN NOMMOY DELTA ;
  389.  DRETRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  390.    NOMDEN NOMRET DELTA ;
  391.  
  392.  DRGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  393.    NOMMOX NOMDEN DELTA ;
  394.  DGXGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  395.    NOMMOX NOMMOX DELTA ;
  396.  DGYGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  397.    NOMMOX NOMMOY DELTA ;
  398.  DRETGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  399.    NOMMOX NOMRET DELTA ;
  400.  
  401.  DRGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  402.    NOMMOY NOMDEN DELTA ;
  403.  DGXGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  404.    NOMMOY NOMMOX DELTA ;
  405.  DGYGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  406.    NOMMOY NOMMOY DELTA ;
  407.  DRETGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  408.    NOMMOY NOMRET DELTA ;
  409.  
  410.  DRRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  411.    NOMRET NOMDEN DELTA ;
  412.  DGXRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  413.    NOMRET NOMMOX DELTA ;
  414.  DGYRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5 PCEN5
  415.    NOMRET NOMMOY DELTA ;
  416.  DRETRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN5
  417.    PCEN5
  418.    NOMRET NOMRET DELTA ;
  419.  
  420. * Test des comparaisons jacobien exact-jacobien numerique
  421.  
  422.  ERR1 = ('ABS' (DRRN '-' DRR)) '*' (ro0 '/' (ro0 '*' cson0)) ;
  423.  'SI' (ERR1 > ERRTOL) ;
  424.       'MESSAGE' 'Problem 1';
  425.       'ERREUR' 5 ;
  426.  'FINSI' ;
  427.  
  428.  ERR1 = ('ABS' (DGXRN '-' DGXR)) '*'
  429.           (ro0 '/' (ro0 '*' cson0 '*' cson0)) ;
  430.  'SI' (ERR1 > ERRTOL) ;
  431.       'MESSAGE' 'Problem 2';
  432.       'ERREUR' 5 ;
  433.  'FINSI' ;
  434.  
  435.  ERR1 = ('ABS' (DGYRN '-' DGYR)) '*'
  436.           (ro0 '/' (ro0 '*' cson0 '*' cson0)) ;
  437.  'SI' (ERR1 > ERRTOL) ;
  438.       'MESSAGE' 'Problem 3';
  439.        'ERREUR' 5 ;
  440.  'FINSI' ;
  441.  
  442.  ERR1 = ('ABS' (DRETRN '-' DRETR)) '*'
  443.         (ro0 '/' (ret0 '*' cson0)) ; 
  444.  'SI' (ERR1 > ERRTOL) ;
  445.       'MESSAGE' 'Problem 4';
  446.       'ERREUR' 5 ;
  447.  'FINSI' ;
  448.  
  449.  
  450.  ERR1 = ('ABS' (DRGXN '-' DRGX))  ;
  451.  'SI' (ERR1 > ERRTOL) ;
  452.       'MESSAGE' 'Problem 5';
  453.       'ERREUR' 5 ;
  454.  'FINSI' ;
  455.  
  456.  ERR1 = ('ABS' (DGXGXN '-' DGXGX)) '/' cson0 ;
  457.  'SI' (ERR1 > ERRTOL) ;
  458.         'MESSAGE' 'Problem 6';
  459.         'ERREUR' 5 ;
  460.  'FINSI' ;
  461.  
  462.  ERR1 = ('ABS' (DGYGXN '-' DGYGX)) '/' cson0 ;
  463.  'SI' (ERR1 > ERRTOL) ;
  464.         'MESSAGE' 'Problem 7';
  465.         'ERREUR' 5 ;
  466.  'FINSI' ;
  467.  
  468.  
  469.  ERR1 = ('ABS' (DRETGXN '-' DRETGX)) '*'
  470.         (ro0  / ret0) ;
  471.  'SI' (ERR1 > ERRTOL) ;
  472.         'MESSAGE' 'Problem 8';
  473.         'ERREUR' 5 ;
  474.  'FINSI' ;
  475.  
  476.  ERR1 = ('ABS' (DRGYN '-' DRGY))  ;
  477.  'SI' (ERR1 > ERRTOL) ;
  478.         'MESSAGE' 'Problem 9';
  479.         'ERREUR' 5 ;
  480.  'FINSI' ;
  481.  
  482.  ERR1 = ('ABS' (DGXGYN '-' DGXGY))  '/' cson0 ;
  483.  'SI' (ERR1 > ERRTOL) ;
  484.         'MESSAGE' 'Problem 10';
  485.         'ERREUR' 5 ;
  486.  'FINSI' ;
  487.  
  488.  ERR1 = ('ABS' (DGYGYN '-' DGYGY)) '/' cson0 ;
  489.  'SI' (ERR1 > ERRTOL) ;
  490.         'MESSAGE' 'Problem 11';
  491.         'ERREUR' 5 ;
  492.  'FINSI' ;
  493.  
  494.  ERR1 = ('ABS' (DRETGYN '-' DRETGY))  '*'
  495.         (ro0 / ret0)  ;
  496.  'SI' (ERR1 > ERRTOL) ;
  497.         'MESSAGE' 'Problem 12';
  498.         'ERREUR' 5 ;
  499.  'FINSI' ;
  500.  
  501.  ERR1 = ('ABS' (DRRETN '-' DRRET)) '*' (ret0 '/' (ro0 '*' cson0))   ;
  502.  'SI' (ERR1 > ERRTOL) ;
  503.         'MESSAGE' 'Problem 13';
  504.         'ERREUR' 5 ;
  505.  'FINSI' ;
  506.  
  507.  ERR1 = ('ABS' (DGXRETN '-' DGXRET)) '*'
  508.           (ret0 '/' (ro0 '*' cson0 * cson0))   ;
  509.  'SI' (ERR1 > ERRTOL) ;
  510.         'MESSAGE' 'Problem 14';
  511.         'ERREUR' 5 ;
  512.  'FINSI' ;
  513.  
  514.  ERR1 = ('ABS' (DGYRETN '-' DGYRET)) '*' 
  515.              (ret0 '/' (ro0 '*' cson0 * cson0)) ;
  516.  'SI' (ERR1 > ERRTOL) ;
  517.         'MESSAGE' 'Problem 15';
  518.         'ERREUR' 5 ;
  519.  'FINSI' ;
  520.  
  521.  
  522.  ERR1 = ('ABS' (DRETRETN '-' DRETRET)) '/' cson0 ;
  523.  'SI' (ERR1 > ERRTOL) ;
  524.         'MESSAGE' 'Problem 16';
  525.         'ERREUR' 5 ;
  526.  'FINSI' ;
  527.  
  528. *****************************************************
  529. *****************************************************
  530. ******* TEST2 ***************************************
  531. *****************************************************
  532. *****************************************************
  533. *
  534. * On compare le jacobien et la variation du residu
  535. * en $DOM5 'CENTRE' par rapport à une variation
  536. * infinitésimal en $DOM6 'CENTRE'
  537. *
  538.  
  539.  PCEN6 = ('DOMA' $DOM6 'CENTRE') 'POIN' 1 ;
  540.  
  541.  ro0 = 'EXTRAIRE' RN0 PCEN6 'SCAL' ;
  542.  cson0 = ('EXTRAIRE' CN20 PCEN6 'SCAL') '**' 0.5 ;
  543.  ret0 = 'EXTRAIRE' RETN0 PCEN6 'SCAL' ;
  544.  
  545. *
  546. * Le jacobien exact. 
  547. *
  548.  
  549.  DRR   = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMDEN NOMDEN ;
  550.  DGXR  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMDEN NOMMOX ;
  551.  DGYR  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMDEN NOMMOY ;
  552.  DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMDEN NOMRET ;
  553.  
  554.  DRGX   = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMMOX NOMDEN ;
  555.  DGXGX  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMMOX NOMMOX ;
  556.  DGYGX  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMMOX NOMMOY ;
  557.  DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMMOX NOMRET ;
  558.  
  559.  DRGY   = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMMOY NOMDEN ;
  560.  DGXGY  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMMOY NOMMOX ;
  561.  DGYGY  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMMOY NOMMOY ;
  562.  DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMMOY NOMRET ;
  563.  
  564.  DRRET   = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMRET NOMDEN ;
  565.  DGXRET  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMRET NOMMOX ;
  566.  DGYRET  = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMRET NOMMOY ;
  567.  DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN5 NOMRET NOMRET ;
  568.  
  569. * Le jacobien numerique
  570.  
  571.  
  572.  DELTA = 1.0D-5 ;
  573.  
  574.  DRRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  575.    NOMDEN NOMDEN DELTA ;
  576.  DGXRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  577.    NOMDEN NOMMOX DELTA ;
  578.  DGYRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  579.    NOMDEN NOMMOY DELTA ;
  580.  DRETRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  581.    NOMDEN NOMRET DELTA ;
  582.  
  583.  DRGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  584.    NOMMOX NOMDEN DELTA ;
  585.  DGXGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  586.    NOMMOX NOMMOX DELTA ;
  587.  DGYGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  588.    NOMMOX NOMMOY DELTA ;
  589.  DRETGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  590.    NOMMOX NOMRET DELTA ;
  591.  
  592.  DRGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  593.    NOMMOY NOMDEN DELTA ;
  594.  DGXGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  595.    NOMMOY NOMMOX DELTA ;
  596.  DGYGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  597.    NOMMOY NOMMOY DELTA ;
  598.  DRETGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  599.    NOMMOY NOMRET DELTA ;
  600.  
  601.  DRRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  602.    NOMRET NOMDEN DELTA ;
  603.  DGXRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  604.    NOMRET NOMMOX DELTA ;
  605.  DGYRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN5
  606.    NOMRET NOMMOY DELTA ;
  607.  DRETRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6
  608.    PCEN5
  609.    NOMRET NOMRET DELTA ;
  610.  
  611. * Test des comparaisons jacobien exact-jacobien numerique
  612.  
  613.  ERR1 = ('ABS' (DRRN '-' DRR)) '*' (ro0 '/' (ro0 '*' cson0)) ;
  614.  'SI' (ERR1 > ERRTOL) ;
  615.       'MESSAGE' 'Problem 1';
  616.       'ERREUR' 5 ;
  617.  'FINSI' ;
  618.  
  619.  ERR1 = ('ABS' (DGXRN '-' DGXR)) '*'
  620.           (ro0 '/' (ro0 '*' cson0 '*' cson0)) ;
  621.  'SI' (ERR1 > ERRTOL) ;
  622.       'MESSAGE' 'Problem 2';
  623.       'ERREUR' 5 ;
  624.  'FINSI' ;
  625.  
  626.  ERR1 = ('ABS' (DGYRN '-' DGYR)) '*'
  627.           (ro0 '/' (ro0 '*' cson0 '*' cson0)) ;
  628.  'SI' (ERR1 > ERRTOL) ;
  629.       'MESSAGE' 'Problem 3';
  630.        'ERREUR' 5 ;
  631.  'FINSI' ;
  632.  
  633.  ERR1 = ('ABS' (DRETRN '-' DRETR)) '*'
  634.         (ro0 '/' (ret0 '*' cson0)) ; 
  635.  'SI' (ERR1 > ERRTOL) ;
  636.       'MESSAGE' 'Problem 4';
  637.       'ERREUR' 5 ;
  638.  'FINSI' ;
  639.  
  640.  
  641.  ERR1 = ('ABS' (DRGXN '-' DRGX))  ;
  642.  'SI' (ERR1 > ERRTOL) ;
  643.       'MESSAGE' 'Problem 5';
  644.       'ERREUR' 5 ;
  645.  'FINSI' ;
  646.  
  647.  ERR1 = ('ABS' (DGXGXN '-' DGXGX)) '/' cson0 ;
  648.  'SI' (ERR1 > ERRTOL) ;
  649.         'MESSAGE' 'Problem 6';
  650.         'ERREUR' 5 ;
  651.  'FINSI' ;
  652.  
  653.  ERR1 = ('ABS' (DGYGXN '-' DGYGX)) '/' cson0 ;
  654.  'SI' (ERR1 > ERRTOL) ;
  655.         'MESSAGE' 'Problem 7';
  656.         'ERREUR' 5 ;
  657.  'FINSI' ;
  658.  
  659.  
  660.  ERR1 = ('ABS' (DRETGXN '-' DRETGX)) '*'
  661.         (ro0  / ret0) ;
  662.  'SI' (ERR1 > ERRTOL) ;
  663.         'MESSAGE' 'Problem 8';
  664.         'ERREUR' 5 ;
  665.  'FINSI' ;
  666.  
  667.  ERR1 = ('ABS' (DRGYN '-' DRGY))  ;
  668.  'SI' (ERR1 > ERRTOL) ;
  669.         'MESSAGE' 'Problem 9';
  670.         'ERREUR' 5 ;
  671.  'FINSI' ;
  672.  
  673.  ERR1 = ('ABS' (DGXGYN '-' DGXGY))  '/' cson0 ;
  674.  'SI' (ERR1 > ERRTOL) ;
  675.         'MESSAGE' 'Problem 10';
  676.         'ERREUR' 5 ;
  677.  'FINSI' ;
  678.  
  679.  ERR1 = ('ABS' (DGYGYN '-' DGYGY)) '/' cson0 ;
  680.  'SI' (ERR1 > ERRTOL) ;
  681.         'MESSAGE' 'Problem 11';
  682.         'ERREUR' 5 ;
  683.  'FINSI' ;
  684.  
  685.  ERR1 = ('ABS' (DRETGYN '-' DRETGY))  '*'
  686.         (ro0 / ret0)  ;
  687.  'SI' (ERR1 > ERRTOL) ;
  688.         'MESSAGE' 'Problem 12';
  689.         'ERREUR' 5 ;
  690.  'FINSI' ;
  691.  
  692.  ERR1 = ('ABS' (DRRETN '-' DRRET)) '*' (ret0 '/' (ro0 '*' cson0))   ;
  693.  'SI' (ERR1 > ERRTOL) ;
  694.         'MESSAGE' 'Problem 13';
  695.         'ERREUR' 5 ;
  696.  'FINSI' ;
  697.  
  698.  ERR1 = ('ABS' (DGXRETN '-' DGXRET)) '*'
  699.           (ret0 '/' (ro0 '*' cson0 * cson0))   ;
  700.  'SI' (ERR1 > ERRTOL) ;
  701.         'MESSAGE' 'Problem 14';
  702.         'ERREUR' 5 ;
  703.  'FINSI' ;
  704.  
  705.  ERR1 = ('ABS' (DGYRETN '-' DGYRET)) '*' 
  706.              (ret0 '/' (ro0 '*' cson0 * cson0)) ;
  707.  'SI' (ERR1 > ERRTOL) ;
  708.         'MESSAGE' 'Problem 15';
  709.         'ERREUR' 5 ;
  710.  'FINSI' ;
  711.  
  712.  
  713.  ERR1 = ('ABS' (DRETRETN '-' DRETRET)) '/' cson0 ;
  714.  'SI' (ERR1 > ERRTOL) ;
  715.         'MESSAGE' 'Problem 16';
  716.         'ERREUR' 5 ;
  717.  'FINSI' ;
  718.  
  719.  
  720. *****************************************************
  721. *****************************************************
  722. ******* TEST3 ***************************************
  723. *****************************************************
  724. *****************************************************
  725. *
  726. * On observe que la variation du residu
  727. * en $DOM5 'CENTRE' par rapport à une variation
  728. * infinitésimal en $DOM7 'CENTRE' doit etre nulle.
  729. * Donc les elements du jacobien relatifs à ça
  730. * doivent etre nuls.
  731. *
  732.  
  733.  PCEN7 = ('DOMA' $DOM7 'CENTRE') 'POIN' 1 ;
  734.  
  735.  ro0 = 'EXTRAIRE' RN0 PCEN7 'SCAL' ;
  736.  cson0 = ('EXTRAIRE' CN20 PCEN7 'SCAL') '**' 0.5 ;
  737.  ret0 = 'EXTRAIRE' RETN0 PCEN7 'SCAL' ;
  738.  
  739. *
  740. * Le jacobien exact. 
  741. *
  742.  
  743.  DRR   = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMDEN NOMDEN ;
  744.  DGXR  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMDEN NOMMOX ;
  745.  DGYR  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMDEN NOMMOY ;
  746.  DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMDEN NOMRET ;
  747.  
  748.  DRGX   = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMMOX NOMDEN ;
  749.  DGXGX  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMMOX NOMMOX ;
  750.  DGYGX  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMMOX NOMMOY ;
  751.  DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMMOX NOMRET ;
  752.  
  753.  DRGY   = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMMOY NOMDEN ;
  754.  DGXGY  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMMOY NOMMOX ;
  755.  DGYGY  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMMOY NOMMOY ;
  756.  DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMMOY NOMRET ;
  757.  
  758.  DRRET   = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMRET NOMDEN ;
  759.  DGXRET  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMRET NOMMOX ;
  760.  DGYRET  = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMRET NOMMOY ;
  761.  DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN5 NOMRET NOMRET ;
  762.  
  763.  
  764.  'SI' (('MAXIMUM'
  765.     ('PROG'  DRR  DGXR  DGYR  DRETR
  766.              DRGX DGXGX DGYGX DRETGX
  767.              DRGY DGXGY DGYGY DRETGY
  768.              DRRET DGXRET DGYRET DRETRET ) 'ABS' )
  769.       > ERRTOL) ;
  770.       'MESSAGE' 'Problem final' ;
  771.       'ERREUR' 5 ;
  772.  'FINSI' ;    
  773.  
  774.  'FIN' ;
  775.  
  776.  
  777.  
  778.  
  779.  
  780.  
  781.  
  782.  
  783.  

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