Télécharger la plaquette Cast3M


fr en

Télécharger fixhc.procedur

Retour à la liste

Numérotation des lignes : 

  1. * FIXHC     PROCEDUR  BR232186  16/12/05    21:39:28     9238           
  2. ****************************************************
  3. * Procedure de tracking par resolution d'un probleme
  4. * de conduction thermique stationnaire
  5. * Application au modele E-FEM
  6. * Version  : 
  7. *     Fixed global tracking
  8. * Creation : 
  9. *     Francesco RICCARDI 22/08/2016
  10. *     Benjamin RICHARD
  11. * Contact :  
  12. *     Francesco.Riccardi[at]cea.fr
  13. *     Benjamin.Richard[at]cea.fr 
  14. * Institution : 
  15. *     CEA\DEN\DANS\DM2S\SEMT\EMSI
  16. ****************************************************
  17. * Voir notice TRACKING pour plus d'informations
  18. ****************************************************
  19. *
  20. *******----DEBUT Procedure FIXHC----*******
  21. *
  22. DEBP FIXHC TAB2*'TABLE';
  23.  
  24.  
  25. * Indice d'evolution
  26. nn = dime tab2.deplacements;
  27. nn = nn-1;
  28.  
  29. * Modeles
  30. MODEL1 = TAB2. 'MODELE';
  31. MODELS = EXTR MODEL1 'ZONE';
  32.  
  33. * Modele et materiau EFEM
  34. MODEL11 = MODELS. 1;
  35. SURFACE1 = MODELS. 2;
  36. MATr1 = TAB2. 'WTABLE'. 'CARACTERISTIQUES';
  37. MATr11 = REDU MATr1 SURFACE1;
  38. E11 = EXCO MATr11 YOUN;
  39. NU11 = EXCO MATr11 NU;
  40. P11 = EXCO MATr11 FT;
  41. P12 = EXCO MATR11 XNX;
  42. P13 = EXCO MATR11 XNY;
  43. P14 = EXCO MATr11 IND1;
  44. P15 = EXCO MATR11 XIE;
  45. P16 = EXCO MATR11 YIE;
  46.  
  47. * Resistance en traction
  48. STMAX = MAXI P11;
  49.  
  50. * Coordonnees de la geometrie
  51. X0 = COOR 1 SURFACE1;
  52. Y0 = COOR 2 SURFACE1;
  53.  
  54. * Maximum des contraintes principales
  55. SIG1 = TAB2. 'CONTINUATION'. 'CONTRAINTES';
  56. PRIN1 = PRIN SIG1 MODEL11;
  57. SIGMA11 = EXCO 'SI11' PRIN1;
  58. SIGMAX = MAXI SIGMA11;
  59.  
  60. * Variables internes : NFLA
  61. FLAG = EXCO (TAB2. 'VARIABLES_INTERNES'.NN) 'NFLA';
  62. VFLAG = MAXI FLAG;
  63. *TRAC FLAG MODEL11;
  64.  
  65. * Nombre de root elements
  66. NROOTS = TAB2. 'NROOTS'; 
  67. *LIST NROOTS;
  68.  
  69. * Nombre maximum de fissures
  70. MAXR = TAB2. 'NRMAX'; 
  71.  
  72. * 
  73. *----- Debut algorithme de tracking --------------------------- 
  74. * 
  75. * Cas1 : on fixe le nombre de fissures
  76. SI ((SIGMAX > 0) ET (NROOTS < MAXR));
  77. * Cas2 : on fixe le pas de temps
  78. *SI ((SIGMAX > 0) ET (NN < (TAB2. 'NPAS_TRACKING')));
  79. * Cas3 : pas de limitations
  80. *SI (SIGMAX > 0);
  81.  
  82.  
  83. * Points ou la contrainte depasse la resistance en traction
  84. SI (VFLAG EGA 0);
  85.         PSLIM = SIGMA11 POIN 'MAXI';    
  86.         NLIM = 1;
  87.         SIGMAX0 = (SIGMA11 MASQ 'EGSUPE' SIGMAX) * SIGMA11;
  88. SINON;
  89.         MTEST = SURFACE1 DIFF (TAB2. 'MESH');
  90.         SIGTEST = REDU SIGMA11 MTEST; 
  91.         EL_LIM = SIGTEST ELEM 'EGSUPE' STMAX;
  92.         PSLIM = SIGMA11 POIN 'EGSUPE' STMAX;
  93.         SIGMAX0 = (SIGTEST MASQ 'EGSUPE' STMAX) * SIGTEST;
  94.         SI (VFLAG EGA 0);
  95.                 NLIM = 0;
  96.         SINON;
  97.                 NLIM = NBEL EL_LIM;
  98.         FINSI;
  99. FINSI;
  100.  
  101. * On ordonne les possibles root elements par contrainte decr
  102. Prr = TABLE;
  103. Prr. 'coor_x' = TABLE;
  104. Prr. 'coor_y' = TABLE;
  105. Srr = TABLE;
  106. Trr = TABLE;
  107. ff = 0;
  108. REPETER ORDINE (NROOTS + NLIM);
  109.         ff = (ff + 1);
  110.         SI ((NROOTS > 0) ET (ff &lt;EG NROOTS));
  111.                 Prr. 'coor_x'. ff = 
  112.                 COOR 1 (TAB2. 'PROOTS'. ff);
  113.                 Prr. 'coor_y'. ff = 
  114.                 COOR 2 (TAB2. 'PROOTS'. ff);
  115.                 Trr. ff = CHAN 'POI1' (TAB2. 'EROOTS'. ff);
  116.                 Srr. ff = 0;
  117.         SINON;
  118.                 SIGMAX1 = MAXI SIGMAX0;         
  119.                 Prr. 'coor_x'. ff = MAXI (COOR 1 (SIGMAX0 POIN 'MAXI'));
  120.                 Prr. 'coor_y'. ff = MAXI (COOR 2 (SIGMAX0 POIN 'MAXI'));
  121.                 Trr. ff = CHAN 'POI1' (SIGMAX0 ELEM 'MAXI');            
  122.                 Srr. ff = SIGMAX1;              
  123.                 XAX = SIGMAX0 MASQ 'EGALE' SIGMAX1;             
  124.                 DIFF1 = XAX * SIGMAX1;          
  125.                 SIGMAX0 = (SIGMAX0 - DIFF1);            
  126.         FINSI;
  127. FIN ORDINE;
  128.  
  129. *
  130. *\\\\ Resolution du probleme de conduction ////*
  131. *
  132. * Cosinus des tangentes
  133. COS2X = EXCO PRIN1 'COX2';
  134. COS2Y = EXCO PRIN1 'COY2';
  135. COS2Z = EXCO PRIN1 'COZ2';
  136.  
  137. * Champ vectoriel des tangentes
  138. *SIGMA22 = MANU 'CHML' MODEL11 'SI22' 1. 'TYPE' 
  139. *'CONTRAINTES PRINCIPALES' 'STRESSES';
  140. *PRIN12 = (SIGMA22 ET COS2X ET COS2Y ET COS2Z);
  141. *VECTAN = VECT PRIN12 MODEL11 0.001 'SI22';  
  142.  
  143. * Modele thermique
  144. MOTH = MODE SURFACE1 THERMIQUE ANISOTROPE;
  145.  
  146. * Perturbation conductivite numerique
  147. TTPETIT = 1.E-6;
  148.  
  149. * Tenseur conductivite anisotrope
  150. KXX=(COS2X**2)+((COS2X**0)*TTPETIT);
  151. KYY=(COS2Y**2)+((COS2Y**0)*TTPETIT);
  152. TTX=NOMC 'UX' COS2X;
  153. TTY=NOMC 'UX' COS2Y;
  154. KXY=TTX*TTY;
  155. PX = 1.0 0.;
  156. MATH= MATE MOTH 'DIRECTION' PX K11 KXX K21 KXY K22 KYY;
  157. * Matrice de rigidite
  158. KTH = COND MOTH MATH;
  159. * Conditions aux limites
  160. PCL1 = CHAN 'POI1' (TAB2. 'LBC'); 
  161. CLT1 = BLOQ PCL1 'T';
  162. CLTH = CLT1;
  163. * Temperature imposee
  164. VBC = TAB2. 'BCSTH';
  165. DE1 = DEPI CLT1 VBC;
  166. DETH = DE1;
  167. * Assemblage
  168. KTHTOT = KTH ET CLTH;
  169. *Resolution
  170. RESTH = RESO KTHTOT DETH;
  171. TAB2. 'RESTH' = RESTH;
  172. *TRAC RESTH SURFACE1;
  173.  
  174. * 
  175. *\\\\ Extraction des isovaleurs ////*   
  176. * 
  177. T0 = CHAN 'CHAM' RESTH MOTH 'NOEUD';
  178. T0 = EXCO T0 'T';
  179. i = 0;
  180. NEW_MESH = VIDE 'MAILLAGE';
  181. EL0 = SURFACE1 ELEM 'CONTENANT' (0. 0.);
  182. P0 = SURFACE1 POIN 'PROC' (0. 0.);
  183. P0 = MANU 'POI1' P0;
  184. REPETER BOUCLEEL (NROOTS + NLIM);
  185.         i = i + 1;      
  186.         PPP = (Prr. 'coor_x'. i) (Prr. 'coor_y'. i);
  187.         ELACTIF = 0;    
  188.         TRACCIA = FAUX;         
  189. * Verification des elements
  190. ***On a trois cas : 
  191. ***ELACTIF = 0 : possibles elements root*** 
  192. ***ELACTIF = 1 : elements appartenant a une isovaleur*** 
  193. ***ELACTIF = 2 : elements a negliger***
  194.         SI (i &lt;EG NROOTS);
  195.                 ELACTIF = 1;
  196.                 TRACCIA = VRAI; 
  197.         SINON;
  198.                 CHECK1 = NEW_MESH ELEM 'CONTENANT' PPP 'NOVERIF';
  199. * CHECK1 : on verifie si PPP appartient a une isovaleur active          
  200.                 SI ((NBEL CHECK1) > 0);
  201.                         ELACTIF = 1;    
  202.                 FINSI;
  203.         FINSI;  
  204. * CHECK2 : les elements root doivent etre bien separes  
  205.         SI (ELACTIF EGA 0);     
  206.                 SI (NROOTS EGA 0);
  207.                         TRACCIA = VRAI;
  208.                 SINON;
  209.                         CHECK2 = NEW_MESH ELEM 'APPUYE' 'LARGEMENT'
  210.                         (Trr. i) 'NOVERIF';
  211.                         SI ((NBEL CHECK2) EGA 0);
  212.                                 TRACCIA = VRAI;
  213.                         SINON;
  214.                                 ELACTIF = 2;
  215.                         FINSI;
  216.                 FINSI;
  217.         FINSI;  
  218. * Isovaleurs    
  219.         SI (TRACCIA EGA VRAI);
  220.                 TREDUC = REDU RESTH (Trr. i);
  221.                 NA = (Trr. i) POIN 1;
  222.                 NB = (Trr. i) POIN 2;
  223.                 NC = (Trr. i) POIN 3;   
  224.                 TA = EXTR TREDUC 'T' NA;
  225.                 TB = EXTR TREDUC 'T' NB;
  226.                 TC = EXTR TREDUC 'T' NC;
  227.                 TG = ((TA + TB + TC) / 3.);
  228.                 ISOVG = ISOV T0 TG;
  229.                 P_ISOVG = CHAN 'POI1' ISOVG; 
  230. * Extraction des elements traverses par l'isovaleur 
  231.                 T0G = MANU 'CHML' MOTH 'T' TG 'NOEUD';
  232.                 T01 = T0 - T0G;
  233.                 T02 = T01 MASQ 'SUPERIEUR' 0.;
  234.                 MAILA = T02 ELEM 'EGALE' 1.;
  235.                 MAILB = T02 ELEM 'EGALE' 0.;
  236.                 MESH_EL = SURFACE1 DIFF (MAILA ET MAILB);       
  237. * CHECK3 : les isovaleurs doivent etre bien separees 
  238. * Si le seul element en commun est EL0 alors on continue
  239.                 SI (i > NROOTS);
  240.                         MM1 = (CHAN 'POI1' MESH_EL) ET P0;
  241.                         MMT = (CHAN 'POI1' NEW_MESH) ET P0; 
  242.                         CHECK3 = MMT INTE MM1;
  243.                         SI ((NBEL CHECK3) > 1);
  244.                                 ELACTIF = 2;
  245.                         FINSI;
  246.                 FINSI;          
  247.                 SI (ELACTIF EGA 1);
  248.             NEW_MESH = NEW_MESH ET (TAB2. 'MESH_ISO'. i);
  249.                 FINSI;  
  250.                 SI (ELACTIF EGA 0);
  251.                         NEW_MESH = NEW_MESH ET MESH_EL;
  252.                         TAB2. 'NROOTS' = (TAB2. 'NROOTS') + 1;     
  253.             TAB2. 'ISOTOT'. (TAB2. 'NROOTS') = ISOVG;
  254.             TAB2. 'MESH_ISO'. (TAB2. 'NROOTS') = MESH_EL;
  255.                         TAB2. 'PROOTS'. (TAB2. 'NROOTS') = PPP;
  256.                         TAB2. 'EROOTS'. (TAB2. 'NROOTS') = SURFACE1 ELEM 
  257.                         'CONTENANT' PPP;
  258.                 FINSI;
  259.         FINSI;
  260.         DROOTS = (TAB2. 'NROOTS') - NROOTS;
  261.         SI (DROOTS EGA 1);
  262.                 QUIT BOUCLEEL;
  263.         FINSI;
  264.         SI (i EGA (NROOTS + NLIM));
  265.                 QUIT BOUCLEEL;
  266.         FINSI;
  267. FIN BOUCLEEL;           
  268. TAB2. 'MESH' = NEW_MESH;
  269.  
  270. *
  271. *\\\\ Creation des chamelem ////*       
  272. *
  273. CHM = MANU 'CHML' model11 'SI11' 0. 
  274. 'TYPE' 'CONTRAINTES PRINCIPALES' 'STRESSES';
  275. X0S = CHAN 'CHAM' X0 MODEL11 'STRESSES';
  276. X0S = CHAN 'TYPE' X0S 'CONTRAINTES PRINCIPALES';
  277. X0S = CHAN 'COMP' 'SI11' X0S;
  278. Y0S = CHAN 'CHAM' Y0 MODEL11 'STRESSES';
  279. Y0S = CHAN 'TYPE' Y0S 'CONTRAINTES PRINCIPALES';
  280. Y0S = CHAN 'COMP' 'SI11' Y0S;
  281. nlm1 = nbel NEW_MESH;
  282. i=0;
  283. repeter BOUCLEB nlm1;
  284.         i = i + 1;
  285.         ELEMi = NEW_MESH ELEM i;
  286.         VALXi = MAXI (REDU X0S ELEMi);
  287.         VALYi = MAXI (REDU Y0S ELEMi);
  288.         CHPX = X0S MASQ 'EGALE' VALXi;
  289.         CHPY = Y0S MASQ 'EGALE' VALYi;
  290.         CHP = CHPX * CHPY;
  291.         CHM = CHM + CHP;
  292. fin BOUCLEB;
  293.  
  294. * Abscisses des points d'entree des fissures
  295. CHIX = MANU 'CHML' MODEL11 'SI11' 0. 'TYPE' 
  296. 'CONTRAINTES PRINCIPALES' 'STRESSES';
  297. * Ordonnees des points d'entree des fissures
  298. CHIY = CHIX;
  299. * Composante x des normales
  300. CHNX = CHIX;
  301. * Composante y des normales
  302. CHNY = CHIX;
  303. j = 0;
  304. * Debut de la boucle sur les root elements
  305. REPETER BOUCLE1 (TAB2. 'NROOTS');
  306.     j = j + 1; 
  307.     ISOVGj = ORDO (TAB2. 'ISOTOT'. j);
  308.     ROOTj = TAB2. 'EROOTS'. j;
  309.     MESHj = TAB2. 'MESH_ISO'. j;
  310.     NUMEL = NBEL MESHj;   
  311.     CHIXj = CHIX * 0.;
  312.     CHIYj = CHIY * 0.;
  313.     CHNXj = CHNX * 0.;
  314.     CHNYj = CHNY * 0.;
  315.     PINI = ISOVGj POIN 'INITIAL';
  316.     PFIN = ISOVGj POIN 'FINAL';
  317.     TEST1 = ROOTj ELEM 'CONTENANT' PINI 'NOVERIF';
  318.     SI ((NBEL TEST1) EGA 1);
  319.         SEG1 = ISOVGj ELEM 'CONTENANT' PINI;
  320.     SINON;
  321.         SEG1 = ISOVGj ELEM 'CONTENANT' PFIN;
  322.         PINI = PFIN;
  323.     FINSI;
  324. * Debut de la boucle sur les elements de l'isovaleur j
  325.     i = 0;
  326.     REPETER BOUCLE2 NUMEL;
  327.         i = i + 1;
  328.         SI (i EGA 1);
  329.             ELi = ROOTj;
  330.             MTEST = ROOTj;
  331.             SEGi = SEG1;
  332.             STEST = SEG1;
  333.             PTEST = PINI;
  334.         SINON;
  335.             MTEST = MTEST ET ELi;
  336.             SI (i > 2);
  337.                 STEST = STEST ET SEGi;
  338.             FINSI;
  339.             PTEST1 = STEST POIN 'INITIAL';
  340.             PTEST2 = STEST POIN 'FINAL';
  341.             SI (PTEST1 EGA PINI);
  342.                 PTEST = PTEST2;
  343.             SINON;
  344.                 PTEST = PTEST1;
  345.             FINSI;
  346.             SDIFF = ISOVGj DIFF STEST;        
  347.             SEGi = SDIFF ELEM 'CONTENANT' PTEST;
  348.         FINSI;
  349.         PP = MANU 'POI1' PTEST;
  350.         XTEST YTEST = COOR PTEST;
  351.         COORXi = REDU X0S ELi;
  352.         COORYi = REDU Y0S ELi;
  353.         VALXi = MAXI COORXi;
  354.         VALYi = MAXI COORYi;    
  355.         PP1 = SEGi POIN 1;
  356.         XP1 YP1 = COOR PP1; 
  357.         PP2 = SEGi POIN 2;
  358.         XP2 YP2 = COOR PP2;
  359. * Calcul de la tangente et de la normale a la fissure 
  360.         LL = (((XP1 - XP2) ** 2) + ((YP1 - YP2) ** 2)) ** 0.5;
  361.         SI (YP1 > YP2);
  362.             TY = (YP1 - YP2) / LL;
  363.             TX = (XP1 - XP2) / LL;
  364.         SINON;
  365.             TY = (YP2 - YP1) / LL;
  366.             TX = (XP2 - XP1) / LL;
  367.         FINSI;
  368.         SI (YP1 EGA YP2);
  369.             TY = 0.;
  370.             TX = -1.;
  371.         FINSI;
  372.         NX = TY;
  373.         NY = -1.*TX;  
  374.         CHPX = X0S MASQ 'EGALE' VALXi;
  375.         CHPY = Y0S MASQ 'EGALE' VALYi;
  376.         CH = CHPX * CHPY;
  377.         CHIXj = CHIXj + (CH * XTEST);
  378.         CHIYj = CHIYj + (CH * YTEST);
  379.         CHNXj = CHNXj + (CH * NX);
  380.         CHNYj = CHNYj + (CH * NY);
  381.         MREDj = REDU MODEL11 MESHj;  
  382. * Pour le pas suivant    
  383.         SI (i < NUMEL);
  384.             CTEST = CONT MTEST;
  385.             MDIFF = MESHj DIFF MTEST;
  386.             CDIFF = CONT MDIFF;
  387.             CINTE = CDIFF INTE CTEST;
  388.             GINTE = BARY CINTE;
  389.             ELi = MDIFF ELEM 'CONTENANT' GINTE;
  390.         FINSI;
  391.     FIN BOUCLE2;
  392.     CHIX = CHIX + CHIXj;
  393.     CHIY = CHIY + CHIYj;
  394.     CHNX = CHNX + CHNXj;
  395.     CHNY = CHNY + CHNYj; 
  396. FIN BOUCLE1;
  397. MATr2 = MATE MODEL11 YOUN E11 NU NU11 FT P11 XNX CHNX
  398. XNY CHNY IND1 CHM XIE CHIX YIE CHIY;
  399. MATr1 = MATr2;
  400. TAB2. 'WTABLE'. 'CARACTERISTIQUES' = MATR1;
  401.  
  402.  
  403. FINSI;
  404.  
  405.  
  406. FINP MATR1;
  407.  
  408.  

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