37
38
39
40#include "implicit_f.inc"
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80#include "mvsiz_p.inc"
81#include "scr17_c.inc"
82#include "units_c.inc"
83#include "comlock.inc"
84#include "param_c.inc"
85#include "impl1_c.inc"
86
87
88
89
90 INTEGER NEL, NUPARAM, NUVAR,NGL(NEL)
92 . time,timestep,uparam(*),
93 . signxx(nel),signyy(nel),signzz(nel),
94 . signxy(nel),signyz(nel),signzx(nel)
95
96
97
98
99
100
101
102
103 my_real uvar(nel,nuvar), off(nel),dfmax(nel),tdele(nel)
104
105
106
107 INTEGER I,J,IDEL,,IFLAG,INDX(MVSIZ),IADBUF,NINDX,
108 . NINDEX,INDEX(MVSIZ),IFAIL,IR,JJ
110 . tba,tbk,sigr,svm,scale,sxx,syy,szz,
111 . e1,e2,e3,e4,e5,e6,
112 . epst2, a, cc1, b, y ,yp, d, e42, e52, c, e62,p,sigmax
113
114
115 tba = uparam(1)
116 tbk = uparam(2)
117 sigr = uparam(3)
118 iflag = int(uparam(5))
119
120
121 idel=0
122 idev=0
123 scale = zero
124 IF(iflag==1)THEN
125 idel=1
126 ELSEIF(iflag==2)THEN
127 idev =1
128 END IF
129
130 IF(idel==1)THEN
131 DO i=1,nel
132 IF(off(i)<0.1) off(i)=0.0
133 IF(off(i)<1.0) off(i)=off(i)*0.8
134 END DO
135 END IF
136
137 IF(idel==1)THEN
138 nindx=0
139
140
141
142
143 DO i=1,nel
144 IF(iflag==1.AND.off(i)==1.)THEN
145
146
147
148 p = third*(signxx(i) + signyy(i) + signzz(i))
149 e1 = signxx(i) - p
150 e2 = signyy(i) - p
151 e3 = signzz(i) - p
152 e4 = signxy(i)
153 e5 = signyz(i)
154 e6 = signzx(i)
155
156
157
158
159
160
161 e42 = e4*e4
162 e52 = e5*e5
163 e62 = e6*e6
164 c = - half * (e1*e1 + e2*e2 + e3*e3) - e42 - e52 - e62
165 d = - e1*e2*e3 + e1*e52 + e2*e62 + e3*e42
166 & - two*e4*e5*e6
167 cc1 = c*third
168 sigmax = sqrt(-cc1)
169 epst2 = sigmax * sigmax
170 y = (epst2 + c)* sigmax + d
171 IF(abs(y)>em8)THEN
172 sigmax = 1.75 * sigmax
173 epst2 = sigmax * sigmax
174 y = (epst2 + c)* sigmax + d
175 yp = three*epst2 + c
176 IF(yp/=zero)sigmax = sigmax - y/yp
177 epst2 = sigmax * sigmax
178 y = (epst2 + c)* sigmax + d
179 yp = three*epst2 + c
180 IF(yp/=zero)sigmax = sigmax - y/yp
181 epst2 = sigmax * sigmax
182 y = (epst2 + c)* sigmax + d
183 yp = three*epst2 + c
184 IF(yp/=zero)sigmax = sigmax - y/yp
185 epst2 = sigmax * sigmax
186 y = (epst2 + c)* sigmax + d
187 yp = three*epst2 + c
188 IF(yp/=zero)sigmax = sigmax - y/yp
189 ENDIF
190 sigmax = sigmax + p
191 IF(sigmax>=sigr)
192 . uvar(i,1)=uvar(i,1) + timestep*(sigmax - sigr)**tba
193
194 IF (uvar(i,1)>tbk) THEN
195 off(i) = four_over_5
196 nindx=nindx+1
197 indx(nindx)=i
198 tdele(i) = time
199 ENDIF
200 ENDIF
201 ENDDO
202 IF(nindx>0.AND.imconv==1)THEN
203 DO j=1,nindx
204 i=indx(j)
205#include "lockon.inc"
206 WRITE(istdo,1000)ngl(i)
207 WRITE(iout,1100)ngl(i),time
208#include "lockoff.inc"
209 END DO
210 END IF
211
212 END IF
213
214
215 IF(idev==1)THEN
216 nindx=0
217 DO i=1,nel
218 IF(iflag==2.AND.off(i)==1.)THEN
219 IF(uvar(i,1)<tbk)THEN
220 p = third*(signxx(i) + signyy(i) + signzz(i))
221 e1 = signxx(i) - p
222 e2 = signyy(i) - p
223 e3 = signzz(i) - p
224 e4 = signxy(i)
225 e5 = signyz(i)
226 e6 = signzx(i)
227
228
229
230
231
232
233 e42 = e4*e4
234 e52 = e5*e5
235 e62 = e6*e6
236 c = - half * (e1*e1 + e2*e2 + e3*e3) - e42 - e52 - e62
237 d = - e1*e2*e3 + e1*e52 + e2*e62 + e3*e42
238 & - two*e4*e5*e6
239 cc1 = c*third
240 sigmax = sqrt(-cc1)
241 epst2 = sigmax * sigmax
242 y = (epst2 + c)* sigmax + d
243 IF(abs(y)>em8)THEN
244 sigmax = 1.75 * sigmax
245 epst2 = sigmax * sigmax
246 y = (epst2 + c)* sigmax + d
247 yp = three*epst2 + c
248 IF(yp/=zero)sigmax = sigmax - y/yp
249 epst2 = sigmax * sigmax
250 y = (epst2 + c)* sigmax + d
251 yp = three*epst2 + c
252 IF(yp/=zero)sigmax = sigmax - y/yp
253 epst2 = sigmax * sigmax
254 y = (epst2 + c)* sigmax + d
255 yp = three*epst2 + c
256 IF(yp/=zero)sigmax = sigmax - y/yp
257 epst2 = sigmax * sigmax
258 y = (epst2 + c)* sigmax + d
259 yp = three*epst2 + c
260 IF(yp/=zero)sigmax = sigmax - y/yp
261 ENDIF
262 sigmax = sigmax + p
263 IF(sigmax>=sigr)
264 . uvar(i,1)=uvar(i,1) + timestep*(sigmax - sigr)**tba
265
266 IF(uvar(i,1)>tbk) THEN
267 nindx=nindx+1
268 indx(nindx)=i
269 signxx(i) = p
270 signyy(i) = p
271 signzz(i) = p
272 signxy(i) = zero
273 signyz(i) = zero
274 signzx(i) = zero
275 ENDIF
276
277 ELSE
278 p= third*(signxx(i) + signyy(i) + signzz(i))
279 signxx(i) = p
280 signyy(i) = p
281 signzz(i) = p
282 signxy(i) = zero
283 signyz(i) = zero
284 signzx(i) = zero
285 ENDIF
286 ENDIF
287 ENDDO
288 IF(nindx>0.AND.imconv==1)THEN
289 DO j=1,nindx
290 i = indx(j)
291#include "lockon.inc"
292 WRITE(iout, 2000) ngl(i)
293 WRITE(istdo,2100) ngl(i),time
294#include "lockoff.inc"
295 END DO
296 END IF
297 ENDIF
298
299 DO i=1,nel
300 dfmax(i)=
min(one,
max(dfmax(i),uvar(i,1)/tbk))
301 ENDDO
302
303 1000 FORMAT(1x,'DELETE SOLID ELEMENT NUMBER ',i10)
304 1100 FORMAT(1x,'DELETE SOLID ELEMENT NUMBER ',i10,
305 . ' AT TIME :',1pe20.13)
306
307 2000 FORMAT(1x,' DEVIATORIC STRESS WILL BE VANISHED',i10)
308 2100 FORMAT(1x,' DEVIATORIC STRESS WILL BE VANISHED',i10,
309 . ' AT TIME :',1pe20.13)
310 RETURN