OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
fail_tab2_c.F File Reference
#include "implicit_f.inc"
#include "units_c.inc"
#include "comlock.inc"
#include "com04_c.inc"
#include "tabsiz_c.inc"
#include "lockon.inc"
#include "lockoff.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine fail_tab2_c (nel, nuparam, nuvar, nfunc, ifunc, npf, table, tf, time, uparam, ngl, aldt, dpla, epsp, uvar, signxx, signyy, signxy, signyz, signzx, temp, foff, dfmax, tdele, ipt, ipg, dmg_flag, dmg_scale, ntablf, itablf)

Function/Subroutine Documentation

◆ fail_tab2_c()

subroutine fail_tab2_c ( integer, intent(in) nel,
integer, intent(in) nuparam,
integer, intent(in) nuvar,
integer, intent(in) nfunc,
integer, dimension(nfunc), intent(in) ifunc,
integer, dimension(snpc), intent(in) npf,
type (ttable), dimension(ntable), intent(in) table,
dimension(stf), intent(in) tf,
intent(in) time,
dimension(nuparam), intent(in) uparam,
integer, dimension(nel), intent(in) ngl,
dimension(nel), intent(in) aldt,
dimension(nel), intent(in) dpla,
dimension(nel), intent(in) epsp,
dimension(nel,nuvar), intent(inout) uvar,
dimension(nel), intent(in) signxx,
dimension(nel), intent(in) signyy,
dimension(nel), intent(in) signxy,
dimension(nel), intent(in) signyz,
dimension(nel), intent(in) signzx,
dimension(nel), intent(in) temp,
integer, dimension(nel), intent(inout) foff,
dimension(nel), intent(inout) dfmax,
dimension(nel), intent(inout) tdele,
integer, intent(in) ipt,
integer, intent(in) ipg,
integer, intent(inout) dmg_flag,
dimension(nel), intent(inout) dmg_scale,
integer, intent(in) ntablf,
integer, dimension(ntablf), intent(in) itablf )

Definition at line 36 of file fail_tab2_c.F.

