Télécharger ktest-calp.dgibi

Retour à la liste

Numérotation des lignes :

  1. * fichier : ktest-calp.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. ***********************************************************
  5. * *
  6. * ktest-calp.dgibi *
  7. * *
  8. * Test de fonctionnement de l'opérateur CALP sur un *
  9. * simple exemple de plaque carrée en flexion pure. *
  10. * Plusieurs types de coques sont testées en même temps. *
  11. * *
  12. * Auteur : Michel Bulik *
  13. * *
  14. * Date : Août 1995 *
  15. * *
  16. ***********************************************************
  17.  
  18. *** Options ...
  19.  
  20. opti dime 3 mode trid elem seg2 echo 1 ;
  21.  
  22. *** Paramètres ...
  23.  
  24. * longueur et largeur de la plaque
  25. L = 1. ;
  26.  
  27. * nombre de divisions le long de ses côtés
  28. nbdiv = 10 ;
  29.  
  30. * module de YOUNG
  31. valyoun = 2.e+11 ;
  32.  
  33. * coeff. de Poisson nul pour qu'il n'y ait pas
  34. * de courbure anticlastique
  35. valnu = 0. ;
  36.  
  37. * épaisseur de la plaque
  38. ep = 1.e-2 ;
  39.  
  40. * calcul linéaire ou non
  41. calclin = vrai ;
  42.  
  43. * affichage graphique ou non
  44. graph = faux ;
  45.  
  46. * types d'éléments testés
  47. typdel = mots COQ4 DKT DST COQ3 ;
  48.  
  49. *** Points ...
  50.  
  51. dens (L / nbdiv) ;
  52.  
  53. p1 = 0 0 0 ;
  54. p2 = L 0 0 ;
  55.  
  56. vectz = 0 0 L ;
  57. vectxz = L 0 L ;
  58.  
  59. *** Lignes ...
  60.  
  61. li1 = p1 d nbdiv p2 ;
  62.  
  63. *** Surface ...
  64.  
  65. opti elem qua4 ;
  66.  
  67. su1 = li1 tran nbdiv vectz ;
  68.  
  69. *** Début de la boucle sur les types d'éléments ...
  70.  
  71. nbfois = dime typdel ;
  72. repeter surelem nbfois ;
  73.  
  74. elemact = extr typdel &surelem ;
  75.  
  76. si(ega &surelem 2) ;
  77. mess 'SU1 passe au type TRI3' ;
  78. su1 = chan tri3 su1 ;
  79. finsi ;
  80.  
  81. mess '--------------------------------------' ;
  82. mess 'On analyse actuellement l élément ' elemact ;
  83.  
  84. *** Modèle ...
  85.  
  86. mo1 = mode su1 mecanique elastique elemact ;
  87.  
  88. *** Matériau ...
  89.  
  90. ma1 = mate mo1 YOUN valyoun NU valnu ;
  91.  
  92. ca1 = cara mo1 EPAI ep ;
  93.  
  94. ma1 = ma1 et ca1 ;
  95.  
  96. *** Rigidité + Blocages ...
  97.  
  98. si(calclin) ;
  99. ri1 = rigi mo1 ma1 ;
  100. finsi ;
  101.  
  102. cl1 = (bloq uy uz li1) et
  103. (bloq rx li1) et
  104. (bloq ux p1 ) ;
  105.  
  106. *** Chargement ...
  107.  
  108. * un champ de forces obtenus par l'intégration avec les
  109. * fonctions de forme - on va le modifier pour obtenir
  110. * le chargement réparti uniformément sur le bord supérieur
  111. * de la plaque
  112.  
  113. MOP = 'MODE' SU1 'CHARGEMENT' 'PRESSION' elemact ;
  114. MAP = 'PRES' MOP 'PRES' 10000. ;
  115. fopres = 'BSIG' MOP MAP ;
  116.  
  117. mafo1 = su1 poin droi vectz vectxz ;
  118. fo1 = exco fy (redu fopres mafo1) mx ;
  119. * moment appliqué par unité de longueur
  120. momtot = (abs (maxi (resu fo1))) / L ;
  121. mess 'Moment appliqué par (par unité de longueur) = '
  122. momtot ;
  123.  
  124. *** L'affichage éventuel du chargement ...
  125.  
  126. vecchar = faux ;
  127. si(exis (extr fo1 comp) FX) ;
  128. titr 'Chargement : forces nodales' ;
  129. vecfo1 = vect fo1 FX FY FZ (1 / (maxi (abs fo1))) ;
  130. vecchar = vrai ;
  131. finsi ;
  132. si (exis (extr fo1 comp) MX) ;
  133. titr 'Chargement : couples nodaux' ;
  134. vecfo1 = vect fo1 MX MY MZ (1 / (maxi (abs fo1))) ;
  135. vecchar = vrai ;
  136. finsi ;
  137. si(graph et vecchar) ;
  138. trac vecfo1 su1 ;
  139. finsi ;
  140.  
  141. *** Résolution ...
  142.  
  143. si(calclin) ;
  144. de1 = reso (ri1 et cl1) fo1 ;
  145. sinon ;
  146. tabnl = table ;
  147. tabnl . MODELE = mo1 ;
  148. tabnl . CARACTERISTIQUES = ma1 ;
  149. tabnl . BLOCAGES_MECANIQUES = cl1 ;
  150. lr1 = prog 0. 1. ;
  151. ev1 = evol manu 't' lr1 'f(t)' lr1 ;
  152. tabnl . TEMPS_CALCULES = prog 0.0 pas 1.e-1 1.0 ;
  153. tabnl . TEMPS_SAUVES = lr1 ;
  154. tabnl . GRANDS_DEPLACEMENTS = vrai ;
  155. tabn1 . 'HYPOTHESE_DEFORMATIONS'='LINEAIRE' ;
  156.  
  157. tabnl . MTOL = 1.e+30 ;
  158.  
  159. pasapas tabnl ;
  160.  
  161. de1 = tabnl . DEPLACEMENTS . 1 ;
  162. finsi ;
  163.  
  164. *** Affichage éventuel de la déformée et des DDL ...
  165.  
  166. defo0 = defo su1 de1 0. ;
  167. titr 'La deformee' ;
  168. defo1 = defo su1 de1 1. ;
  169. si(graph) ;
  170. trac (defo0 et defo1) ;
  171.  
  172. * lescomp = extr (redu de1 su1) comp ;
  173. lescomp = mots UY RX ;
  174. nbcomp = dime lescomp ;
  175. repe surcomp nbcomp ;
  176. lacomp = extr lescomp &surcomp ;
  177. titr 'Element' elemact ' : Degre de liberte'
  178. lacomp ;
  179. trac (exco lacomp de1) su1 ;
  180. fin surcomp ;
  181.  
  182. finsi ;
  183.  
  184. *** Contraintes ...
  185.  
  186. si(calclin) ;
  187. si1 = sigm mo1 ma1 de1 'LINE' ;
  188. sinon ;
  189. si1 = sigm mo1 ma1 de1 'QUAD' ;
  190. finsi ;
  191. si1 = rten si1 mo1 (1 0 0) (0 0 1) ;
  192.  
  193. *** L'affichage éventuel des contraintes ...
  194.  
  195. si(graph) ;
  196.  
  197. lescomp = extr si1 comp ;
  198. nbcomp = dime lescomp ;
  199. repe surcom2 nbcomp ;
  200. lacomp = extr lescomp &surcom2 ;
  201. titr 'Element' elemact ' : Contrainte generalisee'
  202. lacomp ;
  203. trac (exco lacomp si1) mo1 ;
  204. fin surcom2 ;
  205.  
  206. finsi ;
  207.  
  208. *** Vérif de la valeur de la contrainte ...
  209.  
  210. valmax = maxi (exco m22 si1) ;
  211. valmin = mini (exco m22 si1) ;
  212. si( (valmax - valmin) > 1.e-6 ) ;
  213. mess 'Attention ! Résultat non homogène (M22)' ;
  214. * valm22 = (valmin + valmax) / 2 ;
  215. valm22 = (intg (exco m22 si1) mo1) / (L * L) ;
  216. sinon ;
  217. valm22 = valmax ;
  218. finsi ;
  219. ecartrel = 100 * ((valm22 - momtot) / momtot) ;
  220. mess 'Ecart relatif (M22) = ' ecartrel '%' ;
  221. si(ecartrel > 1.e-5) ; erreur 5 ; finsi ;
  222.  
  223. *** Application de CALP ...
  224.  
  225. si1sup = calp si1 ca1 mo1 SUPE ;
  226. si1inf = calp si1 ca1 mo1 INFE ;
  227.  
  228. valtheor = (6 * momtot) / (ep * ep) ;
  229.  
  230. valsup = (intg (exco smyy si1sup) mo1) / (L * L) ;
  231. valinf = (intg (exco smyy si1inf) mo1) / (L * L) ;
  232.  
  233. ecsup = 100 * ((valsup - valtheor) / valtheor) ;
  234. mess 'Ecart relatif (peau supérieure) = ' ecsup '%' ;
  235. si(ecsup > 1.e-5) ; erreur 5 ; finsi ;
  236.  
  237. ecinf = 100 * ((valinf + valtheor) / valtheor) ;
  238. mess 'Ecart relatif (peau inférieure) = ' ecinf '%' ;
  239. si(ecinf > 1.e-5) ; erreur 5 ; finsi ;
  240.  
  241.  
  242. *** Résultats analytiques pour le déplacement
  243. *** et la rotation ...
  244.  
  245. * inertie par unité de longueur ...
  246. I = (ep * ep * ep) / 12 ;
  247.  
  248. depanal = (momtot * L * L) / (2 * valyoun * I) ;
  249. depcal = (maxi (exco uy de1)) ;
  250. ecdep = 100 * ((depcal - depanal) / depanal) ;
  251. mess 'Écart relatif (déplacement maxi) = ' ecdep '%' ;
  252. si(ecdep > 0.5) ; erreur 5 ; finsi ;
  253.  
  254. rotanal = (momtot * L) / (valyoun * I) ;
  255. rotcal = abs (mini (exco rx de1)) ;
  256. ecrot = 100 * ((rotcal - rotanal) / rotanal) ;
  257. mess 'Écart relatif (rotation maxi) = ' ecrot '%' ;
  258. si(ecrot > 2.0) ; erreur 5 ; finsi ;
  259.  
  260. *** Déformations généralisées ...
  261.  
  262. si (calclin) ;
  263. ep1 = epsi mo1 de1 ma1 'LINE' ;
  264. sinon ;
  265. ep1 = epsi mo1 de1 ma1 'QUAD' ;
  266. finsi ;
  267. ep1 = rten ep1 mo1 (1 0 0) (0 0 1) ;
  268.  
  269. *** et sur les peaux ...
  270.  
  271. ep1sup = calp ep1 ca1 mo1 SUPE ;
  272. ep1inf = calp ep1 ca1 mo1 INFE ;
  273.  
  274. *** Comparaison avec les valeurs analytiques ...
  275.  
  276. * RTTT
  277. valmax = maxi (exco rttt ep1) ;
  278. valmin = mini (exco rttt ep1) ;
  279. si( (valmax - valmin) > 1.e-6 ) ;
  280. mess 'Attention ! Résultat non homogène (RTTT)' ;
  281. valrttt = (intg (exco rttt ep1) mo1) / (L * L) ;
  282. sinon ;
  283. valrttt = valmax ;
  284. finsi ;
  285. valanal = momtot / (valyoun * I) ;
  286. ecartrel = 100 * ((valrttt - valanal) / valanal) ;
  287. mess 'Ecart relatif (RTTT) = ' ecartrel '%' ;
  288. si(ecartrel > 1.e-5) ; erreur 5 ; finsi ;
  289.  
  290. * EPYY (peau supérieure)
  291. valanal = (ep * valrttt) / 2 ;
  292.  
  293. valmax = maxi (exco epyy ep1sup) ;
  294. valmin = mini (exco epyy ep1sup) ;
  295. si( (valmax - valmin) > 1.e-6 ) ;
  296. mess 'Attention ! Résultat non homogène (EPYY sup)' ;
  297. valepyys = (intg (exco epyy ep1sup) mo1) / (L * L) ;
  298. sinon ;
  299. valepyys = valmax ;
  300. finsi ;
  301. ecartyys = 100 * ((valepyys - valanal) / valanal) ;
  302. mess 'Ecart relatif (EPYY sup) = ' ecartyys '%' ;
  303. si(ecartyys > 1.e-5) ; erreur 5 ; finsi ;
  304.  
  305. valmax = maxi (exco epyy ep1inf) ;
  306. valmin = mini (exco epyy ep1inf) ;
  307. si( (valmax - valmin) > 1.e-6 ) ;
  308. mess 'Attention ! Résultat non homogène (EPYY inf)' ;
  309. valepyyi = (intg (exco epyy ep1inf) mo1) / (L * L) ;
  310. sinon ;
  311. valepyyi = valmax ;
  312. finsi ;
  313. ecartyyi = 100 * ((valepyyi + valanal) / valanal) ;
  314. mess 'Ecart relatif (EPYY inf) = ' ecartyyi '%' ;
  315. si(ecartyyi > 1.e-5) ; erreur 5 ; finsi ;
  316.  
  317. *** Fin de la boucle sur les types d'éléments ...
  318.  
  319. fin surelem ;
  320.  
  321. *** Bye ...
  322.  
  323. fin ;
  324.  
  325.  
  326.  
  327.  
  328.  
  329.  
  330.  
  331.  
  332.  
  333.  
  334.  
  335.  
  336.  

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