Télécharger konv_impl_centre.dgibi

Retour à la liste

Numérotation des lignes :

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

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