OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
fail_orthstrain_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_orthstrain_c ../engine/source/materials/fail/orthstrain/fail_orthstrain_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!|| vinter2 ../engine/source/tools/curve/vinter.F
31!||====================================================================
32 SUBROUTINE fail_orthstrain_c(
33 1 NEL ,NUPARAM ,NUVAR ,UPARAM ,UVAR ,
34 2 NFUNC ,IFUNC ,NPF ,TF ,NGL ,
35 3 TIME ,TIMESTEP ,IPG ,ILAY ,IPT ,
36 4 EPSXX ,EPSYY ,EPSXY ,DMG_FLAG ,DMG_SCALE,
37 5 EPSPXX ,EPSPYY ,EPSPXY ,ALDT ,ISMSTR ,
38 6 SIGNXX ,SIGNYY ,SIGNXY ,LF_DAMMX ,
39 7 OFF ,OFFLY ,FOFF ,DFMAX ,TDEL )
40C-----------------------------------------------
41C Orthotropic strain failure model
42C-----------------------------------------------
43C I m p l i c i t T y p e s
44C-----------------------------------------------
45#include "implicit_f.inc"
46C-----------------------------------------------
47C G l o b a l P a r a m e t e r s
48C-----------------------------------------------
49#include "units_c.inc"
50#include "comlock.inc"
51C---------+---------+---+---+--------------------------------------------
52C VAR | SIZE |TYP| RW| DEFINITION
53C---------+---------+---+---+--------------------------------------------
54C NEL | 1 | I | R | SIZE OF THE ELEMENT GROUP NEL
55C NUPARAM | 1 | I | R | SIZE OF THE USER PARAMETER ARRAY
56C UPARAM | NUPARAM | F | R | USER MATERIAL PARAMETER ARRAY
57C NUVAR | 1 | I | R | NUMBER OF USER ELEMENT VARIABLES
58C UVAR |NEL*NUVAR| F |R/W| USER ELEMENT VARIABLE ARRAY
59C---------+---------+---+---+--------------------------------------------
60C NFUNC | 1 | I | R | NUMBER FUNCTION USED FOR THIS USER LAW not used
61C IFUNC | NFUNC | I | R | FUNCTION INDEX not used
62C NPF | * | I | R | FUNCTION ARRAY
63C TF | * | F | R | FUNCTION ARRAY
64C---------+---------+---+---+--------------------------------------------
65C EPSXX | NEL | F | R | STRAIN XX
66C EPSYY | NEL | F | R | STRAIN YY
67C ... | | | |
68C EPSPXX | NEL | F | R | STRAIN RATE XX
69C EPSPYY | NEL | F | R | STRAIN RATE YY
70C ... | | | |
71C SIGNXX | NEL | F | W | NEW ELASTO PLASTIC STRESS XX
72C SIGNYY | NEL | F | W | NEW ELASTO PLASTIC STRESS YY
73C ... | | | |
74C---------+---------+---+---+--------------------------------------------
75C OFF | NEL | F | R | DELETED ELEMENT FLAG (=1. ON, =0. OFF)
76C FOFF | NEL | I |R/W| DELETED INTEGRATION POINT FLAG (=1 ON, =0 OFF)
77C DFMAX | NEL | F |R/W| MAX DAMAGE FACTOR
78C TDEL | NEL | F | W | FAILURE TIME
79C---------+---------+---+---+--------------------------------------------
80C NGL ELEMENT ID
81C IPG CURRENT GAUSS POINT (in plane)
82C ILAY CURRENT LAYER
83C IPT CURRENT INTEGRATION POINT IN THE LAYER
84C-----------------------------------------------
85C I N P U T A r g u m e n t s
86C-----------------------------------------------
87 INTEGER ,INTENT(IN) :: NEL,NUPARAM,NUVAR,IPG,ILAY,IPT,ISMSTR
88 INTEGER ,INTENT(IN) :: LF_DAMMX
89 INTEGER ,DIMENSION(NEL) ,INTENT(IN) :: NGL
90 my_real ,INTENT(IN) :: TIME,TIMESTEP
91 my_real ,DIMENSION(NEL) ,INTENT(IN) :: OFF,EPSXX,EPSYY,EPSXY,
92 . EPSPXX,EPSPYY,EPSPXY,ALDT
93 my_real ,DIMENSION(NUPARAM) ,INTENT(IN) :: UPARAM
94 TARGET :: uvar
95C-----------------------------------------------
96C I N P U T O U T P U T A r g u m e n t s
97C-----------------------------------------------
98 INTEGER ,INTENT(OUT) :: DMG_FLAG
99 INTEGER ,DIMENSION(NEL) ,INTENT(INOUT) :: OFFLY,FOFF
100 my_real ,DIMENSION(NEL) ,INTENT(INOUT) :: signxx,signyy,signxy
101 my_real ,DIMENSION(NEL,LF_DAMMX), INTENT(INOUT) :: dfmax
102 my_real ,DIMENSION(NEL) ,INTENT(OUT) :: tdel,dmg_scale
103 my_real ,DIMENSION(NEL,NUVAR) ,INTENT(INOUT) :: uvar
104C-----------------------------------------------
105C VARIABLES FOR FUNCTION INTERPOLATION
106C-----------------------------------------------
107 INTEGER NPF(*), NFUNC, IFUNC(NFUNC)
108 my_real finter ,tf(*)
109 EXTERNAL finter
110C-----------------------------------------------
111C Y = FINTER(IFUNC(J),X,NPF,TF,DYDX)
112C Y : y = f(x)
113C X : x
114C DYDX : f'(x) = dy/dx
115C IFUNC(J): FUNCTION INDEX
116C J : FIRST(J=1), SECOND(J=2) .. FUNCTION USED FOR THIS LAW
117C NPF,TF : FUNCTION PARAMETER
118C-----------------------------------------------
119C L o c a l V a r i a b l e s
120C-----------------------------------------------
121 INTEGER I,J,NINDX,FAILI,MODE,XX,XY,YY,IFUNC_SIZ,STRDEF,STRFLAG
122 INTEGER ,DIMENSION(NEL) :: INDX,IPOSP,IADP,ILENP
123 INTEGER ,DIMENSION(NEL,6) :: FMODE
124 my_real ODAM(NEL,6),EPSP(NEL,6),
125 . SIGOXX(NEL),SIGOYY(NEL),SIGOXY(NEL),ESCAL(NEL),DYDXE(NEL)
126 my_real DYDX,FRATE,ET1,ET2,ET12,EC1,EC2,EC12,ET3,EC3,EC23,ET23,EC31,
127 . ET31,ET1M,ET2M,ET12M,EC1M,EC2M,EC12M,ET3M,EC3M,EC23M,ET23M,EC31M,
128 . et31m,di,alpha,damxx,damyy,damxy,odamxx,odamyy,odamxy,theta,c,s,
129 . epspref,epdot,deps,depsm,fact,fscale_siz,ref_siz,epsd1,epsd2,
130 . epst1,epst2,epst_a,epst_b
131 my_real ,DIMENSION(NEL) :: eps11,eps22,eps12,epsp11,epsp22,epsp12,
132 . epsm1,epsm2,epsm3,epsm4,epsm5,epsm6
133 my_real ,DIMENSION(:), POINTER :: fsize
134C=======================================================================
135C USER VARIABLES INITIALIZATION
136C-----------------------------------------------
137 dmg_flag = int(uparam(1))
138 alpha = min(two*pi*uparam(16)*timestep,one)
139 epspref = uparam(17)
140 fscale_siz = uparam(26)
141 ref_siz = uparam(27)
142 strdef = nint(uparam(28))
143 ifunc_siz = ifunc(13)
144c
145c et1 = mode1 , et2 = mode 2, ec1 = mode 3, ec2 = mode 4,
146c et12 = mode 5, ec12 = mode 6
147c xx = dir 1, yy = dir 2, xy = dir 3
148c xx => et1, ec1
149c yy => et2, ec2
150c xy => et12, ec12
151c
152c----------------------------------------------
153c strain transformation following input definition
154c-------------------
155 strflag = 0
156 IF (strdef == 2) THEN ! failure defined as engineering strain
157 IF (ismstr /= 1 .AND. ismstr /= 3 .AND. ismstr /= 11) THEN
158 strflag = 2
159 END IF
160 ELSE IF (strdef == 3) THEN ! failure defined as true strain
161 IF (ismstr == 1 .OR. ismstr == 3 .OR. ismstr == 11) THEN
162 strflag = 3
163 END IF
164 END IF
165c--------------------------
166 SELECT CASE (strflag)
167c
168 CASE (2) ! transform true strain to engineering
169 DO i=1,nel
170 IF (off(i) == one ) THEN
171 theta = half*atan(epsxy(i) / max((epsxx(i)-epsyy(i)),em20))
172 c = cos(theta)
173 s = sin(theta)
174c strain principal + transformation
175 epst_a = half*(epsxx(i)+epsyy(i))
176 epst_b = sqrt( (half*(epsxx(i)-epsyy(i)))**2 + (half*epsxy(i))**2)
177 epst1 = epst_a + epst_b
178 epst2 = epst_a - epst_b
179 epst1 = exp(epst1) - one
180 epst2 = exp(epst2) - one
181 eps11(i) = c*c*epst1 + s*s*epst2
182 eps22(i) = s*s*epst1 + c*c*epst2
183 eps12(i) = s*c*(epst2 - epst1)
184c strain rate
185 epst_a = half*(epspxx(i)+epspyy(i))
186 epst_b = sqrt( (half*(epspxx(i)-epspyy(i)))**2 + (half*epspxy(i))**2)
187 epsd1 = epst_a + epst_b
188 epsd2 = epst_a - epst_b
189 epsd1 = (one+epst1)*epsd1
190 epsd2 = (one+epst2)*epsd2
191 epsp11(i) = c*c*epsd1 + s*s*epsd2
192 epsp22(i) = s*s*epsd1 + c*c*epsd2
193 epsp12(i) = s*c*(epsd2 - epsd1)
194 END IF
195 ENDDO
196c
197 CASE (3) ! transform engineering to true strain
198 DO i=1,nel
199 IF (off(i) == one ) THEN
200 theta = half*atan(epsxy(i) / max((epsxx(i)-epsyy(i)),em20))
201 c = cos(theta)
202 s = sin(theta)
203c strain principal + transformation
204 epst_a = half*(epsxx(i)+epsyy(i))
205 epst_b = sqrt( (half*(epsxx(i)-epsyy(i)))**2 + (half*epsxy(i))**2)
206 epst1 = epst_a + epst_b
207 epst2 = epst_a - epst_b
208 epst1 = log(one + epst1)
209 epst2 = log(one + epst2)
210 eps11(i) = c*c*epst1 + s*s*epst2
211 eps22(i) = s*s*epst1 + c*c*epst2
212 eps12(i) = s*c*(epst2 - epst1)
213c strain rate
214 epst_a = half*(epspxx(i)+epspyy(i))
215 epst_b = sqrt( (half*(epspxx(i)-epspyy(i)))**2 + (half*epspxy(i))**2)
216 epsd1 = epst_a + epst_b
217 epsd2 = epst_a - epst_b
218 epsd1 = (one/exp(epst1))*epsd1
219 epsd2 = (one/exp(epst2))*epsd2
220 epsp11(i) = c*c*epsd1 + s*s*epsd2
221 epsp22(i) = s*s*epsd1 + c*c*epsd2
222 epsp12(i) = s*c*(epsd2 - epsd1)
223 END IF
224 ENDDO
225c
226 CASE DEFAULT
227 ! no transformation : failure strain measure is defined by Ismstr
228 eps11(1:nel) = epsxx(1:nel)
229 eps22(1:nel) = epsyy(1:nel)
230 eps12(1:nel) = half*epsxy(1:nel)
231 epsp11(1:nel) = epspxx(1:nel)
232 epsp22(1:nel) = epspyy(1:nel)
233 epsp12(1:nel) = half*epspxy(1:nel)
234 END SELECT
235c----------------------------------------------
236c Element size scale factor
237 fsize => uvar(1:nel,13)
238c
239 IF (fsize(1)==zero) THEN
240 IF (ifunc_siz > 0) THEN
241 escal(1:nel) = aldt(1:nel) / ref_siz
242 iposp(1:nel) = 0
243 iadp(1:nel) = npf(ifunc_siz) / 2 + 1
244 ilenp(1:nel) = npf(ifunc_siz + 1) / 2 - iadp(1:nel)
245 CALL vinter2(tf,iadp,iposp,ilenp,nel,escal,dydxe,fsize)
246 fsize(1:nel) = fscale_siz * fsize(1:nel)
247 ELSE
248 fsize(1:nel) = one
249 ENDIF
250 ENDIF
251c
252 xx=1
253 yy=2
254 xy=3
255c variable initialisation, current time
256 DO j=1 ,6
257 DO i=1, nel
258 odam(i,j) = min(dfmax(i,1+j),one-em3)
259 epsp(i,j) = uvar(i,j)
260 ENDDO
261 ENDDO
262 DO i =1, nel
263 sigoxx(i) = uvar(i,6+xx)
264 sigoyy(i) = uvar(i,6+yy)
265 sigoxy(i) = uvar(i,6+xy)
266 ENDDO
267c
268c strain rate filtering
269 DO i=1,nel
270 epsp(i,xx) = alpha * abs(epsp11(i)) + (one-alpha)*epsp(i,xx)
271 epsp(i,yy) = alpha * abs(epsp22(i)) + (one-alpha)*epsp(i,yy)
272 epsp(i,xy) = alpha * abs(epsp12(i)) + (one-alpha)*epsp(i,xy)
273 ENDDO
274c
275c failure and damage calculation in each mode (dir)
276c
277 nindx = 0
278 DO i=1, nel
279 IF (off(i) == one .AND. offly(i) == 1 .AND. foff(i) == 1) THEN
280 et1 = uparam(4)
281 et1m = uparam(5)
282 et2 = uparam(6)
283 et2m = uparam(7)
284 ec1 = uparam(8)
285 ec1m = uparam(9)
286 ec2 = uparam(10)
287 ec2m = uparam(11)
288 et12 = uparam(12)
289 et12m = uparam(13)
290 ec12 = uparam(14)
291 ec12m = uparam(15)
292 et3 = uparam(18)
293 et3m = uparam(19)
294 ec3 = uparam(20)
295 ec3m = uparam(21)
296 et23 = uparam(22)
297 et23m = uparam(23)
298 ec23 = uparam(24)
299 ec23m = uparam(25)
300 et31 = uparam(22)
301 et31m = uparam(23)
302 ec31 = uparam(24)
303 ec31m = uparam(25)
304c
305 fmode(i,1:6) = 0
306 faili = 0
307c
308 mode = 1 ! traction xx
309 epdot = abs(epsp(i,xx))/epspref
310 IF (ifunc(mode) > 0) THEN
311 frate = finter(ifunc(mode),epdot,npf,tf,dydx)
312 ELSE
313 frate = one
314 ENDIF
315 frate = frate*fsize(i)
316 et1m = frate * et1m
317 et1 = frate * et1
318 deps = eps11(i) - et1
319 IF (deps > zero .AND. dfmax(i,1+mode) < one) THEN
320 depsm = et1m - et1
321 fact = et1m / eps11(i)
322 di = fact * deps / depsm
323 dfmax(i,1+mode) = max(dfmax(i,1+mode), di)
324 IF (dfmax(i,1+mode) >= one) THEN
325 faili = 1
326 fmode(i,mode) = 1
327 dfmax(i,1+mode) = one
328 epsm1(i) = et1m
329 ENDIF
330 ENDIF
331c
332 mode = 2 ! traction yy
333 epdot = abs(epsp(i,yy))/epspref
334 IF (ifunc(mode) > 0) THEN
335 frate = finter(ifunc(mode),epdot,npf,tf,dydx)
336 ELSE
337 frate = one
338 ENDIF
339 frate = frate*fsize(i)
340 et2m = frate * et2m
341 et2 = frate * et2
342 deps = max(eps22(i) - et2, zero)
343 IF (deps > zero .AND. dfmax(i,1+mode) < one) THEN
344 depsm = (et2m - et2)
345 fact = et2m / eps22(i)
346 di = fact * deps / depsm
347 dfmax(i,1+mode) = max(dfmax(i,1+mode), di)
348 IF (dfmax(i,1+mode) >= one) THEN
349 faili = 1
350 fmode(i,mode) = 1
351 dfmax(i,1+mode) = one
352 epsm2(i) = et2m
353 ENDIF
354 ENDIF
355c
356 mode = 3 ! traction xy
357 epdot = abs(epsp(i,xy))/epspref
358 IF (ifunc(mode) > 0) THEN
359 frate = finter(ifunc(mode),epdot,npf,tf,dydx)
360 ELSE
361 frate = one
362 ENDIF
363 frate = frate*fsize(i)
364 et12m = frate * et12m
365 et12 = frate * et12
366 deps = max(eps12(i) - et12, zero)
367 IF (deps > zero .AND. dfmax(i,1+mode) < one) THEN
368 depsm = (et12m - et12)
369 fact = et12m / eps12(i)
370 di = fact * deps / depsm
371 dfmax(i,1+mode) = max(dfmax(i,1+mode), di)
372 IF (dfmax(i,1+mode) >= one) THEN
373 faili = 1
374 fmode(i,mode) = 1
375 dfmax(i,1+mode) = one
376 epsm3(i) = et12m
377 ENDIF
378 ENDIF
379c
380 mode = 4 ! compression xx
381 epdot = abs(epsp(i,xx))/epspref
382 IF (ifunc(mode) > 0) THEN
383 frate = finter(ifunc(mode),epdot,npf,tf,dydx)
384 ELSE
385 frate = one
386 ENDIF
387 frate = frate*fsize(i)
388 ec1m = frate * ec1m
389 ec1 = frate * ec1
390 deps = -eps11(i) - ec1
391 IF (deps > zero .AND. dfmax(i,1+mode) < one) THEN
392 depsm = (ec1m - ec1)
393 fact = ec1m / abs(eps11(i))
394 di = fact * deps / depsm
395 dfmax(i,1+mode) = max(dfmax(i,1+mode), di)
396 IF (dfmax(i,1+mode) >= one) THEN
397 faili = 1
398 fmode(i,mode) = 1
399 dfmax(i,1+mode) = one
400 epsm4(i) = -ec1m
401 ENDIF
402 ENDIF
403c
404 mode = 5 ! compression yy
405 epdot = abs(epsp(i,yy))/epspref
406 IF (ifunc(mode) > 0) THEN
407 frate = finter(ifunc(mode),epdot,npf,tf,dydx)
408 ELSE
409 frate = one
410 ENDIF
411 frate = frate*fsize(i)
412 ec2m = frate * ec2m
413 ec2 = frate * ec2
414 deps = max(-eps22(i) - ec2, zero)
415 IF (deps > zero .AND. dfmax(i,1+mode) < one) THEN
416 depsm = (ec2m - ec2)
417 fact = ec2m / abs(eps22(i))
418 di = fact * deps / depsm
419 dfmax(i,1+mode) = max(dfmax(i,1+mode), di)
420 IF (dfmax(i,1+mode) >= one) THEN
421 faili = 1
422 fmode(i,mode) = 1
423 dfmax(i,1+mode) = one
424 epsm5(i) = -ec2m
425 ENDIF
426 ENDIF
427c
428 mode = 6 ! compression xy
429 epdot = abs(epsp(i,xy))/epspref
430 IF (ifunc(mode) > 0) THEN
431 frate = finter(ifunc(mode),epdot,npf,tf,dydx)
432 ELSE
433 frate = one
434 ENDIF
435 frate = frate*fsize(i)
436 ec12m = frate * ec12m
437 ec12 = frate * ec12
438 deps = max(-eps12(i) - ec12, zero)
439 IF (deps > zero .AND. dfmax(i,1+mode) < one) THEN
440 depsm = (ec12m - ec12)
441 fact = ec12m / abs(eps12(i))
442 di = fact * deps / depsm
443 dfmax(i,1+mode) = max(dfmax(i,1+mode), di)
444 IF (dfmax(i,1+mode) >= one) THEN
445 faili = 1
446 fmode(i,mode) = 1
447 dfmax(i,1+mode) = one
448 epsm6(i) = -ec12m
449 ENDIF
450 ENDIF
451c-----
452 IF (faili == 1) THEN
453 foff(i) = 0
454 tdel(i) = time
455 nindx = nindx + 1
456 indx(nindx) = i
457 ENDIF
458 ENDIF
459 ENDDO
460c------------------------
461 IF (nindx > 0) THEN
462 DO j=1,nindx
463 i = indx(j)
464#include "lockon.inc"
465 WRITE(iout, 1000) ngl(i),ipg,ilay,ipt
466 WRITE(istdo,1000) ngl(i),ipg,ilay,ipt
467 IF (fmode(i,1)==1) WRITE(iout, 2000) '1 - TRACTION XX',eps11(i),epsm1(i),epsp(i,xx)
468 IF (fmode(i,2)==1) WRITE(iout, 2000) '2 - TRACTION YY',eps22(i),epsm2(i),epsp(i,yy)
469 IF (fmode(i,3)==1) WRITE(iout, 2000) '3 - POSITIVE SHEAR XY',eps12(i),epsm3(i),epsp(i,xy)
470 IF (fmode(i,4)==1) WRITE(iout, 2000) '4 - COMPRESSION XX',eps11(i),epsm4(i),epsp(i,xx)
471 IF (fmode(i,5)==1) WRITE(iout, 2000) '5 - COMPRESSION YY',eps22(i),epsm5(i),epsp(i,yy)
472 IF (fmode(i,6)==1) WRITE(iout, 2000) '6 - NEGATIVE SHEAR XY',eps12(i),epsm6(i),epsp(i,xy)
473#include "lockoff.inc"
474 END DO
475 END IF
476c-----------------------------------------------------------------------
477c Damage stress reduction
478c damage / direction is equal to the max of damages par mode in the same direction
479c-----------------------------------------------------------------------
480 IF (dmg_flag == 1) THEN ! scalar damage applied to nodal force contribution
481 DO i=1, nel
482 damxx = max(dfmax(i,2),dfmax(i,5))
483 damyy = max(dfmax(i,3),dfmax(i,6))
484 damxy = max(dfmax(i,4),dfmax(i,7))
485 dfmax(i,1) = max(dfmax(i,1), damxx)
486 dfmax(i,1) = max(dfmax(i,1), damyy)
487 dfmax(i,1) = max(dfmax(i,1), damxy)
488 dfmax(i,1) = min(one,dfmax(i,1))
489 dmg_scale(i) = one - dfmax(i,1)
490 ENDDO
491 ELSE ! directional damage applied to stress
492 DO i=1, nel
493 damxx = max(dfmax(i,2),dfmax(i,5))
494 damyy = max(dfmax(i,3),dfmax(i,6))
495 damxy = max(dfmax(i,4),dfmax(i,7))
496 dfmax(i,1) = max(dfmax(i,1), damxx)
497 dfmax(i,1) = max(dfmax(i,1), damyy)
498 dfmax(i,1) = max(dfmax(i,1), damxy)
499 dfmax(i,1) = min(one,dfmax(i,1))
500 odamxx = max(odam(i,1),odam(i,4))
501 odamyy = max(odam(i,2),odam(i,5))
502 odamxy = max(odam(i,3),odam(i,6))
503c Undamaged stress estimate
504 signxx(i) = sigoxx(i)/(one-odamxx) + (signxx(i) - sigoxx(i))
505 signyy(i) = sigoyy(i)/(one-odamyy) + (signyy(i) - sigoyy(i))
506 signxy(i) = sigoxy(i)/(one-odamxy) + (signxy(i) - sigoxy(i))
507c Damaged stress
508 signxx(i) = signxx(i)* (one-damxx)
509 signyy(i) = signyy(i)* (one-damyy)
510 signxy(i) = signxy(i)* (one-damxy)
511 uvar(i,7) = signxx(i)
512 uvar(i,8) = signyy(i)
513 uvar(i,9) = signxy(i)
514 ENDDO
515 ENDIF
516c Update state variables
517 DO j=1 ,6
518 DO i=1, nel
519 uvar(i,j) = epsp(i,j)
520 ENDDO
521 ENDDO
522c-----------------------------------------------------------------------
523 1000 FORMAT(1x,'FAILURE (ORTHSTR) OF SHELL ELEMENT ',i10,1x,',GAUSS PT',
524 . i2,1x,',LAYER',i3,1x,',INTEGRATION PT',i3)
525 2000 FORMAT(1x,'MODE',1x,a,', STRAIN',1pe12.4,', CRITICAL VALUE',1pe12.4,
526 . ', STRAIN RATE',1pe12.4)
527c-----------------------------------------------------------------------
528 RETURN
529 END
alpha
#define alpha
Definition eval.h:35
fail_orthstrain_c
subroutine fail_orthstrain_c(nel, nuparam, nuvar, uparam, uvar, nfunc, ifunc, npf, tf, ngl, time, timestep, ipg, ilay, ipt, epsxx, epsyy, epsxy, dmg_flag, dmg_scale, epspxx, epspyy, epspxy, aldt, ismstr, signxx, signyy, signxy, lf_dammx, off, offly, foff, dfmax, tdel)
Definition fail_orthstrain_c.F:40
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
vinter2
subroutine vinter2(tf, iad, ipos, ilen, nel0, x, dydx, y)
Definition vinter.F:143