Télécharger deriv.eso

Retour à la liste

Numérotation des lignes :

deriv
  1. C DERIV SOURCE CB215821 18/09/21 21:15:34 9930
  2. SUBROUTINE DERIV(TAU,SIG,EPSV,VAR,SIGD,EPSVD,VARD,EPSTHD,
  3. &DSPT,XMAT,FKX,NPTEUR,NSTRS,NVARI,NCOMAT,NKX,NC,DD,DDV,DDINV,
  4. &YOUNG,NYOUNG,XNU,NXNU,MFR,XCAR,ICARA,IFOURB,INDIC,
  5. &TI,TPOINT,ITHHER,TRUC,NCOURB)
  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
  10. C-----------------------------------------------------------------------
  11. C
  12. C-----------------------------------------------------------------------
  13. C Entree: TAU pas de temps
  14. C TI temperature au debut du pas
  15. C TPOINT derivee de la temperature
  16. C SIG(NSTRS) tenseur des contraintes
  17. C EPSV(NSTRS) tenseur des deformations
  18. C VAR(NVARI) tableau contenant les variables internes
  19. C pilotant les equations
  20. C VAR(1)=p ; VAR(2)=r ;VAR(3)=D
  21. C EPSTHD vitesse de dilatation thermique au debut du pas
  22. C DSPT(NSTRS,NSTRS) vitesse de deformation totale
  23. C XMAT(NCOMAT) tableau des parametres scalaires du materiau
  24. C a une temperature TI donnee
  25. C XMAT(1)=YOUNG ;XMAT(2)=XNU ;XMAT(3)=n
  26. C XMAT(4)=m ;XMAT(5)=KK; XMAT(6)=ALPHA
  27. C XMAT(7)=BETA ;XMAT(8)=r; XMAT(9)=A ;
  28. C XMAT(10)=EVOL ;XMAT(11)=ALPHAT
  29. C XMAT(12)=RHO; XMAT(13)=SIGY
  30. C FKX(NKX) tableau de la fonction k(X,T) contenant les courbes
  31. C de k a differentes temperatures T1,T2,T3,...
  32. C NPTEUR(NC) tableau de pointeurs sur les differentes courbes
  33. C de la fonction k
  34. C MFR indice de la formulation mecanique(seulement massif
  35. C ou coque pour les materiaux endommageables
  36. C ICARA nombre de caracteristiques geometriques des elements
  37. C finis
  38. C XCAR(ICARA) tableau des caracteristiques geometriques des
  39. C elements finis
  40. C IFOURB= -2 EN CONTR. PLANES
  41. C -1 EN DEFORM. PLANES
  42. C 0 EN AXISYMETRIE
  43. C 1 EN SERIE DE FOURIER
  44. C 2 EN TRIDIM
  45. C INDIC=0, 1 OU -1 pour plasticite avec endommagement
  46. C =2 OU -2 pour viscoplasticite avec endommagement
  47. C ITHHER = 0 pas de chargement thermique et materiau constant
  48. C = 1 chargement thermique et materiau constant
  49. C = 2 chargement thermique et materiau(T)
  50. C------------------------------------------------------------------------
  51. C
  52. C------------------------------------------------------------------------
  53. C Sortie: EPSVD(NSTRS) derivee du tenseur des deformations
  54. C VARD(NVARI) tableau contenant les derivees des variables
  55. C internes
  56. C SIGD(NSTRS) derivee des contraintes
  57. C DD(NSTRS,NSTRS) matrice de Hooke au debut du pas
  58. C DDV(NSTRS,NSTRS) derivee de DD
  59. C DDINV(NSTRS,NSTRS) inverse de DD
  60. C------------------------------------------------------------------------
  61. C
  62. IMPLICIT INTEGER(I-N)
  63. IMPLICIT REAL*8(A-H,O-Z)
  64. DIMENSION SIG(*),EPSV(*),VAR(*),XCAR(*),YOUNG(*),XNU(*)
  65. DIMENSION SIGD(*),EPSVD(*),VARD(*),DDINV(NSTRS,NSTRS)
  66. DIMENSION FKX(*),NPTEUR(*),XMAT(*),DSPT(*),TRUC(*)
  67. DIMENSION EPSTHD(*),DD(NSTRS,NSTRS),DDV(NSTRS,NSTRS)
  68. DIMENSION XX(6),A(3,3),D(3),S(3,3)
  69. C
  70. C----------------------------------------------------------------------
  71. C CALCULATE J0,J1,J2
  72. C J0: THE BIGEST PRINCIPAL STRESS
  73. C J1: MEAN STRESS
  74. C J2: CONTR. DE VON MISES
  75. C----------------------------------------------------------------------
  76. * cas des massif et coques
  77. * write(6,*) ' entree dans deriv mfr ',mfr
  78. if(mfr.eq.1.or.mfr.eq.3) then
  79. AJ1 = SIG(1)+SIG(2)+SIG(3)
  80. AJ1 = AJ1/3.0D0
  81. DO 30 I=1,3
  82. 30 XX(I) = SIG(I)-AJ1
  83. DO 40 I=4,NSTRS
  84. 40 XX(I) = SIG(I)
  85. AJ2 = PROCON(XX,XX,NSTRS)
  86. AJ2 = SQRT(1.5D0*AJ2)
  87. elseif (mfr.eq.5) then
  88. * cas des coques epaisses <= a priori non utilise dans ce modele !
  89. alfah = 1.D0
  90. ep1 = xcar(1)
  91. AJ1= (SIG(1)+SIG(2)) / 3.D0
  92. AJ2 = vonmis0(sig,nstrs,mfr,ifourb,ep1,alfah)
  93. XX(1)= sig(1) - aj1
  94. xx(2) = sig(2)-aj1
  95. DO 41 I=3,NSTRS
  96. 41 XX(I) = SIG(I)
  97. endif
  98. CALL ZERO(A,3,3)
  99. if(mfr.eq.1.or.mfr.eq.3) then
  100. DO 50 I=1,3
  101. 50 A(I,I) = SIG(I)
  102. A(1,2) = SIG(4)
  103. IF(NSTRS.NE.6) GOTO 54
  104. A(1,3) = SIG(5)
  105. A(2,3) = SIG(6)
  106. 54 CONTINUE
  107. elseif(mfr.eq.5) then
  108. a(1,1)=sig(1)
  109. a(2,2)=sig(2)
  110. a(1,2)= sig(3)
  111. a(1,3)= sig(4)
  112. a(2,3)= sig(5)
  113. endif
  114. A(2,1) = A(1,2)
  115. A(3,1) = A(1,3)
  116. A(3,2) = A(2,3)
  117. IF (NSTRS.EQ.4) CALL JACOB3(A,2,D,S)
  118. * on appele jacod3 car on n'a pas besoin des vecteurs propres
  119. * write(6,*) ' entree de jacob3'
  120. IF (NSTRS.GE.5) CALL JACOD3(A,3,D)
  121. * write(6,*) ' sortie de jacob3'
  122. C
  123. C---------------------------------------------------------------------
  124. C CALCULATE PRINCIPLE STRESS : D(*)
  125. C---------------------------------------------------------------------
  126. AJ0 = MAX(D(1),D(2),D(3))
  127. C
  128. C---------------------------------------------------------------------
  129. C CALCULATE VARIABLE X
  130. C---------------------------------------------------------------------
  131. C
  132. ALFA = XMAT(6)
  133. BETA = XMAT(7)
  134. X = (ALFA*AJ0)+(BETA*AJ1)+((1.0D0-ALFA-BETA)*AJ2)
  135. C
  136. C---------------------------------------------------------------------
  137. C CALCULATE PARAMETER k(X,T)=AK
  138. C---------------------------------------------------------------------
  139. C
  140. CALL VALPAR(FKX,NPTEUR,NKX,NC,X,AK,AKV,TI,TINF,TSUP,TRUC,NCOURB)
  141.  
  142. C
  143. C---------------------------------------------------------------------
  144. C CALCULATE RATE OF EQUAL EFFECTIVE STRAIN VARD(1)
  145. C RATE OF REAL EQUAL EFFECTIVE STRAIN VARD(2)
  146. C WITH VARD(2)=(1-D)*VARD(1)
  147. C---------------------------------------------------------------------
  148. DEN=1.D0-VAR(3)
  149. EM=XMAT(4)
  150. EN=XMAT(3)
  151. EK=XMAT(5)
  152. SIGY=XMAT(13)
  153. IF ( VAR(2).EQ.0.D0 ) VAR(2)=VAR(2)+1.D-12
  154. VARR=AJ2-(SIGY*DEN)
  155. D1=DEN*EK*(VAR(2)**(1.D0/EM))
  156. VARR=VARR/D1
  157. IF ( VARR.GT.0.D0 ) THEN
  158. VARD(2)=VARR**EN
  159. ELSE
  160. VARD(2)=0.D0
  161. ENDIF
  162.  
  163. VARD(1)=VARD(2)/DEN
  164. C
  165. C---------------------------------------------------------------------
  166. C CALCULATE RATE OF DAMAGE
  167. C---------------------------------------------------------------------
  168. ER=XMAT(8)
  169. EA=XMAT(9)
  170. IF (X.LE.0.D0) THEN
  171. VARD(3)=0.D0
  172. ELSE
  173. VARD(3)=((X/EA)**ER)/(DEN**AK)
  174. ENDIF
  175. C
  176. C---------------------------------------------------------------------
  177. C CALCULATE RATE OF PLASTIC STRAIN
  178. C---------------------------------------------------------------------
  179. IF ( AJ2.EQ.0.D0 ) THEN
  180. AA3=0.D0
  181. ELSE
  182. AA3=1.5D0*VARD(1)/AJ2
  183. ENDIF
  184. if(mfr.eq.1.or.mfr.eq.3) then
  185. DO 10 I=1,NSTRS
  186. AGAM=1.D0
  187. IF (I.GE.4) AGAM=2.D0
  188. EPSVD(I)=AA3*XX(I)*AGAM
  189. 10 CONTINUE
  190. elseif(mfr.eq.5) then
  191. DO 11 I=1,NSTRS
  192. AGAM=1.D0
  193. IF (I.GE.3) AGAM=2.D0
  194. EPSVD(I)=AA3*XX(I)*AGAM
  195. 11 CONTINUE
  196. endif
  197. C
  198. C----------------------------------------------------------------------
  199. C CALCULATE RATE OF STRESS
  200. C----------------------------------------------------------------------
  201.  
  202. CALL DERTRA(NYOUNG,YOUNG,TI,YUNG,YUNGV,TINF,TSUP)
  203. CALL DERTRA(NXNU,XNU,TI,ENU,ENUV,TINF,TSUP)
  204. C
  205. C Calcul de l'inverse de la matrice de Hooke au debut du pas
  206. CALL ELAST1(2,IFOURB,VAR,NVARI,XMAT,NCOMAT,YUNGV,ENUV,
  207. &XCAR,ICARA,MFR,NSTRS,DDINV,DDV,KERRE,INDIC,ITHHER)
  208. C
  209. C Calcul de la matrice de Hooke et de sa derivee au debut du pas
  210. CALL ELAST1(1,IFOURB,VAR,NVARI,XMAT,NCOMAT,YUNGV,ENUV,
  211. &XCAR,ICARA,MFR,NSTRS,DD,DDV,KERRE,INDIC,ITHHER)
  212. C
  213. CALL ZDANUL(SIGD,NSTRS)
  214. DO I=1,NSTRS
  215. XFLO=( DSPT(I)-EPSVD(I)-EPSTHD(I) )
  216. DO J=1,NSTRS
  217. SIGD(J)=SIGD(J)+ DD(J,I)*XFLO
  218. ENDDO
  219. ENDDO
  220.  
  221. CALL MULMAT(DD,DDV,DDINV,NSTRS,NSTRS,NSTRS)
  222.  
  223. XFLO1 = VARD(3)/DEN
  224. DO I=1,NSTRS
  225. XFLO=SIG(I)
  226. DO J=1,NSTRS
  227. SIGD(J)=SIGD(J) + TPOINT*DD(J,I)*XFLO
  228. ENDDO
  229. SIGD(I)=SIGD(I) - XFLO1*SIG(I)
  230. ENDDO
  231.  
  232. RETURN
  233. END
  234.  
  235.  
  236.  

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