OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
fail_tab2_s.F File Reference
#include "implicit_f.inc"
#include "scr17_c.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_s (nel, nuparam, nuvar, nfunc, ifunc, npf, table, tf, time, uparam, ngl, aldt, dpla, epsp, uvar, signxx, signyy, signzz, signxy, signyz, signzx, temp, off, dfmax, tdele, dmg_scale, uelr, ipg, npg, loff, ntablf, itablf, noff, voln)

Function/Subroutine Documentation

◆ fail_tab2_s()

subroutine fail_tab2_s ( 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) signzz,
dimension(nel), intent(in) signxy,
dimension(nel), intent(in) signyz,
dimension(nel), intent(in) signzx,
dimension(nel), intent(in) temp,
dimension(nel), intent(inout) off,
dimension(nel), intent(inout) dfmax,
dimension(nel), intent(inout) tdele,
dimension(nel), intent(inout) dmg_scale,
dimension(nel), intent(inout) uelr,
integer, intent(in) ipg,
integer, intent(in) npg,
dimension(nel), intent(inout) loff,
integer, intent(in) ntablf,
integer, dimension(ntablf), intent(in) itablf,
integer, dimension(nel), intent(inout) noff,
dimension(nel), intent(in) voln )

Definition at line 37 of file fail_tab2_s.F.

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