Télécharger konv_impl2.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : konv_impl2.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: AUSMPLUS ****
  15. **** ****
  16. **** A. BECCANTINI DRN/DMT/SEMT/LTMF AOUT 2000 ****
  17. **** A. KUDRIAKOV DRN/DMT/SEMT/LTMF AOUT 2000 ****
  18. ***********************************************************
  19. ***********************************************************
  20.  
  21. 'OPTION' 'DIME' 2 ;
  22. 'OPTION' 'ELEM' QUA4 ;
  23. 'OPTION' 'ECHO' 0 ;
  24. 'OPTION' 'TRAC' 'X' ;
  25.  
  26. *
  27. *** GRAPH
  28. *
  29.  
  30. GRAPH = FAUX ;
  31. * GRAPH = VRAI ;
  32.  
  33. ERRTOL = 1.0D-3 ;
  34.  
  35. ***************************
  36. ***** DOMAINE SPATIAL ****
  37. ***************************
  38.  
  39.  
  40. A0 = 0.0D0 0.0D0;
  41. A1 = 1.0D0 0.0D0;
  42. A2 = 2.0D0 0.0D0;
  43. A3 = 3.0D0 0.0D0;
  44.  
  45. A0A1 = A0 'DROIT' 1 A1;
  46. A1A2 = A1 'DROIT' 1 A2;
  47. A2A3 = A2 'DROIT' 1 A3;
  48.  
  49.  
  50. DOM1 = 'TRANSLATION' A0A1 1 (0.0 1.0) ;
  51. DOM2 = DOM1 'PLUS' (0.0 1.0) ;
  52. DOM3 = DOM2 'PLUS' (0.0 1.0) ;
  53. DOM4 = 'TRANSLATION' A1A2 1 (0.0 1.0) ;
  54. DOM5 = DOM4 'PLUS' (0.0 1.0) ;
  55. DOM6 = DOM5 'PLUS' (0.0 1.0) ;
  56. DOM7 = 'TRANSLATION' A2A3 1 (0.0 1.0) ;
  57. DOM8 = DOM7 'PLUS' (0.0 1.0) ;
  58. DOM9 = DOM8 'PLUS' (0.0 1.0) ;
  59.  
  60.  
  61. DOMTOT = DOM1 'ET' DOM2 'ET' DOM3 'ET' DOM4 'ET' DOM5 'ET'
  62. DOM6 'ET' DOM7 'ET' DOM8 'ET' DOM9 'ELIMINATION' 0.0001 ;
  63.  
  64. $DOMTOT = 'MODELISER' DOMTOT 'EULER';
  65.  
  66. $DOM6 = 'MODELISER' DOM6 'EULER';
  67. $DOM7 = 'MODELISER' DOM7 'EULER';
  68. $DOM9 = 'MODELISER' DOM9 'EULER';
  69.  
  70. TDOMTOT = 'DOMA' $DOMTOT 'VF';
  71.  
  72. TDOM6 = 'DOMA' $DOM6 'VF';
  73. TDOM7 = 'DOMA' $DOM7 'VF';
  74. TDOM9 = 'DOMA' $DOM9 'VF';
  75.  
  76. MDOM6 = TDOM6 . 'QUAF' ;
  77. MDOM7 = TDOM7 . 'QUAF' ;
  78. MDOM9 = TDOM9 . '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 MDOM9) 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 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' 1.4 ;
  95. CSONN = (GAMMAN '*' PN0) '/' RN0 ;
  96. CSONN = 'KOPS' CSONN '**' 0.5 ;
  97. UXN0 = 0.2 * CSONN ;
  98. UYN0 = 0.3 * CSONN ;
  99.  
  100. GN0 = 'KCHT' $DOMTOT 'VECT' 'CENTRE' (
  101. ('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 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE'
  106. ((PN0 '/' (GAMMAN '-' 1.0)) '+' ECIN) ;
  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. 'AUSMPLUS' 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 $DOM9 'CENTRE' par rapport à une variation
  337. * infinitésimal en $DOM9 'CENTRE'
  338. *
  339.  
  340. PCEN9 = 'POIN' 1 ('DOMA' $DOM9 'CENTRE') ;
  341.  
  342. * Les grandeurs pour adimesionner les erreurs
  343.  
  344. ro0 = 'EXTRAIRE' RN0 PCEN9 'SCAL' ;
  345.  
  346. VN0 PN0 = 'PRIM' 'PERFMONO' RN0 GN0 RETN0 GAMMAN ;
  347. CN20 = GAMMAN '*' (PN0 '/' RN0) ;
  348. cson0 = ('EXTRAIRE' CN20 PCEN9 'SCAL') '**' 0.5 ;
  349.  
  350. ret0 = ('EXTRAIRE' RETN0 PCEN9 'SCAL') '**' 0.5 ;
  351.  
  352. *
  353. * Le jacobien exact.
  354. * DRR = d(RES_ro)/dro (variable primale en PCEN9, variable duale en PCEN9) ;
  355. * DGXR = d(RES_gx)/dro (variable primale en PCEN9, variable duale en PCEN9) ;
  356. * ...
  357.  
  358. DRR = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMDEN NOMDEN ;
  359. DGXR = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMDEN NOMMOX ;
  360. DGYR = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMDEN NOMMOY ;
  361. DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMDEN NOMRET ;
  362.  
  363. DRGX = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOX NOMDEN ;
  364. DGXGX = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOX NOMMOX ;
  365. DGYGX = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOX NOMMOY ;
  366. DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOX NOMRET ;
  367.  
  368. DRGY = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOY NOMDEN ;
  369. DGXGY = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOY NOMMOX ;
  370. DGYGY = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOY NOMMOY ;
  371. DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMMOY NOMRET ;
  372.  
  373. DRRET = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMRET NOMDEN ;
  374. DGXRET = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMRET NOMMOX ;
  375. DGYRET = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 NOMRET NOMMOY ;
  376. DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN9 PCEN9 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 PCEN9 PCEN9
  384. NOMDEN NOMDEN DELTA ;
  385. DGXRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  386. NOMDEN NOMMOX DELTA ;
  387. DGYRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  388. NOMDEN NOMMOY DELTA ;
  389. DRETRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  390. NOMDEN NOMRET DELTA ;
  391.  
  392. DRGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  393. NOMMOX NOMDEN DELTA ;
  394. DGXGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  395. NOMMOX NOMMOX DELTA ;
  396. DGYGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  397. NOMMOX NOMMOY DELTA ;
  398. DRETGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  399. NOMMOX NOMRET DELTA ;
  400.  
  401. DRGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  402. NOMMOY NOMDEN DELTA ;
  403. DGXGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  404. NOMMOY NOMMOX DELTA ;
  405. DGYGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  406. NOMMOY NOMMOY DELTA ;
  407. DRETGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  408. NOMMOY NOMRET DELTA ;
  409.  
  410. DRRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  411. NOMRET NOMDEN DELTA ;
  412. DGXRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  413. NOMRET NOMMOX DELTA ;
  414. DGYRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9 PCEN9
  415. NOMRET NOMMOY DELTA ;
  416. DRETRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN9
  417. PCEN9
  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 $DOM9 '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 PCEN9 NOMDEN NOMDEN ;
  550. DGXR = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMDEN NOMMOX ;
  551. DGYR = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMDEN NOMMOY ;
  552. DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMDEN NOMRET ;
  553.  
  554. DRGX = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOX NOMDEN ;
  555. DGXGX = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOX NOMMOX ;
  556. DGYGX = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOX NOMMOY ;
  557. DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOX NOMRET ;
  558.  
  559. DRGY = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOY NOMDEN ;
  560. DGXGY = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOY NOMMOX ;
  561. DGYGY = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOY NOMMOY ;
  562. DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMMOY NOMRET ;
  563.  
  564. DRRET = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMRET NOMDEN ;
  565. DGXRET = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMRET NOMMOX ;
  566. DGYRET = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 NOMRET NOMMOY ;
  567. DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN6 PCEN9 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 PCEN9
  575. NOMDEN NOMDEN DELTA ;
  576. DGXRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  577. NOMDEN NOMMOX DELTA ;
  578. DGYRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  579. NOMDEN NOMMOY DELTA ;
  580. DRETRN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  581. NOMDEN NOMRET DELTA ;
  582.  
  583. DRGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  584. NOMMOX NOMDEN DELTA ;
  585. DGXGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  586. NOMMOX NOMMOX DELTA ;
  587. DGYGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  588. NOMMOX NOMMOY DELTA ;
  589. DRETGXN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  590. NOMMOX NOMRET DELTA ;
  591.  
  592. DRGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  593. NOMMOY NOMDEN DELTA ;
  594. DGXGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  595. NOMMOY NOMMOX DELTA ;
  596. DGYGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  597. NOMMOY NOMMOY DELTA ;
  598. DRETGYN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  599. NOMMOY NOMRET DELTA ;
  600.  
  601. DRRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  602. NOMRET NOMDEN DELTA ;
  603. DGXRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  604. NOMRET NOMMOX DELTA ;
  605. DGYRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6 PCEN9
  606. NOMRET NOMMOY DELTA ;
  607. DRETRETN = JACNUM $DOMTOT METO RN0 GN0 RETN0 GAMMAN LISTINCO PCEN6
  608. PCEN9
  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 $DOM9 'CENTRE' par rapport à une variation
  728. * infinitésimal en $DOM7 'CENTRE' doit etre nulle.
  729. *
  730.  
  731. PCEN7 = ('DOMA' $DOM7 'CENTRE') 'POIN' 1 ;
  732.  
  733. ro0 = 'EXTRAIRE' RN0 PCEN7 'SCAL' ;
  734. cson0 = ('EXTRAIRE' CN20 PCEN7 'SCAL') '**' 0.5 ;
  735. ret0 = 'EXTRAIRE' RETN0 PCEN7 'SCAL' ;
  736.  
  737. *
  738. * Le jacobien exact.
  739. *
  740.  
  741. DRR = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMDEN NOMDEN ;
  742. DGXR = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMDEN NOMMOX ;
  743. DGYR = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMDEN NOMMOY ;
  744. DRETR = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMDEN NOMRET ;
  745.  
  746. DRGX = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOX NOMDEN ;
  747. DGXGX = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOX NOMMOX ;
  748. DGYGX = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOX NOMMOY ;
  749. DRETGX = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOX NOMRET ;
  750.  
  751. DRGY = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOY NOMDEN ;
  752. DGXGY = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOY NOMMOX ;
  753. DGYGY = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOY NOMMOY ;
  754. DRETGY = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMMOY NOMRET ;
  755.  
  756. DRRET = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMRET NOMDEN ;
  757. DGXRET = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMRET NOMMOX ;
  758. DGYRET = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMRET NOMMOY ;
  759. DRETRET = JACOVA JACO $DOMTOT LISTINCO PCEN7 PCEN9 NOMRET NOMRET ;
  760.  
  761.  
  762. 'SI' (('MAXIMUM'
  763. ('PROG' DRR DGXR DGYR DRETR
  764. DRGX DGXGX DGYGX DRETGX
  765. DRGY DGXGY DGYGY DRETGY
  766. DRRET DGXRET DGYRET DRETRET ) 'ABS' )
  767. > ERRTOL) ;
  768. 'MESSAGE' 'Problem final' ;
  769. 'ERREUR' 5 ;
  770. 'FINSI' ;
  771.  
  772. *****************************************************
  773. *****************************************************
  774. ******* TEST4 ***************************************
  775. *****************************************************
  776. *****************************************************
  777. * Primitive variables
  778. * Redefinition of VN0
  779. *
  780.  
  781. RN0 = 'BRUI' 'BLAN' 'UNIF' 1.11 0.5 ('DOMA' $DOMTOT 'CENTRE') ;
  782. PN0 = 'BRUI' 'BLAN' 'UNIF' 1234.1 800 ('DOMA' $DOMTOT 'CENTRE') ;
  783. GAMMAN = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 1 'SCAL' 1.4
  784. 'NATU' 'DISCRET' ;
  785. CN0 = (GAMMAN '*' PN0) '/' RN0 ;
  786. CN0 = 'KOPS' CN0 '**' 0.5 ;
  787. VN0 = ('NOMC' (CN0 * 0.11) 'UX' 'NATU' 'DISCRET') 'ET'
  788. ('NOMC' (CN0 * 2.11) 'UY' 'NATU' 'DISCRET') ;
  789.  
  790. LISTCONS = 'MOTS' 'RN' 'RUX' 'RUY' 'RETN' ;
  791.  
  792. JACO = 'KONV' 'VF' 'PERFMONO' 'JACOCONS'
  793. $DOMTOT LISTCONS 'AUSMPLUS' RN0 VN0 PN0 GAMMAN ;
  794.  
  795. LISTPRIM = 'MOTS' 'RN' 'UX' 'UY' 'PN' ;
  796.  
  797. DUDV = 'PRIM' 'CONSPRIM' ('DOMA' $DOMTOT 'CENTRE')
  798. LISTCONS LISTPRIM RN0 VN0 PN0 GAMMAN ;
  799.  
  800. JACOP = 'KONV' 'VF' 'PERFMONO' 'JACOPRIM'
  801. $DOMTOT LISTCONS LISTPRIM 'AUSMPLUS' RN0 VN0 PN0 GAMMAN ;
  802.  
  803. UNP = 'MANUEL' 'CHPO' ('DOMA' $DOMTOT 'CENTRE') 4 'RN' 0.11
  804. 'UX' 0.15 'UY' 0.19 'PN' 0.21 ;
  805.  
  806. * JACOP * UNP = JACO * DUDV * UNP, \forall UNP
  807.  
  808. ERRO = 'MAXIMUM' (
  809. ('KOPS' 'MULT' JACOP UNP) '-'
  810. ('KOPS' 'MULT' JACO ('KOPS' 'MULT' DUDV UNP))
  811. ) 'ABS' ;
  812.  
  813. 'SI' (ERRO > ERRTOL) ;
  814. 'MESSAGE' 'Problem final' ;
  815. 'ERREUR' 5 ;
  816. 'FINSI' ;
  817.  
  818. *****************************************************
  819. *****************************************************
  820. ******* TEST5 ***************************************
  821. *****************************************************
  822. *****************************************************
  823. *
  824. * We compute the flux contribution nowhere
  825. *
  826.  
  827. JACO = 'KONV' 'VF' 'PERFMONO' 'JACOCONS'
  828. $DOMTOT LISTCONS ('DOMA' $DOMTOT 'FACE')
  829. 'AUSMPLUS' RN0 VN0 PN0 GAMMAN ;
  830.  
  831. JACOP = 'KONV' 'VF' 'PERFMONO' 'JACOPRIM'
  832. $DOMTOT LISTCONS LISTPRIM ('DOMA' $DOMTOT 'FACE')
  833. 'AUSMPLUS' RN0 VN0 PN0 GAMMAN ;
  834.  
  835. ERRO = 'MAXIMUM'
  836. ('KOPS' 'MULT' JACO ('KOPS' 'MULT' DUDV UNP))
  837. 'ABS' ;
  838.  
  839. 'SI' (ERRO > ERRTOL) ;
  840. 'MESSAGE' 'Problem final, mais tres final' ;
  841. 'ERREUR' 5 ;
  842. 'FINSI' ;
  843.  
  844. ERRO = 'MAXIMUM'
  845. (('KOPS' 'MULT' JACOP UNP))
  846. 'ABS' ;
  847.  
  848. 'SI' (ERRO > ERRTOL) ;
  849. 'MESSAGE' 'Problem final, mais tres final' ;
  850. 'ERREUR' 5 ;
  851. 'FINSI' ;
  852.  
  853. 'FIN' ;
  854.  
  855.  
  856.  
  857.  
  858.  
  859.  
  860.  
  861.  
  862.  

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