Télécharger la plaquette Cast3M


fr en

Télécharger rupt29.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier :  rupt29.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. opti trac psc eptr 5 ;
  5.  
  6. GRAPH = 'O' ;
  7.  
  8. *******************************************************
  9. *           Test rupt29.dgibi: Jeux de données        *
  10. *           ---------------------------------         *
  11. *******************************************************
  12.  
  13. *Cas test de validation pour le calcul de J sous plusieurs chargement 
  14. *avec les procedures g_theta.procedur et g_calcul.procedur
  15. *
  16. *- chargement avec pression sur levres avec plusieurs types d elements
  17. *  sur le maillage de la fissure
  18.  
  19. *Calcul en dimension 3 avec des elements CU20 et PRI15 sur un maillage 
  20. *complet non symetrique
  21.  
  22.  
  23. * DONNEES GEOMETRIQUES
  24. *-------------------- 
  25.  opti dime 3 elem cu20 echo 1;
  26.  zz = 1e-6 ;
  27.  
  28.  aa = 2. ; Rsa = 5. ;
  29. *
  30.  R = Rsa * aa ;
  31.  dda = aa / 40. ;
  32. *
  33. * ----- MAILLAGE
  34. *
  35.  pfond = aa 0. 0. ;
  36.  np = 4 ; nt = 6 ;
  37.  p1 = pfond plus ((np*dda) 0. 0.) ;
  38.  p2 = p1 plus (0. 0. (np*dda)) ;
  39.  p4 = pfond moin ((np*dda) 0. 0.) ;
  40.  p3 = p4 plus (0. 0. (np*dda)) ;
  41.  ss1 = (p1 d np pfond d np p4) tran np (0. 0. (np*dda)) ;
  42.  p5 = pfond plus ((aa / 2.) 0. 0.) ;
  43.  p6 = p5 plus (0. 0. (aa / 2.)) ;
  44.  p8 = pfond moin ((aa / 2.) 0. 0.) ;
  45.  p7 = p8 plus (0. 0. (aa / 2.)) ;
  46.  ll1 = p1 d np p2 d (2*np) p3 d np p4 ;
  47.  ll2 = p5 d np p6 d (2*np) p7 d np p8 ;
  48.  ss2 = ll1 regl ll2 dini dda dfin (aa/8.) ;
  49.  
  50.  dx = R - (1.5 * aa) ;
  51.  ss3 = (p5 d np p6) tran (dx 0. 0.) dini (aa/4.) dfin aa ;
  52.  ss4 = ss1 et ss2 et ss3 ; elim zz ss4 ;
  53.  vv1 = ss4 volu nt rota 90. (0. 0. 0.) (0. 0. 1.) ;
  54.  
  55.  pc = 0. 0. 0. ;
  56.  p9 = 0. (aa / 2.) 0. ;
  57.  nq = 4 ;
  58.  ll3 = (p8 c nt pc p9) et (p9 d nq pc d nq p8) ;
  59.  ss5 = surf ll3 plane ;
  60.  vv2 = ss5 volu tran np (0. 0. (aa / 2.)) ;
  61.  elim zz (vv1 et vv2) ;
  62.  
  63.  ss6 = (vv1 et vv2) poin plan p6 p7 (0. 1. (aa / 2.)) zz ;
  64.  ss7 = (enve (vv1 et vv2)) elem appuye strictement ss6 ;
  65.  dz = R - (0.5 * aa) ;
  66.  vv3 = ss7 volu tran (0. 0. dz) dini (aa/4.) dfin aa ;
  67.  
  68.  vtot1 = vv1 et vv2 et vv3 ;
  69.  
  70. *---------- Surface de la fissure
  71.  ss0 = vtot1 poin plan pc p8 p9 zz ;
  72.  ss1 = (enve vtot1) elem appuye strictement ss0 ;
  73.  chx = ss1 coor 1 ; chy = ss1 coor 2 ;
  74.  chr = ((chx ** 2) + (chy ** 2)) ** 0.5 ;
  75.  chrr = chan cham chr ss1 ;
  76.  sulev1 = chrr elem infe (aa - zz) ;
  77.  sulig = diff ss1 sulev1 ;
  78.  ffis = (cont sulev1) elem appuye strictement sulig ;
  79.  
  80. *---------- Fin du modèle
  81.  
  82.  vtot2 = vtot1 syme plan pc p8 p9 ;
  83.  sulev2 = sulev1 syme plan pc p8 p9 ;
  84.  sulev2 = INVE sulev2 ;
  85.  elim zz (vtot2 et sulev2 et sulig) ;
  86.  
  87.  vtot = vtot1 et vtot2 ;
  88.  
  89. *---------- Lignes et surfaces nommées
  90.  sulev = sulev1 ET sulev2 ;
  91.  
  92.  ss0 = vtot poin plan pc p8 p7 zz ;
  93.  symy = (enve vtot) elem appuye strictement ss0 ;
  94.  
  95.  ss0 = vtot poin plan pc p9 (0. 0. 1.) zz ;
  96.  symx = (enve vtot) elem appuye strictement ss0 ;
  97.  
  98.  ss0 = vtot point plan (0. 0. R) (1. 0. R) (0. 1. R) zz ;
  99.  ssup = (enve vtot) elem appuye strictement ss0 ;
  100.  Rp = -1. * R ;
  101.  ss0 = vtot point plan (0. 0. Rp) (1. 0. Rp) (0. 1. Rp) zz ;
  102.  sinf = (enve vtot) elem appuye strictement ss0 ;
  103.  
  104. si (ega GRAPH 'O') ;
  105.   trac vtot qual cach ;
  106. fins ;
  107.  
  108. ******************************************************************
  109.  
  110. *-- Modele mecanique et parametres materiau
  111. EI20 = 197000. ; ALI20 = 15.54E-6 ;
  112. mo2 = MODE vtot MECANIQUE ELASTIQUE ISOTROPE ;
  113. ma2 = MATE mo2 'YOUN' EI20 'NU' 0.3 'ALPH' AlI20 ;
  114. ri1 = RIGI mo2 ma2 ;
  115.  
  116. *--- Conditions aux limites
  117. * Blocages
  118. cd1 = BLOQ symx  'UX' ;
  119. cd2 = BLOQ symy  'UY' ;
  120. cd3 = BLOQ sinf  'UZ' ;
  121. cdl = cd1 ET cd2 ET cd3 ;
  122. *Traction uniaxiale (via un modele de pression)
  123. pp = 100. ;
  124. moph = MODE ssup 'CHARGEMENT' 'PRESSION' 'CONS' 'HAUT' ;
  125. *Pression sur les levres (via un modele de pression)
  126. mopl = MODE sulev 'CHARGEMENT' 'PRESSION' 'CONS' 'LEVRES' ;
  127.  
  128. * Construction des second membres
  129. maph  = PRES moph 'PRES' pp ;
  130. effN  = BSIG moph maph ;
  131. mapl1 = PRES mopl 'PRES' pp ;
  132. ff0   = BSIG mopl mapl1 ;
  133. zer1  = ZERO mopl 'RIGIDITE' ;
  134. chx chy chz = COOR zer1 ;
  135. chr = (((chx ** 2) + (chy ** 2)) ** 0.5) / aa ;
  136. *ch2 = pp * ((1. - (chr ** 2)) ** 0.5) ;
  137. ch2 = pp * (chr ** 2) ;
  138. mapl2 = PRES mopl ch2 ;
  139. ff2  = BSIG mopl mapl2 ;
  140.  
  141.  
  142. *--- 11 - RESOLUTION
  143. depN dep0 dep2 = RESO (ri1 ET cdl) effN ff0 ff2 ;
  144.  
  145.  
  146. *--- 12 - Calcul de J pour les trois conditions aux limites
  147. *----------------------------------------------------------
  148. EE = EI20 / 910. ; rpa = (3.1416 * aa / 1000.) ** 0.5 ;
  149. * cas 1 : traction seule
  150. SUPTAB = TABLE;
  151. SUPTAB.'MODELE' = mo2 ET moph ;
  152. SUPTAB.'CARACTERISTIQUES' = ma2 ;
  153. SUPTAB.'PRESSION'         = maph ;
  154. SUPTAB.'BLOCAGES_MECANIQUES' = cdl ;
  155. SUPTAB.'OBJECTIF' = 'J' ;
  156. SUPTAB.'LEVRE_SUPERIEURE' = sulev1 ;
  157. SUPTAB.'LEVRE_INFERIEURE' = sulev2 ;
  158. SUPTAB.'COUCHE' = 4 ;
  159. SUPTAB.'FRONT_FISSURE' = ffis ;
  160. SUPTAB.'SOLUTION_RESO' = depN ;
  161. G_THETA SUPTAB ;
  162. Jc0 = suptab.resultats.global ;
  163. kc0 = (Jc0 * EE) ** 0.5 ;
  164. ic0 = kc0 / rpa / pp ;
  165.  
  166. * cas 2 : pression uniforme sur les levres
  167. SUPTAB.'MODELE' = mo2 ET mopl ;
  168. SUPTAB.'CARACTERISTIQUES' = ma2 ;
  169. SUPTAB.'PRESSION'         = mapl1 ;
  170. SUPTAB.'SOLUTION_RESO' = dep0 ;
  171. G_THETA SUPTAB ;
  172. Jc1 = suptab.resultats.global ;
  173. kc1 = (Jc1 * EE) ** 0.5 ;
  174. ic1 = kc1 / rpa / pp ;
  175.  
  176. * cas 3 : pression parabolique sur les levres
  177. SUPTAB.'MODELE' = mo2 ET mopl ;
  178. SUPTAB.'CARACTERISTIQUES' = ma2 ;
  179. SUPTAB.'PRESSION'         = mapl2 ;
  180. SUPTAB.'SOLUTION_RESO' = dep2 ;
  181. G_THETA SUPTAB ;
  182. Jc2 = suptab.resultats.global ;
  183. kc2 = (Jc2 * EE) ** 0.5 ;
  184. ic2 = kc2 / rpa / pp ;
  185.  
  186.  
  187.  
  188.  
  189. ********************
  190. * CALCUL DE L'ERREUR
  191. ********************
  192. **---> Les 2 résultats doivent être identiques
  193. **---> et converger vers la solution d'Irwin : 2/Pi
  194. **---> Le résultat ic2 doit lui converger vers 4/(3.Pi)
  195. err0 = ABS (100. * ((ic0 * pi / 2.) - 1.)) ;
  196. err1 = ABS (100. * ((ic1 * pi / 2.) - 1.)) ;
  197. err2 = ABS (100. * ((ic2 * 3. * pi / 4.) - 1.)) ;
  198.  
  199.  
  200. ****************************************
  201. * AFFICHAGE DES RESULTATS ET DES ERREURS
  202. ****************************************
  203. SAUT 5 'LIGNE' ;
  204.  
  205. mess 'Erreurs en % : ' err0 err1 err2 ;
  207.  
  208. SI (err0 > 0.4) ;
  209.   ERRE 'Erreur sur le calcul de Jc0' ;
  210. FINS ;
  211. SI (err1 > 0.4) ;
  212.   ERRE 'Erreur sur le calcul de Jc1' ;
  213. FINS ;
  214. SI (err2 > 0.5) ;
  215.   ERRE 'Erreur sur le calcul de Jc2' ;
  216. FINS ;
  217.  
  218.  
  219. FIN;
  220.  
  221.  
  222.  
  223.  
  224.  
  225.  
  226.  
  227.  
  228.  
  229.  
  230.  

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