Télécharger bc30.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : bc30.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. *$$$ BC30
  5. *X BC30 (Jeux de donnees)
  6.  
  7. **
  8. ** BC30
  9. **
  10. ** CANAL CHAUFFE INCLINE ANGL = 30. LONGUEUR H0=10 LARGEUR 1.
  11. ** K-E PLAN --- RE=2.E+5
  12. **
  13. ** teste NSKE(SUPG) TSCAL(PSI) LAPN FPU ECHI et VNIMP
  14. ** Formulation EFM1 algorithme semi implicite
  15. ** teste aussi les méthodes d'inversion de la pression :
  16. ** Choleski et gradient conjugué.
  17. **
  18. **
  19. ** --- 26 JANVIER 1987 --- modifie --- 26 MAI 1993 ---
  20. ** modifie --- 17 DECEMBRE 1997 ---
  21. **
  22. ** A noter :
  23. **
  24. ** 1/ la manière dont sont donnés les différents coefficients
  25. ** du fait que df/dt a comme coefficient 1. pour toutes les
  26. ** équations.
  27. ** 2/ La manière dont est réalisé le couplage : construction
  28. ** du maillage de raccord (uniquement valable en 2D)
  29. ** et procédure CALCUL .
  30. ** 3/ Le bilan thermique dans ce cas est assez mal vérifié .
  31. ** Les causes : maillage grossier, flux thermique tres faible
  32. ** par rapport au debit enthalpique, terme d'accumulation
  33. ** d teta / dt non pris en compte dans le bilan.
  34. **
  35. **
  36.  
  37. GRAPH = 'N' ;
  38. COMPLET = FAUX ;
  39.  
  40. NITMA=200;
  41.  
  42. SI (NON COMPLET) ;
  43. err1=5.E-4;
  44. NITMA = 40 ;
  45. FINSI ;
  46.  
  47. option dime 2 elem tri6 ;
  48. option dime 2 elem qua8 ;
  49. p1=0 0 ;
  50. p2=1. 0.;
  51.  
  52. nbpe=7 ; nbpp=10 ;
  53. *nbpe=40 ; nbpp=80 ;
  54. *nbpe=50 ; nbpp=150 ;
  55. tole=1.e-5 ;
  56. H0=10.;
  57.  
  58. entree= p1 d nbpe p2 ;
  59. *entree= p1 d dini 0.05 dfin 0.2 p2 ;
  60. c1= 5. 0. ;
  61.  
  62. q1=p1 tour c1 -90. ;
  63. q2=p2 tour c1 -90. ;
  64. pp1=p1 c c1 q1 nbpp ;
  65. pp2=p2 c c1 q2 nbpp ;
  66.  
  67. q1=p1 plus (0 H0) ;
  68. q2=p2 plus (0 H0) ;
  69. pp1=p1 d q1 nbpp ;
  70. eps=0.1 ;
  71. pri=pp1 plus ((-1.*eps) 0) ;
  72. pp2=p2 d q2 nbpp ;
  73.  
  74. sortie=entree tour c1 -90 ;
  75. sortie=entree plus (0 H0) ;
  76. pp1=inve pp1 ;
  77. sortie=inve sortie;
  78. elim tole (sortie et pp1 et pp2 );
  79. cont=entree et pp2 et sortie et pp1 ;
  80. *bell=surf cont ;
  81. bell=daller entree pp2 sortie pp1 ;
  82. paroi=pri trans 4 (-0.1 0) ;
  83. pr= cote 3 paroi ;
  84. entref=entree elem seg3 (lect 2 pas 1 6) ;
  85.  
  86. bell=bell tour c1 -20. ;
  87. entree=entree tour c1 -20. ;
  88. entref=entref tour c1 -20. ;
  89. sortie=sortie tour c1 -20. ;
  90. pp1=pp1 tour c1 -20. ;
  91. pp2=pp2 tour c1 -20. ;
  92.  
  93. paroi=paroi tour c1 -20. ;
  94. paroi=paroi coul rouge ;
  95. pr= pr tour c1 -20. ;
  96. pri= pri tour c1 -20. ;
  97. MT=bell et paroi ;
  98. pp1=inve pp1 ;
  99.  
  100. elim (pp1 et pp2 et entree et entref et sortie et bell) tole ;
  101. elim ( pri et pr et paroi ) tole ;
  102.  
  103. titred= 'PLAQUES // COUDEES : TEST TRIO-EF NSKE FPU VNIMP' ;
  104. titre= 'PLAQUES // COUDEES : TEST TRIO-EF NSKE FPU VNIMP' ;
  105. Mmt = CHAN mt QUAF ;
  106. Mbell = CHAN bell QUAF ;
  107. Mparoi = CHAN paroi QUAF ;
  108. Mentree= CHAN entree QUAF ;
  109. Mentref= CHAN entref QUAF ;
  110. Msortie= CHAN sortie QUAF ;
  111. Mpri = CHAN pri QUAF ;
  112. Mpp1 = CHAN pp1 QUAF ;
  113. Mpp2 = CHAN pp2 QUAF ;
  114.  
  115. elim (Mmt et Mbell et Mparoi et Mentree et Mentref et Msortie
  116. et Mpri et Mpp1 et Mpp2) tole ;
  117.  
  118. $mt = MODE MMT 'NAVIER_STOKES' 'MACRO' ;
  119. $bell = MODE Mbell 'NAVIER_STOKES' 'MACRO' ;
  120. $paroi= MODE Mparoi 'NAVIER_STOKES' 'MACRO' ;
  121. paroi = doma $paroi 'MAILLAGE' ;
  122. DOMA $mt 'IMPR' ; DOMA $bell 'IMPR' ; DOMA $paroi 'IMPR' ;
  123.  
  124. $pri = MODE Mpri 'NAVIER_STOKES' 'MACRO' ;
  125. $pp1 = MODE Mpp1 'NAVIER_STOKES' 'MACRO' ;
  126. $pp2 = MODE Mpp2 'NAVIER_STOKES' 'MACRO' ;
  127. $entree = MODE Mentree 'NAVIER_STOKES' 'MACRO' ;
  128. $entref = MODE Mentref 'NAVIER_STOKES' 'MACRO' ;
  129. $sortie = MODE Msortie 'NAVIER_STOKES' 'MACRO' ;
  130.  
  131. pp1c=pp1 elem (lect 2 pas 1 (nbel pp1));
  132. pric=pri elem (lect 2 pas 1 (nbel pri));
  133.  
  134. Mpp1c = CHAN pp1c QUAF ;
  135. Mpric = CHAN pric QUAF ;
  136. elim (Mpric et Mpp1c et Mmt ) tole ;
  137.  
  138. $pp1c=MODE Mpp1c 'NAVIER_STOKES' 'MACRO' ;
  139. $pric=MODE Mpric 'NAVIER_STOKES' 'MACRO' ;
  140. pp1l=doma $pp1c MAILLAGE ;
  141. pril=doma $pric MAILLAGE ;
  142.  
  143. rfp=racc eps pp1l pril ;
  144. rfpl= rfp chan ligne ;
  145. rfp= diff rfpl (pp1l et pril) ;
  146. rfpi= inve rfp ;
  147. ************************* PROCEDURE CALCUL *****************
  148. DEBP CALCUL ;
  149. ARGU RX*TABLE ;
  150. rv=rx.'EQEX' ;
  151.  
  152. tn=rv.'INCO'.'TN' ;
  153. un=rv.'INCO'.'UN' ;
  154.  
  155. TF=redu tn pp1c ;
  156. rv.'INCO'.'TF'= kcht $pric scal sommet (tf kpro rfp) ;
  157. TP=redu tn pric ;
  158. rv.'INCO'.'TP'= kcht $pp1c scal sommet (tp kpro rfpi) ;
  159. rv.'INCO'.'UF'= kcht $bell vect sommet (redu un bell) ;
  160.  
  161. as2 ama1 = 'KOPS' 'MATRIK' ;
  162. finproc as2 ama1 ;
  163.  
  164. ************************* PROCEDURE CALCUL *****************
  165.  
  166. *option donn 5 ;
  167. NU =7.3E-7 ;
  168. LAMBDAP=15.2 ;
  169. ROCPP=3.7E6 ;
  170. ALFP=LAMBDAP/ROCPP ;
  171. *ALFP=4.06E-6 ;
  172. LAMBDAE=0.62 ;
  173. ROCPE=4.15E6 ;
  174. ALFE=LAMBDAE/ROCPE ;
  175. *ALFE=1.5E-7 ;
  176. mess ' ALFE=' ALFE ' ALFP=' ALFP ;
  177. YP =5.E-2 ;
  178. R0 =0.5 ;
  179. TPSC =1540/span>.;
  180. TREF=0. ;
  181. ZERO =1.E-20 ;
  182. gb = 0. 0. ;
  183. H=1.E4 ;
  184. HF=H/ROCPE ;
  185. HP=H/ROCPP ;
  186. HP=H/ROCPE ;
  187. CMU = 0.09 ;
  188. UREF = 4. ;
  189. KE = (0.05 * UREF)*(0.05 * UREF) ;
  190. EE = (KE**1.5) / R0 ;
  191. nut = 0.09 * KE * KE / EE ;
  192. mess ' NUT entree ' nut ;
  193. UET1 = KCHT $PP1 SCAL CENTRE 5.E-2 ;
  194. UET2 = KCHT $PP2 SCAL CENTRE 5.E-2 ;
  195. NUT = KCHT $BELL SCAL CENTRE nut ;
  196.  
  197. RV= EQEX 'DUMP' 'ITMA' NITMA 'ALFA' 0.9 OPTI 'SUPG'
  198. 'ZONE' $MT 'OPER' 'CALCUL'
  199. 'ZONE' $BELL 'OPER' 'NSKE' NU NUT 'INCO' 'UN' 'KN' 'EN'
  200. 'ZONE' $PP1 'OPER' 'FPU' NU UET1 YP 'INCO' 'UN' 'KN' 'EN'
  201. 'ZONE' $PP2 'OPER' 'FPU' NU UET2 YP 'INCO' 'UN' 'KN' 'EN'
  202. ;
  203. RV=EQEX RV OPTI EFM1 'PSI'
  204. 'ZONE' $BELL 'OPER' 'TSCA' ALFE 'UF' 0. NUT 0.7 'INCO' 'TN'
  205. 'ZONE' $PAROI 'OPER' 'LAPN' ALFP 'INCO' 'TN'
  206. 'ZONE' $PRIC 'OPER' 'ECHI' HP 'TF' 'INCO' 'TN'
  207. 'ZONE' $PP1C 'OPER' 'ECHI' HF 'TP' 'INCO' 'TN'
  208. ;
  209. RV=EQEX RV OPTI EFM1 'CENTREE'
  210. 'ZONE' $BELL 'OPER' 'DFDT' 1. 'UN' 'DELTAT' 'INCO' 'UN'
  211. 'ZONE' $MT 'OPER' 'DFDT' 1. 'TN' 'DELTAT' 'INCO' 'TN'
  212. 'ZONE' $BELL 'OPER' 'DFDT' 1. 'KN' 'DELTAT' 'INCO' 'KN'
  213. 'ZONE' $BELL 'OPER' 'DFDT' 1. 'EN' 'DELTAT' 'INCO' 'EN'
  214. ;
  215. RV=EQEX RV CLIM
  216. 'KN' TIMP ENTREE KE 'EN' TIMP ENTREE EE
  217. 'TN' TIMP PR 100. 'TN' TIMP ENTREE 20.
  218. ;
  219. RV.INCO=TABLE 'INCO' ;
  220. UF = kcht $BELL VECT SOMMET (UREF*(0.5 0.9));
  221. RV.INCO.'UN'=kcht $BELL VECT SOMMET (0. 0.) UF ;
  222. RV.INCO.'KN'=kcht $BELL SCAL SOMMET KE ;
  223. RV.INCO.'EN'=kcht $BELL SCAL SOMMET EE ;
  224. tp=kcht $paroi scal sommet 100. ;
  225. RV.INCO.'TN'=kcht $MT SCAL SOMMET 20. tp ;
  226.  
  227. * Suivi temporel de quelques points
  228. lh= elem (doma $bell sommet) 'POI1' (lect 1 pas 20 101) ;
  229. lc1= elem (doma $pp1 sommet) 'POI1' (lect 4 pas 4 20) ;
  230. lc2= elem (doma $pri sommet) 'POI1' (lect 4 pas 4 20) ;
  231. lc=lc1 et lc2 ;
  232. his = khis 'UN' 1 lh
  233. 'UN' 2 lc1
  234. 'KN' lh
  235. 'EN' lh
  236. 'TN' lc1 ;
  237. rv.'HIST'=his ;
  238.  
  239.  
  240.  
  241.  
  242. DEBPROC TEST ityp*entier ksto*entier rv*table ;
  243.  
  244. * 1 : Choleski (valeur par defaut)
  245.  
  246. **** Le reste est non maintenu et ne fonctionne pas sur machine 64bytes**
  247. * 2 : gradient conjugue sans preconditionnement
  248. * 3 : gradient conjugue avec preconditionnement diagonal
  249. * 4 : gradient conjugue avec preconditionnement par Choleski incomplet
  250. * 7 : gradient conjugue avec preconditionnement diagonal (variante)
  251. *
  252. *Pour 2,3,4, la matrice de pression est calculee et stockee une fois
  253. *pour toutes. Pour 7, elle n'est jamais stockee explicitement et permet
  254. *de ce fait d'atteindre les plus grosses tailles de problemes.
  255. *
  256. *ityp=1 ;
  257. * obtenir ityp*entier ;
  258. *
  259. *Pour 2,3,4, on propose deux types de stockage : morse et compresse. Le
  260. *second type n'est interessant que pour les machines vectorielles.
  261. * ksto=0 morse ksto=1 compresse
  262. * ksto=0 Valeur par defaut
  263. *
  264. *ksto=0 ;
  265. *obtenir ksto*entier ;
  266.  
  267. RVP = EQPR $BELL KTYPI ITYP 'ZONE' $BELL 'OPER' 'PRESSION' 0.
  268. 'ZONE' $PP1 'OPER' 'VNIMP' 0. 'ZONE' $PP2 'OPER' 'VNIMP' 0.
  269. 'ZONE' $ENTREE 'OPER' 'VNIMP' UREF 'ZONE' $ENTREF 'OPER' 'VTIMP' 0.;
  270.  
  271. si (non (ega ityp 1));
  272. rvp.METHODE.KSTOCK=ksto ;
  273. rvp.METHODE.JTEPS=5.E-4;
  274. rvp.METHODE.'IPAT'=1;
  275. finsi ;
  276.  
  277. RV.'PRESSION' =RVP ;
  278.  
  279. EXEC RV ;
  280.  
  281. mess ' Fin PROC ' ;
  282. FINPROC RV ;
  283.  
  284. DEBPROC POST RV*TABLE GRAPH*MOT ;
  285.  
  286. un=rv.'INCO'.'UN' ;
  287. tn=rv.'INCO'.'TN' ;
  288.  
  289. qe=dbit un $entree ;
  290. qs=dbit un $sortie ;
  291. dq=(abs qs) - (abs qe) ;
  292. dqr=dq / (abs qe) ;
  293. mess ' Erreur sur la masse ' dqr ;
  294. ER1= abs (abs(dqr) - 1.E-10) ;
  295.  
  296. srti=inve (doma $sortie MAILLAGE) ;
  297. lcu = extr un 'COMP';
  298. mu2=un lcu 'PSCA' un lcu;
  299. mu= mu2 ** 0.5 ;
  300. evolV = EVOL 'CHPO' mu SCAL srti ;
  301. evy=extr evolV 'ORDO' ;
  302. list evy ;
  303. ley=prog 3.4573 3.8455 4.0470 4.0806 4.0848 4.0855 4.0853 4.0852
  304. 4.0853 4.0855 4.0848 4.0806 4.0470 3.8455 3.4573 ;
  305.  
  306.  
  307. ER2=SOMM( (abs (ley - evy)) / 4. ) ;
  308. mess ' Ecart sur profil de vitesse en sortie ' er2 ;
  309.  
  310. evolT = EVOL 'CHPO' tn SCAL srti ;
  311. evt=extr evolT 'ORDO' ;
  312. list evt ;
  313. let=prog 20.985 20.501 20.179 20.078 20.045 20.026 20.015 20.008
  314. 20.005 20.002 20.001 20.001 20.000 20.000 20.000 ;
  315.  
  316.  
  317. ER3=SOMM( (abs (let - evt)) / 20. ) ;
  318. mess ' Ecart sur profil de température en sortie ' er3 ;
  319.  
  320. D=doma $pp1c 'XXDIAGSI' ;
  321. TF= kcht $pp1c scal sommet (redu tn pp1c) ;
  322. TP=rv.inco.'TP' ;
  323. DTFP=kops tp '-' tf ;
  324. FP=somt ( hf* (kops d '*' (kops tp '-' tf) ) ) ;
  325. mess 'FP=' fp ;
  326. uteta= kops ('KCHT' $bell 'SCAL' 'SOMMET' tn) '*' un ;
  327. ute=dbit uteta $entree ;
  328. uts=dbit uteta $sortie ;
  329. dut=(abs uts) - (abs ute) - fp ;
  330. dutr=dut/FP ;
  331. mess (' Bilan : erreur relative = ') dutr ;
  332. ER4= abs (abs(dutr) - .58735) ;
  333.  
  334. si (non COMPLET ) ;
  335. si ( ER1 > 1.E-9 ) ; erreur 5 ; finsi ;
  336. si ( ER2 > err1 ) ; erreur 5 ; finsi ;
  337. si ( ER3 > err1 ) ; erreur 5 ; finsi ;
  338. si ( ER4 > err1 ) ; erreur 5 ; finsi ;
  339. finsi ;
  340.  
  341. SI ('EGA' GRAPH 'O' );
  342.  
  343. TAB1=TABLE;
  344. TAB1.'TITRE'=TABLE ;
  345. TAB1 . 1 ='MARQ REGU ' ;
  346. TAB1.'TITRE' . 1 = mot 'Module de U ' ;
  347. DESS evolV 'TITX' 'R (m)' 'TITY' 'V (m/s)' LEGE TAB1 ;
  348.  
  349. TAB1.'TITRE' . 1 = mot 'Temperature sortie ' ;
  350. DESS evolT 'TITX' 'R (m)' 'TITY' 'T °C' LEGE TAB1 ;
  351.  
  352.  
  353. ung1=vect un 0.1 ux uy jaune;
  354. pn=elno (kcht $bell scal centre (rv.'PRESSION'.pression)) $bell ;
  355.  
  356. his=rv.'HIST' ;
  357. dessin his.'TABD' his.'1UN' ;
  358. dessin his.'TABD' his.'2UN' ;
  359. dessin his.'TABD' his.'TN' ;
  360. dessin his.'TABD' his.'KN' ;
  361. dessin his.'TABD' his.'EN' ;
  362.  
  363. trace ung1 bell ;
  364. trace pn bell ;
  365. trace tn bell ;
  366. trace tn paroi;
  367. kn=rv.inco.'KN' ;
  368. en=rv.inco.'EN' ;
  369. nut=0.09*kn*kn/en ;
  370. trace nut bell ;
  371. trace kn bell ;
  372. trace en bell ;
  373. FINSI ;
  374.  
  375. FINPROC ;
  376.  
  377. ******************** TESTS **************************
  378.  
  379. ityp=1 ;
  380. ksto = 0 ;
  381. rv = test ityp ksto rv ;
  382.  
  383. * ityp=4 ;non maintenu, ne marche pas sur machine 64bytes
  384. ityp=1 ;
  385. ksto = 0 ;
  386. rv = test ityp ksto rv ;
  387.  
  388. * ityp=4 ;non maintenu, ne marche pas sur machine 64bytes
  389. ityp=1 ;
  390. ksto = 1 ;
  391. rv = test ityp ksto rv ;
  392.  
  393. * ityp=3 ;non maintenu, ne marche pas sur machine 64bytes
  394. ityp=1 ;
  395. ksto = 1 ;
  396. rv = test ityp ksto rv ;
  397.  
  398. * ityp=2 ;non maintenu, ne marche pas sur machine 64bytes
  399. ityp=1 ;
  400. ksto = 1 ;
  401. rv = test ityp ksto rv ;
  402.  
  403. post rv graph ;
  404. FIN ;
  405.  
  406.  
  407.  
  408.  
  409.  
  410.  
  411.  
  412.  
  413.  
  414.  
  415.  
  416.  
  417.  
  418.  
  419.  
  420.  
  421.  

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