Télécharger xxcre1.eso

Retour à la liste

Numérotation des lignes :

  1. C XXCRE1 SOURCE CB215821 17/11/30 21:17:47 9639
  2. SUBROUTINE XXCRE1(TAU,SIG,EPSV,VAR,SIGD,EPSVD,VARD,EPSTHD,
  3. &DSPT,XMAT,XMAINF,XMASUP,NSTRS,NVARI,NCOMAT,DD,DDV,DDINV,
  4. &DDINVp,YOUNG,NYOUNG,XNU,NXNU,MFR,XCAR,ICARA,IFOURB,INDIC,TI,
  5. &TPOINT,TINF,TSUP,ITEST,ITHHER,VARF,IB,IGAU,KERRE)
  6. C-----------------------------------------------------------------------
  7. C Objet: Cette subroutine calcule les derivees des variables internes
  8. C pour un materiau viscoplastique a endommagement et ecrouissage
  9. C isotrope en regime anisotherme de 2 manieres au choix suivant
  10. C la valeur de ITEST,ITHHER
  11. C-----------------------------------------------------------------------
  12. C
  13. C-----------------------------------------------------------------------
  14. C Entree: TAU pas de temps
  15. C TI temperature au debut du pas
  16. C TPOINT derivee de la temperature
  17. C TINF et TSUP bornes entre lesquelles TI est comprise
  18. C SIG(NSTRS) tenseur des contraintes
  19. C EPSV(NSTRS) tenseur des deformations
  20. C VAR(NVARI) tableau contenant les variables internes
  21. C pilotant les equations
  22. C VAR(1)=p ; VAR(2)=D11 ; VAR(3)=D22; VAR(4)=D33;
  23. C VAR(5)=D12; VAR(6)=D13; VAR(7)=D23; VAR(8)=DMAX
  24. C EPSTHD vitesse de dilatation thermique au debut du pas
  25. C DSPT(NSTRS,NSTRS) vitesse de deformation totale
  26. C XMAT(NCOMAT) tableau des parametres scalaires du materiau
  27. C a une temperature T donnee
  28. C XMAINF(NCOMAT) tableau des parametres scalaires du materiau
  29. C a la temperature TINF
  30. C XMASUP(NCOMAT) tableau des parametres scalaires du materiau
  31. C a la temperature TSUP
  32. C XMAT(1)=YOUNG ;XMAT(2)=XNU ;XMAT(3)=N
  33. C XMAT(4)=M ;XMAT(5)=KK; XMAT(6)=ALF2;
  34. C XMAT(7)=BET2 ;XMAT(8)=r; XMAT(9)=A ;
  35. C XMAT(10)=q ;XMAT(11)=ALPHAT
  36. C XMAT(12)=RHO; XMAT(13)=SIGY
  37. C MFR indice de la formulation mecanique(seulement massif
  38. C ou coque pour les materiaux endommageables
  39. C ICARA nombre de caracteristiques geometriques des elements
  40. C finis
  41. C XCAR(ICARA) tableau des caracteristiques geometriques des
  42. C elements finis
  43. C IFOURB= -2 EN CONTR. PLANES
  44. C -1 EN DEFORM. PLANES
  45. C 0 EN AXISYMETRIE
  46. C 1 EN SERIE DE FOURIER
  47. C 2 EN TRIDIM
  48. C INDIC=0, 1 OU -1 pour plasticite avec endommagement
  49. C =2 OU -2 pour viscoplasticite avec endommagement
  50. C ITEST = 0 pas uniformite des listes de temperatures pour tous les
  51. C coefficients non lineaires
  52. C d'ou interpolation sur les coefficients du materiau
  53. C = 1 uniformite des listes de temperatures pour tous les
  54. C coefficients non lineaires
  55. C d'ou moyenne ponderee sur la variable elle-meme
  56. C ex: D(T)=A(T)*K**B(T)
  57. C T=TETA*TINF+(1-TETA)*TSUP
  58. C si ITEST=0 on remplace A(T) par (TETA*A(TINF)+(1-TETA))*A(TSUP)
  59. C et B(T) par (TETA*B(TINF)+(1-TETA))*B(TSUP)
  60. C si ITEST=1 D(T)=TETA*D(TINF)+(1-TETA)*D(TSUP)
  61. C ITHHER = 0 pas de chargement thermique et materiau constant
  62. C = 1 chargement thermique et materiau constant
  63. C = 2 chargement thermique et materiau(T)
  64. C------------------------------------------------------------------------
  65. C
  66. C------------------------------------------------------------------------
  67. C Sortie: EPSVD(NSTRS) derivee du tenseur des deformations
  68. C VARD(NVARI) tableau contenant les derivees des variables
  69. C internes
  70. C VARF(NVARI) état des variables internes à la fin du sous programme
  71. C SIGD(NSTRS) derivee des contraintes
  72. C DD(NSTRS,NSTRS) matrice de Hooke au debut du pas
  73. C DDV(NSTRS,NSTRS) derivee de DD
  74. C DDINV(NSTRS,NSTRS) inverse de DD
  75. C------------------------------------------------------------------------
  76. IMPLICIT INTEGER(I-N)
  77. IMPLICIT REAL*8 (A-H,O-Z)
  78. DIMENSION XMAINF(*),XMASUP(*),XCAR(*)
  79. DIMENSION SIG(*),EPSV(*),VAR(*)
  80. DIMENSION SIGD(*),EPSVD(*),VARD(*),VARF(NVARI)
  81. DIMENSION XMAT(*),DSPT(*),YOUNG(*),XNU(*),EPSTHD(*)
  82. DIMENSION DDV(NSTRS,NSTRS),DDINVp(NSTRS,NSTRS)
  83. DIMENSION DD(NSTRS,NSTRS),DDINV(NSTRS,NSTRS)
  84. DIMENSION VARD1(3),EPSVD1(6),SIGD1(6),VARD2(3),EPSVD2(6)
  85. DIMENSION EPSED(6),XX(6),Di(3),VaP(3),AA(3),BB(6)
  86. DIMENSION ROT(3,3),ROTT(3,3),CO(3,3),COM(3,3)
  87. DIMENSION D1(3,3),D2(3,3),DV1(3,3),DI1(3,3),DI2(3,3)
  88. DIMENSION AAA(3,3),BBB(3,3),CCC(3,3),VcP(3,3),FTH(3,3)
  89. PARAMETER (UN=1.D0)
  90. LOGICAL endonul
  91. C
  92. C----- le materiau est constant ou pas, on derive normalement les contraintes
  93. C
  94. IFLAG=0
  95. IF ( (ITEST.EQ.1).AND.(ITHHER.EQ.2) ) THEN
  96. C----- si le materiau depend de la temperature et si les listes sont identiques
  97. IF (TINF.NE.TSUP) THEN
  98. C-------------- Moyenne ponderee sur les vitesses des variables
  99. IFLAG=1
  100. TETA=(TI-TSUP)/(TINF-TSUP)
  101. CALL DERIV1(TAU,SIG,EPSV,VAR,SIGD1,EPSVD1,VARD1,EPSTHD,
  102. & DSPT,XMAINF,NSTRS,NVARI,NCOMAT,DD,DDV,DDINV,DDINVp,
  103. & YOUNG,NYOUNG,XNU,NXNU,MFR,XCAR,ICARA,IFOURB,INDIC,
  104. & TINF,TPOINT,ITHHER,VARF,IB,IGAU,KERRRE)
  105. IF(KERRE.NE.0) GO TO 100
  106.  
  107. DO 10 I=1,NVARI
  108. 10 VARD2(I)=TETA*VARD1(I)
  109. DO 20 I=1,NSTRS
  110. 20 EPSVD2(I)=TETA*EPSVD1(I)
  111. CALL DERIV1(TAU,SIG,EPSV,VAR,SIGD1,EPSVD1,VARD1,EPSTHD,
  112. & DSPT,XMASUP,NSTRS,NVARI,NCOMAT,DD,DDV,DDINV,DDINVp,
  113. & YOUNG,NYOUNG,XNU,NXNU,MFR,XCAR,ICARA,IFOURB,INDIC,
  114. & TSUP,TPOINT,ITHHER,VARF,IB,IGAU,KERRRE)
  115. IF(KERRE.NE.0) GO TO 100
  116.  
  117. DO 30 I=1,NVARI
  118. 30 VARD2(I)=VARD2(I)+(1.D0-TETA)*VARD1(I)
  119. DO 40 I=1,NSTRS
  120. 40 EPSVD2(I)=EPSVD2(I)+(1.D0-TETA)*EPSVD1(I)
  121. ENDIF
  122. ENDIF
  123. C
  124. CALL DERIV1(TAU,SIG,EPSV,VAR,SIGD,EPSVD,VARD,EPSTHD,
  125. &DSPT,XMAT,NSTRS,NVARI,NCOMAT,DD,DDV,DDINV,DDINVp,YOUNG,
  126. &NYOUNG,XNU,NXNU,MFR,XCAR,ICARA,IFOURB,INDIC,TI,TPOINT,
  127. &ITHHER,VARF,IB,IGAU,KERRE)
  128. IF(KERRE.NE.0) GO TO 100
  129. C
  130. IF (IFLAG.EQ.0) GOTO 100
  131. C
  132. C --- Cas ou ITEST = 1 et ITHHER = 2
  133. C
  134. DO I=1,NVARI
  135. VARD(I)=VARD2(I)
  136. ENDDO
  137. DO I=1,NSTRS
  138. EPSVD(I)=EPSVD2(I)
  139. ENDDO
  140. C
  141. endonul = .FALSE.
  142. IF ((VAR(2).EQ.0).AND.(VAR(3).EQ.0).AND.(VAR(4).EQ.0).AND.
  143. & (VAR(5).EQ.0).AND.(VAR(6).EQ.0).AND.(VAR(7).EQ.0).AND.
  144. & (VAR(8).EQ.0)) THEN
  145. endonul = .TRUE.
  146. ENDIF
  147. C
  148. C----------------------------------------------
  149. C TENSEUR DES CONTRAINTES (dans repere general)
  150. C----------------------------------------------
  151. C
  152. CALL ZERO(CO,3,3)
  153. IF (mfr.eq.1) THEN
  154. DO I=1,3
  155. CO(I,I) = SIG(I)
  156. ENDDO
  157. CO(1,2) = SIG(4)
  158. IF(NSTRS.NE.6) GOTO 50
  159. CO(1,3) = SIG(5)
  160. CO(2,3) = SIG(6)
  161. 50 CONTINUE
  162. CO(2,1) = CO(1,2)
  163. CO(3,1) = CO(1,3)
  164. CO(3,2) = CO(2,3)
  165. ENDIF
  166. C
  167. C------------------------------------
  168. C TENSEUR D ENDOMMAGEMENT ET
  169. C TENSEUR D ENDOMMAGEMENT DIAGONALISE
  170. C------------------------------------
  171. C
  172. XD11=VAR(2)
  173. XD22=VAR(3)
  174. XD33=VAR(4)
  175. XD12=VAR(5)
  176. XD13=VAR(6)
  177. XD23=VAR(7)
  178. CALL ZERO(D1,3,3)
  179. D1(1,1) = XD11
  180. D1(2,2) = XD22
  181. D1(3,3) = XD33
  182. D1(1,2) = XD12
  183. D1(1,3) = XD13
  184. D1(2,3) = XD23
  185. D1(2,1) = D1(1,2)
  186. D1(3,1) = D1(1,3)
  187. D1(3,2) = D1(2,3)
  188. C
  189. CALL JACOB3(D1,3,Di,ROT)
  190. C
  191. CALL ZERO(D2,3,3)
  192. D2(1,1) = Di(1)
  193. D2(2,2) = Di(2)
  194. D2(3,3) = Di(3)
  195. w1 = 1.D0 - Di(1)
  196. w2 = 1.D0 - Di(2)
  197. w3 = 1.D0 - Di(3)
  198. IF (D2(1,1).EQ.1.D0) w1=1.D-6
  199. C
  200. ADMAX = MAX(Di(1),Di(2),Di(3))
  201. C
  202. C-------------------------------------------------------
  203. C TENSEUR DES CONTRAINTES MOTRICES (dans repere general)
  204. C-------------------------------------------------------
  205. C
  206. EQ=XMAT(10)
  207. DO I=1,3
  208. DO J=1,3
  209. ROTT(I,J)=ROT(I,J)
  210. ENDDO
  211. ENDDO
  212. CALL INVALM (ROTT, 3 , 3 , KER, 1.D-12)
  213. C
  214. IF (endonul) THEN
  215. DO I=1,3
  216. DO J=1,3
  217. COM(I,J) = CO(I,J)
  218. ENDDO
  219. ENDDO
  220. ELSE
  221. C --- Matrice de [(I-D)/det(I-D)]**(q/2) dans le repère principal de D ---
  222. CALL ZERO(DI1,3,3)
  223. DO I=1,3
  224. XNB=(1.D0-Di(I))/(w1*w2*w3)
  225. IF (XNB.GE.0.D0) THEN
  226. DI1(I,I)=(XNB)**(EQ/2.D0)
  227. ELSE
  228. DI1(I,I)=0.D0
  229. WRITE(6,*) 'erreur'
  230. ENDIF
  231. ENDDO
  232. C --- passage dans le repère de général ---
  233. CALL MULMAT(AAA,ROT,DI1,3,3,3)
  234. CALL MULMAT(DI2,AAA,ROTT,3,3,3)
  235. C
  236. CALL MULMAT(AAA,DI2,CO,3,3,3)
  237. CALL MULMAT(COM,AAA,DI2,3,3,3)
  238. ENDIF
  239. C
  240. C----------------------------------------------------------------------
  241. C
  242. C CALCUL DE J1,J2
  243. C J1: CONTRAINTE MOYENNE
  244. C J2: CONTRAINTE DE VON MISES
  245. C
  246. C----------------------------------------------------------------------
  247. C
  248. AJ1 = (SIG(1)+SIG(2)+SIG(3))/(3.D0)
  249. DO I=1,3
  250. XX(I) = SIG(I)-AJ1
  251. ENDDO
  252. DO I=4,NSTRS
  253. XX(I) = SIG(I)
  254. ENDDO
  255. AJ2 = PROCON(XX,XX,NSTRS)
  256. AJ2 = SQRT(1.5D0*AJ2)
  257. C
  258. C----------------------
  259. C FORCE THERMODYNAMIQUE
  260. C----------------------
  261. C
  262. C --- Etat de triaxialité de la contrainte ---
  263. tri = (AJ2)/(SIG(1)+SIG(2)+SIG(3))
  264. C --- Valeurs propres et vecteurs propres du tenseur des contraintes motrices ---
  265. CALL JACOB3(COM,3,VaP,VcP)
  266. ER=XMAT(8)
  267. EA=XMAT(9)
  268. DO I=1,3
  269. AA(I) = ((tri*VaP(I)/(EA))+ABS(tri*VaP(I)/(EA)))/(2.D0)
  270. AA(I) = (AA(I))**(ER)
  271. ENDDO
  272. C -- Dans le repère général ---
  273. CALL ZERO(FTH,3,3)
  274. FTH(1,1) = AA(1)*(VcP(1,1))**2+AA(2)*(VcP(1,2))**2
  275. FTH(1,1) = FTH(1,1)+AA(3)*(VcP(1,3))**2
  276. FTH(2,2) = AA(1)*(VcP(2,1))**2+AA(2)*(VcP(2,2))**2
  277. FTH(2,2) = FTH(2,2)+AA(3)*(VcP(2,3))**2
  278. FTH(3,3) = AA(1)*(VcP(3,1))**2+AA(2)*(VcP(3,2))**2
  279. FTH(3,3) = FTH(3,3)+AA(3)*(VcP(3,3))**2
  280. FTH(1,2) = AA(1)*VcP(1,1)*VcP(2,1)+AA(2)*VcP(1,2)*VcP(2,2)
  281. FTH(1,2) = FTH(1,2)+AA(3)*VcP(1,3)*VcP(2,3)
  282. FTH(1,3) = AA(1)*VcP(1,1)*VcP(3,1)+AA(2)*VcP(1,2)*VcP(3,2)
  283. FTH(1,3) = FTH(1,3)+AA(3)*VcP(1,3)*VcP(3,3)
  284. FTH(2,3) = AA(1)*VcP(2,1)*VcP(3,1)+AA(2)*VcP(2,2)*VcP(3,2)
  285. FTH(2,3) = FTH(1,3)+AA(3)*VcP(2,3)*VcP(3,3)
  286. FTH(3,1) = FTH(1,3)
  287. FTH(2,1) = FTH(1,2)
  288. FTH(3,2) = FTH(2,3)
  289. C
  290. BET2=XMAT(7)
  291. DV1(1,1) = (BET2*FTH(1,1))+FTH(2,2)+FTH(3,3)
  292. DV1(2,2) = (BET2*FTH(2,2))+FTH(1,1)+FTH(3,3)
  293. DV1(3,3) = (BET2*FTH(3,3))+FTH(2,2)+FTH(1,1)
  294. DV1(1,2) = (BET2-1.D0)*FTH(1,2)
  295. DV1(1,3) = (BET2-1.D0)*FTH(1,3)
  296. DV1(2,3) = (BET2-1.D0)*FTH(2,3)
  297. DV1(2,1) = DV1(1,2)
  298. DV1(3,1) = DV1(1,3)
  299. DV1(3,2) = DV1(2,3)
  300. C
  301. C----------------------------------------------------------------------
  302. C
  303. C CALCUL DE LA VARIATION DES CONTRAINTES
  304. C
  305. C----------------------------------------------------------------------
  306. C
  307. CALL DERTRA(NYOUNG,YOUNG,TI,YUNG,YUNGV,TINF,TSUP)
  308. CALL DERTRA(NXNU,XNU,TI,ENU,ENUV,TINF,TSUP)
  309. C
  310. C --- Calcul de l inverse de la matrice de Hooke au debut du pas ---
  311. CALL ELAST4(2,IFOURB,VAR,NVARI,XMAT,NCOMAT,YUNGV,ENUV,
  312. &XCAR,ICARA,MFR,NSTRS,DDINV,DDV,KERRE,INDIC,ITHHER)
  313. IF(KERRE.NE.0) GO TO 100
  314. C
  315. C --- Calcul de la matrice de Hooke et de sa derivee au debut du pas ---
  316. CALL ELAST4(1,IFOURB,VAR,NVARI,XMAT,NCOMAT,YUNGV,ENUV,
  317. &XCAR,ICARA,MFR,NSTRS,DD,DDV,KERRE,INDIC,ITHHER)
  318. IF(KERRE.NE.0) GO TO 100
  319. C
  320. C --- Détermination de la variation de la matrive de Hooke inversée ---
  321. CALL HOOKP(VAR,VARD,NVARI,XMAT,NCOMAT,MFR,NSTRS,DDINVp)
  322. C
  323. C --- Calcul de la variation de contrainte ---
  324. CALL ZDANUL(SIGD,NSTRS)
  325. DO I=1,3
  326. XX(I) = CO(I,I)
  327. ENDDO
  328. XX(4) = CO(1,2)
  329. XX(5) = CO(1,3)
  330. XX(6) = CO(2,3)
  331. C
  332. CALL ZDANUL(EPSED,NSTRS)
  333. DO I=1,NSTRS
  334. EPSED(I)=DSPT(I)-EPSVD(I)-EPSTHD(I)
  335. ENDDO
  336. C
  337. CALL MULMAT(BB,DDINVp,XX,NSTRS,1,NSTRS)
  338. C
  339. DO I=1,NSTRS
  340. XX(I)=EPSED(I)-BB(I)
  341. ENDDO
  342. C
  343. CALL MULMAT(SIGD,DD,XX,NSTRS,1,NSTRS)
  344. C
  345. 100 CONTINUE
  346. RETURN
  347. END
  348.  
  349.  
  350.  
  351.  
  352.  
  353.  
  354.  
  355.  
  356.  
  357.  

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