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