Télécharger source3.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : SOURCE3.dgibi
  2. * section : thermique
  3. *----------------------------------------------------------------------*
  4. * SOURCE3.DGIBI *
  5. *----------------------------------------------------------------------*
  6. *
  7. * Objet :
  8. * -------
  9. *
  10. * Verfication / validation d'un modele de source de chaleur.
  11. * Cas d'une SOURCE GAUSSIENNE ELLIPTIQUE.
  12. *
  13. * Description :
  14. * -------------
  15. * Comparaison des flux nodaux equivalents obtenus avec le modele a
  16. * ceux obtenus en construisant le champ "a la main", puis en l'integrant
  17. * avec l'operateur SOURCE.
  18. *
  19. * Type de calcul : Aucun
  20. * Mode de calcul : 2D PLAN, AXIS et 3D
  21. * Type d'element : TRI3 TRI6 QUA4 QUA8 TET4 TE10 PYR5 PY13 PRI6 PR15
  22. * CUB8 CU20
  23. * Objectifs : Ecart relatif entre flux integres < 1.e-12
  24. * Ecart relatif entre la chaleur totale fournie (QTOT)
  25. * et la resultante sur le champ integre < 1.e-12
  26. *
  27. *----------------------------------------------------------------------*
  28. *
  29. * IG1 vrai : traces actives
  30. IG1 = faux ;
  31. *
  32. *------------------------ 2D ELEMENTS LINEAIRES -----------------------*
  33. *
  34. opti dime 2 elem qua4 ;
  35.  
  36. * Parametres du maillage :
  37. lo1 = 50.e-3 ;
  38. ep1 = 20.e-3 ;
  39.  
  40. * Parametres de la source :
  41. QT1 = 1.e3 ;
  42. OR1 = 30.e-3 ep1 ;
  43. DR1 = 0. 2. ;
  44. RG1 = 5.e-3 ;
  45. ZG1 = 9.e-3 ;
  46.  
  47. * Maillage :
  48. l1 = (0 0) droi 10 (0 ep1) ;
  49. s1 = l1 tran 50 (lo1 0) ;
  50. l2 = s1 cote 3 ;
  51. s2 = l2 tran 50 (lo1 ep1) ;
  52. s2 = s2 chan tri3 ;
  53. s0 = s1 et s2 ;
  54.  
  55. si IG1 ;
  56. ltyp1 = s0 elem type ;
  57. mot1 = mot 'ELEM =' ;
  58. repe bm1 (dime ltyp1) ;
  59. mot1 = chai mot1 ' ' (ltyp1 extr &bm1) ;
  60. fin bm1 ;
  61. mot1 = chai ' Maillage test en 2D (' mot1 ')' ;
  62. trac s0 titr mot1 ;
  63. fins;
  64.  
  65. * MODE PLAN :
  66. * -----------
  67. *
  68. opti mode plan ;
  69. *
  70. moq1 = mode s0 thermique source gaussienne elliptique ;
  71. maq1 = mate moq1 'QTOT' QT1 'ORIG' OR1 'RGAU' RG1 'DIRE' DR1 'ZGAU' ZG1 ;
  72. sq1 = sour moq1 maq1 ;
  73.  
  74. qtmoq1 = sq1 resu maxi ;
  75.  
  76. mo1 = mode s0 thermique ;
  77. ma1 = mate mo1 k 25. ;
  78.  
  79. x1 y1 = s0 coor ;
  80. x1 y1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) ;
  81. x1 y1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) ;
  82. r1 = (x1 * x1) + (y1 * y1) ;
  83. z1 = (x1 * (DR1 coor 1)) + (y1 * (DR1 coor 2)) / (norm DR1) ;
  84. r1 = r1 - (z1 * z1) ;
  85. xqx1 = 4. / PI / RG1 / ZG1 ;
  86. sqref1 = exp (-2. * ((r1 / RG1 / RG1) + (z1 * z1 / ZG1 / ZG1))) * QT1 * xqx1 ;
  87. sqref1 = sour mo1 sqref1 ;
  88.  
  89. si IG1 ;
  90. trac sq1 S0 titr ' Flux nodaux du champ de source Gaussienne du MODELE ';
  91. trac sqref1 s0 titr ' Flux nodaux du champ de source Gaussienne REFERENCE' ;
  92. trac (sq1 - sqref1) s0 titr ' Difference des 2 champs' ;
  93. fins ;
  94.  
  95. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  96. list err1 ;
  97. err2 = (abs (qtmoq1 - qt1)) / (abs qt1) ;
  98. list err2 ;
  99. si ((err1 > 1.e-12) ou (err2 > 1.e-4)) ;
  100. ERRE ' Probleme MODE SOURCE' ;
  101. fins ;
  102.  
  103. * MODE AXIS :
  104. * -----------
  105. *
  106. opti mode axis ;
  107. *
  108. moq1 = mode s0 thermique source gaussienne elliptique ;
  109. maq1 = mate moq1 'QTOT' QT1 'ORIG' OR1 'RGAU' RG1 'DIRE' DR1 'ZGAU' ZG1 ;
  110. sq1 = sour moq1 maq1 ;
  111.  
  112. qtmoq1 = sq1 resu maxi ;
  113.  
  114. mo1 = mode s0 thermique ;
  115. ma1 = mate mo1 k 25. ;
  116.  
  117. x0 y0 = (or1 coor 1) (or1 coor 2) ;
  118. x1 y1 = s0 coor ;
  119. x1 y1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) ;
  120. x1 y1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) ;
  121. r1 = (x1 * x1) + (y1 * y1) ;
  122. z1 = (x1 * (DR1 coor 1)) + (y1 * (DR1 coor 2)) / (norm DR1) ;
  123. r1 = r1 - (z1 * z1) ;
  124. rc2 = 2. ** 0.5 ;
  125. rcpi1 = pi ** 0.5 ;
  126. xik1 = 0.25 * RG1 * RG1 * (exp (-2. * x0 * x0 / RG1 / RG1)) ;
  127. xik2 = 0.5 * RG1 * x0 * rcpi1 / rc2 ;
  128. xik3 = RG1 * x0 * (erf (x0 * rc2 / RG1)) * rcpi1 / (rc2 * rc2 * rc2) ;
  129. xqx1 = PI * ZG1 * rcpi1 / rc2 ;
  130. xqx1 = xqx1 * (xik1 + xik2 + xik3) ;
  131. sqref1 = exp (-2. * ((r1 / RG1 / RG1) + (z1 * z1 / ZG1 / ZG1))) * QT1 / xqx1 ;
  132. sqref1 = sour mo1 sqref1 ;
  133.  
  134. si IG1 ;
  135. trac sq1 S0 titr ' Flux nodaux du champ de source Gaussienne du MODELE ';
  136. trac sqref1 s0 titr ' Flux nodaux du champ de source Gaussienne REFERENCE' ;
  137. trac (sq1 - sqref1) s0 titr ' Difference des 2 champs' ;
  138. fins ;
  139.  
  140. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  141. list err1 ;
  142. err2 = (abs (qtmoq1 - qt1)) / (abs qt1) ;
  143. list err2 ;
  144. si ((err1 > 1.e-12) ou (err2 > 1.e-4)) ;
  145. ERRE ' Probleme MODE SOURCE' ;
  146. fins ;
  147.  
  148. *---------------------- 2D ELEMENTS QUADRATIQUES ----------------------*
  149.  
  150. * Passage en EF quadratiques :
  151. s0 = s0 chan quad ;
  152.  
  153. si IG1 ;
  154. ltyp1 = s0 elem type ;
  155. mot1 = mot 'ELEM =' ;
  156. repe bm1 (dime ltyp1) ;
  157. mot1 = chai mot1 ' ' (ltyp1 extr &bm1) ;
  158. fin bm1 ;
  159. mot1 = chai ' Maillage test en 2D (' mot1 ')' ;
  160. trac s0 titr mot1 ;
  161. fins;
  162.  
  163. * MODE PLAN :
  164. * -----------
  165. *
  166. opti mode plan ;
  167. *
  168. moq1 = mode s0 thermique source gaussienne elliptique ;
  169. maq1 = mate moq1 'QTOT' QT1 'ORIG' OR1 'RGAU' RG1 'DIRE' DR1 'ZGAU' ZG1 ;
  170. sq1 = sour moq1 maq1 ;
  171.  
  172. qtmoq1 = sq1 resu maxi ;
  173.  
  174. mo1 = mode s0 thermique ;
  175. ma1 = mate mo1 k 25. ;
  176.  
  177. x1 y1 = s0 coor ;
  178. x1 y1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) ;
  179. x1 y1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) ;
  180. r1 = (x1 * x1) + (y1 * y1) ;
  181. z1 = (x1 * (DR1 coor 1)) + (y1 * (DR1 coor 2)) / (norm DR1) ;
  182. r1 = r1 - (z1 * z1) ;
  183. xqx1 = 4. / PI / RG1 / ZG1 ;
  184. sqref1 = exp (-2. * ((r1 / RG1 / RG1) + (z1 * z1 / ZG1 / ZG1))) * QT1 * xqx1 ;
  185. sqref1 = sour mo1 sqref1 ;
  186.  
  187. si IG1 ;
  188. trac sq1 S0 titr ' Flux nodaux du champ de source Gaussienne du MODELE ';
  189. trac sqref1 s0 titr ' Flux nodaux du champ de source Gaussienne REFERENCE' ;
  190. trac (sq1 - sqref1) s0 titr ' Difference des 2 champs' ;
  191. fins ;
  192.  
  193. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  194. list err1 ;
  195. err2 = (abs (qtmoq1 - qt1)) / (abs qt1) ;
  196. list err2 ;
  197. si ((err1 > 1.e-12) ou (err2 > 1.e-4)) ;
  198. ERRE ' Probleme MODE SOURCE' ;
  199. fins ;
  200.  
  201. * MODE AXIS :
  202. * -----------
  203. *
  204. opti mode axis ;
  205. *
  206. moq1 = mode s0 thermique source gaussienne elliptique ;
  207. maq1 = mate moq1 'QTOT' QT1 'ORIG' OR1 'RGAU' RG1 'DIRE' DR1 'ZGAU' ZG1 ;
  208. sq1 = sour moq1 maq1 ;
  209.  
  210. qtmoq1 = sq1 resu maxi ;
  211.  
  212. mo1 = mode s0 thermique ;
  213. ma1 = mate mo1 k 25. ;
  214.  
  215. x0 y0 = (or1 coor 1) (or1 coor 2) ;
  216. x1 y1 = s0 coor ;
  217. x1 y1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) ;
  218. x1 y1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) ;
  219. r1 = (x1 * x1) + (y1 * y1) ;
  220. z1 = (x1 * (DR1 coor 1)) + (y1 * (DR1 coor 2)) / (norm DR1) ;
  221. r1 = r1 - (z1 * z1) ;
  222. rc2 = 2. ** 0.5 ;
  223. rcpi1 = pi ** 0.5 ;
  224. xik1 = 0.25 * RG1 * RG1 * (exp (-2. * x0 * x0 / RG1 / RG1)) ;
  225. xik2 = 0.5 * RG1 * x0 * rcpi1 / rc2 ;
  226. xik3 = RG1 * x0 * (erf (x0 * rc2 / RG1)) * rcpi1 / (rc2 * rc2 * rc2) ;
  227. xqx1 = PI * ZG1 * rcpi1 / rc2 ;
  228. xqx1 = xqx1 * (xik1 + xik2 + xik3) ;
  229. sqref1 = exp (-2. * ((r1 / RG1 / RG1) + (z1 * z1 / ZG1 / ZG1))) * QT1 / xqx1 ;
  230. sqref1 = sour mo1 sqref1 ;
  231.  
  232. si IG1 ;
  233. trac sq1 S0 titr ' Flux nodaux du champ de source Gaussienne du MODELE ';
  234. trac sqref1 s0 titr ' Flux nodaux du champ de source Gaussienne REFERENCE' ;
  235. trac (sq1 - sqref1) s0 titr ' Difference des 2 champs' ;
  236. fins ;
  237.  
  238. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  239. list err1 ;
  240. err2 = (abs (qtmoq1 - qt1)) / (abs qt1) ;
  241. list err2 ;
  242. si ((err1 > 1.e-12) ou (err2 > 1.e-4)) ;
  243. ERRE ' Probleme MODE SOURCE' ;
  244. fins ;
  245.  
  246. *------------------------ 3D ELEMENTS LINEAIRES -----------------------*
  247. *
  248. opti dime 3 elem cub8 ;
  249. s0 = s0 chan line ;
  250.  
  251. * Maillage du volume :
  252. v0 = s0 volu tran 5 (0 0 ep1) ;
  253. x1 = v0 coor 1 ;
  254.  
  255. * Ajour de tetraedres :
  256. v3 = v0 elem appu stri (x1 poin supe (1.5 * lo1)) ;
  257. v0 = v0 diff v3 ;
  258. v0 = v0 et (v3 chan tet4) ;
  259.  
  260. * Ajout d'un element pyramide :
  261. el1 = (v0 elem cub8) elem 1 ;
  262. pe1 = el1 chan poi1 ;
  263. ps1 = pe1 poin 8 ;
  264. py1 = manu pyr5 (pe1 poin 1) (pe1 poin 2) (pe1 poin 3) (pe1 poin 4) ps1 ;
  265. py2 = manu pyr5 (pe1 poin 1) (pe1 poin 2) (pe1 poin 6) (pe1 poin 5) ps1 ;
  266. py3 = manu pyr5 (pe1 poin 2) (pe1 poin 3) (pe1 poin 7) (pe1 poin 6) ps1 ;
  267. py0 = py1 et py2 et py3 ;
  268. v0 = (v0 diff el1) et py0 ;
  269.  
  270. si IG1 ;
  271. ltyp1 = v0 elem type ;
  272. mot1 = mot 'ELEM =' ;
  273. repe bm1 (dime ltyp1) ;
  274. mot1 = chai mot1 ' ' (ltyp1 extr &bm1) ;
  275. fin bm1 ;
  276. mot1 = chai ' Maillage test en 3D (' mot1 ')' ;
  277. trac v0 titr mot1 ;
  278. fins;
  279.  
  280. * Parametres de la source en 3D
  281. OR1 = 30.e-3 (0.5 * ep1) ep1 ;
  282. DR1 = 0. 0. 1. ;
  283.  
  284. * Modeles :
  285. moq1 = mode V0 thermique source gaussienne elliptique ;
  286. maq1 = mate moq1 'QTOT' QT1 'ORIG' OR1 'RGAU' RG1 'DIRE' DR1 'ZGAU' ZG1 ;
  287. sq1 = sour moq1 maq1 ;
  288.  
  289. qtmoq1 = sq1 resu maxi ;
  290.  
  291. mo1 = mode V0 thermique ;
  292. ma1 = mate mo1 k 25. ;
  293.  
  294. x1 y1 z1 = V0 coor ;
  295. x1 y1 z1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) (z1 - (or1 coor 3)) ;
  296. x1 y1 z1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) (chan cham z1 mo1 stresses) ;
  297. r1 = (x1 * x1) + (y1 * y1) + (z1 * z1) ;
  298. z1 = (x1 * (DR1 coor 1)) + (y1 * (DR1 coor 2)) + (z1 * (DR1 coor 3)) / (norm DR1) ;
  299. r1 = r1 - (z1 * z1) ;
  300. xqx1 = (32. / (PI * PI * PI)) ** 0.5 / RG1 / RG1 / ZG1 ;
  301. sqref1 = exp (-2. * ((r1 / RG1 / RG1) + (z1 * z1 / ZG1 / ZG1))) * QT1 * xqx1 ;
  302. sqref1 = sour mo1 sqref1 ;
  303.  
  304. si IG1 ;
  305. trac sq1 V0 titr ' Flux nodaux du champ de source Gaussienne du MODELE ' coup or1 (0 1 0) (0 0 1) ;
  306. trac sqref1 V0 titr ' Flux nodaux du champ de source Gaussienne REFERENCE' coup or1 (0 1 0) (0 0 1) ;
  307. trac (sq1 - sqref1) V0 titr ' Difference des 2 champs' coup or1 (0 1 0) (0 0 1) ;
  308. fins ;
  309.  
  310. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  311. list err1 ;
  312. err2 = (abs (qtmoq1 - qt1)) / (abs qt1) ;
  313. list err2 ;
  314. si ((err1 > 1.e-12) ou (err2 > 1.e-4)) ;
  315. ERRE ' Probleme MODE SOURCE' ;
  316. fins ;
  317.  
  318. *---------------------- 3D ELEMENTS QUADRATIQUES ----------------------*
  319. *
  320. v0 = v0 chan quad ;
  321.  
  322. si IG1 ;
  323. ltyp1 = v0 elem type ;
  324. mot1 = mot 'ELEM =' ;
  325. repe bm1 (dime ltyp1) ;
  326. mot1 = chai mot1 ' ' (ltyp1 extr &bm1) ;
  327. fin bm1 ;
  328. mot1 = chai ' Maillage test en 3D (' mot1 ')' ;
  329. trac v0 titr mot1 ;
  330. fins;
  331.  
  332. * Modeles :
  333. moq1 = mode V0 thermique source gaussienne elliptique ;
  334. maq1 = mate moq1 'QTOT' QT1 'ORIG' OR1 'RGAU' RG1 'DIRE' DR1 'ZGAU' ZG1 ;
  335. sq1 = sour moq1 maq1 ;
  336.  
  337. qtmoq1 = sq1 resu maxi ;
  338.  
  339. mo1 = mode V0 thermique ;
  340. ma1 = mate mo1 k 25. ;
  341.  
  342. x1 y1 z1 = V0 coor ;
  343. x1 y1 z1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) (z1 - (or1 coor 3)) ;
  344. x1 y1 z1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) (chan cham z1 mo1 stresses) ;
  345. r1 = (x1 * x1) + (y1 * y1) + (z1 * z1) ;
  346. z1 = (x1 * (DR1 coor 1)) + (y1 * (DR1 coor 2)) + (z1 * (DR1 coor 3)) / (norm DR1) ;
  347. r1 = r1 - (z1 * z1) ;
  348. xqx1 = (32. / (PI * PI * PI)) ** 0.5 / RG1 / RG1 / ZG1 ;
  349. sqref1 = exp (-2. * ((r1 / RG1 / RG1) + (z1 * z1 / ZG1 / ZG1))) * QT1 * xqx1 ;
  350. sqref1 = sour mo1 sqref1 ;
  351.  
  352. si IG1 ;
  353. trac sq1 V0 titr ' Flux nodaux du champ de source Gaussienne du MODELE ' ;
  354. trac sqref1 V0 titr ' Flux nodaux du champ de source Gaussienne REFERENCE' ;
  355. trac (sq1 - sqref1) V0 titr ' Difference des 2 champs' ;
  356. fins ;
  357.  
  358. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  359. list err1 ;
  360. err2 = (abs (qtmoq1 - qt1)) / (abs qt1) ;
  361. list err2 ;
  362. si ((err1 > 1.e-12) ou (err2 > 1.e-4)) ;
  363. ERRE ' Probleme MODE SOURCE' ;
  364. fins ;
  365.  
  366. *---------------/!\ SOURCES SUR PLUSIEURS SOUS-ZONES /!\---------------*
  367. *
  368. * Maillage :
  369. v0 = v0 chan line ;
  370. ptsv1 = (v0 coor 1) poin egin lo1 ;
  371. v1 = v0 elem appu stri ptsv1 ;
  372. v1 = v1 coul vert ;
  373. v2 = v0 diff v1;
  374. v2 = v2 coul oran ;
  375.  
  376. si IG1 ;
  377. trac (v1 et v2) titr ' Maillage a 2 sous-zones' ;
  378. fins ;
  379.  
  380. * DEUX SOURCES : GAUSSIENNE & UNIFORME :
  381. * --------------------------------------
  382. *
  383. * Sur 2 sous-zones disctinctes :
  384. * ------------------------------
  385. *
  386. * Modeles / caracteristiques / sources :
  387. moq1 = mode v1 thermique source gaussienne ;
  388. moq2 = mode v2 thermique source ;
  389. maq1 = mate moq1 qtot qt1 orig or1 rgau rg1 ;
  390. maq2 = mate moq2 qvol 1.e9 ;
  391. moq0 = moq1 et moq2 ;
  392. maq0 = maq1 et maq2 ;
  393. sq1 = sour moq1 maq1 ;
  394. sq2 = sour moq2 maq2 ;
  395. sq0 = sour moq0 maq0 ;
  396.  
  397. mo1 = mode v1 thermique ;
  398. ma1 = mate mo1 k 25. ;
  399. x1 y1 z1 = V1 coor ;
  400. x1 y1 z1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) (z1 - (or1 coor 3)) ;
  401. x1 y1 z1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) (chan cham z1 mo1 stresses) ;
  402. r1 = (x1 * x1) + (y1 * y1) + (z1 * z1) ;
  403. xqx1 = (32. / (PI * PI * PI)) ** 0.5 / RG1 / RG1 / RG1 ;
  404. sqref1 = exp (-2. * (r1 / RG1 / RG1)) * QT1 * xqx1 ;
  405. sqref1 = sour mo1 sqref1 ;
  406.  
  407. mo2 = mode v2 thermique ;
  408. ma2 = mate mo2 k 25. ;
  409. sqref2 = sour mo2 v2 1.e9 ;
  410.  
  411. sqref0 = sqref1 et sqref2 ;
  412.  
  413. si IG1 ;
  414. trac sq0 v0 titr ' Flux nodaux des deux modeles de source ' ;
  415. trac sqref0 v0 titr ' Flux nodaux REFERENCE' ;
  416. trac (sq0 - sqref0) v0 titr ' Difference des deux champs' ;
  417. fins ;
  418.  
  419. err0 = (maxi abs (sq0 - sqref0)) / (maxi abs sqref0) ;
  420. list err0 ;
  421. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  422. list err1 ;
  423. err2 = (maxi abs (sq2 - sqref2)) / (maxi abs sqref2) ;
  424. list err2 ;
  425. err1 = maxi (prog err0 err1 err2) ;
  426. si (err1 > 1.e-12) ;
  427. ERRE ' Probleme MODE SOURCE' ;
  428. fins ;
  429.  
  430. * Sur 2 sous-zones "a cheval" :
  431. * -----------------------------
  432. *
  433. * Modeles / caracteristiques / sources :
  434. moq1 = mode v0 thermique source gaussienne ;
  435. moq2 = mode v2 thermique source ;
  436. maq1 = mate moq1 qtot qt1 orig or1 rgau rg1 ;
  437. maq2 = mate moq2 qvol 1.e9 ;
  438. moq0 = moq1 et moq2 ;
  439. maq0 = maq1 et maq2 ;
  440. sq1 = sour moq1 maq1 ;
  441. sq2 = sour moq2 maq2 ;
  442. sq0 = sour moq0 maq0 ;
  443.  
  444. mo1 = mode v0 thermique ;
  445. x1 y1 z1 = V0 coor ;
  446. x1 y1 z1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) (z1 - (or1 coor 3)) ;
  447. x1 y1 z1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) (chan cham z1 mo1 stresses) ;
  448. r1 = (x1 * x1) + (y1 * y1) + (z1 * z1) ;
  449. xqx1 = (32. / (PI * PI * PI)) ** 0.5 / RG1 / RG1 / RG1 ;
  450. sqref1 = exp (-2. * (r1 / RG1 / RG1)) * QT1 * xqx1 ;
  451. sqref1 = sour mo1 sqref1 ;
  452.  
  453. mo2 = mode v2 thermique ;
  454. sqref2 = sour mo2 v2 1.e9 ;
  455.  
  456. sqref0 = sqref1 et sqref2 ;
  457.  
  458. si IG1 ;
  459. trac sq0 V0 titr ' Flux nodaux des deux modeles de source ' ;
  460. trac sqref0 V0 titr ' Flux nodaux REFERENCE' ;
  461. trac (sq0 - sqref0) V0 titr ' Difference des deux champs' ;
  462. fins ;
  463.  
  464. err0 = (maxi abs (sq0 - sqref0)) / (maxi abs sqref0) ;
  465. list err0 ;
  466. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  467. list err1 ;
  468. err2 = (maxi abs (sq2 - sqref2)) / (maxi abs sqref2) ;
  469. list err2 ;
  470. err1 = maxi (prog err0 err1 err2) ;
  471. si (err1 > 1.e-12) ;
  472. ERRE ' Probleme MODE SOURCE' ;
  473. fins ;
  474.  
  475.  
  476. * DEUX SOURCES GAUSSIENNES : SPHERIQUE & ELLIPTIQUE :
  477. * ---------------------------------------------------
  478.  
  479. * Sur 2 sous-zones disctinctes :
  480. * ------------------------------
  481. *
  482. * Parametres de la source en 3D
  483. qt2 = 0.5 * qt1 ;
  484. OR2 = 60.e-3 (0.6*ep1) ep1 ;
  485. DR2 = 0. 0. 2. ;
  486.  
  487. * Modeles / caracteristiques / sources :
  488. moq1 = mode v1 thermique source gaussienne ;
  489. moq2 = mode v2 thermique source gaussienne elliptique ;
  490. maq1 = mate moq1 qtot qt1 orig or1 rgau rg1 ;
  491. maq2 = mate moq2 qtot qt2 orig or2 rgau rg1 dire dr2 zgau zg1 ;
  492. moq0 = moq1 et moq2 ;
  493. maq0 = maq1 et maq2 ;
  494. sq1 = sour moq1 maq1 ;
  495. sq2 = sour moq2 maq2 ;
  496. sq0 = sour moq0 maq0 ;
  497.  
  498. mo1 = mode v1 thermique ;
  499. ma1 = mate mo1 k 25. ;
  500. x1 y1 z1 = V1 coor ;
  501. x1 y1 z1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) (z1 - (or1 coor 3)) ;
  502. x1 y1 z1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) (chan cham z1 mo1 stresses) ;
  503. r1 = (x1 * x1) + (y1 * y1) + (z1 * z1) ;
  504. xqx1 = (32. / (PI * PI * PI)) ** 0.5 / RG1 / RG1 / RG1 ;
  505. sqref1 = exp (-2. * (r1 / RG1 / RG1)) * QT1 * xqx1 ;
  506. sqref1 = sour mo1 sqref1 ;
  507.  
  508. mo2 = mode v2 thermique ;
  509. ma2 = mate mo2 k 25. ;
  510. x2 y2 z2 = V2 coor ;
  511. x2 y2 z2 = (x2 - (or2 coor 1)) (y2 - (or2 coor 2)) (z2 - (or2 coor 3)) ;
  512. x2 y2 z2 = (chan cham x2 mo2 stresses) (chan cham y2 mo2 stresses) (chan cham z2 mo2 stresses) ;
  513. r2 = (x2 * x2) + (y2 * y2) + (z2 * z2) ;
  514. z2 = (x2 * (dr2 coor 1)) + (y2 * (dr2 coor 2)) + (z2 * (dr2 coor 3)) / (norm dr2) ;
  515. r2 = r2 - (z2 * z2) ;
  516. xqx2 = (32. / (PI * PI * PI)) ** 0.5 / RG1 / RG1 / ZG1 ;
  517. sqref2 = exp (-2. * ((r2 / RG1 / RG1) + (z2 * z2 / ZG1 / ZG1))) * QT2 * xqx2 ;
  518. sqref2 = sour mo2 sqref2 ;
  519.  
  520. sqref0 = sqref1 et sqref2 ;
  521.  
  522. si IG1 ;
  523. trac sq0 V0 titr ' Flux nodaux des deux modeles de source ' ;
  524. trac sqref0 V0 titr ' Flux nodaux REFERENCE' ;
  525. trac (sq0 - sqref0) V0 titr ' Difference des deux champs' ;
  526. fins ;
  527.  
  528. err0 = (maxi abs (sq0 - sqref0)) / (maxi abs sqref0) ;
  529. list err0 ;
  530. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  531. list err1 ;
  532. err2 = (maxi abs (sq2 - sqref2)) / (maxi abs sqref2) ;
  533. list err2 ;
  534. err1 = maxi (prog err0 err1 err2) ;
  535. si (err1 > 1.e-12) ;
  536. ERRE ' Probleme MODE SOURCE' ;
  537. fins ;
  538.  
  539. * Sur 2 sous-zones "a cheval" :
  540. * -----------------------------
  541. *
  542. * Modeles / caracteristiques / sources :
  543. moq1 = mode v0 thermique source gaussienne ;
  544. moq2 = mode v0 thermique source gaussienne elliptique ;
  545. maq1 = mate moq1 qtot qt1 orig or1 rgau rg1 ;
  546. maq2 = mate moq2 qtot qt2 orig or2 rgau rg1 dire dr2 zgau zg1 ;
  547. moq0 = moq1 et moq2 ;
  548. maq0 = maq1 et maq2 ;
  549. sq1 = sour moq1 maq1 ;
  550. sq2 = sour moq2 maq2 ;
  551. sq0 = sour moq0 maq0 ;
  552.  
  553. mo1 = mode v0 thermique ;
  554. ma1 = mate mo1 k 25. ;
  555. x1 y1 z1 = v0 coor ;
  556. x1 y1 z1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) (z1 - (or1 coor 3)) ;
  557. x1 y1 z1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) (chan cham z1 mo1 stresses) ;
  558. r1 = (x1 * x1) + (y1 * y1) + (z1 * z1) ;
  559. xqx1 = (32. / (PI * PI * PI)) ** 0.5 / RG1 / RG1 / RG1 ;
  560. sqref1 = exp (-2. * (r1 / RG1 / RG1)) * QT1 * xqx1 ;
  561. sqref1 = sour mo1 sqref1 ;
  562.  
  563. mo2 = mode v0 thermique ;
  564. ma2 = mate mo2 k 25. ;
  565. x2 y2 z2 = v0 coor ;
  566. x2 y2 z2 = (x2 - (or2 coor 1)) (y2 - (or2 coor 2)) (z2 - (or2 coor 3)) ;
  567. x2 y2 z2 = (chan cham x2 mo2 stresses) (chan cham y2 mo2 stresses) (chan cham z2 mo2 stresses) ;
  568. r2 = (x2 * x2) + (y2 * y2) + (z2 * z2) ;
  569. z2 = (x2 * (dr2 coor 1)) + (y2 * (dr2 coor 2)) + (z2 * (dr2 coor 3)) / (norm dr2) ;
  570. r2 = r2 - (z2 * z2) ;
  571. xqx2 = (32. / (PI * PI * PI)) ** 0.5 / RG1 / RG1 / ZG1 ;
  572. sqref2 = exp (-2. * ((r2 / RG1 / RG1) + (z2 * z2 / ZG1 / ZG1))) * QT2 * xqx2 ;
  573. sqref2 = sour mo2 sqref2 ;
  574.  
  575. sqref0 = sqref1 et sqref2 ;
  576.  
  577. si IG1 ;
  578. trac sq0 V0 titr ' Flux nodaux des deux modeles de source ' ;
  579. trac sqref0 V0 titr ' Flux nodaux REFERENCE' ;
  580. trac (sq0 - sqref0) V0 titr ' Difference des deux champs' ;
  581. fins ;
  582.  
  583. err0 = (maxi abs (sq0 - sqref0)) / (maxi abs sqref0) ;
  584. list err0 ;
  585. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  586. list err1 ;
  587. err2 = (maxi abs (sq2 - sqref2)) / (maxi abs sqref2) ;
  588. list err2 ;
  589. err1 = maxi (prog err0 err1 err2) ;
  590. si (err1 > 1.e-12) ;
  591. ERRE ' Probleme MODE SOURCE' ;
  592. fins ;
  593.  
  594. *-------------------------------- FIN ---------------------------------*
  595. fin ;
  596.  
  597.  
  598.  
  599.  
  600.  
  601.  
  602.  

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