Télécharger source2.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : source2.dgibi
  2. * section : thermique
  3. *----------------------------------------------------------------------*
  4. * SOURCE2.DGIBI *
  5. *----------------------------------------------------------------------*
  6. *
  7. * Objet :
  8. * -------
  9. *
  10. * Verfication / validation d'un modele de source de chaleur.
  11. * Cas d'une SOURCE GAUSSIENNE ISOTROPE.
  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. RG1 = 5.e-3 ;
  44.  
  45. * Maillage :
  46. l1 = (0 0) droi 10 (0 ep1) ;
  47. s1 = l1 tran 50 (lo1 0) ;
  48. l2 = s1 cote 3 ;
  49. s2 = l2 tran 50 (lo1 ep1) ;
  50. s2 = s2 chan tri3 ;
  51. s0 = s1 et s2 ;
  52.  
  53. si IG1 ;
  54. ltyp1 = s0 elem type ;
  55. mot1 = mot 'ELEM =' ;
  56. repe bm1 (dime ltyp1) ;
  57. mot1 = chai mot1 ' ' (ltyp1 extr &bm1) ;
  58. fin bm1 ;
  59. mot1 = chai ' Maillage test en 2D (' mot1 ')' ;
  60. trac s0 titr mot1 ;
  61. fins;
  62.  
  63. * MODE PLAN :
  64. * -----------
  65. *
  66. opti mode plan ;
  67. *
  68. moq1 = mode s0 thermique source gaussienne ;
  69. maq1 = mate moq1 'QTOT' QT1 'ORIG' OR1 'RGAU' RG1 ;
  70. sq1 = sour moq1 maq1 ;
  71.  
  72. qtmoq1 = sq1 resu maxi ;
  73.  
  74. mo1 = mode s0 thermique ;
  75. ma1 = mate mo1 k 25. ;
  76.  
  77. x1 y1 = s0 coor ;
  78. x1 y1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) ;
  79. x1 y1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) ;
  80. xqx1 = 4. / PI / RG1 / RG1 ;
  81. sqref1 = exp (-2. * ((x1 * x1) + (y1 * y1)) / RG1 / RG1) * QT1 * xqx1 ;
  82. sqref1 = sour mo1 sqref1 ;
  83.  
  84. si IG1 ;
  85. trac sq1 S0 titr ' Flux nodaux du champ de source Gaussienne du MODELE ';
  86. trac sqref1 s0 titr ' Flux nodaux du champ de source Gaussienne REFERENCE' ;
  87. trac (sq1 - sqref1) s0 titr ' Difference des 2 champs' ;
  88. fins ;
  89.  
  90. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  91. list err1 ;
  92. err2 = (abs (qtmoq1 - qt1)) / (abs qt1) ;
  93. list err2 ;
  94. si ((err1 > 1.e-12) ou (err2 > 1.e-4)) ;
  95. ERRE ' Probleme MODE SOURCE' ;
  96. fins ;
  97.  
  98. * MODE AXIS :
  99. * -----------
  100. *
  101. opti mode axis ;
  102. *
  103. moq1 = mode s0 thermique source gaussienne isotrope ;
  104. maq1 = mate moq1 'QTOT' QT1 'ORIG' OR1 'RGAU' RG1 ;
  105. sq1 = sour moq1 maq1 ;
  106.  
  107. qtmoq1 = sq1 resu maxi ;
  108.  
  109. mo1 = mode s0 thermique ;
  110. ma1 = mate mo1 k 25. ;
  111.  
  112. x0 y0 = (or1 coor 1) (or1 coor 2) ;
  113. x1 y1 = s0 coor ;
  114. x1 y1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) ;
  115. x1 y1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) ;
  116. rc2 = 2. ** 0.5 ;
  117. rcpi1 = pi ** 0.5 ;
  118. xik1 = 0.25 * RG1 * RG1 * (exp (-2. * x0 * x0 / RG1 / RG1)) ;
  119. xik2 = 0.5 * RG1 * x0 * rcpi1 / rc2 ;
  120. xik3 = RG1 * x0 * (erf (x0 * rc2 / RG1)) * rcpi1 / (rc2 * rc2 * rc2) ;
  121. xqx1 = PI * RG1 * rcpi1 / rc2 ;
  122. xqx1 = xqx1 * (xik1 + xik2 + xik3) ;
  123. sqref1 = exp (-2. * ((x1 * x1) + (y1 * y1)) / RG1 / RG1) * QT1 / xqx1 ;
  124. sqref1 = sour mo1 sqref1 ;
  125.  
  126. si IG1 ;
  127. trac sq1 S0 titr ' Flux nodaux du champ de source Gaussienne du MODELE ';
  128. trac sqref1 s0 titr ' Flux nodaux du champ de source Gaussienne REFERENCE' ;
  129. trac (sq1 - sqref1) s0 titr ' Difference des 2 champs' ;
  130. fins ;
  131.  
  132. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  133. list err1 ;
  134. err2 = (abs (qtmoq1 - qt1)) / (abs qt1) ;
  135. list err2 ;
  136. si ((err1 > 1.e-12) ou (err2 > 1.e-4)) ;
  137. ERRE ' Probleme MODE SOURCE' ;
  138. fins ;
  139.  
  140. *---------------------- 2D ELEMENTS QUADRATIQUES ----------------------*
  141.  
  142. * Passage en EF quadratiques :
  143. s0 = s0 chan quad ;
  144.  
  145. si IG1 ;
  146. ltyp1 = s0 elem type ;
  147. mot1 = mot 'ELEM =' ;
  148. repe bm1 (dime ltyp1) ;
  149. mot1 = chai mot1 ' ' (ltyp1 extr &bm1) ;
  150. fin bm1 ;
  151. mot1 = chai ' Maillage test en 2D (' mot1 ')' ;
  152. trac s0 titr mot1 ;
  153. fins;
  154.  
  155. * MODE PLAN :
  156. * -----------
  157. *
  158. opti mode plan ;
  159. *
  160. moq1 = mode s0 thermique source gaussienne isotrope ;
  161. maq1 = mate moq1 'QTOT' QT1 'ORIG' OR1 'RGAU' RG1 ;
  162. sq1 = sour moq1 maq1 ;
  163.  
  164. qtmoq1 = sq1 resu maxi ;
  165.  
  166. mo1 = mode s0 thermique ;
  167. ma1 = mate mo1 k 25. ;
  168.  
  169. x1 y1 = s0 coor ;
  170. x1 y1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) ;
  171. x1 y1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) ;
  172. xqx1 = 4. / PI / RG1 / RG1 ;
  173. sqref1 = exp (-2. * ((x1 * x1) + (y1 * y1)) / RG1 / RG1) * QT1 * xqx1 ;
  174. sqref1 = sour mo1 sqref1 ;
  175.  
  176. si IG1 ;
  177. trac sq1 S0 titr ' Flux nodaux du champ de source Gaussienne du MODELE ';
  178. trac sqref1 s0 titr ' Flux nodaux du champ de source Gaussienne REFERENCE' ;
  179. trac (sq1 - sqref1) s0 titr ' Difference des 2 champs' ;
  180. fins ;
  181.  
  182. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  183. list err1 ;
  184. err2 = (abs (qtmoq1 - qt1)) / (abs qt1) ;
  185. list err2 ;
  186. si ((err1 > 1.e-12) ou (err2 > 1.e-4)) ;
  187. ERRE ' Probleme MODE SOURCE' ;
  188. fins ;
  189.  
  190. * MODE AXIS :
  191. * -----------
  192. *
  193. opti mode axis ;
  194. *
  195. moq1 = mode s0 thermique source gaussienne ;
  196. maq1 = mate moq1 'QTOT' QT1 'ORIG' OR1 'RGAU' RG1 ;
  197. sq1 = sour moq1 maq1 ;
  198.  
  199. qtmoq1 = sq1 resu maxi ;
  200.  
  201. mo1 = mode s0 thermique ;
  202. ma1 = mate mo1 k 25. ;
  203.  
  204. x0 y0 = (or1 coor 1) (or1 coor 2) ;
  205. x1 y1 = s0 coor ;
  206. x1 y1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) ;
  207. x1 y1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) ;
  208. rc2 = 2. ** 0.5 ;
  209. rcpi1 = pi ** 0.5 ;
  210. xik1 = 0.25 * RG1 * RG1 * (exp (-2. * x0 * x0 / RG1 / RG1)) ;
  211. xik2 = 0.5 * RG1 * x0 * rcpi1 / rc2 ;
  212. xik3 = RG1 * x0 * (erf (x0 * rc2 / RG1)) * rcpi1 / (rc2 * rc2 * rc2) ;
  213. xqx1 = PI * RG1 * rcpi1 / rc2 ;
  214. xqx1 = xqx1 * (xik1 + xik2 + xik3) ;
  215. sqref1 = exp (-2. * ((x1 * x1) + (y1 * y1)) / RG1 / RG1) * QT1 / xqx1 ;
  216. sqref1 = sour mo1 sqref1 ;
  217.  
  218. si IG1 ;
  219. trac sq1 S0 titr ' Flux nodaux du champ de source Gaussienne du MODELE ';
  220. trac sqref1 s0 titr ' Flux nodaux du champ de source Gaussienne REFERENCE' ;
  221. trac (sq1 - sqref1) s0 titr ' Difference des 2 champs' ;
  222. fins ;
  223.  
  224. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  225. list err1 ;
  226. err2 = (abs (qtmoq1 - qt1)) / (abs qt1) ;
  227. list err2 ;
  228. si ((err1 > 1.e-12) ou (err2 > 1.e-4)) ;
  229. ERRE ' Probleme MODE SOURCE' ;
  230. fins ;
  231.  
  232. *------------------------ 3D ELEMENTS LINEAIRES -----------------------*
  233. *
  234. opti dime 3 elem cub8 ;
  235. s0 = s0 chan line ;
  236.  
  237. * Maillage du volume :
  238. v0 = s0 volu tran 5 (0 0 ep1) ;
  239. x1 = v0 coor 1 ;
  240.  
  241. * Ajour de tetraedres :
  242. v3 = v0 elem appu stri (x1 poin supe (1.5 * lo1)) ;
  243. v0 = v0 diff v3 ;
  244. v0 = v0 et (v3 chan tet4) ;
  245.  
  246. * Ajout d'un element pyramide :
  247. el1 = (v0 elem cub8) elem 1 ;
  248. pe1 = el1 chan poi1 ;
  249. ps1 = pe1 poin 8 ;
  250. py1 = manu pyr5 (pe1 poin 1) (pe1 poin 2) (pe1 poin 3) (pe1 poin 4) ps1 ;
  251. py2 = manu pyr5 (pe1 poin 1) (pe1 poin 2) (pe1 poin 6) (pe1 poin 5) ps1 ;
  252. py3 = manu pyr5 (pe1 poin 2) (pe1 poin 3) (pe1 poin 7) (pe1 poin 6) ps1 ;
  253. py0 = py1 et py2 et py3 ;
  254. v0 = (v0 diff el1) et py0 ;
  255.  
  256. si IG1 ;
  257. ltyp1 = v0 elem type ;
  258. mot1 = mot 'ELEM =' ;
  259. repe bm1 (dime ltyp1) ;
  260. mot1 = chai mot1 ' ' (ltyp1 extr &bm1) ;
  261. fin bm1 ;
  262. mot1 = chai ' Maillage test en 3D (' mot1 ')' ;
  263. trac v0 titr mot1 ;
  264. fins;
  265.  
  266. * Parametres de la source en 3D
  267. OR1 = 30.e-3 (0.5 * ep1) ep1 ;
  268.  
  269. * Modeles :
  270. moq1 = mode V0 thermique source gaussienne ;
  271. maq1 = mate moq1 'QTOT' QT1 'ORIG' OR1 'RGAU' RG1 ;
  272. sq1 = sour moq1 maq1 ;
  273.  
  274. qtmoq1 = sq1 resu maxi ;
  275.  
  276. mo1 = mode V0 thermique ;
  277. ma1 = mate mo1 k 25. ;
  278.  
  279. x1 y1 z1 = V0 coor ;
  280. x1 y1 z1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) (z1 - (or1 coor 3)) ;
  281. x1 y1 z1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) (chan cham z1 mo1 stresses) ;
  282. xqx1 = (32. / (PI * PI * PI)) ** 0.5 / RG1 / RG1 / RG1 ;
  283. sqref1 = exp (-2. *((x1 * x1) + (y1 * y1) + (z1 * z1)) / RG1 / RG1) * QT1 * xqx1 ;
  284. sqref1 = sour mo1 sqref1 ;
  285.  
  286. si IG1 ;
  287. trac sq1 V0 titr ' Flux nodaux du champ de source Gaussienne du MODELE ' coup or1 (0 1 0) (0 0 1) ;
  288. trac sqref1 V0 titr ' Flux nodaux du champ de source Gaussienne REFERENCE' coup or1 (0 1 0) (0 0 1) ;
  289. trac (sq1 - sqref1) V0 titr ' Difference des 2 champs' coup or1 (0 1 0) (0 0 1) ;
  290. fins ;
  291.  
  292. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  293. list err1 ;
  294. err2 = (abs (qtmoq1 - qt1)) / (abs qt1) ;
  295. list err2 ;
  296. si ((err1 > 1.e-12) ou (err2 > 2.e-4)) ;
  297. ERRE ' Probleme MODE SOURCE' ;
  298. fins ;
  299.  
  300. *---------------------- 3D ELEMENTS QUADRATIQUES ----------------------*
  301. *
  302. v0 = v0 chan quad ;
  303.  
  304. si IG1 ;
  305. ltyp1 = v0 elem type ;
  306. mot1 = mot 'ELEM =' ;
  307. repe bm1 (dime ltyp1) ;
  308. mot1 = chai mot1 ' ' (ltyp1 extr &bm1) ;
  309. fin bm1 ;
  310. mot1 = chai ' Maillage test en 3D (' mot1 ')' ;
  311. trac v0 titr mot1 ;
  312. fins;
  313.  
  314. * Modeles :
  315. moq1 = mode V0 thermique source gaussienne isotrope ;
  316. maq1 = mate moq1 'QTOT' QT1 'ORIG' OR1 'RGAU' RG1 ;
  317. sq1 = sour moq1 maq1 ;
  318.  
  319. qtmoq1 = sq1 resu maxi ;
  320.  
  321. mo1 = mode V0 thermique ;
  322. ma1 = mate mo1 k 25. ;
  323.  
  324. x1 y1 z1 = V0 coor ;
  325. x1 y1 z1 = (x1 - (or1 coor 1)) (y1 - (or1 coor 2)) (z1 - (or1 coor 3)) ;
  326. x1 y1 z1 = (chan cham x1 mo1 stresses) (chan cham y1 mo1 stresses) (chan cham z1 mo1 stresses) ;
  327. xqx1 = (32. / (PI * PI * PI)) ** 0.5 / RG1 / RG1 / RG1 ;
  328. sqref1 = exp (-2. *((x1 * x1) + (y1 * y1) + (z1 * z1)) / RG1 / RG1) * QT1 * xqx1 ;
  329. sqref1 = sour mo1 sqref1 ;
  330.  
  331. si IG1 ;
  332. trac sq1 V0 titr ' Flux nodaux du champ de source Gaussienne du MODELE ' ;
  333. trac sqref1 V0 titr ' Flux nodaux du champ de source Gaussienne REFERENCE' ;
  334. trac (sq1 - sqref1) V0 titr ' Difference des 2 champs' ;
  335. fins ;
  336.  
  337. err1 = (maxi abs (sq1 - sqref1)) / (maxi abs sqref1) ;
  338. list err1 ;
  339. err2 = (abs (qtmoq1 - qt1)) / (abs qt1) ;
  340. list err2 ;
  341. si ((err1 > 1.e-12) ou (err2 > 1.e-4)) ;
  342. ERRE ' Probleme MODE SOURCE' ;
  343. fins ;
  344.  
  345. *-------------------------------- FIN ---------------------------------*
  346. fin ;
  347.  
  348.  
  349.  
  350.  
  351.  

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