43C---------+---------+---+---+--------------------------------------------
44 USE table_mod
46C-----------------------------------------------
47C I m p l i c i t T y p e s
48C-----------------------------------------------
49#include "implicit_f.inc"
50C-----------------------------------------------
51C C o m m o n B l o c k s
52C-----------------------------------------------
53#include "units_c.inc"
54#include "comlock.inc"
55#include "com04_c.inc"
56#include "tabsiz_c.inc"
57C-----------------------------------------------
58C I N P U T A r g u m e n t s
59C-----------------------------------------------
60 INTEGER, INTENT(IN) ::
61 . NEL,NUPARAM,NUVAR,NTABLF,NGL(NEL),IPT,IPG
62 INTEGER, DIMENSION(NTABLF) ,INTENT(IN) :: ITABLF
63 INTEGER, INTENT(INOUT) ::
64 . DMG_FLAG,FOFF(NEL)
65 my_real, INTENT(IN) ::
66 . time,uparam(nuparam),aldt(nel),
67 . dpla(nel),epsp(nel),temp(nel),
68 . signxx(nel),signyy(nel),signxy(nel),
69 . signyz(nel),signzx(nel)
70 my_real, INTENT(INOUT) ::
71 . uvar(nel,nuvar),dfmax(nel),tdele(nel),
72 . dmg_scale(nel)
73 TYPE (TTABLE), INTENT(IN), DIMENSION(NTABLE) :: TABLE
74C!-----------------------------------------------
75C! VARIABLES FOR FUNCTION INTERPOLATION
76C!-----------------------------------------------
77 INTEGER, INTENT(IN) :: NPF(SNPC),NFUNC,IFUNC(NFUNC)
78 my_real, INTENT(IN) :: tf(stf)
80 . finter
81 EXTERNAL finter
82C-----------------------------------------------
83C L o c a l V a r i a b l e s
84C-----------------------------------------------
85 INTEGER I,J,INDX(NEL),NINDX,ITAB_EPSF,
86 . ITAB_INST,ITAB_SIZE,IREG,
87 . IPOS(NEL,3),NDIM,IPOS2(NEL),
88 . IAD(NEL),ILEN(NEL),LOG_SCALE1,
89 . LOG_SCALE2
90 my_real
91 . fcrit ,dn,dcrit,ecrit,exp_ref,expo,el_ref,
92 . sr_ref1,fscale_el,shrf,biaxf ,sr_ref2,
93 . fscale_sr,cjc,fscale_dlim,temp_ref, fscale_temp
95 . lambda,fac,df,inst(nel) ,dc(nel) ,l0(nel) ,
96 . triax(nel) ,xi(nel) ,epsf(nel) ,epsl(nel) ,
97 . depsf(nel) ,depsl(nel) ,xvec(nel,3) ,dpl_def ,
98 . cos3theta ,p ,svm ,sizefac(nel),
99 . ratefac(nel),dsize(nel) ,softexp(nel),dlim(nel) ,
100 . tempfac(nel),tempfac2(nel),dft(nel) ,var(nel)
101C!--------------------------------------------------------------
102 !=======================================================================
103 ! - INITIALISATION OF COMPUTATION ON TIME STEP
104 !=======================================================================
105 ! Recovering failure criterion parameters
106 fcrit = uparam(1)
107 dn = uparam(4)
108 dcrit = uparam(5)
109 ecrit = uparam(6)
110 exp_ref = uparam(7)
111 expo = uparam(8)
112 ireg = nint(uparam(9))
113 el_ref = uparam(10)
114 sr_ref1 = uparam(11)
115 fscale_el = uparam(12)
116 shrf = uparam(13)
117 biaxf = uparam(14)
118 sr_ref2 = uparam(15)
119 fscale_sr = uparam(16)
120 cjc = uparam(17)
121 fscale_dlim = uparam(18)
122 temp_ref = uparam(19)
123 fscale_temp = uparam(20)
124 log_scale1 = nint(uparam(21))
125 log_scale2 = nint(uparam(22))
126c
127 itab_epsf = itablf(1)
128 itab_inst = itablf(2)
129 itab_size = itablf(3)
130c
131 ! Set flag for stress softening
132 dmg_flag = 1
133c
134 ! Checking element failure and recovering user variable
135 DO i=1,nel
136 ! If necking control is activated
137 IF ((itab_inst > 0).OR.(ecrit > zero)) THEN
138 ! Instability damage
139 inst(i) = uvar(i,1)
140 ! Necking critical damage
141 IF (uvar(i,2) == zero) uvar(i,2) = one
142 dc(i) = uvar(i,2)
143 ENDIF
144 END DO
145c
146 !====================================================================
147 ! - LOOP OVER THE ELEMENT TO COMPUTE THE STRESS STATE QUANTITIES
148 !====================================================================
149 DO i=1,nel
150c
151 ! Computation of hydrostatic stress, Von Mises stress, and stress triaxiality
152 p = third*(signxx(i) + signyy(i))
153 ! von mises equivalent stress
154 svm = signxx(i)**2 + signyy(i)**2 - signxx(i)*signyy(i) +
155 . three*signxy(i)**2
156 svm = sqrt(max(svm,zero))
157 triax(i) = p/max(em20,svm)
158 IF (triax(i) < -two_third) triax(i) = -two_third
159 IF (triax(i) > two_third) triax(i) = two_third
160c
161 ! Computation of Lode parameter
162 cos3theta = -half*twenty7*triax(i)*(triax(i)**2 - third)
163 IF (cos3theta < -one ) cos3theta = -one
164 IF (cos3theta > one ) cos3theta = one
165 xi(i) = one - two*acos(cos3theta)/pi
166c
167 ENDDO
168c
169 !====================================================================
170 ! - COMPUTE FACTORS FOR ELEMENT SIZE, STRAIN RATE AND TEMPERATURE
171 !====================================================================
172 ! At initial time, save the element size
173 IF (uvar(1,3) == zero) uvar(1:nel,3) = aldt(1:nel)
174 l0(1:nel) = uvar(1:nel,3)
175c
176 ! Compute the softening exponent
177 IF (ifunc(1) > 0) THEN
178 DO i=1,nel
179 lambda = l0(i)/exp_ref
180 softexp(i) = finter(ifunc(1),lambda,npf,tf,df)
181 softexp(i) = expo*softexp(i)
182 ENDDO
183 ELSE
184 softexp(1:nel) = expo
185 ENDIF
186c
187 ! Compute the temperature dependency factor
188 IF (ifunc(4) > 0) THEN
189 var(1:nel) = temp(1:nel)/temp_ref
190 ipos2(1:nel) = 1
191 iad(1:nel) = npf(ifunc(4)) / 2 + 1
192 ilen(1:nel) = npf(ifunc(4)+1) / 2 - iad(1:nel) - ipos2(1:nel)
193 CALL vinter2(tf,iad,ipos2,ilen,nel,var,dft,tempfac)
194 tempfac(1:nel) = fscale_temp*tempfac(1:nel)
195 tempfac2(1:nel) = tempfac(1:nel)
196 ELSE
197 tempfac(1:nel) = one
198 tempfac2(1:nel) = one
199 ENDIF
200c
201 ! Compute the element size regularization factor
202 IF (itab_size > 0) THEN
203 ! Element size scaling dependency
204 ndim = table(itab_size)%NDIM
205 IF (ireg == 1) THEN
206 SELECT CASE (ndim)
207 ! Scale factor vs element size
208 CASE(1)
209 xvec(1:nel,1) = l0(1:nel)/el_ref
210 xvec(1:nel,2:3) = zero
211 ipos(1:nel,1:3) = 1
212 ! Scale factor vs element size vs strain rate
213 CASE(2)
214 xvec(1:nel,1) = l0(1:nel)/el_ref
215 IF (log_scale1 > 0) THEN
216 DO i = 1,nel
217 xvec(i,2) = log(max(epsp(i),em20)/sr_ref1)
218 ENDDO
219 ELSE
220 xvec(1:nel,2) = epsp(1:nel)/sr_ref1
221 ENDIF
222 xvec(1:nel,3) = zero
223 ipos(1:nel,1:3) = 1
224 END SELECT
225 ELSEIF (ireg == 2) THEN
226 SELECT CASE (ndim)
227 ! Scale factor vs element size
228 CASE(1)
229 xvec(1:nel,1) = l0(1:nel)/el_ref
230 xvec(1:nel,2:3) = zero
231 ipos(1:nel,1:3) = 1
232 ! Scale factor vs element size vs triaxiality
233 CASE(2)
234 xvec(1:nel,1) = l0(1:nel)/el_ref
235 xvec(1:nel,2) = triax(1:nel)
236 xvec(1:nel,3) = zero
237 ipos(1:nel,1:3) = 1
238 ! Scale factor vs element size vs triaxiality vs Lode parameter
239 CASE(3)
240 xvec(1:nel,1) = l0(1:nel)/el_ref
241 xvec(1:nel,2) = triax(1:nel)
242 xvec(1:nel,3) = xi(1:nel)
243 ipos(1:nel,1:3) = 1
244 END SELECT
245 ENDIF
246 CALL table_vinterp(table(itab_size),nel,nel,ipos,xvec,sizefac,dsize)
247 sizefac(1:nel) = sizefac(1:nel)*fscale_el
248 IF (ireg == 1) THEN
249 DO i = 1,nel
250 IF (triax(i) < shrf) THEN
251 sizefac(i) = one
252 ELSEIF (triax(i) > biaxf) THEN
253 sizefac(i) = one
254 ENDIF
255 ENDDO
256 ENDIF
257 ELSE
258 sizefac(1:nel) = one
259 ENDIF
260c
261 ! Compute the strain rate dependency factor
262 IF (ifunc(2) > 0) THEN
263 IF (log_scale2 > 0) THEN
264 DO i = 1,nel
265 var(i) = log(max(epsp(i),em20)/sr_ref2)
266 ENDDO
267 ELSE
268 var(1:nel) = epsp(1:nel)/sr_ref2
269 ENDIF
270 ipos2(1:nel) = 1
271 iad(1:nel) = npf(ifunc(2)) / 2 + 1
272 ilen(1:nel) = npf(ifunc(2)+1) / 2 - iad(1:nel) - ipos2(1:nel)
273 CALL vinter2(tf,iad,ipos2,ilen,nel,var,dft,ratefac)
274 ratefac(1:nel) = fscale_sr*ratefac(1:nel)
275 ELSEIF (cjc > zero) THEN
276 DO i=1,nel
277 IF (epsp(i) > sr_ref2) THEN
278 ratefac(i) = one + cjc*log(epsp(i)/sr_ref2)
279 ELSE
280 ratefac(i) = one
281 ENDIF
282 ENDDO
283 ELSE
284 ratefac(1:nel) = one
285 ENDIF
286c
287 ! Compute the damage limit value
288 IF (ifunc(3) > 0) THEN
289 DO i = 1,nel
290 lambda = triax(i)
291 dlim(i) = finter(ifunc(3),lambda,npf,tf,df)
292 dlim(i) = fscale_dlim*dlim(i)
293 dlim(i) = min(dlim(i),one)
294 dlim(i) = max(dlim(i),zero)
295 ENDDO
296 ELSE
297 dlim(1:nel) = one
298 ENDIF
299c
300 !====================================================================
301 ! - COMPUTATION OF PLASTIC STRAIN AT FAILURE
302 !====================================================================
303 IF (itab_epsf > 0) THEN
304 ! Failure plastic strain map dependency
305 ndim = table(itab_epsf)%NDIM
306 SELECT CASE (ndim)
307 ! Failure plastic strain vs triaxiality
308 CASE (1)
309 xvec(1:nel,1) = triax(1:nel)
310 xvec(1:nel,2:3) = zero
311 ipos(1:nel,1:3) = 1
312 ! Failure plastic strain vs triaxiality vs Lode parameter
313 CASE (2)
314 xvec(1:nel,1) = triax(1:nel)
315 xvec(1:nel,2) = xi(1:nel)
316 xvec(1:nel,3) = zero
317 ipos(1:nel,1:3) = 1
318 ! Failure plastic strain vs triaxiality vs Lode parameter vs temperature
319 CASE (3)
320 xvec(1:nel,1) = triax(1:nel)
321 xvec(1:nel,2) = xi(1:nel)
322 xvec(1:nel,3) = temp(1:nel)/temp_ref
323 ipos(1:nel,1:3) = 1
324 tempfac(1:nel) = one
325 END SELECT
326 CALL table_vinterp(table(itab_epsf),nel,nel,ipos,xvec,epsf,depsf)
327 epsf(1:nel) = epsf(1:nel)*fcrit
328 ELSE
329 epsf(1:nel) = fcrit
330 ENDIF
331c
332 !====================================================================
333 ! - COMPUTATION OF PLASTIC STRAIN AT NECKING
334 !====================================================================
335 IF (itab_inst > 0) THEN
336 ! Instability plastic strain map dependency
337 ndim = table(itab_inst)%NDIM
338 SELECT CASE (ndim)
339 ! Instability plastic strain vs triaxiality
340 CASE(1)
341 xvec(1:nel,1) = triax(1:nel)
342 xvec(1:nel,2:3) = zero
343 ipos(1:nel,1:3) = 1
344 ! Instability plastic strain vs triaxiality vs Lode
345 CASE(2)
346 xvec(1:nel,1) = triax(1:nel)
347 xvec(1:nel,2) = xi(1:nel)
348 xvec(1:nel,3) = zero
349 ipos(1:nel,1:3) = 1
350 ! Instability plastic strain vs triaxiality vs Lode vs temperature
351 CASE(3)
352 xvec(1:nel,1) = triax(1:nel)
353 xvec(1:nel,2) = xi(1:nel)
354 xvec(1:nel,3) = temp(1:nel)/temp_ref
355 ipos(1:nel,1:3) = 1
356 tempfac2(1:nel) = one
357 END SELECT
358 CALL table_vinterp(table(itab_inst),nel,nel,ipos,xvec,epsl,depsl)
359 epsl(1:nel) = epsl(1:nel)*ecrit
360 ELSEIF (ecrit > zero) THEN
361 epsl(1:nel) = ecrit
362 ENDIF
363c
364 !====================================================================
365 ! - COMPUTATION OF THE DAMAGE VARIABLE EVOLUTION
366 !====================================================================
367 ! Initialization of element failure index
368 nindx = 0
369 indx(1:nel) = 0
370c
371 ! Loop over the elements
372 DO i=1,nel
373c
374 ! If the element is not broken
375 IF (foff(i) /= 0 .AND. dpla(i) > zero) THEN
376c
377 ! Needs to initialize damage at a very small value the first time
378 IF (dfmax(i) == zero) dfmax(i) = em20
379 IF (inst(i) == zero) inst(i) = em20
380c
381 ! Compute failure strain damage variable
382 dpl_def = dpla(i)/max(epsf(i)*ratefac(i)*sizefac(i)*tempfac(i),em20)
383 dfmax(i) = dfmax(i) + dpl_def*dn*(dfmax(i)**(one-(one/dn)))
384 dfmax(i) = min(dfmax(i),dlim(i))
385 IF (dfmax(i) >= one) THEN
386 nindx = nindx + 1
387 indx(nindx) = i
388 foff(i) = 0
389 tdele(i) = time
390 ENDIF
391c
392 ! Compute the control necking instability damage
393 IF ((itab_inst > 0).OR.(ecrit > zero)) THEN
394 dpl_def = dpla(i)/max(epsl(i)*ratefac(i)*sizefac(i)*tempfac2(i),em20)
395 inst(i) = inst(i) + dpl_def*dn*(inst(i)**(one-(one/dn)))
396 inst(i) = min(inst(i),one)
397 IF ((inst(i) >= one).AND.(dc(i) == one)) THEN
398 dc(i) = dfmax(i)
399 ENDIF
400 ENDIF
401c
402 ENDIF
403 ENDDO
404c
405 !====================================================================
406 ! - UPDATE UVAR AND THE STRESS TENSOR
407 !====================================================================
408 DO i = 1,nel
409 IF ((itab_inst > 0).OR.(ecrit > zero)) THEN
410 uvar(i,1) = inst(i)
411 uvar(i,2) = dc(i)
412 IF (dfmax(i) >= dc(i)) THEN
413 IF (dc(i) < one) THEN
414 dmg_scale(i) = one - ((dfmax(i)-dc(i))/max(one-dc(i),em20))**softexp(i)
415 ELSE
416 dmg_scale(i) = zero
417 ENDIF
418 ELSE
419 dmg_scale(i) = one
420 ENDIF
421 ELSE
422 IF (dfmax(i) >= dcrit) THEN
423 IF (dcrit < one) THEN
424 dmg_scale(i) = one - ((dfmax(i)-dcrit)/max(one-dcrit,em20))**softexp(i)
425 ELSE
426 dmg_scale(i) = zero
427 ENDIF
428 ELSE
429 dmg_scale(i) = one
430 ENDIF
431 ENDIF
432 ENDDO
433c
434 !====================================================================
435 ! - PRINTOUT DATA ABOUT FAILED ELEMENTS
436 !====================================================================
437 IF (nindx > 0) THEN
438 DO j=1,nindx
439 i = indx(j)
440#include "lockon.inc"
441 WRITE(iout, 1000) ngl(i),ipg,ipt,time
442 WRITE(istdo,1000) ngl(i),ipg,ipt,time
443#include "lockoff.inc"
444 END DO
445 END IF
446c------------------------
447 1000 FORMAT(1x,'FOR SHELL ELEMENT NUMBER el#',i10,
448 . ' FAILURE (TAB2) AT GAUSS POINT ',i3,' LAYER ',i3,
449 . ' AT TIME :',1pe12.4)
450c
#define my_real
Definition cppsort.cpp:32
subroutine ecrit(timers, partsav, ms, v, in, r, dmas, weight, enintot, ekintot, a, ar, fxbipm, fxbrpm, monvol, xmom_sms, sensors, qfricint, ipari, weight_md, wfexth, iflag, ms_2d, multi_fvm, mas_nd, kend, h3d_data, dynain_data, usreint, output)
Definition ecrit.F:52
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
subroutine vinter2(tf, iad, ipos, ilen, nel0, x, dydx, y)
Definition vinter.F:143