Télécharger rupthe.eso

Retour à la liste

Numérotation des lignes :

rupthe
  1. C RUPTHE SOURCE PV 17/12/08 21:17:45 9660
  2. C--------------------------------------------------------------------
  3. subroutine RUPTHE(WRK52,WRK53,WRK54,WRK27,IB,IGAU,NBPGAU)
  4.  
  5. C ***************************************************
  6. C Elément macro de rupteur (derniÚre version)
  7. C crée par Huyen Nguyen - décembre 2010
  8. C *************************************
  9. IMPLICIT INTEGER(I-N)
  10. IMPLICIT REAL*8(A-H,O-Z)
  11.  
  12. -INC PPARAM
  13. -INC CCOPTIO
  14. -INC DECHE
  15. C-----------------------------------------------------------------------
  16. C GLOBAL MODEL
  17. DIMENSION XMKBel(6,6)
  18. DIMENSION XMKAel(6,6)
  19. * real*8 XKBel
  20. * real*8 XKAel
  21.  
  22. C Endommagement du béton
  23. * real*8 Zplus, Zmoin
  24. * real*8 Dplus, Dmoin, Dh
  25. * real*8 Ym, Ymp, Fonc1
  26.  
  27. C Frottement du béton
  28. logical conv
  29. * integer iter
  30.  
  31. * real*8 up, Ecrb, Fp
  32. * real*8 Fonc2, Xmul1, Fonc3, Xlamda1
  33. * real*8 dfdfb, dpidfb, dfdxb, dpidxb
  34.  
  35. C Frottement du béton
  36.  
  37. * real*8 ua, Ecra, Fa
  38. * real*8 Fonc4, Xmul2, Fonc5, Xlamda2
  39. * real*8 dfdfa, dpidfa, dfdxa, dpidxa
  40.  
  41. C
  42. DIMENSION uo(6)
  43. DIMENSION duo(6)
  44. DIMENSION force(6)
  45. DIMENSION da(10)
  46.  
  47.  
  48. C----------------------------------------------------------------------
  49. C Affectation des données du modÚle
  50.  
  51.  
  52. da(1) =XMAT(13)
  53. da(2) =XMAT(14)
  54. da(3) =XMAT(15)
  55. da(4) =XMAT(16)
  56. da(5) =XMAT(17)
  57. da(6) =XMAT(18)
  58. da(7) =XMAT(19)
  59. da(8) =XMAT(20)
  60. da(9) =XMAT(21)
  61. da(10) =XMAT(22)
  62.  
  63. XKBel=da(1)
  64. Ybo=da(2)
  65. C1=da(3)
  66. D1=da(4)
  67. alphab=da(5)
  68. betab=da(6)
  69. XKAel=da(7)
  70. Yao=da(8)
  71. alphaa=da(9)
  72. betaa=da(10)
  73. C---------------------------------------------------------------------
  74. C Affectation des données historiques courantes
  75.  
  76. force(1)=SIG0(1)
  77. force(2)=SIG0(2)
  78. force(3)=SIG0(3)
  79. force(4)=SIG0(4)
  80. force(5)=SIG0(5)
  81. force(6)=SIG0(6)
  82.  
  83. Dplus=VAR0(2)
  84. Dmoin=VAR0(3)
  85. Dh=VAR0(4)
  86. Zplus=VAR0(5)
  87. Zmoin=VAR0(6)
  88. up=VAR0(7)
  89. Ecrb=VAR0(8)
  90. Fp=VAR0(9)
  91. ua=VAR0(10)
  92. Ecra=VAR0(11)
  93. Fa=VAR0(12)
  94. C---------------------------------------------------------------------
  95.  
  96. uo(1)=epst0(1)
  97. uo(2)=epst0(2)
  98. uo(3)=epst0(3)
  99. uo(4)=epst0(4)
  100. uo(5)=epst0(5)
  101. uo(6)=epst0(6)
  102.  
  103.  
  104. duo(1)=DEPST(1)
  105. duo(2)=DEPST(2)
  106. duo(3)=DEPST(3)
  107. duo(4)=DEPST(4)
  108. duo(5)=DEPST(5)
  109. duo(6)=DEPST(6)
  110.  
  111. C Correction données d'entrées Castem
  112.  
  113. DO I=1,6
  114. uo(I)=uo(I)+duo(I)
  115. ENDDO
  116.  
  117. C PRINT*, uo
  118. C PRINT*, duo
  119.  
  120.  
  121.  
  122. C *****************************************************************
  123.  
  124.  
  125. do i=1,6
  126. do j=1,6
  127. XMKBel(i,j)=0.
  128. enddo
  129. enddo
  130. XMKBel(2,2)=XKBel
  131. C XMKBel(4,4)=XKBel
  132.  
  133. do i=1,6
  134. do j=1,6
  135. XMKAel(i,j)=0.
  136. enddo
  137. enddo
  138. do i=1,6
  139. XMKAel(i,i)=XMAT(i+6)
  140. enddo
  141. XMKAel(2,2)=XKAel
  142. C XMKAel(4,4)=XKAel
  143.  
  144. C ************************Endomagement du béton
  145.  
  146. Ym=5.d-1*XMKBel(2,2)*(uo(2)**2.d0)
  147. Ymp=Ym
  148.  
  149. IF (uo(2) .ge. 0.d0) THEN
  150.  
  151. Fonc1=(Ymp-Zplus-Ybo)
  152.  
  153. IF (Fonc1 .gt. 0.d0) THEN
  154.  
  155. Dplus=1.d0-(1.d0/(1.d0+C1*(Ym**D1)))
  156. Dmoin=Dmoin
  157. Dh=max(Dplus,Dmoin)
  158. Zplus=Zplus+(5.d-1*XMKBel(2,2)*duo(2)*(2.d0*uo(2)-duo(2)))
  159. Zmoin=Zmoin
  160.  
  161. ELSE
  162. Dh=Dh
  163. Dplus=Dplus
  164. Dmoin=Dmoin
  165. Zplus=Zplus
  166. Zmoin=Zmoin
  167.  
  168. ENDIF
  169.  
  170. ELSE
  171.  
  172. Fonc1=(Ymp-Zmoin-Ybo)
  173.  
  174. C print*, 'Déplacement négatif'
  175.  
  176. IF (Fonc1 .gt. 0.d0) THEN
  177.  
  178. Dplus=Dplus
  179. Dmoin=1.d0-(1.d0/(1.d0+C1*(Ym**D1)))
  180. Dh=max(Dplus,Dmoin)
  181. Zmoin=Zmoin+(5.d-1*XMKBel(2,2)*duo(2)*(2.d0*uo(2)-duo(2)))
  182. Zplus=Zplus
  183.  
  184. ELSE
  185.  
  186. Dh=Dh
  187. Dplus=Dplus
  188. Dmoin=Dmoin
  189. Zplus=Zplus
  190. Zmoin=Zmoin
  191.  
  192. ENDIF
  193. ENDIF
  194.  
  195. C***************Frottement du beton
  196.  
  197. Fp=Dh*XMKBel(2,2)*(uo(2)-up)
  198.  
  199. Fonc2=ABS(Fp-Ecrb)
  200.  
  201. IF (Fonc2 .gt. 0.d0) THEN
  202. conv=.false.
  203. iter=0
  204. DO WHILE (.not. conv .and. iter .le. 500)
  205.  
  206. iter=iter+1
  207.  
  208. IF ( (Fp-Ecrb) .GE. 0.d0) THEN
  209.  
  210. Xmul1=1.d0
  211. ELSE
  212. Xmul1=-1.d0
  213. ENDIF
  214.  
  215. Fonc3=ABS(Fp-Ecrb)
  216.  
  217. dfdfb=Xmul1
  218. dpidfb=Xmul1
  219. dfdxb=-1.0*Xmul1
  220. dpidxb=-1.0*Xmul1+alphab*Ecrb
  221.  
  222. Xlamda1=Fonc3/(XMKBel(2,2)*Dh*dfdfb*dpidfb+betab*dfdxb*dpidxb)
  223.  
  224. up=up+Xlamda1*dpidfb
  225.  
  226. Ecrb=Ecrb-betab*Xlamda1*dpidxb
  227. Fp=Fp-XMKBel(2,2)*Dh*Xlamda1*dpidfb
  228.  
  229. Fonc3=ABS(Fp-Ecrb)
  230.  
  231. conv=((ABS(Fonc3/Fonc2) .LE. 1.d-7). or.
  232. & (Xlamda1 .LE. 1.d-7))
  233.  
  234. ENDDO
  235.  
  236. IF (.not. conv) THEN
  237. KERRE = 22
  238. RETURN
  239. ENDIF
  240. ENDIF
  241.  
  242.  
  243. C**************Acier plastique
  244.  
  245. Fa=XMKAel(2,2)*(uo(2)-ua)
  246. Fonc4=ABS(Fa-Ecra)-Yao
  247.  
  248. IF (Fonc4 .gt. 0.d0) THEN
  249. conv=.false.
  250. iter=0
  251. DO WHILE (.not. conv .and. iter .le. 500)
  252. iter=iter+1
  253.  
  254. IF ( (Fa-Ecra) .GE. 0.d0) THEN
  255.  
  256. Xmul2=1.d0
  257. ELSE
  258. Xmul2=-1.d0
  259. ENDIF
  260.  
  261. Fonc5=ABS(Fa-Ecra)-Yao
  262.  
  263. dfdfa=Xmul2
  264. dpidfa=Xmul2
  265. dfdxa=-1.0*Xmul2
  266. dpidxa=-1.0*Xmul2+alphaa*Ecra
  267.  
  268. Xlamda2=Fonc5/(XMKAel(2,2)*dfdfa*dpidfa+betaa*dfdxa*dpidxa)
  269.  
  270. ua=ua+Xlamda2*dpidfa
  271. Ecra=Ecra-betaa*Xlamda2*dpidxa
  272. Fa=Fa-XMKAel(2,2)*Xlamda2*dpidfa
  273.  
  274. Fonc5=ABS(Fa-Ecra)-Yao
  275.  
  276. conv=((ABS(Fonc5/Fonc4) .LE. 1.d-7). or.
  277. & (Xlamda2 .LE. 1.D-7))
  278.  
  279. ENDDO
  280.  
  281. IF (.not. conv) THEN
  282. KERRE = 22
  283. RETURN
  284. ENDIF
  285. ENDIF
  286.  
  287. Fa=XMKAel(2,2)*(uo(2)-ua)
  288.  
  289.  
  290. C*****************************************************************
  291. IF (uo(2) .ge. 0.d0) THEN
  292. force(2)=XMKBel(2,2)*uo(2)*(1.d0-Dplus)
  293. & +Dh*XMKBel(2,2)*(uo(2)-up)+Fa
  294. ELSE
  295. force(2)=XMKBel(2,2)*uo(2)*(1.d0-Dmoin)
  296. & +Dh*XMKBel(2,2)*(uo(2)-up)+Fa
  297. ENDIF
  298.  
  299. force(1)=XMAT(7)*uo(1)
  300. force(3)=XMAT(9)*uo(3)
  301. force(4)=XMAT(10)*uo(4)
  302. force(5)=XMAT(11)*uo(5)
  303. force(6)=XMAT(12)*uo(6)
  304.  
  305. C*****************************************************************
  306. EPSE=sqrt(EPSE)
  307. VARF(1)=EPSE
  308. VARF(2)=Dplus
  309. VARF(3)=Dmoin
  310. VARF(4)=Dh
  311. VARF(5)=Zplus
  312. VARF(6)=Zmoin
  313. VARF(7)=up
  314. VARF(8)=Ecrb
  315. VARF(9)=Fp
  316. VARF(10)=ua
  317. VARF(11)=Ecra
  318. VARF(12)=Fa
  319.  
  320. SIGF(1)=force(1)
  321. SIGF(2)=force(2)
  322. SIGF(3)=force(3)
  323. SIGF(4)=force(4)
  324. SIGF(5)=force(5)
  325. SIGF(6)=force(6)
  326.  
  327. return
  328.  
  329. C#####################################################################
  330. C#####################################################################
  331.  
  332. end
  333.  
  334.  
  335.  
  336.  
  337.  
  338.  
  339.  
  340.  
  341.  
  342.  

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