38
39
40
41
42
44
45
46
47#include "implicit_f.inc"
48
49
50
51#include "units_c.inc"
52#include "comlock.inc"
53#include "com_xfem1.inc"
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
81
82
83
84
85
86 INTEGER NEL, NUPARAM, ,IXFEM,IXEL,ILAY,NPTT,
87 . IPT,IPTT,NPTOT
88 INTEGER NGL(NEL),ELCRKINI(NXLAYMAX,*),NOFF(NEL)
89 my_real time,uparam(*),dpla(nel),npttf(nel),
90 . epsp(nel),tstar(nel),tens(nel,5)
91
92
93
94 my_real uvar(nel,nuvar),off(nel),offl(nel),dfmax(nel),tdel(nel),
95 . signxx(nel),signyy(nel),signxy(nel),signyz(nel),signzx(nel)
96
97
98
99 INTEGER I,J,K,L,JJ,ISHELL,NINDX,NINDXP,IADR,LAYXFEM
100 INTEGER INDX(NEL),INDXP(NEL),RFLAG(),RFLAGP(NEL)
101 my_real d1,d2,d3,d4,d5,epsp0,p,epsf,svm,thk_tmp,
102 . dadv,cc,bb,cr,orm,ss1,ss2,epsf_min
103 CHARACTER (LEN=3) :: XCHAR
104
105 iadr = (ipt-1)*nel
106 xchar = ' '
107
108 d1 = uparam(1)
109 d2 = uparam(2)
110 d3 = uparam(3)
111 d4 = uparam(4)
112 d5 = uparam(5)
113 epsp0 = uparam(6)
114 dadv = uparam(10)
115 epsf_min = uparam(12)
116 ishell = int(uparam(7))
117
118 layxfem = ixfem
119 IF (layxfem == 1 .and. ishell == 1) ishell=2
120
121 nindx = 0
122 nindxp = 0
123 DO i=1,nel
124 rflag(i) = 0
125 rflagp(i)= 0
126 indxp(i) = 0
127 indx(i) = 0
128 ENDDO
129
130 DO i=1,nel
131 tens(i,1) = signxx(i)
132 tens(i,2) = signyy(i)
133 tens(i,3) = signxy(i)
134 tens(i,4) = signyz(i)
135 tens(i,5) = signzx(i)
136 END DO
137
138 IF (ixel > 0) THEN
139 IF (ixel == 1) THEN
140 xchar = '1st'
141 ELSEIF (ixel == 2) THEN
142 xchar = '2nd'
143 ELSEIF (ixel == 3) THEN
144 xchar = '3rd'
145 ENDIF
146 ENDIF
147
148 SELECT CASE (layxfem)
149
150 CASE (1)
151
152 DO i=1,nel
153 IF (off(i) == one) THEN
154 IF (dfmaxTHEN
155 p = third*(signxx(i) + signyy(i))
156 svm = sqrt(signxx(i)*signxx(i)
157 . + signyy(i)*signyy(i)
158 . - signxx(i)*signyy(i)
159 . + three*signxy(i)*signxy(i))
160 epsf = d3*p/
max(em20,svm)
161 epsf = (d1 + d2*exp(epsf))*(one
162 . + d4*log(
max(one,epsp(i)/epsp0)))
163 . * (one + d5*tstar(i))
164 epsf =
max(epsf,epsf_min)
165 IF (epsf > zero) dfmax(i) = dfmax(i) + dpla(i)/epsf
166
167 IF (ixel == 0) THEN
168 IF (elcrkini(ilay,i)==0 .AND. dfmax(i) >= oneTHEN
169 nindxp = nindxp+1
170 indxp(nindxp) = i
171 offl(i) = zero
172 noff(i) = noff(i) + 1
173 npttf(i) = npttf(i) + one
174 rflagp(i) = 1
175 IF (int(npttf(i)) == nptt) THEN
176 nindx = nindx+1
177 indx(nindx) = i
178 elcrkini(ilay,i) = -1
179 rflag(i) = 1
180 ENDIF
181 IF (noff(i) == nptot) THEN
182 off(i) = four_over_5
183 tdel(i)= time
184 ENDIF
185 ELSEIF (elcrkini(ilay,i) == 2 .AND. dfmax(i) >= dadv) THEN
186 nindxp = nindxp + 1
187 indxp(nindxp) = i
188 offl(i) = zero
189 noff(i) = noff(i) + 1
190 npttf(i) = npttf(i)+ one
191 rflagp(i) = 1
192 IF (int(npttf(i)) == nptt) THEN
193 nindx = nindx+1
194 indx(nindx) = i
195 elcrkini(ilay,i) = 1
196 rflag(i) = -1
197 ENDIF
198 IF (noff(i) == nptot) THEN
199 off(i) = four_over_5
200 tdel(i)= time
201 ENDIF
202 ENDIF
203 ELSEIF (dfmax(i) >= one) THEN
204 nindxp = nindxp + 1
205 indxp(nindxp) = i
206 offl(i) = zero
207 npttf(i) = npttf(i) + one
208 rflagp(i) = 1
209 IF (int(npttf(i)) == nptt) THEN
210 nindx = nindx+1
211 indx(nindx) = i
212 off(i) = four_over_5
213 tdel(i) = time
214 rflag(i) = 3
215 ENDIF
216 ENDIF
217 ELSE
218
219 signyy(i) = zero
220 signxy(i) = zero
221 signyz(i) = zero
222 signzx(i) = zero
223 ENDIF
224 ENDIF
225 ENDDO
226
227 IF (nindxp > 0) THEN
228 DO j=1,nindxp
229 i = indxp(j)
230 signxx(i) = zero
231 signyy(i) = zero
232 signxy(i) = zero
233 signyz(i) = zero
234 signzx(i) = zero
235#include "lockon.inc"
236 IF (rflagp(i) == 1) WRITE(iout, 3800)ngl(i),ilay,iptt
237 IF (rflagp(i) == 1) WRITE(istdo,3900)ngl(i),ilay,iptt,time
238#include "lockoff.inc"
239 ENDDO
240 ENDIF
241 IF (nindx > 0) THEN
242 DO j=1,nindx
243 i = indx(j)
244#include "lockon.inc"
245
246 IF (rflag(i) == 1) WRITE(iout ,3000) ngl(i),ilay
247 IF (rflag(i) == 1) WRITE(istdo,3100) ngl(i),ilay,time
248
249 IF (rflag(i) == -1) WRITE(iout ,3200) ngl(i),ilay
250 IF (rflag(i) == -1) WRITE(istdo,3300) ngl(i),ilay,time
251
252 IF (rflag(i) == 3) WRITE(iout, 3400)xchar,ngl(i),ilay
253 IF (rflag(i) == 3) WRITE(istdo,3500)xchar,ngl(i),ilay,time
254#include "lockoff.inc"
255 ENDDO
256 ENDIF
257
258
259 CASE (2)
260
261 IF (ishell == 1) THEN
262
263 DO i=1,nel
264 IF (off(i) == one) THEN
265 p = third*(signxx(i) + signyy(i))
266 svm = sqrt(signxx(i)*signxx(i)
267 . + signyy(i)*signyy(i)
268 . - signxx(i)*signyy(i)
269 . + three*signxy(i)*signxy(i))
270 epsf = d3*p/
max(em20,svm)
271 epsf = (d1 + d2*exp(epsf))*(one
272 . + d4*log(
max(one,epsp(i)/epsp0)))
273 . * (one + d5*tstar(i))
274 epsf =
max(epsf,epsf_min)
275 IF (epsf > zero) dfmax(i) = dfmax(i) + dpla(i)/epsf
276
277 IF (ixel == 0) THEN
278 IF (elcrkini(ilay,i)==0 .AND. dfmax(i)>=one) THEN
279 elcrkini(ilay,i) = -1
280 nindx = nindx + 1
281 indx(nindx) = i
282 rflag(i) = 1
283 tdel(i)= time
284 off(i) = four_over_5
285 ELSEIF (elcrkini(ilay,i) == 2 .AND.
286 . dfmax(i) >= dadv) THEN
287 elcrkini(ilay,i) = 1
288 nindx = nindx + 1
289 indx(nindx) = i
290 rflag(i) =-1
291 tdel(i)= time
292 off(i) = four_over_5
293 ENDIF
294 ELSEIF (dfmax(i) >= one) THEN
295 nindx = nindx+1
296 indx(nindx) = i
297 rflag(i) = 3
298 off(i) = four_over_5
299 ENDIF
300 ENDIF
301 ENDDO
302
303 IF (nindx > 0) THEN
304 DO j=1,nindx
305 i = indx(j)
306#include "lockon.inc"
307 WRITE(iout, 4800)ngl(i),iptt
308 WRITE(istdo,4900)ngl(i),iptt,time
309
310 IF (rflag(i) == 1) WRITE(iout, 4000) ngl(i)
311 IF (rflag(i) == 1) WRITE(istdo,4100) ngl(i),time
312
313 IF (rflag(i) ==-1) WRITE(iout, 4200) ngl(i)
314 IF (rflag(i) ==-1) WRITE(istdo,4300) ngl(i),time
315
316 IF (rflag(i) == 3) WRITE(iout, 4400) xchar,ngl(i)
317 IF (rflag(i) == 3) WRITE(istdo,4500) xchar,ngl(i),time
318#include "lockoff.inc"
319 ENDDO
320 ENDIF
321
322 ELSEIF (ishell == 2) THEN
323
324 DO i=1,nel
325 IF (off(i) == one) THEN
326 IF (dfmax(i) < one) THEN
327 p = third*(signxx(i) + signyy(i))
328 svm = sqrt(signxx(i)*signxx(i)
329 . + signyy(i)*signyy(i)
330 . - signxx(i)*signyy(i)
331 . + three*signxy(i)*signxy(i))
332 epsf = d3*p/
max(em20,svm)
333 epsf = (d1 + d2*exp(epsf))*(one
335 . * (one + d5*tstar(i))
336 epsf =
max(epsf,epsf_min)
337 IF (epsf > zero) dfmax(i) = dfmax(i) + dpla(i)/epsf
338
339 IF (ixel == 0) THEN
340 IF (elcrkini(ilay,i)==0 .and. dfmax(i) >= one) THEN
341 nindxp = nindxp + 1
342 indxp(nindxp) = i
343 noff(i) = noff(i) + 1
344 offl(i) = zero
345 rflagp(i) = 1
346 IF (noff(i) == nptot) THEN
347 nindx = nindx + 1
348 indx(nindx) = i
349 elcrkini(ilay,i) = -1
350 rflag(i) = 1
351 tdel(i)= time
352 off(i) = four_over_5
353 ENDIF
354 ELSEIF (elcrkini(ilay,i) == 2 .AND. dfmax(i) >= dadv) THEN
355 nindxp = nindxp + 1
356 indxp(nindxp) = i
357 noff(i) = noff(i) + 1
358 offl(i) = zero
359 dfmax(i) = one
360 rflagp(i) = 1
361 IF (noff(i) == nptot) THEN
362 nindx = nindx + 1
363 indx(nindx) = i
364 elcrkini(ilay,i) = 1
365 rflag(i) = -1
366 tdel(i)= time
367 off(i) = four_over_5
368 ENDIF
369 ENDIF
370 ELSEIF (dfmax(i) >= one) THEN
371 nindxp = nindxp + 1
372 indxp(nindxp) = i
373 noff(i) = noff(i) + 1
374 rflagp(i) = 1
375 IF (noff(i) == nptot) THEN
376 nindx = nindx + 1
377 indx(nindx) = i
378 off(i) = four_over_5
379 rflag(i) = 3
380 ENDIF
381 ENDIF
382 ELSE
383 signxx(i) = zero
384 signyy(i) = zero
385 signxy(i) = zero
386 signyz(i) = zero
387 signzx(i) = zero
388 ENDIF
389 ENDIF
390 ENDDO
391
392 IF (nindxp > 0) THEN
393 DO j=1,nindxp
394 i = indxp(j)
395 signxx(i) = zero
396 signyy(i) = zero
397 signxy(i) = zero
398 signyz(i) = zero
399 signzx(i) = zero
400#include "lockon.inc"
401 IF (rflagp(i) == 1) WRITE(iout, 4800)ngl(i),iptt
402 IF (rflagp(i) == 1) WRITE(istdo,4900)ngl(i),iptt,time
403#include "lockoff.inc"
404 ENDDO
405 ENDIF
406 IF (nindx > 0) THEN
407 DO j=1,nindx
408 i = indx(j)
409#include "lockon.inc"
410
411 IF (rflag(i) == 1) WRITE(iout ,4000) ngl(i)
412 IF (rflag(i) == 1) WRITE(istdo,4100) ngl(i),time
413
414 IF (rflag(i) == -1) WRITE(iout, 4200) ngl(i)
415 IF (rflag(i) == -1) WRITE(istdo,4300) ngl(i),time
416
417 IF (rflag(i) == 3) WRITE(iout, 4400) xchar,ngl(i)
418 IF (rflag(i) == 3) WRITE(istdo,4500) xchar,ngl(i),time
419#include "lockoff.inc"
420 ENDDO
421 ENDIF
422
423 ENDIF
424
425 END SELECT
426
427
428 DO i=1,nel
429 dfmax(i) =
min(one,dfmax(i))
430 ENDDO
431
432 2000 FORMAT(1x,'FOR SHELL ELEMENT (JOHNS)',i10,1x,'LAYER',i3,':',/
433 . 1x, 'STRESS TENSOR SET TO ZERO')
434 2100 FORMAT(1x,'FOR SHELL ELEMENT (JOHNS)',i10,1x,'LAYER',i3,':',/,
435 . 1x,'STRESS TENSOR SET TO ZERO',1x,'AT TIME :',1pe12.4)
436 2400 FORMAT(1x,1pg20.13,' % OF THICKNESS OF SHELL BROKEN ')
437 2410 FORMAT(1x,1pg20.13,' % OF THICKNESS OF SHELL ',i10,' BROKEN ')
438 2500 FORMAT(1x,' LOWER SKIN -> UPPER SKIN ')
439 2600 FORMAT(1x,' UPPER SKIN -> LOWER SKIN ')
440
441 3000 FORMAT(1x,'CRACK INITIALIZATION IN SHELL ELEMENT',i10,1x,'LAYER',i3)
442 3100 FORMAT(1x,'CRACK INITIALIZATION IN SHELL ELEMENT',i10,1x,'LAYER',i3,/
443 . 1x,'AT TIME :',1pe12.4)
444 3200 FORMAT(1x,'CRACK ADVANCEMENT IN SHELL ELEMENT',i10,' LAYER',i3)
445 3300 FORMAT(1x,'CRACK ADVANCEMENT IN SHELL ELEMENT',i10,' LAYER',i3/
446 . 1x,'AT TIME :',1pe12.4)
447 3400 FORMAT(1x,'DELETE ',a4,' PHANTOM ELEMENT, SHELL ID=',i10,' LAYER',i3)
448 3500 FORMAT(1x,'DELETE ',a4,' PHANTOM ELEMENT, SHELL ID=',i10,' LAYER',i3,/
449 . 1x,'AT TIME :',1pe12.4)
450 3800 FORMAT(1x,'JC FAILURE IN SHELL',i10,1x,'LAYER',i3,1x,'INT POINT',i2)
451 3900 FORMAT(1x,'JC FAILURE IN SHELL',i10,1x,'LAYER',i3,1x,'INT POINT',i2,/
452 . 1x,'AT TIME :',1pe12.4)
453
454 4000 FORMAT(1x,'CRACK INITIALIZATION IN SHELL ELEMENT',i10)
455 4100 FORMAT(1x,'CRACK INITIALIZATION IN SHELL ELEMENT',i10,1x,'AT TIME :',1pe12.4)
456 4200 FORMAT(1x,'CRACK ADVANCEMENT IN SHELL ELEMENT',i10)
457 4300 FORMAT(1x,'CRACK ADVANCEMENT IN SHELL ELEMENT',i10,1x,'AT TIME :',1pe12.4)
458 4400 FORMAT(1x,'DELETE ',a4,' PHANTOM ELEMENT, SHELL ID=',i10)
459 4500 FORMAT(1x,'DELETE ',a4,' PHANTOM ELEMENT, SHELL ID=',i10,1x,'AT TIME :',1pe12.4)
460 4800 FORMAT(1x,'JC FAILURE IN SHELL',i10,1x,'INT POINT',i2)
461 4900 FORMAT(1x,'JC FAILURE IN SHELL',i10,1x,'INT POINT',i2,1x,'AT TIME :',1pe12.4)
462
463 RETURN