Télécharger konv_resi_ther_cons.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : konv_resi_ther_cons.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. ***********************************************************
  5. **** APPROCHE VF "Cell-Centred Formulation" pour la ****
  6. **** solution des ****
  7. **** Equations d'Euler pour un gaz parfait ****
  8. **** OPERATEURS PRET, KONV ****
  9. **** ****
  10. **** Cas gaz monoespece "thermally perfect" ****
  11. **** Consistence ****
  12. **** ****
  13. **** Methodes: VLH/SS ****
  14. **** ****
  15. **** A. BECCANTINI DM2S/SEMT/LTMF AVRIL 2001 ****
  16. ***********************************************************
  17.  
  18.  
  19. 'OPTION' 'DIME' 2 ;
  20. 'OPTION' 'ELEM' QUA4 ;
  21. 'OPTION' 'ECHO' 0 ;
  22. 'OPTION' 'TRAC' 'X' ;
  23.  
  24. *
  25. *** GRAPH
  26. *
  27.  
  28. GRAPH = FAUX ;
  29. * GRAPH = VRAI ;
  30.  
  31. **** Cas monoespece : la table proprieté de gaz
  32.  
  33. PGAZ = 'TABLE' ;
  34. PGAZ2 = 'TABLE' ;
  35.  
  36. *
  37. **** Ordre des polynoms cv_i
  38. *
  39.  
  40. PGAZ . 'NORD' = 4 ;
  41. PGAZ2 . 'NORD' = 4 ;
  42.  
  43. *
  44. **** Especes qui sont dans les equations d'Euler (rien)
  45. *
  46.  
  47. *
  48. **** Espece qui n'y est pas
  49. *
  50.  
  51.  
  52. PGAZ . 'ESPNEULE' = 'H2 ';
  53. PGAZ2 . 'ESPNEULE' = 'H2 ';
  54.  
  55. *
  56.  
  57. PGAZ . 'H2 ' = 'TABLE' ;
  58. PGAZ2 . 'H2 ' = 'TABLE' ;
  59.  
  60. *
  61. **** R (J/Kg/K)
  62. *
  63.  
  64. PGAZ . 'H2 ' . 'R' = 4130.0 ;
  65. PGAZ2 . 'H2 ' . 'R' = 4130.0 ;
  66.  
  67. *
  68. **** Regressions polynomials
  69. *
  70.  
  71. PGAZ . 'H2 ' . 'A' = 'PROG' 9834.91866 0.54273926 0.000862203836
  72. -2.37281455E-07 1.84701105E-11 ;
  73. PGAZ2 . 'H2 ' . 'A' = 'PROG' 9834.91866 0.54273926 0.000862203836
  74. -2.37281455E-07 1.84701105E-11 ;
  75.  
  76. *
  77. **** "Enthalpies" (ou energies) de formations a OK (J/Kg)
  78. * Note: ce sont des entites numeriques
  79. *
  80.  
  81. PGAZ . 'H2 ' . 'H0K' = -4.195D6 ;
  82. PGAZ2 . 'H2 ' . 'H0K' = -4.195D6 ;
  83.  
  84.  
  85. * PGAZ2: transport de scalaires passifs
  86.  
  87. PGAZ2 . 'SCALPASS' = 'MOTS' 'TUTU' 'TATA' ;
  88.  
  89. ***************************
  90. ***** DOMAINE SPATIAL ****
  91. ***************************
  92.  
  93.  
  94. A1 = 0.0D0 0.0D0;
  95. A2 = 1.0D0 0.0D0;
  96. A3 = 2.0D0 0.0D0;
  97. A4 = 2.0D0 1.0D0;
  98. A5 = 1.0D0 1.0D0;
  99. A6 = 0.0D0 1.0D0;
  100.  
  101. L12 = A1 'DROIT' 1 A2;
  102. L23 = A2 'DROIT' 1 A3;
  103. L34 = A3 'DROIT' 1 A4;
  104. L45 = A4 'DROIT' 1 A5;
  105. L56 = A5 'DROIT' 1 A6;
  106. L61 = A6 'DROIT' 1 A1;
  107. L25 = A2 'DROIT' 1 A5;
  108.  
  109.  
  110. DOM10 = 'DALL' L12 L25 L56 L61
  111. 'PLANE';
  112. DOM20 = 'DALL' L23 L34 L45 ('INVERSE' L25)
  113. 'PLANE';
  114. *
  115. *** Point ou on controlle la consistence
  116. *
  117.  
  118. P10 = 1.0 0.5;
  119.  
  120. *
  121. *** Etats gauche et droite
  122. *
  123.  
  124. rogd = 1.412 ;
  125. pgd = 1.17D6 ;
  126. ungd = 8.34 ;
  127. utgd = 18.2 ;
  128. tutugd = 14.9 ;
  129. tatagd = 23.8 ;
  130.  
  131. *
  132. **** tg
  133. *
  134.  
  135. tgd = pgd '/' (rogd * (PGAZ . 'H2 ' . 'R')) ;
  136.  
  137.  
  138. *
  139. **** gamg, gamd , htg , htd
  140. *
  141. *
  142.  
  143. A0H2 = 'EXTRAIRE' (PGAZ . 'H2 ' . 'A') 1 ;
  144. A1H2 = 'EXTRAIRE' (PGAZ . 'H2 ' . 'A') 2 ;
  145. A2H2 = 'EXTRAIRE' (PGAZ . 'H2 ' . 'A') 3 ;
  146. A3H2 = 'EXTRAIRE' (PGAZ . 'H2 ' . 'A') 4 ;
  147. A4H2 = 'EXTRAIRE' (PGAZ . 'H2 ' . 'A') 5 ;
  148.  
  149. T = tgd ;
  150. T2 = T * T ;
  151. T3 = T2 * T ;
  152. T4 = T3 * T ;
  153. T5 = T4 * T ;
  154.  
  155. cvgd = A0H2 '+' (A1H2 * T) '+' (A2H2 * T2) '+' (A3H2 * T3) '+'
  156. (A4H2 * T4) ;
  157.  
  158. ETGd = (A0H2 * T) '+' (A1H2 * T2 '/' 2) '+' (A2H2 * T3 '/' 3) '+'
  159. (A3H2 * T4 '/' 4) '+' (A4H2 * T5 '/' 5) ;
  160.  
  161. ecingd = 0.5D0 '*' ((ungd '*' ungd) '+' (utgd '*' utgd));
  162.  
  163. retgd = rogd '*' (ETGd '+' ecingd) ;
  164.  
  165. gamgd = (cvgd '+' (PGAZ . 'H2 ' . 'R')) '/' cvgd ;
  166.  
  167. *
  168. *** flux en (n,t)
  169. *
  170.  
  171. f1gd = ungd '*' rogd ;
  172. f2gd = (rogd '*' (ungd '*' ungd)) '+' pgd ;
  173. f3gd = rogd '*' ungd '*' utgd ;
  174. f4gd = ungd '*' (retgd '+' pgd);
  175. ftutugd = ungd '*' rogd '*' tutugd ;
  176. ftatagd = ungd '*' rogd '*' tatagd ;
  177.  
  178.  
  179.  
  180. ****************************************************
  181. ****************************************************
  182. ******** Boucle sur les angles *********
  183. ****************************************************
  184. ****************************************************
  185.  
  186. DANGLE = 360 '/' 7.15;
  187. ANGLE = 11.3 ;
  188.  
  189. 'REPETER' BLOC 8;
  190.  
  191. *
  192. *** Rotation
  193. *
  194.  
  195.  
  196. ANGLE = ANGLE '+' DANGLE;
  197. ORIG = 0.0D0 0.0D0;
  198.  
  199. 'MESSAGE' ;
  200. 'MESSAGE' (CHAIN 'Angle de rotation= ' ANGLE);
  201. 'MESSAGE' ;
  202.  
  203. DOM1 = DOM10 'TOURNER' ANGLE ORIG;
  204. DOM2 = DOM20 'TOURNER' ANGLE ORIG;
  205. P1 = P10 'TOURNER' ANGLE ORIG;
  206.  
  207. DOMTOT = DOM1 ET DOM2;
  208. 'ELIMINATION' DOMTOT 1D-6;
  209.  
  210. $DOMTOT = 'MODELISER' DOMTOT 'EULER';
  211. $DOM1 = 'MODELISER' DOM1 'EULER';
  212. $DOM2 = 'MODELISER' DOM2 'EULER';
  213. TDOMTOT = 'DOMA' $DOMTOT 'VF';
  214. TDOM1 = 'DOMA' $DOM1 'VF';
  215. TDOM2 = 'DOMA' $DOM2 'VF';
  216. MDOM1 = TDOM1 . 'QUAF' ;
  217. MDOM2 = TDOM2 . 'QUAF' ;
  218. MDOMTOT = TDOMTOT . 'QUAF' ;
  219.  
  220. 'ELIM' (MDOMTOT 'ET' MDOM1 'ET' MDOM2) 1.E-6 ;
  221.  
  222. 'SI' GRAPH;
  223. 'TRACER' (('DOMA' $DOMTOT 'MAILLAGE') 'ET'
  224. ('DOMA' $DOMTOT 'FACEL') 'ET' P1) 'TITRE' 'Domaine et FACEL';
  225. 'FINSI' ;
  226.  
  227. *
  228. **** Redefinition de P1 dans $DOMTOT 'FACE'
  229. *
  230.  
  231. P1 = ('DOMA' $DOMTOT 'FACE') 'POIN' 'PROC' P1;
  232.  
  233. ***********************
  234. **** Les CHPOINTs ****
  235. ***********************
  236.  
  237. uxgd = (ungd '*' ('COS' ANGLE)) '-' (utgd '*' ('SIN' ANGLE));
  238. uygd = (ungd '*' ('SIN' ANGLE)) '+' (utgd '*' ('COS' ANGLE));
  239.  
  240. RN1 = 'KCHT' $DOM1 'SCAL' 'CENTRE' rogd;
  241. RN2 = 'KCHT' $DOM2 'SCAL' 'CENTRE' rogd;;
  242. RN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' (RN1 'ET' RN2);
  243.  
  244. VN1 = 'KCHT' $DOM1 'VECT' 'CENTRE' (uxgd uygd);
  245. VN2 = 'KCHT' $DOM2 'VECT' 'CENTRE' (uxgd uygd);
  246. VN = 'KCHT' $DOMTOT 'VECT' 'CENTRE' (VN1 'ET' VN2);
  247.  
  248. PN1 = 'KCHT' $DOM1 'SCAL' 'CENTRE' pgd;
  249. PN2 = 'KCHT' $DOM2 'SCAL' 'CENTRE' pgd;
  250. PN = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' (PN1 'ET' PN2);
  251.  
  252. SCAN1 = 'KCHT' $DOM1 'VECT' 'CENTRE'
  253. 'COMP' (PGAZ2 . 'SCALPASS') (tutugd tatagd) ;
  254. SCAN2 = 'KCHT' $DOM2 'VECT' 'CENTRE'
  255. 'COMP' (PGAZ2 . 'SCALPASS') (tutugd tatagd) ;
  256. SCAN = 'KCHT' $DOMTOT 'VECT' 'CENTRE'
  257. 'COMP' (PGAZ2 . 'SCALPASS') (SCAN1 'ET' SCAN2);
  258.  
  259. ***************************
  260. **** L'operateur PRET****
  261. ***************************
  262.  
  263. ORDESP = 1;
  264. ORDTEM = 1;
  265.  
  266.  
  267. ROF VITF PF = 'PRET' 'PERFTEMP' ORDESP ORDTEM
  268. $DOMTOT PGAZ RN VN PN ;
  269.  
  270. ROF2 VITF2 PF2 SCAF = 'PRET' 'PERFTEMP' ORDESP ORDTEM
  271. $DOMTOT PGAZ2 RN VN PN SCAN ;
  272.  
  273. GEOP1 = ('DOMA' $DOMTOT 'FACEL') 'ELEM' 'APPUYE' 'LARGEMENT' P1;
  274. GEOP2 = ('DOMA' $DOMTOT 'FACE') 'ELEM' 'APPUYE' 'LARGEMENT' P1;
  275.  
  276.  
  277. REFGEOP1 = 'REDU' VITF GEOP2;
  278.  
  279. *
  280. **** Orientation de la normal n de castem par raport a la
  281. * notre; t est par consequence
  282. *
  283.  
  284. NCOS = 'EXTRAIRE' REFGEOP1 'NX' 1 1 1;
  285. NSIN = 'EXTRAIRE' REFGEOP1 'NY' 1 1 1;
  286.  
  287. 'SI' (('ABS' NCOS) > ('ABS' NSIN));
  288. ORIENT = ('COS' ANGLE) '/' NCOS;
  289. 'SINON';
  290. ORIENT = ('SIN' ANGLE) '/' NSIN;
  291. 'FINSI' ;
  292.  
  293. ORIENT = 'SIGNE' 'FLOTTANT' ORIENT;
  294.  
  295.  
  296.  
  297. ********************************
  298. ********************************
  299. **** 'KONV' ****
  300. ********************************
  301. ********************************
  302.  
  303.  
  304. 'REPETER' BLMETO 2 ;
  305.  
  306. 'SI' ('EGA' &BLMETO 1) ;
  307. METO = 'VLH' ;
  308. 'FINSI' ;
  309.  
  310. 'SI' ('EGA' &BLMETO 2) ;
  311. METO = 'SS' ;
  312. 'FINSI' ;
  313.  
  314. 'MESSAGE' ('CHAINE' 'METO = ' METO) ;
  315. LISTINC1 = ('MOTS' 'RN' 'RUNX' 'RUNY' 'RETN') ;
  316. CHPFLU DT = 'KONV' 'VF' 'PERFTEMP' 'FLUX' METO
  317. $DOMTOT PGAZ LISTINC1 ROF VITF PF ;
  318.  
  319. FLUX1 = 'EXCO' ('EXTRAIRE' LISTINC1 1) CHPFLU ;
  320. FLUX2X = 'EXCO' ('EXTRAIRE' LISTINC1 2) CHPFLU ;
  321. FLUX2Y = 'EXCO' ('EXTRAIRE' LISTINC1 3) CHPFLU ;
  322. FLUX3 = 'EXCO' ('EXTRAIRE' LISTINC1 4) CHPFLU ;
  323.  
  324. FLUX2N = (FLUX2X '*' ('COS' ANGLE)) '+' (FLUX2Y * ('SIN' ANGLE));
  325. FLUX2T = (FLUX2Y '*' ('COS' ANGLE)) '-' (FLUX2X * ('SIN' ANGLE));
  326.  
  327.  
  328. f1 = 'EXTRAIRE' FLUX1 'SCAL' P1;
  329. ERRO = 1D-8 '*' f1gd;
  330. LOGI1 = ('ABS' ((f1 '*' ORIENT) '-' f1gd)) < ERRO ;
  331.  
  332. f2 = 'EXTRAIRE' FLUX2N 'SCAL' P1;
  333. ERRO = 1D-8 '*' f2gd;
  334. LOGI2 = ('ABS' ((f2 '*' ORIENT) '-' f2gd)) < ERRO;
  335. LOGI1 = LOGI1 'ET' LOGI2;
  336.  
  337. f3 = 'EXTRAIRE' FLUX2T 'SCAL' P1;
  338. ERRO = 1D-8 '*' f3gd;
  339. LOGI2 = ('ABS' ((f3 '*' ORIENT) '-' f3gd)) < ERRO;
  340. LOGI1 = LOGI1 'ET' LOGI2;
  341.  
  342. f4 = 'EXTRAIRE' FLUX3 'SCAL' P1;
  343. ERRO = 1D-8 '*' f4gd;
  344. LOGI2 = ('ABS' ((f4 '*' ORIENT) '-' f4gd)) < ERRO;
  345. LOGI1 = LOGI1 'ET' LOGI2;
  346.  
  347. 'SI' ('NON' LOGI1);
  348. 'MESSAGE' ;
  349. 'MESSAGE' 'OPERATEUR KONV';
  350. 'MESSAGE' ('CHAINE' METO);
  351. 'MESSAGE' ;
  352. 'MESSAGE' ('CHAINE' 'df1 ' ('ABS' (f1gd '-' (ORIENT '*'f1)))
  353. ' f1 ' f1);
  354. 'MESSAGE' ('CHAINE' 'df2 ' ('ABS' (f2gd '-' (ORIENT '*'f2)))
  355. ' f2 ' f2);
  356. 'MESSAGE' ('CHAINE' 'df3 ' ('ABS' (f3gd '-' (ORIENT '*'f3)))
  357. ' f3 ' f3);
  358. 'MESSAGE' ('CHAINE' 'df4 ' ('ABS' (f4gd '-' (ORIENT '*'f4)))
  359. ' f4 ' f4);
  360. 'ERREUR' 5;
  361. 'FINSI' ;
  362.  
  363. *
  364. **** Scalaires passifs
  365. *
  366.  
  367. 'MESSAGE' ('CHAINE' 'METO = ' METO ' (scalaires passifs)') ;
  368. LISTINC2 = ('MOTS' 'RN' 'RUNX' 'RUNY' 'RETN' 'TUTU' 'TATA') ;
  369. CHPFLU2 DT2 = 'KONV' 'VF' 'PERFTEMP' 'FLUX' METO
  370. $DOMTOT PGAZ2 LISTINC2 ROF VITF PF SCAF ;
  371.  
  372. CHPFLU3 = 'EXCO' LISTINC1 CHPFLU2 LISTINC1 ;
  373.  
  374. 'SI' (('MAXIMUM' (CHPFLU3 '-' CHPFLU) 'ABS') > 1.0D-15) ;
  375. 'MESSAGE' 'Probleme dans les transport de scalaires passifs' ;
  376. 'ERREUR' 5;
  377. 'FINSI' ;
  378.  
  379. fceltutu = 'EXTRAIRE' CHPFLU2 'TUTU' P1;
  380. LOGI1 = ('ABS' ((fceltutu '*' ORIENT) '-' ftutugd)) <
  381. ('ABS' (1.0D-8 * fceltutu)) ;
  382.  
  383. fceltata = 'EXTRAIRE' CHPFLU2 'TATA' P1;
  384. LOGI2 = ('ABS' ((fceltata '*' ORIENT) '-' ftatagd)) <
  385. ('ABS' (1.0D-8 * fceltata)) ;
  386. LOGI1 = LOGI1 'ET' LOGI2;
  387.  
  388. 'SI' ('NON' LOGI1);
  389. 'MESSAGE' ;
  390. 'MESSAGE' 'OPERATEUR KONV';
  391. 'MESSAGE' 'Transport de scalaires passifs' ;
  392. 'ERREUR' 5;
  393. 'FINSI' ;
  394.  
  395.  
  396. *
  397. ****** Le residu
  398. *
  399. *
  400. ****** CHPDIR = chpoint sur les face de chaque elts
  401. *
  402. *
  403. * NB Dans ce cas particulier, DOMTOT n'a que de QUA4,
  404. * donc ($DOMTOT . 'XXNORMAE') contient une seule
  405. * zona elementaire!
  406.  
  407.  
  408. RESI MATCEL = 'KOPS' 'MATRIK' ;
  409. MAICHA = 'EXTRAIRE' ('DOMA' $DOMTOT 'XXNORMAE') 'MAILLAGE' ;
  410. NBELEM = 'NBEL' MAICHA ;
  411.  
  412. 'REPETER' BLELE NBELEM ;
  413. CHPDIR MATCEL = 'KOPS' 'MATRIK' ;
  414. ELEM1 = &BLELE 'ELEM' MAICHA ;
  415. ELEMC1 = 'CHANGER' ELEM1 'POI1' ;
  416. NBPOIN = 'NBNO' ELEMC1 ;
  417. PCEN = ('DOMA' $DOMTOT 'CENTRE') 'POIN' &BLELE ;
  418. XVOL = ('DOMA' $DOMTOT 'XXVOLUM') 'EXTRAIRE' PCEN 'SCAL' ;
  419. MAICEN = 'MANU' 'POI1'PCEN ;
  420. 'REPETER' BLEPO NBPOIN ;
  421. P1 = &BLEPO 'POIN' ELEMC1 ;
  422. ELEP1 = 'MANUEL' 'POI1' P1 ;
  423. COEFNO = 'EXTRAIRE' ('DOMA' $DOMTOT 'XXNORMAE')
  424. 'SCAL' 1 &BLELE &BLEPO ;
  425. CHPDIR = CHPDIR 'ET' ('MANUEL' 'CHPO' ELEP1 1 'SCAL'
  426. COEFNO 'NATU' 'DISCRET') ;
  427. 'FIN' BLEPO ;
  428. MAIDIR = 'EXTRAIRE' CHPDIR 'MAILLAGE' ;
  429. *
  430. ********* Le flux qui sort
  431. *
  432. FLUSOR = CHPDIR '*' ('REDU' CHPFLU2 MAIDIR)
  433. ('MOTS' 'SCAL' 'SCAL' 'SCAL' 'SCAL' 'SCAL' 'SCAL')
  434. LISTINC2 LISTINC2 ;
  435. CHPUN = 'MANUEL' 'CHPO' MAIDIR 1 'SCAL' (-1./XVOL) ;
  436. NBCOMP = 'DIME' LISTINC2 ;
  437. 'REPETER' BLCOM NBCOMP ;
  438. NOMCOM = 'EXTRAIRE' LISTINC2 &BLCOM ;
  439. CELSCA = 'XTY' FLUSOR CHPUN ('MOTS' NOMCOM)
  440. ('MOTS' 'SCAL') ;
  441. RESI = RESI 'ET' ('MANUEL' 'CHPO' MAICEN 1
  442. NOMCOM CELSCA 'NATU' 'DISCRET') ;
  443. 'FIN' BLCOM ;
  444. 'FIN' BLELE ;
  445.  
  446. CHPRES DT1 = 'KONV' 'VF' 'PERFTEMP' 'RESI' METO
  447. $DOMTOT PGAZ2 LISTINC2 ROF VITF PF SCAF ;
  448.  
  449. ERRO = 'MAXIMUM' (RESI '-' CHPRES) 'ABS' ;
  450. 'SI' (ERRO > 1.0D-8) ;
  451. 'MESSAGE' ;
  452. 'MESSAGE' 'Probleme dans KONV' ;
  453. 'ERREUR' 5 ;
  454. 'FINSI' ;
  455.  
  456. ERRO = 'ABS' (DT '-' DT1) ;
  457. 'SI' (ERRO > 1.0D-15) ;
  458. 'MESSAGE' ;
  459. 'MESSAGE' 'Probleme dans KONV' ;
  460. 'ERREUR' 5 ;
  461. 'FINSI' ;
  462.  
  463.  
  464.  
  465. 'FIN' BLMETO ;
  466.  
  467. ****************************************************
  468. ****************************************************
  469. ******** Fin boucle sur les angles *********
  470. ****************************************************
  471. ****************************************************
  472.  
  473. 'FIN' BLOC;
  474.  
  475. 'FIN' ;
  476.  
  477.  
  478.  
  479.  
  480.  
  481.  
  482.  
  483.  
  484.  
  485.  

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