Télécharger microp.eso

Retour à la liste

Numérotation des lignes :

  1. C MICROP SOURCE CB215821 16/04/21 21:17:45 8920
  2. SUBROUTINE MICROP (WRK0,WRK1,NSTRS,NVARI,NMATT,
  3. & ICARA,KERRE,MFR1,IFOURB)
  4. *
  5. * modele d'endommagement microplan couple a la plasticite
  6. * C. La Borderie + S. Fichant Oct. 95
  7. * routines utilisees:
  8. * micro1: plasticite nadai
  9. * IDECAL=3 DANS LE CAS ISO IDECAL=8 DANS LE CAS ANISO
  10. * jacob3: diagonalisation:
  11. * attention jacob3 modifie la matrice a diagonaliser!!
  12. * prodt et prodt2
  13. * attention prodt2 ne fonctionne qu'avec la matrice des V. P. !!
  14. *
  15. IMPLICIT INTEGER(I-N)
  16. C IMPLICIT REAL*8(A-H,O-Z)
  17. INTEGER NSTRS,NVARI,NMATT,ICARA,KERRE,MFR1,IFOURB
  18. INTEGER ISTRS,I,J
  19. REAL*8 YOUNG,XNU,EPSD0,BT,LAMB,DEUXMU,ALFA
  20. REAL*8 DOM(6),SIGAN(6),TRSIG,DEF33(3,3),EPSIPP(3),VECP(3,3)
  21. REAL*8 VECPT(3,3),DOM33(3,3),DEFRPDE(3,3),DOMRPDE(3,3)
  22. REAL*8 D1,D2,D3,DOM3(3)
  23. REAL*8 SIGPP(3),SIGPM(3),SIG33(3,3),SIG33P(3,3),SIG33M(3,3)
  24. REAL*8 S33PRD(3,3),S33MRD(3,3)
  25. REAL*8 LAMBDAP(6), LAMBDAM(6)
  26. REAL*8 DSIGT(6)
  27. LOGICAL COMP
  28. INTEGER IDECAL
  29. *
  30. -INC CCOPTIO
  31. *
  32. SEGMENT WRK0
  33. REAL*8 XMAT(NMATT)
  34. ENDSEGMENT
  35. *
  36. SEGMENT WRK1
  37. REAL*8 DDHOOK(LHOOK,LHOOK),SIG0(NSTRS),DEPST(NSTRS)
  38. REAL*8 SIGF(NSTRS),VAR0(NVARI),VARF(NVARI)
  39. REAL*8 DEFP(NSTRS),XCAR(ICARA)
  40. ENDSEGMENT
  41. SEGMENT WRKK1
  42. REAL*8 DEFELA(NSTRS)
  43. REAL*8 DDEFEL(NSTRS)
  44. ENDSEGMENT
  45. *
  46. * on recupere les variable materielles
  47. *
  48. YOUNG=XMAT(1)
  49. XNU=XMAT(2)
  50. EPSD0=XMAT(5)
  51. BT=XMAT(6)
  52. ALFA=XMAT(9)
  53. DEUXMU=YOUNG/(1.D0+XNU)
  54. LAMB=XNU*DEUXMU/(1.D0-2.D0*XNU)
  55. * print*,'dans microp'
  56. *
  57. * recuperation des variables internes d'endommagement
  58. *
  59. DO ISTRS=1,NSTRS
  60. DOM(ISTRS)=VAR0(2+ISTRS)
  61. END DO
  62. SEGINI WRKK1
  63. *
  64. * on ecoule plastiquement sur la contrainte effective
  65. *
  66. IDECAL=8
  67. * print*,'appel a micro1'
  68. CALL MICRO1(WRK0,NMATT,WRK1,NSTRS,DSIGT,
  69. 1 NVARI,IDECAL,KERRE,.false.)
  70. IF (KERRE .NE. 0) THEN
  71. print*,'on n''a pas converge dans micro1'
  72. CALL MICRO1(WRK0,NMATT,WRK1,NSTRS,DSIGT,
  73. 1 NVARI,IDECAL,KERRE,.true.)
  74. SEGSUP WRKK1
  75. RETURN
  76. ENDIF
  77. *
  78. * on ecoule en endommagement sur les deformations elastiques
  79. *
  80. * print*,'apres micro1'
  81. *
  82. * calcul de l'increment deformations elastiques DDEFEL
  83. * 1) on calcule l'increment de deformations totales avec
  84. * l'increment de contraintes elastique DSIGT
  85. * 2) on retranche l'increment de deformations plastiques DEFP
  86. *
  87. TRSIG=0.D0
  88. DO ISTRS=1,3
  89. TRSIG=TRSIG+DSIGT(ISTRS)
  90. END DO
  91. DO ISTRS=1,3
  92. DDEFEL(ISTRS)=( (1.D0+XNU)*DSIGT(ISTRS)-XNU*TRSIG)/YOUNG
  93. 1 - DEFP(ISTRS)
  94. END DO
  95. DO ISTRS=4,NSTRS
  96. DDEFEL(ISTRS)= (1.D0+XNU)*DSIGT(ISTRS)/YOUNG
  97. 1 - 0.5d0*DEFP(ISTRS)
  98. END DO
  99. *
  100. * on diagonalise l'increment de deformations elastiques:
  101. * 1) on met sous forme 3x3 avec endoca
  102. * 2) on diagonalise avec jacob3
  103. *
  104. * print*,'incr de def el dans rpg'
  105. * print*,ddefel
  106. * print*,'avant endoca ddefl'
  107. CALL ENDOCA (DDEFEL,DEF33,1)
  108. * print*,'apres endoca ddefl'
  109. CALL JACOB3 (DEF33,IDIM,EPSIPP,VECP)
  110. * print*,'apres jacob3 def33'
  111. * print*,'incr deformations principales'
  112. * print*,epsipp
  113. *
  114. * calcul des deformations elastiques
  115. *
  116. TRSIG=0.D0
  117. * print*,'-----sigf------'
  118. * print*,sigf(1),sigf(2),sigf(3)
  119. DO ISTRS=1,3
  120. TRSIG=TRSIG+SIGF(ISTRS)
  121. END DO
  122. DO ISTRS=1,3
  123. DEFELA(ISTRS)=( (1.D0+XNU)*SIGF(ISTRS)-XNU*TRSIG)/YOUNG
  124. END DO
  125. * print*,'-----defela------'
  126. * print*,defela(1),defela(2),defela(3)
  127. DO ISTRS=4,NSTRS
  128. DEFELA(ISTRS)= (1.D0+XNU)*SIGF(ISTRS)/YOUNG
  129. END DO
  130. * print*,'deformations elastiques dans rpg'
  131. * print*,defela
  132. *
  133. * on met les deformations sous forme de matrice 3x3
  134. * puis on ecrit la matrice dans le repere de depsilon:DEFRPDE
  135. *
  136. CALL ENDOCA (DEFELA,DEF33,1)
  137. * print*,'deformations elastiques dans rpg(3x3)'
  138. * print*,def33
  139. CALL PRODT (DEFRPDE,DEF33,VECP,3,3)
  140. * print*,'def elast dans rpddeps'
  141. * print*,defrpde
  142. *
  143. * on calcule l'endommagement resultant
  144. *
  145. IF (DEFRPDE(1,1) .GT. EPSD0) THEN
  146. D1=1.D0-EPSD0/DEFRPDE(1,1)*EXP(BT*( EPSD0 - DEFRPDE(1,1)))
  147. ELSE
  148. D1=0.D0
  149. END IF
  150. IF (DEFRPDE(2,2) .GT. EPSD0) THEN
  151. D2=1.D0-EPSD0/DEFRPDE(2,2)*EXP(BT*( EPSD0 - DEFRPDE(2,2)))
  152. ELSE
  153. D2=0.D0
  154. END IF
  155. IF (DEFRPDE(3,3) .GT. EPSD0) THEN
  156. D3=1.D0-EPSD0/DEFRPDE(3,3)*EXP(BT*( EPSD0 - DEFRPDE(3,3)))
  157. ELSE
  158. D3=0.D0
  159. END IF
  160. * print*,defrpde(1,1),epsd0,D1
  161. * print*,defrpde(2,2),epsd0,D2
  162. * print*,defrpde(3,3),epsd0,D3
  163. *
  164. * on met l'endommagement initial dans le meme repere
  165. *
  166. CALL ENDOCA(DOM,DOM33,1)
  167. * print*,'prodt DOMRPDE'
  168. CALL PRODT (DOMRPDE,DOM33,VECP,3,3)
  169. * print*,'endommagement initial dans rpddeps'
  170. * print*,domrpde
  171. *
  172. * et on en deduit l'increment d'endommagement dans RPDE
  173. *
  174. IF(d1.gt.domrpde(1,1))then
  175. domrpde(1,1)=d1
  176. endif
  177. IF(d2.gt.domrpde(2,2))then
  178. domrpde(2,2)=d2
  179. endif
  180. IF(d3.gt.domrpde(3,3))then
  181. domrpde(3,3)=d3
  182. endif
  183. * print*,'endommagement final dans rpddeps'
  184. * print*,domrpde
  185. *
  186. * on remet D dans RPG
  187. *
  188. DO I=1,3
  189. DO J=1,3
  190. VECPT(I,J)=VECP(J,I)
  191. END DO
  192. END DO
  193. call prodt(DOM33,DOMRPDE,VECPT,3,3)
  194. * print*,'endommagement final dans rpg'
  195. * print*,dom33
  196.  
  197. *
  198. * on separe les contraintes effectives en + et - dans rpsigma
  199. *
  200. CALL ENDOCA (SIGF,SIG33,1)
  201. CALL JACOB3 (SIG33,3,SIGPP,VECP)
  202. * print*,'contraintes ppales'
  203. * print*,sigpp
  204. DO I=1,3
  205. IF (SIGPP(I) .LT. 0.D0)THEN
  206. SIGPM(I)=SIGPP(I)
  207. SIGPP(I)=0.D0
  208. ELSE
  209. SIGPM(I)=0.D0
  210. END IF
  211. END DO
  212. CALL PRODT2(SIG33P,SIGPP,VECP,3)
  213. CALL PRODT2(SIG33M,SIGPM,VECP,3)
  214. * print*,'contraintes dans rpg'
  215. * print*,sig33p
  216. * print*,sig33m
  217. *
  218. * on met le tout dans le repere ppal d'endo
  219. * attention jacob3 modifie la matrice fournie
  220. * --> on passe une copie
  221. DO I=1,3
  222. DO J=1,3
  223. SIG33(I,J)=DOM33(I,J)
  224. END DO
  225. END DO
  226. CALL JACOB3(SIG33,3,DOM3,VECP)
  227. * print*,'endom dans rpd'
  228. * print*,dom3
  229. **************
  230. * A REVOIR EN ATENDANT MIEUX ON BORNE LES VALEURS PROPRES DE D ENTRE 0. ET 1.*
  231. **************
  232. DO I=1,3
  233. DOM3(I)=DMAX1 (DOM3(I),0.D0)
  234. DOM3(I)=DMIN1 (DOM3(I),1.D0-1.d-6)
  235. END DO
  236. CALL PRODT2 (DOM33,DOM3,VECP,3)
  237. **************
  238. * FIN A REVOIR
  239. **************
  240. * print*,'dom3=',dom3
  241. * print*,'dom33=',dom33
  242. CALL PRODT (S33PRD,SIG33P,VECP,3,3)
  243. CALL PRODT (S33MRD,SIG33M,VECP,3,3)
  244.  
  245. * print*,'s33prd',s33prd(3,3)
  246. * print*,'s33mrd',s33mrd(3,3)
  247. * print*,'callambdap'
  248. COMP=.FALSE.
  249. CALL CLMBDA(DOM3(1),DOM3(2),DOM3(3),LAMBDAP,YOUNG,XNU,ALFA,
  250. 1 COMP)
  251. COMP=.TRUE.
  252. CALL CLMBDA(DOM3(1),DOM3(2),DOM3(3),LAMBDAM,YOUNG,XNU,ALFA,
  253. 1 COMP)
  254. * print*,'apres callambdap'
  255. SIG33P(1,1)=LAMBDAP(1)*S33PRD(1,1)+LAMBDAM(1)*S33MRD(1,1)
  256. SIG33P(1,2)=LAMBDAP(4)*S33PRD(1,2)+LAMBDAM(4)*S33MRD(1,2)
  257. SIG33P(1,3)=LAMBDAP(5)*S33PRD(1,3)+LAMBDAM(5)*S33MRD(1,3)
  258. SIG33P(2,1)=SIG33P(1,2)
  259. SIG33P(2,2)=LAMBDAP(2)*S33PRD(2,2)+LAMBDAM(2)*S33MRD(2,2)
  260. SIG33P(2,3)=LAMBDAP(6)*S33PRD(2,3)+LAMBDAM(6)*S33MRD(2,3)
  261. SIG33P(3,1)=SIG33P(1,3)
  262. SIG33P(3,2)=SIG33P(2,3)
  263. SIG33P(3,3)=LAMBDAP(3)*S33PRD(3,3)+LAMBDAM(3)*S33MRD(3,3)
  264. * print*,'contraintes dans RPD'
  265. * print*,sig33p
  266. * print*,'sig33p(1,2)',sig33p(1,2)
  267. * print*,'dsigt',dsigt(1),dsigt(2),dsigt(3)
  268. * print*,'dsigt',dsigt(4),dsigt(5),dsigt(6)
  269. * print*,'sigf',sigf(1),sigf(2),sigf(3)
  270. * print*,'sigf',sigf(4),sigf(5),sigf(6)
  271. *
  272. * on remet le tout dans le repere global
  273. *
  274. * print*,'matrice des vect proprs'
  275. * print*,vecp
  276. DO I=1,3
  277. DO J=1,3
  278. VECPT(I,J)=VECP(J,I)
  279. END DO
  280. END DO
  281. CALL PRODT (SIG33,SIG33P,VECPT,3,3)
  282. * print*,'contraintes dans rpg'
  283. * print*,sig33
  284. *
  285. * on rend les contraintes et les variables internes finales
  286. *
  287. SIGAN(1)=SIGF(1)-SIG33(1,1)
  288. SIGF(1)=SIG33(1,1)
  289. VARF(3)=DMAX1 (DOM33(1,1),0.d0)
  290. SIGAN(2)=SIGF(2)-SIG33(2,2)
  291. SIGF(2)=SIG33(2,2)
  292. VARF(4)=DMAX1 (DOM33(2,2),0.d0)
  293. SIGAN(3)=SIGF(3)-SIG33(3,3)
  294. SIGF(3)=SIG33(3,3)
  295. VARF(5)=DMAX1 (DOM33(3,3),0.d0)
  296. SIGAN(4)=SIGF(4)-SIG33(1,2)
  297. SIGF(4)=SIG33(1,2)
  298. VARF(6)=DMAX1 (DOM33(1,2),0.d0)
  299. VARF(7)=DMAX1 (DOM33(1,3),0.d0)
  300. VARF(8)=DMAX1 (DOM33(2,3),0.d0)
  301. IF(IFOUR.GE.1.OR.IFOUR.LE.-3) THEN
  302. SIGAN(5)=SIGF(5)-SIG33(1,3)
  303. SIGF(5)=SIG33(1,3)
  304. SIGAN(6)=SIGF(6)-SIG33(2,3)
  305. SIGF(6)=SIG33(2,3)
  306. ELSE
  307. SIGAN(5)=0.D0
  308. SIGAN(6)=0.D0
  309. END IF
  310. DO ISTRS=1,6
  311. VARF(ISTRS+8)=SIGAN(ISTRS)
  312. END DO
  313. * print*,'sigf',sigf(1),sigf(2),sigf(3)
  314. * print*,'sigf',sigf(4),sigf(5),sigf(6)
  315. * print*,'sigan',sigan(1),sigan(2),sigan(3)
  316. * print*,'sigan',sigan(4),sigan(5),sigan(6)
  317. SEGSUP WRKK1
  318. RETURN
  319. END
  320.  
  321.  
  322.  
  323.  
  324.  

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