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#include "com01_c.inc"
76#include "units_c.inc"
77#include "comlock.inc"
78
79
80
81 INTEGER NEL,NUPARAM,NUVAR,NFUNC,IPG,ISOLID
82 INTEGER NGL(*),NPF(*),IFUNC(*)
83 INTEGER, INTENT(IN) :: LF_DAMMX
85 my_real uparam(nuparam),uvar(nel,nuvar),tf(*)
86 my_real,
DIMENSION(NEL) :: offg,offl,eps1,eps2,eps3,epst,epsp,
87 . signzz,signyz,signzx,dein,deit,tdele,
area,soft
88 my_real ,
DIMENSION(NEL,LF_DAMMX),
INTENT(INOUT) :: dfmax
89
90
91
92 INTEGER I,J,IDEL,NINDX,NINDXA,NINDXE,NINDXS,FUNN,FUNT,IFAILS,IFAILE,
93 . ISYM
94 INTEGER, DIMENSION(NEL) :: INDX,INDXA,INDXE,INDXS
95 my_real c1,c2,c3,c4,cm,dydx,damt,iflag,defo,dtime,dsoft,areascale,
96 . maxdn,maxdt,maxei,maxen,maxet,fac1,fac2,ascale,expn,expt,nn,nt
97 my_real,
DIMENSION(NEL) :: ei,en,et,dam1,dam2,dam3,dam4
98
99
100
102 EXTERNAL finter
103
104 maxdn = uparam(1)
105 maxdt = uparam(2)
106 expn = uparam(3)
107 expt = uparam(4)
108 fac1 = uparam(5)
109 fac2 = uparam(6)
110 ascale = uparam(7)
111 ifails = nint(uparam(8))
112 ifaile = nint(uparam(9))
113 isolid = nint(uparam(10))
114 maxei = uparam(11)
115 maxen = uparam(12)
116 maxet = uparam(13)
117 nn = uparam(14)
118 nt = uparam(15)
119 dtime = uparam(16)
120 dsoft = uparam(17)
121 isym = nint(uparam(18))
122 areascale = uparam(19)
123
124
125
126 IF (isigi == zero) THEN
127 IF ((uvar(1,1) == zero).AND.(offl(1) == one)) THEN
128 DO i=1,nel
129 uvar(i,1)= one
130 ENDDO
131 ENDIF
132 IF (uvar(1,8) == zero) THEN
133 DO i=1,nel
135 ENDDO
136 ENDIF
137 ENDIF
138
139 DO i=1,nel
140 dam1(i) = uvar(i,2)
141 dam2(i) = uvar(i,3)
142 dam3(i) = uvar(i,4)
143 dam4(i) = uvar(i,5)
144 en(i) = uvar(i,6)
145 et(i) = uvar(i,7)
146 ENDDO
147
148 nindx = 0
149 nindxe = 0
150 nindxs = 0
151 nindxa = 0
152
153 DO i=1,nel
154 idel = 0
155 et(i) = et(i) + deit(i)
156 IF (isym == 1 .AND. signzz(i) < zero) THEN
157 fac1 = zero
158 ELSE
159 en(i) = en(i) + dein(i)
160 ENDIF
161 ei(i) = en(i) + et(i)
162 c1 = zero
163 c2 = zero
164 c3 = zero
165 c4 = zero
166 cm = zero
167
168
169
170 IF (ifaile == 1) THEN
171 c1 = ei(i)/maxei
172 dfmax(i,7) =
max(dfmax(i,7),c1)
173 dfmax(i,7) =
min(dfmax(i,7),one)
174 IF (c1 > one) THEN
175 IF (dam1(i) == zero) THEN
176 nindxe = nindxe+1
177 indxe(nindxe) = i
178 ENDIF
179 dam1(i) = dam1(i) + c1 * timestep
180 ENDIF
181
182 ELSEIF (ifaile == 2) THEN
183 c2 = (et(i)/maxet)**nt+ (en(i)/maxen)**nn
184 dfmax(i,5) =
max(dfmax(i,5),en(i)/maxen)
185 dfmax(i,5) =
min(dfmax(i,5),one)
186 dfmax(i,6) =
max(dfmax(i,6),et(i)/maxet)
187 dfmax(i,6) =
min(dfmax(i,6),one)
188 IF (c2 > one) THEN
189 IF (dam1(i) == zero) THEN
190 nindxe = nindxe+1
191 indxe(nindxe) = i
192 ENDIF
193 dam1(i) = dam1(i) + c2*timestep
194 ENDIF
195 ELSEIF (ifaile == 3) THEN
196 c1 = ei(i)/maxei
197 dfmax(i,7) =
max(dfmax(i,7),c1)
198 dfmax(i,7) =
min(dfmax(i,7),one)
199 c2 = (et(i)/maxet)**nt+ (en(i)/maxen)**nn
200 dfmax(i,5) =
max(dfmax(i,5),en(i)/maxen)
201 dfmax(i,5) =
min(dfmax(i,5),one)
202 dfmax(i,6) =
max(dfmax(i,6),et(i)/maxet)
203 dfmax(i,6) =
min(dfmax(i,6),one)
205 dfmax(i,8) =
max(dfmax(i,8),c1)
206 dfmax(i,8) =
min(dfmax(i,8),one)
207 IF (c1 > one) THEN
208 IF (dam1(i) == zero) THEN
209 nindxe = nindxe+1
210 indxe(nindxe) = i
211 ENDIF
212 dam1(i) = dam1(i) + c1*timestep
213 ENDIF
214 ENDIF
215
216
217
218 IF (ifails > 0) THEN
219 funn = ifunc(1)
220 funt = ifunc(2)
221
222 epst(i) = sqrt(eps2(i)**2 + eps3(i)**2)
223 c3 = fac1*eps1(i)/maxdn
224 c4 = fac2*epst(i)/maxdt
225 IF (funn > 0) c3 = c3*finter(funn,epsp(i)*ascale,npf,tf,dydx)
226 IF (funt > 0) c4 = c4*finter(funt,epsp(i)*ascale,npf,tf,dydx)
227 dfmax(i,2) =
max(dfmax(i,2),c3)
228 dfmax(i,2) =
min(dfmax(i,2),one)
229 dfmax(i,3) =
max(dfmax(i,3),c4)
230 dfmax(i,3) =
min(dfmax(i,3),one)
231 IF (ifails == 1) THEN
232
234 dfmax(i,4) =
max(dfmax(i,4),c3)
235 dfmax(i,4) =
min(dfmax(i,4),one)
236 IF (c3 > one) THEN
237 IF (dam3(i) == zero) THEN
238 nindxs = nindxs+1
239 indxs(nindxs) = i
240 ENDIF
241 dam3(i) = dam3(i) + c3*timestep
242 ENDIF
243 ELSEIF (ifails == 2) THEN
244
245 cm = abs(c3)**expn + abs(c4)**expt
246 dfmax(i,4) =
max(dfmax(i,4),cm)
247 dfmax(i,4) =
min(dfmax(i,4),one)
248 IF (cm > one) THEN
249 IF (dam3(i) == zero) THEN
250 nindxs = nindxs+1
251 indxs(nindxs) = i
252 ENDIF
253 dam3(i) = dam3(i) + cm* timestep
254 ENDIF
255 ENDIF
256 ENDIF
257
258 damt =
max(dam1(i),dam3(i))
259 IF (damt > dtime) THEN
260 idel =1
261 soft(i) = zero
262 ELSEIF (dsoft == one) THEN
263 soft(i) = one - damt/
max(dtime,em20)
264 ELSEIF (dsoft > zero) THEN
265 soft(i) = (one - damt/
max(dtime,em20))**dsoft
266 ENDIF
267
268 IF (idel == 1 .AND. offl(i) == one) THEN
269 offl(i) = zero
270 nindx = nindx+1
271 indx(nindx) = i
272 tdele(i) = time
273 ENDIF
274
275
276 IF (areascale > zero .AND. offg(i) == one ) THEN
277 defo = uvar(i,8) * areascale
278 IF (
area(i) > defo)
THEN
279 offl(i) = zero
280 nindxa = nindxa+1
281 indxa(nindxa) = i
282 tdele(i) = time
283 isolid = 1
284 ENDIF
285 ENDIF
286
287 uvar(i,1) = soft(i)
288
289 uvar(i,2) = dam1(i)
290 uvar(i,3) = dam2(i)
291 uvar(i,4) = dam3(i)
292 uvar(i,5) = dam4(i)
293 uvar(i,6) = en(i)
294 uvar(i,7) = et(i)
295
296 dfmax(i,1) =
min(one,
max(dfmax(i,1),damt/
max(dtime,em20)))
297
298 ENDDO
299
300 IF (nindxe > 0) THEN
301 DO j=1,nindxe
302 i = indxe(j)
303#include "lockon.inc"
304 WRITE(iout ,1000) ngl(i),ipg,ei(i)
305 WRITE(istdo,1100) ngl(i),ipg,ei(i),time
306#include "lockoff.inc"
307 END DO
308 ELSEIF (nindxs > 0) THEN
309 DO j=1,nindxs
310 i = indxs(j)
311#include "lockon.inc"
312 WRITE(iout ,1600) ngl(i),ipg,eps1(i),epst(i)
313 WRITE(istdo,1700) ngl(i),ipg,eps1(i),epst(i),time
314#include "lockoff.inc"
315 END DO
316 ELSEIF (nindx > 0) THEN
317 DO j=1,nindx
318 i = indx(j)
319#include "lockon.inc"
320 WRITE(iout ,1200) ngl(i),ipg
321 WRITE(istdo,1300) ngl(i),ipg,time
322#include "lockoff.inc"
323 END DO
324 ELSEIF (nindxa > 0) THEN
325 DO j=1,nindxa
326 i = indxa(j)
327#include "lockon.inc"
328 WRITE(iout ,1400) ngl(i),
area(i)
329 WRITE(istdo,1500) ngl(i),
area(i),time
330#include "lockoff.inc"
331 END DO
332 ENDIF
333
334 1000 FORMAT(5x,'START DAMAGE CONNECTION ELEMENT ',i10,
335 . ' INTEGRATION POINT',i2,', ENERGY=',1pe16.9)
336 1100 FORMAT(5x,'START DAMAGE CONNECTION ELEMENT ',i10,
337 . ' INTEGRATION POINT',i2,', ENERGY=',1pe16.9,
338 . ' AT TIME ',1pe16.9)
339 1200 FORMAT(5x,'FAILURE CONNECTION SOLID ELEMENT ',i10,
340 . ' INTEGRATION POINT',i2)
341 1300 FORMAT(5x,'FAILURE CONNECTION SOLID ELEMENT ',i10,
342 . ' INTEGRATION POINT',i2,' AT TIME ',1pe16.9)
343 1400 FORMAT(5x,'FAILURE CONNECTION SOLID ELEMENT ',i10,
344 . ', AREA(LIMIT REACHED) =',1pe16.9)
345 1500 FORMAT(5x,'FAILURE CONNECTION SOLID ELEMENT ',i10,
346 . ', AREA(LIMIT REACHED) :',1pe16.9,' AT TIME ',1pe16.9)
347 1600 FORMAT(5x,'START DAMAGE CONNECTION ELEMENT ',i10,
348 . ' INTEGRATION POINT',i2,', EPSN=',1pe16.9,', EPST=',1pe16.9)
349 1700 FORMAT(5x,'START DAMAGE CONNECTION ELEMENT ',i10,
350 . ' INTEGRATION POINT',i2,', EPSN=',1pe16.9,', EPST=',1pe16.9,
351 . ' AT TIME ',1pe16.9)
352
353 RETURN
subroutine area(d1, x, x2, y, y2, eint, stif0)