OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
fail_tab2_c.F
Go to the documentation of this file.
1Copyright> OpenRadioss
2Copyright> Copyright (C) 1986-2025 Altair Engineering Inc.
3Copyright>
4Copyright> This program is free software: you can redistribute it and/or modify
5Copyright> it under the terms of the GNU Affero General Public License as published by
6Copyright> the Free Software Foundation, either version 3 of the License, or
7Copyright> (at your option) any later version.
8Copyright>
9Copyright> This program is distributed in the hope that it will be useful,
10Copyright> but WITHOUT ANY WARRANTY; without even the implied warranty of
11Copyright> MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12Copyright> GNU Affero General Public License for more details.
13Copyright>
14Copyright> You should have received a copy of the GNU Affero General Public License
15Copyright> along with this program. If not, see <https://www.gnu.org/licenses/>.
16Copyright>
17Copyright>
18Copyright> Commercial Alternative: Altair Radioss Software
19Copyright>
20Copyright> As an alternative to this open-source version, Altair also offers Altair Radioss
21Copyright> software under a commercial license. Contact Altair to discuss further if the
22Copyright> commercial version may interest you: https://www.altair.com/radioss/.
23!||====================================================================
24!|| fail_tab2_c ../engine/source/materials/fail/tabulated/fail_tab2_c.F
25!||--- called by ------------------------------------------------------
26!|| mulawc ../engine/source/materials/mat_share/mulawc.F90
27!|| usermat_shell ../engine/source/materials/mat_share/usermat_shell.f
28!||--- calls -----------------------------------------------------
29!|| finter ../engine/source/tools/curve/finter.F
30!|| table_vinterp ../engine/source/tools/curve/table_tools.F
31!|| vinter2 ../engine/source/tools/curve/vinter.F
32!||--- uses -----------------------------------------------------
33!|| interface_table_mod ../engine/share/modules/table_mod.F
34!|| table_mod ../engine/share/modules/table_mod.F
35!||====================================================================
36 SUBROUTINE fail_tab2_c (
37 1 NEL ,NUPARAM ,NUVAR ,NFUNC ,IFUNC ,
38 2 NPF ,TABLE ,TF ,TIME ,UPARAM ,
39 3 NGL ,ALDT ,DPLA ,EPSP ,UVAR ,
40 4 SIGNXX ,SIGNYY ,SIGNXY ,SIGNYZ ,SIGNZX ,
41 5 TEMP ,FOFF ,DFMAX ,TDELE ,IPT ,
42 6 IPG ,DMG_FLAG ,DMG_SCALE,NTABLF ,ITABLF )
43C---------+---------+---+---+--------------------------------------------
44 USE table_mod
45 USE interface_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)
79 my_real
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
94 my_real
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
451 END
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
fail_tab2_c
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)
Definition fail_tab2_c.F:43
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
usermat_shell
subroutine usermat_shell(timers, elbuf_str, mat_elem, jft, jlt, nel, pm, for, mom, gstr, thk, eint, off, dir_a, dir_b, mat, area, exx, eyy, exy, exz, eyz, kxx, kyy, kxy, geo, thk_ly, pid, tf, npf, mtn, dt1c, dm, bufmat, ssp, rho, viscmx, ipla, iofc, indx, ngl, thkly, matly, zcfac, ng, shf, gs, sigy, thk0, epsd_pg, posly, igeo, ipm, failwave, fwave_el, ifailure, aldt, tempel, die, r11, r12, r13, r21, r22, r23, r31, r32, r33, table, ixfem, elcrkini, dir1_crk, dir2_crk, iparg, jhbe, ismstr, jthe, tensx, ir, is, nlay, npt, ixlay, ixel, ithk, f_def, ishplyxfem, itask, pm_stack, isubstack, stack, alpe, ply_exx, ply_eyy, ply_exy, ply_exz, ply_eyz, ply_f, varnl, nloc_dmg, nlay_max, laynpt_max, dt)
Definition usermat_shell.F:116
vinter2
subroutine vinter2(tf, iad, ipos, ilen, nel0, x, dydx, y)
Definition vinter.F:143