Télécharger ricj3.eso

Retour à la liste

Numérotation des lignes :

  1. C RICJ3 SOURCE CB215821 16/04/21 21:18:17 8920
  2. C sub ricjoi3d
  3. C RICJOI3D SOURCE CHAT 05/01/12 22:24:09 5004
  4. SUBROUTINE RICJ3(IB,IGAU,NSTRS,SIG0,EPIN0,VAR0,NVARI,
  5. & DEPST,IFOUR,TETA1,TETA2,XMAT,
  6. & NMATT,IVAL,DD,SIGF,DEFP,VARF,KERRE)
  7. C-----------------------------------------------------------------------
  8. C Modele RICJOI en version element joint 3D
  9. C-----------------------------------------------------------------------
  10. C Développé par : Benjamin Richard
  11. C
  12. C Supervisé par : Frederic Ragueneau
  13. C Lucas Hector Adelaide
  14. C Christian Cremona
  15. C Jean Louis Tailhan
  16. C-----------------------------------------------------------------------
  17.  
  18. C
  19. C-----DECLARATION GENERALE----------------------------------------------
  20. C
  21. IMPLICIT INTEGER(I-N)
  22. IMPLICIT REAL*8(A-H,O-Z)
  23. C
  24. C-----APPEL AUX INCLUDES------------------------------------------------
  25. C
  26. -INC CCREEL
  27. C
  28. C-----PARAMETRES
  29. C
  30. PARAMETER(XDILAT = 80.0D0)
  31. C
  32. C-----DIMENSION---------------------------------------------------------
  33. C
  34. DIMENSION SIG0(*),EPIN0(*),VAR0(*),DEPST(*),XMAT(*)
  35. DIMENSION DD(NSTRS,NSTRS)
  36. DIMENSION SIGF(*),DEFP(*),VARF(*)
  37.  
  38. FVARLIM = 2.0D-5
  39. C
  40. C-----COMPTEUR POUR LE NOMBRE DE PASSAGE DAS LA BOUCLE
  41. C
  42. COUN = VAR0(17)
  43. C
  44. C-----ETAT PRECEDENT----------------------------------------------------
  45. C
  46. D = VAR0(1 )
  47. EPSF1 = VAR0(2 )
  48. EPSF2 = VAR0(3 )
  49. XEC1 = VAR0(4 )
  50. XEC2 = VAR0(5 )
  51. EPST1 = VAR0(6 )
  52. EPST2 = VAR0(7 )
  53. EPSN = VAR0(8 )
  54. CIN1 = VAR0(9 )
  55. CIN2 = VAR0(10)
  56. EPSR = VAR0(11)
  57. IFLAG = nint(VAR0(12))
  58. FVAR = VAR0(13)
  59. SIGM = VAR0(14)
  60. XKS = VAR0(15)
  61. XKN = VAR0(16)
  62. crit2tol=1.001d0
  63. C
  64. C-----INITIALISATION SPECIFIQUE POUR LES VARIABLES INTERNES DIFF. DE ZERO
  65. C
  66. IF (COUN.EQ.0.0D0) THEN
  67. FVAR = 0.2D0
  68. SIGM = XMAT(13)
  69. XKS = XMAT(1)
  70. XKN = XMAT(2)
  71. ENDIF
  72. C
  73. C-----CONSTANTES MATERIAUX----------------------------------------------
  74. C
  75. XKS0 = XMAT(1)
  76. XKN0 = XMAT(2)
  77.  
  78. XAD = XMAT(5)
  79. XY0 = XMAT(6)
  80. XAL = XMAT(7)
  81. XGA = XMAT(8)
  82. XAA = XMAT(9)
  83.  
  84. Q1COE = XMAT(10)
  85. Q2COE = XMAT(11)
  86. Q3COE = XMAT(12)
  87. SYCOE = XMAT(13)
  88. NCOEF = nint(XMAT(14) )
  89. KCOEF = nint(XMAT(15))
  90.  
  91. TCCOR = XMAT(16)
  92. NGONF = nint(XMAT(17))
  93.  
  94. CALL APP3(TCCOR,YR,15)
  95.  
  96.  
  97. C
  98. C-----ACTUALISATION DE LA DEFORMATION-----------------------------------
  99. C
  100. EPST1 = EPST1 + DEPST(1)
  101. EPST2 = EPST2 + DEPST(2)
  102. EPSN = EPSN + DEPST(3)
  103.  
  104. IF (NGONF.EQ.1) THEN
  105. GOTO 123
  106. ENDIF
  107.  
  108. XSTE = XKN*EPSN
  109. XSAD = XDILAT*0.5D0*(-ABS(XSTE)+XSTE)
  110.  
  111. C
  112. C-----CALCUL DE L ENERGIE-----------------------------------------------
  113. C
  114. EPSNP = 0.5D0*(ABS(EPSN)+EPSN)
  115. YN = 0.5D0*XKN*EPSNP*EPSNP
  116. YD1 = 0.5D0*XKS*EPST1*EPST1
  117. YD2 = 0.5D0*XKS*EPST2*EPST2
  118.  
  119. C
  120. C-----CALCUL DU SEUIL D ENDOMMAGEMENT 1---------------------------------
  121. C
  122. XFD1 = XAL*YN+YD1-(XY0+XEC1+YR)
  123.  
  124. C
  125. C-----TEST SUR LE SEUIL D ENDOMMAGEMENT---------------------------------
  126. C
  127. IF (XFD1.GT.0.0D0) THEN
  128. Dn = 1.0D0-1.0D0/(1+XAD*(XAL*YN+YD1-XY0-YR))
  129. XEC1 = XAL*YN+YD1-XY0-YR
  130. IF (Dn.GT.D) THEN
  131. D = Dn
  132. ENDIF
  133. ENDIF
  134. IF (D.GE.1.0D0) THEN
  135. D = 0.999D0
  136. ENDIF
  137.  
  138. C
  139. C-----CALCUL DU SEUIL D ENDOMMAGEMENT 2---------------------------------
  140. C
  141. XFD2 = XAL*YN+YD2-(XY0+XEC2+YR)
  142.  
  143. C
  144. C-----TEST SUR LE SEUIL D ENDOMMAGEMENT---------------------------------
  145. C
  146. IF (XFD2.GT.0.0D0) THEN
  147. Dn = 1.0D0-1.0D0/(1+XAD*(XAL*YN+YD2-XY0-YR))
  148. XEC2 = XAL*YN+YD2-XY0-YR
  149. IF (Dn.GT.D) THEN
  150. D = Dn
  151. ENDIF
  152. ENDIF
  153. IF (D.GE.1.0D0) THEN
  154. D = 0.999D0
  155. ENDIF
  156.  
  157. C
  158. C-----FROTTEMENT 1-------------------------------------------------------
  159. C
  160. SPI1 = D*XKS*(EPST1-EPSF1)
  161.  
  162. CRITF1= ABS(SPI1-CIN1)+XSAD
  163.  
  164. IF (CRITF1.GT.(0.0D0)) THEN
  165. XF0 = CRITF1
  166. XDLAM = 1.0D0
  167.  
  168. DO I=1,500
  169.  
  170. SEUIL = (ABS(SPI1-CIN1)+XSAD)/XF0
  171.  
  172. IF ((SEUIL.LE.1.0D-5).OR.(XDLAM.LE.1.0D-10)) THEN
  173. GOTO 900
  174. ELSE
  175. XSG = ABS(SPI1-CIN1)/(SPI1-CIN1)
  176. XDLAM = (ABS(SPI1-CIN1)+XSAD)/
  177. & (XKS*D-XSG*XGA*(-XSG+XAA*CIN1))
  178. CIN1 = CIN1-XGA*XDLAM*(-XSG+XAA*CIN1)
  179. SPI1 = SPI1-D*XKS*XDLAM*XSG
  180. ENDIF
  181.  
  182. EPSF1 = -SPI1/(D*XKS)+EPST1
  183.  
  184. ENDDO
  185. ENDIF
  186.  
  187. 900 CONTINUE
  188.  
  189. C
  190. C-----FROTTEMENT 2-------------------------------------------------------
  191. C
  192. SPI2 = D*XKS*(EPST2-EPSF2)
  193.  
  194. CRITF2= ABS(SPI2-CIN2)+XSAD
  195.  
  196. IF (CRITF2.GT.(0.0D0)) THEN
  197. XF0 = CRITF2
  198. XDLAM = 1.0D0
  199.  
  200. DO I=1,500
  201.  
  202. SEUIL = (ABS(SPI2-CIN2)+XSAD)/XF0
  203.  
  204. IF ((SEUIL.LE.1.0D-5).OR.(XDLAM.LE.1.0D-10)) THEN
  205. GOTO 901
  206. ELSE
  207. XSG = ABS(SPI2-CIN2)/(SPI2-CIN2)
  208. XDLAM = (ABS(SPI2-CIN2)+XSAD)/
  209. & (XKS*D-XSG*XGA*(-XSG+XAA*CIN2))
  210. CIN2 = CIN2-XGA*XDLAM*(-XSG+XAA*CIN2)
  211. SPI2 = SPI2-D*XKS*XDLAM*XSG
  212. ENDIF
  213.  
  214. EPSF2 = -SPI2/(D*XKS)+EPST2
  215.  
  216. ENDDO
  217. ENDIF
  218.  
  219. 901 CONTINUE
  220.  
  221. C
  222. C-----ANELASTICITE--------------------------------------------------------
  223. C
  224. IF (EPSN.GT.0.0D0) THEN
  225. SIGN = (1-D)*XKN*EPSN
  226. ELSE
  227. SIGN = XKN*(EPSN-EPSR)
  228. ENDIF
  229.  
  230. C
  231. C-----CONTRAINTE DUE AUX EFFETS INELASTIQUE DE LA ROUILLE-----------------
  232. C
  233. sigr=sign
  234. IF (EPSN.LE.0.0D0) THEN
  235. SIGR = -1.0D0*SIGN
  236. ENDIF
  237.  
  238. C
  239. C-----CRITERE D ANELASTICITE DE LA ROUILLE--------------------------------
  240. C
  241. IF (IFLAG.EQ.1) THEN
  242. GOTO 20
  243. ELSE
  244. CRIT2 = 2.0D0*Q1COE*FVAR*DCOSH(Q2COE*SIGR/2.0D0/SYCOE)-
  245. & (1.0D0+(Q3COE*FVAR)**2.0D0)
  246. ENDIF
  247.  
  248. * print*,'crit2=',crit2
  249.  
  250. C
  251. C SI LE CRITERE EST DEPASSE ON DEMARRE LE RETURN MAPPING
  252. C
  253. IF ((CRIT2.GT.0.0D0).AND.
  254. 1 (EPSN .LE.0.0D0).AND.
  255. 2 (TCCOR.GT.0.0D0).AND.
  256. 3 (FVAR.GT.FVARLIM).AND.
  257. 4 (IFLAG.EQ.0).AND.
  258. 5 ((EPSN*DEPST(3)).GT.0.0D0)) THEN
  259. CRIT20 = CRIT2
  260. DLAM2 = 1.0D0
  261. GOTO 10
  262. ELSE IF (FVAR.EQ.FVARLIM) THEN
  263. IFLAG = 1
  264. GOTO 20
  265. ELSE
  266. GOTO 30
  267. END IF
  268. C
  269. C-----DEBUT DES ITERATIONS INTERNES
  270. C
  271. 10 DO I=1,10000
  272. C
  273. C EVALUATION DU CRITERE
  274. C
  275. CRIT2 = 2.0D0*Q1COE*FVAR*DCOSH(Q2COE*SIGR/2.0D0/SYCOE)-
  276. & (1.0D0+(Q3COE*FVAR)**2.0D0)
  277. C
  278. C TEST DE CONVERGENCE SUR LA VALEUR DU CRITERE RELATIVE
  279. C
  280. IF ((ABS(CRIT2/CRIT20).LE.CRIT2TOL).OR.
  281. & (DLAM2.LT.0.0D0).OR.
  282. & (CRIT2.LT.0.0D0)) THEN
  283. GOTO 20
  284. END IF
  285. C
  286. C CALCUL DES DERIVEES DU CRITERES
  287. C
  288. C dFr/dSIGR : DERIV21
  289. C dFr/dSIGM : DERIV22
  290. C dFr/dFVAR : DERIV23
  291. C
  292. DERIV21 = DSINH(Q2COE*SIGR/2.0D0/SIGM)*
  293. & Q1COE*Q2COE*FVAR/SIGM
  294.  
  295. DERIV22 = DSINH(Q2COE*SIGR/2.0D0/SIGM)*
  296. & (-1.0D0)*Q1COE*Q2COE*FVAR*SIGR/SIGM**2.0D0
  297.  
  298. DERIV23 = 2.0D0*Q1COE*DCOSH(Q2COE*SIGR/2.0D0/SIGM)-
  299. & 2.0D0*(Q3COE**2.0D0)*FVAR
  300. C
  301. C CALUL DU MODULE TANGENT COHERENTS ENTRE DSIGM ET DEPSM
  302. C
  303. ETAN = XKN/NCOEF*(SYCOE/SIGM)**(NCOEF-1.0D0)
  304. ESTA = 1.0D0/ETAN - 1.0D0/XKN
  305. ESTA = 1.0D0/ESTA
  306. C
  307. C CALCUL DU MULTIPLICATEUR DE LAGRANGE
  308. C
  309. DLAM2 = -CRIT2/(XKN*(DERIV21)**2.0D0-
  310. & DERIV23*KCOEF*FVAR*(1.0D0-FVAR)*DERIV21-
  311. & ESTA/((1.0D0-FVAR)*SIGM)*DERIV21*SIGR*DERIV22)
  312. C
  313. C ACTUALISATION DES VARIABLES FORCES
  314. C
  315. SIGR = SIGR + DLAM2 * XKN * DERIV21
  316. FVAR = FVAR - DLAM2 * KCOEF * FVAR * (1.0D0-FVAR) * DERIV21
  317. SIGM = SIGM - DLAM2 * SIGR * DERIV21 *
  318. & ESTA/(1.0D0-FVAR)/SIGM
  319. C
  320. C SECURITE SUR LA VALEUR DE FVAR (CAS OU ELLE DEVIENT TROP FAIBLE)
  321. C
  322. IF (FVAR.LE.FVARLIM) THEN
  323. FVAR = FVARLIM
  324. MODR = nint(DLAM2 * XKN * DERIV21/DEPST(2) )
  325. END IF
  326. C
  327. C --------------------------------------------
  328. C RETURN MAPPING : FIN DES ITERATIONS INTERNES
  329. C --------------------------------------------
  330. C
  331. END DO
  332. C
  333. C BALISE DE SORTIE DU RETURN MAPPING AVANT ACUALISATION
  334. C
  335. 20 CONTINUE
  336. C
  337. C CAS OU LA FRACTION VOLUMIQUE DES PORES EST EGALE A 0.02
  338. C
  339. IF ((FVAR.EQ.FVARLIM).AND.
  340. & (IFLAG.EQ.1)) THEN
  341. SIGR = SIGR - MODR * DEPST(3)
  342. END IF
  343. C
  344. C CALCUL DES QUANTITES RESULTATS (EPSR, VKACT ET XKN)
  345. C
  346. SIGN = -SIGR
  347. EPSR = -(EPSN - SIGR/XKN)
  348. C
  349. C BALISE DE SORTIE DU RETURN MAPPING APRES ACUALISATION
  350. C
  351. 30 CONTINUE
  352. C
  353. C ACTUALISATION DES PARAMETRES ELASTIQUES
  354. C
  355. XKNACT = 4.0D0*XKN0*XKS0*(1.0D0-FVAR)/
  356. & (4.0D0*XKS0+3.0D0*XKN0*FVAR)
  357. XKSACT = XKS0*(1.0D0-FVAR)/
  358. & (1.0D0+((6.0D0*XKN0+12.0D0*XKS0)/
  359. & (9.0D0*XKN0+8.0D0*XKS0))*FVAR)
  360.  
  361. XKN = XKNACT
  362. XKS = XKSACT
  363.  
  364. C
  365. C-----CALCUL DES CONTRAINTES FINALES------------------------------------
  366. C
  367. SIGF(1) = XKS*(1-D)*EPST1 + SPI1
  368. SIGF(2) = XKS*(1-D)*EPST2 + SPI2
  369.  
  370. IF (EPSN.GT.0.0D0) THEN
  371. SIGF(3) = XKN*(1-D)*EPSN
  372. ELSE
  373. SIGF(3) = XKN*(EPSN-EPSR)
  374. ENDIF
  375. C
  376. C-----STOCKAGE EN SORTIE DES VARIABLES INTERNES-------------------------
  377. C
  378. VARF(1 ) = D
  379. VARF(2 ) = EPSF1
  380. VARF(3 ) = EPSF2
  381. VARF(4 ) = XEC1
  382. VARF(5 ) = XEC2
  383. VARF(9 ) = CIN1
  384. VARF(10) = CIN2
  385. VARF(11) = EPSR
  386. VARF(12) = IFLAG
  387. VARF(13) = FVAR
  388. VARF(14) = SIGM
  389. VARF(15) = XKS
  390. VARF(16) = XKN
  391. VARF(17) = COUN + 1
  392.  
  393. 123 CONTINUE
  394.  
  395. IF (NGONF.EQ.1) THEN
  396. SIGF(1) = XMAT(1)*EPST1
  397. SIGF(2) = XMAT(1)*EPST2
  398. SIGF(3) = XMAT(2)*EPSN
  399. ENDIF
  400.  
  401. VARF(6 ) = EPST1
  402. VARF(7 ) = EPST2
  403. VARF(8 ) = EPSN
  404.  
  405. RETURN
  406. END
  407.  
  408.  
  409.  
  410.  
  411.  
  412.  
  413.  
  414.  
  415.  
  416.  
  417.  
  418.  
  419.  

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