Télécharger la plaquette Cast3M


fr en

Télécharger cmct1.dgibi

Retour à la liste

Numérotation des lignes : 

  1. * fichier :  cmct1.dgibi
  2. ************************************************************************
  3. ************************************************************************
  4. * 
  5. *     CAST TEST PORTANT SUR L'OPÉRATEUR CMCT et l'opérateur *
  6. *
  7. *      Principe :
  8. *        On considère une matrice formée de l'assemblage
  9. *           de matrice de blocage nuls cl1 et cl2
  10. *           d'une matrie de masse diagonale massc1 (LUMP)
  11. *           et un chargement formé 
  12. *           de forces reparties : chpfi
  13. *           de second membre sur les relations 
  14. *
  15. *       On resoud indirectement par condensation sur les
  16. *          multiplicateurs de Lagrange de cl3 puis par
  17. *          une redescente sur tous les noeuds.
  18. *
  19. *       On résoud directement
  20. *       
  21. *       On compare   
  22. *
  23. *
  24. *     degay 22 04 97
  25. *  
  26.  
  27.  
  28. *==========  taille du maillage ======================
  29. m = 30 ;
  30. n =  m ;
  31. E1 = 1.d15 ;
  32. r0 = 7.d35 ;
  33. graph = faux ; 
  34. * coefficient qui sert à abimer le conditionnement du système
  35. alpha = 1.d6 ;
  36. *==========  maillage ================================
  37. 'OPTI' 'DIME' 2  'MODE' 'PLAN' 'DEFO';
  38. 'OPTI' 'ELEM' 'SEG2';
  39. p1 = 0. 0. ;
  40. p2 = 10. 0. ;
  41. li1 = d p1 n p2 ;
  42. 'OPTI' 'ELEM' 'QUA4' ;
  43. su1 = 'TRAN' li1 m ( 0. 1. ) ;
  44. ls1 = 'COTE' 3 su1 ;
  45. P3 = 'POINT' su1 'PROC'  ( 10. 1. ) ;
  46. P4 = 'POINT' su1 'PROC'  ( 0.  1. ) ;
  47. lr0 = 'COTE' 4 su1 ;
  48. lr10 = 'COTE' 2 su1 ;
  49.  
  50. * rotation du maillage de 45  pour avoir des conditions 
  51. * aux limites composées
  52. 'DEPL' su1 'TOUR' 45.  ( 0. 5. ) ; 
  53.  
  54. *============== modele et matériau  ===================
  55. mod1 = 'MODE' su1 mecanique elastique ; 
  56. mod2 = 'MODE' ( li1 et ls1 ) mecanique elastique coq2 ;
  57. mat1 = 'MATE' mod1 'YOUN' E1 'NU' 0.3 'RHO' r0 ;
  58. mat2 = 'MATE' mod2 'YOUN' E1 'NU' 0.3 'RHO' r0 'EPAI' 0.05 ;
  59.  
  60. *===============  matrices classique ==================
  61. mas1 = 'LUMP' mod1 mat1 ;
  62. mas2 = 'LUMP' mod2 mat2 ;
  63. massc1 = mas1 'ET' mas2 ;
  64.  
  65. * stockage sous la forme d'un champ point
  66. chpu1 = 'MANU' 'CHPO'  su1 4 'UX' 1. 'UY' 1.  'UZ' 1.
  67.                              'RZ' 1. ;
  68. chpm1 = massc1  '*'  chpu1 ;
  69. chpm1 = 'NOMC'  chpm1  
  70.               ( 'MOTS' 'FX' 'FY' 'MZ' ) 
  71.               ( 'MOTS' 'UX' 'UY' 'RZ' ) 
  72.                 'NATU' 'DISCRET'  ;
  73. * inversion du champ point
  74. chpm1 = 'INVE' chpm1 ;
  75.  
  76. *================= conditions aux limites =============
  77. cl0 = 'SYMT' lr10 'DEPL' P2  P3 0.01 ;
  78. *  on multiplie la matrice par alpha pour tester le conditionnement
  79. cl1 = ( 'BLOQ' p1 'RZ' ) * alpha ;
  80. cl2 = 'SYMT' lr0 'DEPL' P1  P4 0.01 ; 
  81. cl3 = cl0 'ET' cl1 'ET' cl2 ;
  82.  
  83. * matrice globale du systeme
  84. rigt1 =  cl3 'ET' massc1 ;
  85.  
  86. * maillage support des multiplicateurs
  87. mailmult = 'EXTR' cl3 'MAIL' 'MULT' ;
  88.  
  89.  
  90. *================= Chargement ===========================
  91. * on impose un deplacement unité de lr10 selon l'axe de la plaque
  92. chpd0 = ( 'DEPI' cl0 1.d0 ) 'ET' ( 'DEPI' cl1 ALPHA ) ;
  93. *  on applique des forces sur le maillage
  94. chpfi =  'CHPOINT' 'ALEATOIRE'  ( mas1 et mas2 )  ; 
  95. chpfi = 'NOMC' chpfi ( 'MOTS' 'UX  ' 'UY  ' 'RZ  ' ) 
  96.                     ( 'MOTS' 'FX  ' 'FY  ' 'MZ  ' ) 
  97.          'NATU' 'DISCRET'  ;
  98.  
  99.  
  100. *============ résolution indirecte   ====================
  101. temps zero ;
  102.  
  103. *-- avec SUPER
  104. sup1 = 'SUPE' 'RIGI'  rigt1  mailmult ;
  105. rigc1 = 'EXTR' sup1 'RIGI' ;
  106. **chpfci1 = 'SUPE'  sup1 'CHAR' chpfi ;
  107. **chpd01 = 'SUPE'  sup1 'CHAR' chpd0;
  108. **depfl1 = 'RESO' rigc1 ( chpd01 'ET' chpfci1 ) ;
  109. * on redescend
  110. **chpf1 = 'REAC' cl3 depfl1 ;
  111. **dep1 = 'RESO'  massc1 ( chpf1 'ET' chpfi );
  112. chpfci1 = 'SUPE'  sup1 'CHAR' (chpd0 'ET' chpfi);
  113. depfl1 = 'RESO' rigc1 chpfci1 ;
  114. dep1 = super depl sup1 depfl1 (chpd0 'ET' chpfi);
  115. chpf1 = 'REAC' cl3 dep1 ;
  116. dep1 = 'ENLE' dep1 'LX';
  117.  
  118. TABTPS = TEMP 'NOEC';
  119. t0     = TABTPS.'TEMPS_CPU'.'INITIAL' ;
  120. temps zero ;
  121.  
  122. *-- avec CMCT
  123.  
  124. rigc2 =  'CMCT' cl3 chpm1 ;
  125. mac0 sgac0 = temps sgac ;
  126. depsl2 = '*'
  127.             chpm1 (MOTS 'UX' 'UY' 'RZ') 
  128.             chpfi (MOTS 'FX' 'FY' 'MZ') (MOTS 'UX' 'UY' 'RZ') ;
  129. mac1 sgac1  = temps sgac  ;
  130. difsg = sgac1 '-' sgac0 ;
  131. difmsg = mac1 '-' mac0 ;
  132.  
  133. * second membre du systeme sur les multiplicateurs
  134. fcl2 =  cl3 '*' depsl2 '-' chpd0  ;
  135. depfl2 =  'RESO'  rigc2 fcl2 ;
  136.  
  137. chpf2 = ( 'REAC' cl3 depfl2 ) ;
  138.  
  139. depcor2 =    '*' 
  140.                  chpm1  (MOTS 'UX' 'UY' 'RZ')
  141.                  chpf2  (MOTS 'FX' 'FY' 'MZ')
  142.                         (MOTS 'UX' 'UY' 'RZ')
  143.               'NATURE' 'DIFFUS'  ;
  144.  
  145. dep2 = depsl2 + depcor2 ;
  146.  
  147.  
  148. TABTPS = TEMP 'NOEC';
  149. t01 = TABTPS.'TEMPS_CPU'.'INITIAL' ;
  150. t01 = t01 + 1;
  151.  temps zero ;
  152.  'MESS' 'Nombre de segment Active' difsg ;
  153.  'MESS' 'Memoire correspondante' difmsg ;
  154. errsg = >eg difsg 17;
  155. * 2019/12/11 SG Le nouveau paradigme (fiche anomalie 10273) étant de ne
  156. * pas désactiver les segments, on commente ce test
  157. *new-paradigm 'SI' errsg ;
  158. *new-paradigm      'MESS' 'Erreur dans la desactivation des segments' ; 
  159. *new-paradigm      'ERREUR' 5 ;
  160. *new-paradigm 'FINSI' ;
  161.  
  162. *============== resolution directe =======================
  163. *
  164. dep3 = 'RESO'   rigt1  ( chpd0 'ET' chpfi ) ;
  165. chpf3 = 'REAC' cl3 dep3 ;
  166. dep3 = 'ENLE'  dep3  'LX' ; 
  167.  
  168. TABTPS = TEMP 'NOEC';
  169. t1 = TABTPS.'TEMPS_CPU'.'INITIAL' ;
  170. t1 = t1 + 1;
  171.  
  172. 'OPTI' 'ECHO' 0 ;
  173. 'SAUTER' 1 'LIGNE' ;
  174. 'MESS' 'Par rapport à la méthode directe : ' ;
  175. 'MESS' 'SUPER' ( 100. * t0 / t1 - 100.) '% de temps' ; 
  176. 'MESS' 'CMCT ' ( 100. *  t01 / t1  - 100.)'% de temps - facteur'
  177.  (1. * t1 / t01 ) ; 
  178. 'SAUTER' 1 'LIGNE' ; 
  179. 'OPTI' 'ECHO' 1 ;
  180.  
  181.  
  182.  
  183. *=============== Comparaison  ============================
  184. * 
  185.  
  186. 'OPTI' 'ECHO' 0 ;
  187. * Erreur sur le deplacement 
  188. deperr =  dep3 * 2.  - dep2 - dep1 / 2. ;
  189. erx = 'EXCO' deperr 'UX' 'SCAL' / ( 'MAXI' 'ABS' ( 'EXCO'  dep3 'UX'));
  190. errx = ( 'MAXI' 'ABS' erx) ;
  191. MESS 'Erreur sur X' ( errx  * 100.) '%'  ;
  192. ery = 'EXCO' deperr 'UY' 'SCAL' / ( 'MAXI' 'ABS' ( 'EXCO'  dep3 'UY'));
  193. erry = ( 'MAXI' 'ABS' ery) ;
  194. MESS 'Erreur sur Y' ( erry * 100.) '%'  ;
  195. erz = 'EXCO' deperr 'RZ' 'SCAL' ;
  196. errz = ( 'MAXI' 'ABS' erz) ;
  197. MESS 'Erreur sur RZ' (errz  * 1.) 'difference absolue'  ;
  198. errtot = ( erx ** 2 ) + ( ery ** 2 ) ; 
  199. 'OPTI' 'ECHO' 1 ;
  200.  
  201. 'SI' graph ;
  202. 'TITR' 'Niveau d erreur' ;
  203. 'TRAC' errtot su1 ; 
  204. 'FINSI' ;
  205.  
  206. * erreur sur les forces 
  207. ferr = chpf3 * 2. '-' chpf2 - chpf1 ;
  208. errf = 'MAXI' 'ABS' ferr / ( 'MAXI' 'ABS' chpf3 ) ;
  209. 'MESS' 'Erreur sur les forces' ( errf '*' 100.) '%' ;
  210.  
  211.  
  212. *=====================================
  213. *      Code de fonctionnement
  214.  
  215. 'SI' ((errx > 1.D-9) 'OU' (erry > 1.D-9) 'OU' (errf > 1.D-9 ));
  216.    'ERREUR' 5 ;
  217. 'SINON' 
  218.    'ERREUR' 0 ;
  219. 'FINSI' ;
  220.  
  221. 'FIN' ;
  222.  
  223.  
  224.  
  225.  
  226.  
  227.  
  228.  
  229.  
  230.  
  231.  
  232.  
  233.  
  234.  
  235.  
  236.  
  237.  
  238.  
  239.  
  240.  
  241.  
  242.  

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