Télécharger shock3d.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : shock3d.dgibi
  2. ************************************************************************
  3. * NAME : main3d.dgibi
  4. * DESCRIPTION : file to solve 3D tests:
  5. * 2. shock tube
  6. *
  7. * LANGUAGE : GIBIANE-CAST3M
  8. * AUTHOR : José R. Garcia Cascales
  9. * Universidad Politecnica de Cartagena,
  10. * jr.garcia@upct.es
  11. ************************************************************************
  12. * VERSION : v1, 04/04/2002, version initiale
  13. * HISTORIQUE : vfinal, 22/11/2005
  14. ************************************************************************
  15. *
  16. * MESH DEFINITION
  17. *
  18. ************************************************************************
  19. 'OPTION' 'DIME' 3 'ELEM' 'CUB8' 'ISOV' 'SULI'
  20. 'ECHO' 1 'TRAC' 'X';
  21.  
  22. GRAPH = VRAI ;
  23. GRAPH = FAUX ;
  24. *
  25. ** number of elements in the z axis (after rotation) = NZ
  26. *
  27. NZ = 100 ;
  28.  
  29. NX = 2;
  30. NY = 2;
  31.  
  32. H = 10.D0 ;
  33. DZ = H/NZ ;
  34.  
  35. L = NX*DZ ;
  36. P = NY*DZ ;
  37. xR = 0.5D0*L ;
  38. xL = -1.0D0*xR ;
  39. yR = 0.5D0*P ;
  40. yL = -1.0D0*yR ;
  41.  
  42. P1 = XR YR 0. ;
  43. P2 = XL YR 0. ;
  44. P3 = XL YL 0. ;
  45. P4 = XR YL 0. ;
  46.  
  47. P1P2 = 'DROIT' NX P1 P2 ;
  48. P2P3 = 'DROIT' NY P2 P3 ;
  49. P3P4 = 'DROIT' NX P3 P4 ;
  50. P4P1 = 'DROIT' NY P4 P1 ;
  51.  
  52. BASEJR = 'DALLER' P1P2 P2P3 P3P4 P4P1 ;
  53.  
  54. POSH = H/2. ;
  55. NEGH = (-1.)*POSH ;
  56. MIDDIV = NZ/2 ;
  57.  
  58. V1 = 0. 0. POSH ;
  59. V2 = 0. 0. NEGH ;
  60.  
  61. 'OPTI' 'ELEM' 'CUB8';
  62. DOM1 = BASEJR 'VOLUME' MIDDIV 'TRAN' V1 ;
  63. DOM2 = BASEJR 'VOLUME' MIDDIV 'TRAN' V2 ;
  64. DOMTOT = DOM1 'ET' DOM2 ;
  65.  
  66. MDOMTOT = 'MODELISER' DOMTOT 'EULER';
  67. MDOM1 = 'MODELISER' DOM1 'EULER';
  68. MDOM2 = 'MODELISER' DOM2 'EULER';
  69.  
  70. $DOMTOT = 'DOMA' MDOMTOT 'VF' ;
  71. $DOM1 = 'DOMA' MDOM1 'VF' ;
  72. $DOM2 = 'DOMA' MDOM2 'VF' ;
  73.  
  74. TDOMTOT = $DOMTOT. 'QUAF';
  75. TDOM1 = $DOM1. 'QUAF';
  76. TDOM2 = $DOM2. 'QUAF';
  77. 'ELIMINATION' (TDOMTOT 'ET' TDOM1 'ET' TDOM2) 1D-6;
  78. *
  79. ******* Creation de la ligne Utilisée pour le Post-Traitement
  80. * reliant les points centres
  81. *
  82. XINIT = 0.5D0*DZ ;
  83. YINIT = 0.5D0*DZ ;
  84. ZINIT = NEGH '+' (0.5D0*DZ) ;
  85. XFIN = XINIT ;
  86. YFIN = YINIT ;
  87. ZFIN = POSH '-' (0.5D0*DZ) ;
  88. PINI = XINIT YINIT ZINIT;
  89. PFIN = XFIN YFIN ZFIN;
  90. DIV = NZ '-' 1 ;
  91. COURB = PINI 'DROIT' DIV PFIN;
  92. COURB = COURB 'COULEUR' 'VERT';
  93. 'ELIMINATION' (1.0D-6) COURB $DOMTOT.CENTRE;
  94.  
  95. 'SI' GRAPH ;
  96. 'TRACER' (DOMTOT 'ET' COURB) 'TITRE' 'Domaine total' ;
  97. 'FINSI' ;
  98. ************************************************************************
  99. *
  100. * INITIAL CONDITIONS
  101. *
  102. ************************************************************************
  103. CP = 2.d0 ;
  104. CVM = 0.d0 ;
  105. *
  106. *
  107. VCH1 = ('MANU' 'CHPO' ($DOM1.'CENTRE') 1 'SCAL' 0.1d0 'NATURE'
  108. 'DISCRET') 'ET'
  109. ('MANU' 'CHPO' ($DOM2.'CENTRE') 1 'SCAL' 0.25d0 'NATURE'
  110. 'DISCRET') ;
  111. VCH2 = ('MANU' 'CHPO' ($DOMTOT.'CENTRE') 3 'UVX' 0.d0 'UVY' 0.d0
  112. 'UVZ' 0.d0 'NATURE' 'DISCRET') ;
  113. VCH3 = ('MANU' 'CHPO' ($DOMTOT.'CENTRE') 3 'ULX' 0.d0 'ULY' 0.d0
  114. 'ULZ' 0.d0 'NATURE' 'DISCRET') ;
  115. VCH4 = ('MANU' 'CHPO' ($DOM1.'CENTRE') 1 'SCAL' 10.d6 'NATURE'
  116. 'DISCRET') 'ET'
  117. ('MANU' 'CHPO' ($DOM2.'CENTRE') 1 'SCAL' 20.d6 'NATURE'
  118. 'DISCRET') ;
  119. VCH5 = ('MANU' 'CHPO' ($DOMTOT.'CENTRE') 1 'SCAL' 308.d0 'NATURE'
  120. 'DISCRET') ;
  121. VCH6 = ('MANU' 'CHPO' ($DOMTOT.'CENTRE') 1 'SCAL' 308.d0 'NATURE'
  122. 'DISCRET');
  123. GAMV = ('MANU' 'CHPO' ($DOMTOT.'CENTRE') 1 'SCAL' 1.4d0 'NATURE'
  124. 'DISCRET');
  125. CPV = ('MANU' 'CHPO' ($DOMTOT.'CENTRE') 1 'SCAL' 1008.d0 'NATURE'
  126. 'DISCRET');
  127. GAML = ('MANU' 'CHPO' ($DOMTOT.'CENTRE') 1 'SCAL' 2.8d0 'NATURE'
  128. 'DISCRET');
  129. CPL = ('MANU' 'CHPO' ($DOMTOT.'CENTRE') 1 'SCAL' 4186.d0 'NATURE'
  130. 'DISCRET');
  131. UNO = ('MANU' 'CHPO' ($DOMTOT.'CENTRE') 1 'SCAL' 1.d0 'NATURE'
  132. 'DISCRET');
  133. PIL = ('MANU' 'CHPO' ($DOMTOT.'CENTRE') 1 'SCAL' 8.5d8 'NATURE'
  134. 'DISCRET');
  135.  
  136. * vch7 (rhov) = gammav*p/((gammav - 1.D0)*Cpv*T)
  137. * vch8 (rhol) = gammal*(p + pil)/((gammal - 1.D0)*Cpl*T)
  138. * IMPROVE, VCH7 & VCH8 ARE NOT ASSOCIATED TO ANY MESH!!!!!!!!!
  139.  
  140. NUM1 = 'PSCAL' GAMV VCH4 ('MOTS' 'SCAL') ('MOTS' 'SCAL') ;
  141. NUM2 = 0.4d0 * ('PSCAL' CPV VCH5 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ;
  142. VCH7 = NUM1 '/' NUM2 ;
  143. VCH7 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' VCH7 ;
  144. NUM3 = 2.8d0 * (VCH4 '+' 8.5d8) ;
  145. NUM4 = 1.8d0 * ('PSCAL' CPL VCH6 ('MOTS' 'SCAL') ('MOTS' 'SCAL')) ;
  146. VCH8 = NUM3 '/' NUM4 ;
  147. VCH8 = 'KCHT' $DOMTOT 'SCAL' 'CENTRE' VCH8 ;
  148. ***************************************************************
  149. ***** PROCEDURE POUR CALCULER LES VARIABLES CONSERVATIVES *****
  150. ***************************************************************
  151. 'DEBPROC' CONS ;
  152. 'ARGUMENT'VCH1*'CHPOINT' VCH2*'CHPOINT' VCH3*'CHPOINT' VCH4*'CHPOINT'
  153. VCH5*'CHPOINT' VCH6*'CHPOINT' VCH7*'CHPOINT' VCH8*'CHPOINT'
  154. ;
  155.  
  156. WCH1 = 'PSCAL' VCH1 VCH7 ('MOTS' 'SCAL') ('MOTS' 'SCAL') ;
  157. WCH1 = 'EXCO' 'SCAL' WCH1 'SCAL' 'NATURE' 'DISCRET' ;
  158. AL = VCH1 '-' 1.0 ;
  159. AL = AL * (-1.0) ;
  160. WCH2 = 'PSCAL' AL VCH8 ('MOTS' 'SCAL') ('MOTS' 'SCAL') ;
  161. WCH2 = 'EXCO' 'SCAL' WCH2 'SCAL' 'NATURE' 'DISCRET' ;
  162. WCH3 = WCH1 '*' VCH2 ('MOTS' 'SCAL' 'SCAL' 'SCAL') ('MOTS' 'UVX'
  163. 'UVY' 'UVZ') ('MOTS' 'UVX' 'UVY' 'UVZ') ;
  164. WCH3 = 'EXCO' ('MOTS' 'UVX' 'UVY' 'UVZ') WCH3 ('MOTS' 'UVX'
  165. 'UVY' 'UVZ') 'NATURE' 'DISCRET' ;
  166. WCH4 = WCH2 '*' VCH3 ('MOTS' 'SCAL' 'SCAL' 'SCAL') ('MOTS' 'ULX'
  167. 'ULY' 'ULZ') ('MOTS' 'ULX' 'ULY' 'ULZ') ;
  168. WCH4 = 'EXCO' ('MOTS' 'ULX' 'ULY' 'ULZ') WCH4 ('MOTS' 'ULX' 'ULY'
  169. 'ULZ')'NATURE' 'DISCRET' ;
  170. UVX2 = 'PSCAL' VCH2 VCH2 ('MOTS' 'UVX') ('MOTS' 'UVX') ;
  171. UVY2 = 'PSCAL' VCH2 VCH2 ('MOTS' 'UVY') ('MOTS' 'UVY') ;
  172. UVZ2 = 'PSCAL' VCH2 VCH2 ('MOTS' 'UVZ') ('MOTS' 'UVZ') ;
  173. ECV = 0.5d0 * (UVX2 '+' UVY2 '+' UVY2) ;
  174. * ev = Cpv*v(5)/gammav
  175. * IEV = ('PSCAL' CPV VCH5 ('MOTS' 'SCAL') ('MOTS' 'SCAL'))'/' 1.4d0 ;
  176. IEV = (1008.d0 * VCH5)'/' 1.4d0 ;
  177. TIEV = IEV '+' ECV ;
  178. WCH5 = WCH1 * TIEV ;
  179. WCH5 = 'EXCO' 'SCAL' WCH5 'SCAL' 'NATURE' 'DISCRET' ;
  180. ULX2 = 'PSCAL' VCH3 VCH3 ('MOTS' 'ULX') ('MOTS' 'ULX') ;
  181. ULY2 = 'PSCAL' VCH3 VCH3 ('MOTS' 'ULY') ('MOTS' 'ULY') ;
  182. ULZ2 = 'PSCAL' VCH3 VCH3 ('MOTS' 'ULZ') ('MOTS' 'ULZ') ;
  183. ECL = 0.5d0 * (ULX2 '+' ULY2 '+' ULZ2) ;
  184. * el = Cpl*v(6)/gammal + pil/rl
  185. * NUM9 = ('PSCAL' CPL VCH5 ('MOTS' 'SCAL') ('MOTS' 'SCAL'))'/' 2.8d0 ;
  186. NUM9 = (4186.d0 * VCH6)'/' 2.8d0 ;
  187. NUM0 = 8.5d8 '*'(VCH8 '**' -1.) ;
  188. NUM0 = 'NOMC' 'SCAL' NUM0 ;
  189. IEL = NUM9 '+' NUM0 ;
  190. TIEL = IEL '+' ECL ;
  191. WCH6 = WCH2 * TIEL ;
  192. WCH6 = 'EXCO' 'SCAL' WCH6 'SCAL' 'NATURE' 'DISCRET' ;
  193.  
  194. 'RESPRO' WCH1 WCH2 WCH3 WCH4 WCH5 WCH6 ;
  195.  
  196. 'FINPROC' ;
  197. ***************************************************************
  198. WCH1 WCH2 WCH3 WCH4 WCH5 WCH6 = CONS
  199. VCH1 VCH2 VCH3 VCH4 VCH5 VCH6 VCH7 VCH8 ;
  200. *
  201. * We create old alpha, to satisfy the prim operator
  202. * necessities
  203. *
  204. OALP = VCH1 ;
  205. OTV = VCH5 ;
  206. OTL = VCH6 ;
  207.  
  208. PINT RL RV TL TV P UL UV ALP = 'PRIM' 'TWOFLUID'
  209. WCH1 WCH2 WCH3 WCH4 WCH5 WCH6 OALP OTV OTL CP CVM;
  210.  
  211. ERRO = 'MAXIMUM' (VCH4 '-' P) 'ABS' ;
  212.  
  213. 'SI' (ERRO > 1.0D-6) ;
  214. 'MESSAGE' 'Problem in the stic file!!!'
  215. 'ERREUR' 5 ;
  216. 'FINSI' ;
  217. *
  218. **** Plot of IC
  219. *
  220. 'SI' GRAPH ;
  221.  
  222. MOD1 = 'MODELISER' ($DOMTOT . 'MAILLAGE' ) 'THERMIQUE';
  223.  
  224. CHM_V1 = 'KCHA' $DOMTOT 'CHAM' VCH1 ;
  225. CHM_V2 = 'KCHA' $DOMTOT 'CHAM' VCH2 ;
  226. CHM_V3 = 'KCHA' $DOMTOT 'CHAM' VCH3 ;
  227. CHM_V4 = 'KCHA' $DOMTOT 'CHAM' VCH4 ;
  228. CHM_V5 = 'KCHA' $DOMTOT 'CHAM' VCH5 ;
  229. CHM_V6 = 'KCHA' $DOMTOT 'CHAM' VCH6 ;
  230. CHM_V7 = 'KCHA' $DOMTOT 'CHAM' VCH7 ;
  231. CHM_V8 = 'KCHA' $DOMTOT 'CHAM' VCH8 ;
  232.  
  233. 'TRACER' CHM_V1 MOD1
  234. 'TITR' ('CHAINE' 'alpha at t=' 0.0);
  235. 'TRACER' CHM_V2 MOD1
  236. 'TITR' ('CHAINE' 'uv at t=' 0.0);
  237. 'TRACER' CHM_V3 MOD1
  238. 'TITR' ('CHAINE' 'ul at t=' 0.0);
  239. 'TRACER' CHM_V4 MOD1
  240. 'TITR' ('CHAINE' 'p at t=' 0.0);
  241. 'TRACER' CHM_V5 MOD1
  242. 'TITR' ('CHAINE' 'Tv at t=' 0.0);
  243. 'TRACER' CHM_V6 MOD1
  244. 'TITR' ('CHAINE' 'Tl at t=' 0.0);
  245. 'TRACER' CHM_V7 MOD1
  246. 'TITR' ('CHAINE' 'rv at t=' 0.0);
  247. 'TRACER' CHM_V8 MOD1
  248. 'TITR' ('CHAINE' 'rl at t=' 0.0);
  249.  
  250.  
  251. 'FINSI' ;
  252. *************************************************************************
  253. *
  254. * MAIN PROGRAM
  255. *
  256. *************************************************************************
  257. * Different upwind scheme available
  258. * AUSM +
  259. * METO = 'AUSMP1' ;
  260. * Preconditioned AUSM +
  261. METO = 'AUSMP2';
  262. * AUSMDV
  263. * METO = 'AUSMDV1';
  264. * Preconditioned AUSMDV
  265. * METO = 'AUSMDV2';
  266. *
  267. * Iterations
  268. *
  269. NITER = 10000000 ;
  270. *
  271. * Gravity, Final time, CFL-like parameter
  272. *
  273. g = 0.d0 ;
  274. TFINAL = 0.005d0 ;
  275. MCFL = 0.1D0 ;
  276. *
  277. * Order in space, order in time and a useful counter
  278. *
  279. ORDESP = 2;
  280. ORDTEM = 2;
  281. COUNT = 1;
  282. *
  283. * TIME IS RESETED TO 0 seconds
  284. *
  285. TPS = 0. ;
  286. *
  287. * Names of the residuum components
  288. *
  289. LISTINCO = 'MOTS' 'RVN' 'RLN' 'RUVX' 'RUVY' 'RUVZ' 'RULX' 'RULY' 'RULZ' 'RETV' 'RETL' ;
  290.  
  291. 'MESSAGE' ;
  292. 'MESSAGE' ('CHAINE' 'Methode = ' METO) ;
  293. 'MESSAGE' ;
  294. *
  295. * Some preliminary calculation for the evaluation of gradients
  296. * 1st. for the gradient of void fraction
  297. * 2nd. for the gradients in 2nd order
  298. *
  299. DALP LIMCH GRG = 'PENT' MDOMTOT 'CENTRE'
  300. 'EULESCAL' 'NOLIMITE' ('MOTS' 'SCAL') ALP ;
  301.  
  302. 'SI' (ORDESP 'EGA' 2) ;
  303.  
  304. GRADAL CACA COEFSCAL = 'PENT' MDOMTOT 'CENTRE' 'EULESCAL'
  305. 'NOLIMITE' ('MOTS' 'SCAL') ALP ;
  306. GRADUV CACA COEFVECT = 'PENT' MDOMTOT 'CENTRE' 'EULEVECT'
  307. 'NOLIMITE' ('MOTS' 'UX' 'UY' 'UZ')
  308. ('EXCO' ('MOTS' 'UVX' 'UVY' 'UVZ')
  309. UV ('MOTS' 'UX' 'UY' 'UZ')) ;
  310. 'FINSI' ;
  311.  
  312. *
  313. **** Temporal loop
  314. *
  315. 'REPETER' BL1 NITER ;
  316. *
  317. **** Evaluation of variables at each side of the interface (PRET operator)
  318. *
  319. 'SI' ((ORDESP 'EGA' 1) 'OU' (&BL1 'EGA' 1)) ;
  320.  
  321. RLF RVF TLF TVF PF ULF UVF ALPF = 'PRET'
  322. 'TWOFLUID' 1 1 $DOMTOT
  323. ALP UV UL P TV TL RV RL ;
  324.  
  325. 'SINON' ;
  326.  
  327. GRADAL LIMALP = 'PENT' MDOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR'
  328. ('MOTS' 'SCAL') ALP 'GRADGEO' COEFSCAL ;
  329. GRADUV LIMUV = 'PENT' MDOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR'
  330. ('MOTS' 'UVX' 'UVY' 'UVZ') UV 'GRADGEO' COEFVECT ;
  331. GRADUL LIMUL = 'PENT' MDOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR'
  332. ('MOTS' 'ULX' 'ULY' 'ULZ') UL 'GRADGEO' COEFVECT ;
  333. GRADP LIMP = 'PENT' MDOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR'
  334. ('MOTS' 'SCAL') P 'GRADGEO' COEFSCAL ;
  335. GRADTV LIMTV = 'PENT' MDOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR'
  336. ('MOTS' 'SCAL') TV 'GRADGEO' COEFSCAL ;
  337. GRADTL LIMTL = 'PENT' MDOMTOT 'CENTRE' 'EULESCAL' 'LIMITEUR'
  338. ('MOTS' 'SCAL') TL 'GRADGEO' COEFSCAL ;
  339.  
  340. RLF RVF TLF TVF PF ULF UVF ALPF = 'PRET'
  341. 'TWOFLUID' ORDESP ORDTEM $DOMTOT
  342. ALP GRADAL LIMALP
  343. UV GRADUV LIMUV
  344. UL GRADUL LIMUL
  345. P GRADP LIMP
  346. TV GRADTV LIMTV
  347. TL GRADTL LIMTL
  348. RV RL DELTAT ;
  349.  
  350. 'FINSI' ;
  351. *
  352. **** Evaluation of Residual and Dt
  353. *
  354. RESIDU DELTAT = 'KONV' 'VF' 'TWOFLUID' 'RESI' METO
  355. $DOMTOT LISTINCO ALPF UVF ULF PF TVF TLF RVF RLF ;
  356.  
  357. DT = MCFL '*' DELTAT ;
  358.  
  359. 'SI' (COUNT 'EGA' 1 ) ;
  360. COUNT = 0 ;
  361. OALP = ALP ;
  362. OLDDT = DT ;
  363. 'FINSI' ;
  364. *
  365. **** Increment of the variables (convection)
  366. *
  367. RESIDU = DT '*' RESIDU ;
  368.  
  369. DRVN = 'EXCO' 'RVN' RESIDU 'SCAL' ;
  370. DRLN = 'EXCO' 'RLN' RESIDU 'SCAL' ;
  371. DRUVN = 'EXCO' ('MOTS' 'RUVX' 'RUVY' 'RUVZ') RESIDU ('MOTS' 'UVX'
  372. 'UVY' 'UVZ') ;
  373. DRULN = 'EXCO' ('MOTS' 'RULX' 'RULY' 'RULZ') RESIDU ('MOTS' 'ULX'
  374. 'ULY' 'ULZ') ;
  375. DRETVN = 'EXCO' 'RETV' RESIDU 'SCAL' ;
  376. DRETLN = 'EXCO' 'RETL' RESIDU 'SCAL' ;
  377. *
  378. * Analizamos la conservacion de la masa
  379. *
  380. MASAV = 'XTY' WCH1 UNO ('MOTS' 'SCAL') ('MOTS' 'SCAL') ;
  381. MASAL = 'XTY' WCH2 UNO ('MOTS' 'SCAL') ('MOTS' 'SCAL') ;
  382. MOMV = 'XTY' WCH3 UNO ('MOTS' 'UVY') ('MOTS' 'SCAL') ;
  383. MOML = 'XTY' WCH4 UNO ('MOTS' 'ULY') ('MOTS' 'SCAL') ;
  384. ENERV = 'XTY' WCH5 UNO ('MOTS' 'SCAL') ('MOTS' 'SCAL') ;
  385. ENERL = 'XTY' WCH6 UNO ('MOTS' 'SCAL') ('MOTS' 'SCAL') ;
  386. *
  387. **** Previous update of conserved variables
  388. *
  389. TPS = TPS '+' DT ;
  390. WCH1 = WCH1 '+' DRVN ;
  391. WCH2 = WCH2 '+' DRLN ;
  392. WCH3 = WCH3 '+' DRUVN ;
  393. WCH4 = WCH4 '+' DRULN ;
  394. WCH5 = WCH5 '+' DRETVN ;
  395. WCH6 = WCH6 '+' DRETLN ;
  396. ************************************************************************
  397. **** Adittion of the source terms, final update of conserved variables
  398. ************************************************************************
  399. DALP1 LIMCH = 'PENT' MDOMTOT 'CENTRE'
  400. 'EULESCAL' 'LIMITEUR' ('MOTS' 'SCAL') ALP 'GRADGEO' GRG ;
  401.  
  402. DADP = DT * PINT * DALP1 ;
  403. *
  404. * Source term corresponding to the vapour momentum equation
  405. *
  406. S31 = ('MANU' 'CHPO' ($DOMTOT.'CENTRE') 2 'UVX' 0.d0 'UVY' 0.d0
  407. 'NATURE' 'DISCRET');
  408. S32 = DT * ALP * RV * g ;
  409. S321 = 'EXCO' 'SCAL' S32 'UVZ' 'NATURE' 'DISCRET' ;
  410. S3 = S31 'ET' S321 ;
  411. DADP = 'EXCO' ('MOTS' 'P1DX' 'P1DY' 'P1DZ') DADP ('MOTS' 'UVX' 'UVY'
  412. 'UVZ') 'NATURE' 'DISCRET' ;
  413. SO3 = (S3 '+' DADP) ;
  414. SO3 = 'NOMC' ('MOTS' 'UVX' 'UVY' 'UVZ') SO3 ('MOTS' 'UVX' 'UVY'
  415. 'UVZ') 'NATURE' 'DISCRET' ;
  416. WCH3 = WCH3 '+' SO3 ;
  417.  
  418. * Source term corresponding to the liquid momentum equation
  419.  
  420. S41 = ('MANU' 'CHPO' ($DOMTOT.'CENTRE') 2 'ULX' 0.d0 'ULY' 0.d0
  421. 'NATURE' 'DISCRET');
  422. S42 = DT *(ALP '-' 1.d0) * RL * g ;
  423. S42 = S42 * (-1.) ;
  424. S421 = 'EXCO' 'SCAL' S42 'ULZ' 'NATURE' 'DISCRET' ;
  425. S4 = S41 'ET' S421 ;
  426. DADP = 'EXCO' ('MOTS' 'UVX' 'UVY' 'UVZ') DADP ('MOTS' 'ULX' 'ULY'
  427. 'ULZ') 'NATURE' 'DISCRET' ;
  428. SO4 = (S4 '-' DADP) ;
  429. SO4 = 'NOMC' ('MOTS' 'ULX' 'ULY' 'ULZ') SO4 ('MOTS' 'ULX' 'ULY'
  430. 'ULZ') 'NATURE' 'DISCRET' ;
  431. WCH4 = WCH4 '+' SO4 ;
  432.  
  433. * Source term corresponding to the vapour total energy equation
  434.  
  435. VY = 'EXCO' 'UVZ' UV ;
  436. S5 = S32 * VY ;
  437. DADT = DT * PINT * (ALP '-' OALP) '/' OLDDT ;
  438. SO5 = (S5 '-' DADT) ;
  439. SO5 = 'EXCO' 'SCAL' SO5 'SCAL' 'NATURE' 'DISCRET' ;
  440. WCH5 = WCH5 '+' SO5 ;
  441.  
  442. * Source term corresponding to the liquid total equation
  443.  
  444. LY = 'EXCO' 'ULZ' UL ;
  445. S6 = S42 * LY ;
  446. SO6 = (S6 '+' DADT) ;
  447. SO6 = 'EXCO' 'SCAL' SO6 'SCAL' 'NATURE' 'DISCRET' ;
  448. WCH6 = WCH6 '+' SO6 ;
  449.  
  450. OALP = ALP ;
  451. OLDDT = DT ;
  452. OTV = TV ;
  453. OTL = TL ;
  454. ************************************************************************
  455. **** Updating primitive variables
  456. ************************************************************************
  457. PINT RL RV TL TV P UL UV ALP = 'PRIM' 'TWOFLUID'
  458. WCH1 WCH2 WCH3 WCH4 WCH5 WCH6 OALP OTV OTL CP CVM;
  459. * ANADIR CVM TRAS COMPILAR
  460. ************************************************************************
  461. **** Graphical outputs
  462. ************************************************************************
  463. 'SI' (((&BL1 '/' 1000) '*' 1000) 'EGA' &BL1) ;
  464. 'MESSAGE' &BL1 TPS ;
  465. 'FINSI' ;
  466.  
  467. 'SI' (TPS > TFINAL) ;
  468.  
  469. EVALP = 'EVOL' 'CHPO' ALP COURB ;
  470. EVUVX = 'EVOL' 'CHPO' UV 'UVX' COURB ;
  471. EVUVY = 'EVOL' 'CHPO' UV 'UVY' COURB ;
  472. EVUVZ = 'EVOL' 'CHPO' UV 'UVZ' COURB ;
  473. EVULX = 'EVOL' 'CHPO' UL 'ULX' COURB ;
  474. EVULY = 'EVOL' 'CHPO' UL 'ULY' COURB ;
  475. EVULZ = 'EVOL' 'CHPO' UL 'ULZ' COURB ;
  476. EVP = 'EVOL' 'CHPO' P COURB ;
  477. EVTV = 'EVOL' 'CHPO' TV COURB ;
  478. EVTL = 'EVOL' 'CHPO' TL COURB ;
  479. EVRV = 'EVOL' 'CHPO' RV COURB ;
  480. EVRL = 'EVOL' 'CHPO' RL COURB ;
  481.  
  482. 'MESSAGE' &BL1 TPS ;
  483.  
  484. 'SI' GRAPH ;
  485. 'DESSIN' EVALP 'TITR' 'Shock tube test' 'TITX' 'x (m)' 'TITY'
  486. 'alpha' ;
  487. 'DESSIN' EVUVX 'TITR' 'Shock tube test' 'TITX' 'x (m)' 'TITY'
  488. 'uvx (m/s)' ;
  489. 'DESSIN' EVUVY 'TITR' 'Shock tube test' 'TITX' 'x (m)' 'TITY'
  490. 'uvy (m/s)' ;
  491. 'DESSIN' EVUVZ 'TITR' 'Shock tube test' 'TITX' 'x (m)' 'TITY'
  492. 'uvz (m/s)' ;
  493. 'DESSIN' EVULX 'TITR' 'Shock tube test' 'TITX' 'x (m)' 'TITY'
  494. 'ulx (m/s)' ;
  495. 'DESSIN' EVULY 'TITR' 'Shock tube test' 'TITX' 'x (m)' 'TITY'
  496. 'uly (m/s)' ;
  497. 'DESSIN' EVULZ 'TITR' 'Shock tube test' 'TITX' 'x (m)' 'TITY'
  498. 'ulz (m/s)' ;
  499. 'DESSIN' EVP 'TITR' 'Shock tube test' 'TITX' 'x (m)' 'TITY'
  500. 'p (Pa)' ;
  501. 'DESSIN' EVTV 'TITR' 'Shock tube test' 'TITX' 'x (m)' 'TITY'
  502. 'Tv (K)' ;
  503. 'DESSIN' EVTL 'TITR' 'Shock tube test' 'TITX' 'x (m)' 'TITY'
  504. 'Tl (K)' ;
  505. 'DESSIN' EVRV 'TITR' 'Shock tube test' 'TITX' 'x (m)' 'TITY'
  506. 'rv (kg/m3)' ;
  507. 'DESSIN' EVRL 'TITR' 'Shock tube test' 'TITX' 'x (m)' 'TITY'
  508. 'rl (kg/m3)' ;
  509. 'FINSI' ;
  510.  
  511. 'QUITTER' BL1 ;
  512.  
  513. 'FINSI' ;
  514.  
  515.  
  516. 'FIN' BL1 ;
  517. ************************************************************************
  518. **** End of loop
  519. ************************************************************************
  520.  
  521. *
  522. 'SI' GRAPH ;
  523. MOD1 = 'MODELISER' ($DOMTOT . 'MAILLAGE' ) 'THERMIQUE';
  524.  
  525. CHM_V1 = 'KCHA' $DOMTOT 'CHAM' alp ;
  526. CHM_V2 = 'KCHA' $DOMTOT 'CHAM' uv ;
  527. CHM_V3 = 'KCHA' $DOMTOT 'CHAM' ul ;
  528. CHM_V4 = 'KCHA' $DOMTOT 'CHAM' p ;
  529. CHM_V5 = 'KCHA' $DOMTOT 'CHAM' Tv ;
  530. CHM_V6 = 'KCHA' $DOMTOT 'CHAM' TL ;
  531. CHM_V7 = 'KCHA' $DOMTOT 'CHAM' RV ;
  532. CHM_V8 = 'KCHA' $DOMTOT 'CHAM' RL ;
  533.  
  534. 'TRACER' CHM_V1 MOD1
  535. 'TITR' ('CHAINE' 'alpha at t=' TPS);
  536. 'TRACER' CHM_V2 MOD1
  537. 'TITR' ('CHAINE' 'uv at t=' TPS);
  538. 'TRACER' CHM_V3 MOD1
  539. 'TITR' ('CHAINE' 'ul at t=' TPS);
  540. 'TRACER' CHM_V4 MOD1
  541. 'TITR' ('CHAINE' 'p at t=' TPS);
  542. 'TRACER' CHM_V5 MOD1
  543. 'TITR' ('CHAINE' 'Tv at t=' TPS);
  544. 'TRACER' CHM_V6 MOD1
  545. 'TITR' ('CHAINE' 'Tl at t=' TPS);
  546. 'TRACER' CHM_V7 MOD1
  547. 'TITR' ('CHAINE' 'rv at t=' TPS);
  548. 'TRACER' CHM_V8 MOD1
  549. 'TITR' ('CHAINE' 'rl at t=' TPS);
  550.  
  551. 'FINSI' ;
  552. ************************************************************************
  553. *** saving the results in a file
  554. ************************************************************************
  555.  
  556. * 'OPTION' 'SAUVER' ('CHAINE' 'shock3d.sauv') ;
  557. * 'SAUVER' ;
  558.  
  559. 'FIN' ;
  560.  
  561.  
  562.  
  563.  
  564.  
  565.  
  566.  

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