Télécharger la plaquette Cast3M


fr en

Télécharger konv_impl3.dgibi

Retour à la liste

Numérotation des lignes : 

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

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