OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
sigeps53.F File Reference
#include "implicit_f.inc"
#include "mvsiz_p.inc"
#include "param_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine sigeps53 (nel, nuparam, nuvar, mfunc, kfunc, npf, tf, time, timestep, uparam, rho0, rho, volume, eint, epspxx, epspyy, epspzz, epspxy, epspyz, epspzx, depsxx, depsyy, depszz, depsxy, depsyz, depszx, epsxx, epsyy, epszz, epsxy, epsyz, epszx, sigoxx, sigoyy, sigozz, sigoxy, sigoyz, sigozx, signxx, signyy, signzz, signxy, signyz, signzx, sigvxx, sigvyy, sigvzz, sigvxy, sigvyz, sigvzx, soundsp, viscmax, uvar, off, ngl, ipt, ipm, mat, epsp, ipla, seq_output)

Function/Subroutine Documentation

◆ sigeps53()

subroutine sigeps53 ( integer nel,
integer nuparam,
integer nuvar,
integer mfunc,
integer, dimension(mfunc) kfunc,
integer, dimension(*) npf,
tf,
time,
timestep,
uparam,
rho0,
rho,
volume,
eint,
epspxx,
epspyy,
epspzz,
epspxy,
epspyz,
epspzx,
depsxx,
depsyy,
depszz,
depsxy,
depsyz,
depszx,
epsxx,
epsyy,
epszz,
epsxy,
epsyz,
epszx,
sigoxx,
sigoyy,
sigozz,
sigoxy,
sigoyz,
sigozx,
signxx,
signyy,
signzz,
signxy,
signyz,
signzx,
sigvxx,
sigvyy,
sigvzz,
sigvxy,
sigvyz,
sigvzx,
soundsp,
viscmax,
uvar,
off,
integer, dimension(nel) ngl,
integer ipt,
integer, dimension(npropmi,*) ipm,
integer, dimension(nel) mat,
epsp,
integer ipla,
seq_output )

Definition at line 36 of file sigeps53.F.

48C-----------------------------------------------
49C M o d u l e s
50C-----------------------------------------------
51 USE message_mod
52C-----------------------------------------------
53C I m p l i c i t T y p e s
54C-----------------------------------------------
55#include "implicit_f.inc"
56C-----------------------------------------------
57C G l o b a l P a r a m e t e r s
58C-----------------------------------------------
59#include "mvsiz_p.inc"
60C---------+---------+---+---+--------------------------------------------
61C VAR | SIZE |TYP| RW| DEFINITION
62C---------+---------+---+---+--------------------------------------------
63C NEL | 1 | I | R | SIZE OF THE ELEMENT GROUP NEL
64C NUPARAM | 1 | I | R | SIZE OF THE USER PARAMETER ARRAY
65C NUVAR | 1 | I | R | NUMBER OF USER ELEMENT VARIABLES
66C---------+---------+---+---+--------------------------------------------
67C MFUNC | 1 | I | R | NUMBER FUNCTION USED FOR THIS USER LAW not used
68C KFUNC | NFUNC | I | R | FUNCTION INDEX not used
69C NPF | * | I | R | FUNCTION ARRAY
70C TF | * | F | R | FUNCTION ARRAY
71C---------+---------+---+---+--------------------------------------------
72C TIME | 1 | F | R | CURRENT TIME
73C TIMESTEP| 1 | F | R | CURRENT TIME STEP
74C UPARAM | NUPARAM | F | R | USER MATERIAL PARAMETER ARRAY
75C RHO0 | NEL | F | R | INITIAL DENSITY
76C RHO | NEL | F | R | DENSITY
77C VOLUME | NEL | F | R | VOLUME
78C EINT | NEL | F | R | TOTAL INTERNAL ENERGY
79C EPSPXX | NEL | F | R | STRAIN RATE XX
80C EPSPYY | NEL | F | R | STRAIN RATE YY
81C ... | | | |
82C DEPSXX | NEL | F | R | STRAIN INCREMENT XX
83C DEPSYY | NEL | F | R | STRAIN INCREMENT YY
84C ... | | | |
85C EPSXX | NEL | F | R | STRAIN XX
86C EPSYY | NEL | F | R | STRAIN YY
87C ... | | | |
88C SIGOXX | NEL | F | R | OLD ELASTO PLASTIC STRESS XX
89C SIGOYY | NEL | F | R | OLD ELASTO PLASTIC STRESS YY
90C ... | | | |
91C---------+---------+---+---+--------------------------------------------
92C SIGNXX | NEL | F | W | NEW ELASTO PLASTIC STRESS XX
93C SIGNYY | NEL | F | W | NEW ELASTO PLASTIC STRESS YY
94C ... | | | |
95C SIGVXX | NEL | F | W | VISCOUS STRESS XX
96C SIGVYY | NEL | F | W | VISCOUS STRESS YY
97C ... | | | |
98C SOUNDSP | NEL | F | W | SOUND SPEED (NEEDED FOR TIME STEP)
99C VISCMAX | NEL | F | W | MAXIMUM DAMPING MODULUS(NEEDED FOR TIME STEP)
100C---------+---------+---+---+--------------------------------------------
101C UVAR |NEL*NUVAR| F |R/W| USER ELEMENT VARIABLE ARRAY
102C OFF | NEL | F |R/W| DELETED ELEMENT FLAG (=1. ON, =0. OFF)
103C---------+---------+---+---+--------------------------------------------
104#include "param_c.inc"
105 INTEGER NEL, NUPARAM, NUVAR,IPT,IPLA
106 INTEGER NGL(NEL),MAT(NEL),IPM(NPROPMI,*)
107 my_real
108 . time,timestep,uparam(*),
109 . rho(nel),rho0(nel),volume(nel),eint(nel),
110 . epspxx(nel),epspyy(nel),epspzz(nel),
111 . epspxy(nel),epspyz(nel),epspzx(nel),
112 . depsxx(nel),depsyy(nel),depszz(nel),
113 . depsxy(nel),depsyz(nel),depszx(nel),
114 . epsxx(nel) ,epsyy(nel) ,epszz(nel) ,
115 . epsxy(nel) ,epsyz(nel) ,epszx(nel) ,
116 . sigoxx(nel),sigoyy(nel),sigozz(nel),
117 . sigoxy(nel),sigoyz(nel),sigozx(nel),
118 . epsp(nel),seq_output(nel)
119C-----------------------------------------------
120C O U T P U T A r g u m e n t s
121C-----------------------------------------------
122 my_real
123 . signxx(nel),signyy(nel),signzz(nel),
124 . signxy(nel),signyz(nel),signzx(nel),
125 . sigvxx(nel),sigvyy(nel),sigvzz(nel),
126 . sigvxy(nel),sigvyz(nel),sigvzx(nel),
127 . soundsp(nel),viscmax(nel)
128C-----------------------------------------------
129C I N P U T O U T P U T A r g u m e n t s
130C-----------------------------------------------
131 my_real
132 . uvar(nel,nuvar), off(nel)
133C-----------------------------------------------
134C VARIABLES FOR FUNCTION INTERPOLATION
135C-----------------------------------------------
136 INTEGER NPF(*), MFUNC, KFUNC(MFUNC)
137 my_real
138 . finter2, tf(*)
139 EXTERNAL finter2
140C
141C Y = FINTER(IFUNC(J),X,NPF,TF,DYDX)
142C Y : y = f(x)
143C X : x
144C DYDX : f'(x) = dy/dx
145C IFUNC(J): FUNCTION INDEX
146C J : FIRST(J=1), SECOND(J=2) .. FUNCTION USED FOR THIS LAW
147C NPF,TF : FUNCTION PARAMETER
148C-----------------------------------------------
149C L o c a l V a r i a b l e s
150C-----------------------------------------------
151 INTEGER I,J,IADBUFV(MVSIZ),J1,J2,JJ(MVSIZ),NFUNC,
152 . NRATE(MVSIZ),IPOS1(MVSIZ),IPOS2(MVSIZ),IAD1(MVSIZ),
153 . ILEN1(MVSIZ),IAD2(MVSIZ),ILEN2(MVSIZ),
154 . IFUNC(MVSIZ,100),NFUNCV(MVSIZ),PFUN(MVSIZ),
155 . IPOSP(MVSIZ),IADP(MVSIZ),ILENP(MVSIZ),NFUNCM,NRATEM,
156 . IPFLAG,IPARAM,NPAR, IFLAG(MVSIZ)
157 my_real
158 . rate(mvsiz,2),yfac(mvsiz,2),
159 . y1(mvsiz),y2(mvsiz),h(mvsiz),dydx1(mvsiz),
160 . dydx2(mvsiz),pla(mvsiz),fail(mvsiz),epsr1(mvsiz),
161 . p0(mvsiz),pfac(mvsiz),
162 . dfdp(mvsiz),dpla,
163 . evl(mvsiz),
164 . e11(mvsiz), e22(mvsiz), g12(mvsiz),g23(mvsiz),
165 . f1(mvsiz), f2(mvsiz), f11(mvsiz), f22(mvsiz), f12(mvsiz),
166 . f44(mvsiz), f55(mvsiz), f23(mvsiz),f(mvsiz), s1c(mvsiz),
167 . s2c(mvsiz),s2t(mvsiz), s3c(mvsiz),s3t(mvsiz), s4t(mvsiz),
168 . s4c(mvsiz), s45c(mvsiz),s45t(mvsiz),scale,bb,aa,dd,cc,ss1,
169 . ss2,s1t(mvsiz), evol(mvsiz)
170
171C-----------------------------------------------
172C USER VARIABLES INITIALIZATION
173C-----------------------------------------------
174C
175 DO i=1,nel
176 iadbufv(i) = ipm(7,mat(i))-1
177 e11(i) = uparam(iadbufv(i)+1)
178 e22(i) = uparam(iadbufv(i)+2)
179 g12(i) = uparam(iadbufv(i)+3)
180 g23(i) = uparam(iadbufv(i)+4)
181 iflag(i) = nint(uparam(iadbufv(i)+5))
182C
183 nfunc = ipm(10,mat(i))
184 DO j=1,nfunc
185 ifunc(i,j) = ipm(10+j,mat(i))
186 ENDDO
187 ENDDO
188C
189 IF(time==zero)THEN
190 DO i=1,nel
191 uvar(i,1)=zero
192 DO j=1,nfunc
193 uvar(i,j+1)=zero
194 ENDDO
195 ENDDO
196 ENDIF
197C-----------------------------------------------
198C
199 DO i=1,nel
200 pla(i) = uvar(i,1)
201 pla(i) = pla(i) + depsxx(i) + depsyy(i) + depszz(i)
202 signxx(i)=sigoxx(i)+ e11(i) * depsxx(i)
203 signyy(i)=sigoyy(i)+ e22(i) * depsyy(i)
204 signzz(i)=sigozz(i)+ e22(i) * depszz(i)
205 signxy(i)=sigoxy(i)+ g12(i) * depsxy(i)
206 signyz(i)=sigoyz(i)+ g23(i) * depsyz(i)
207 signzx(i)=sigozx(i)+ g12(i) * depszx(i)
208C
209 soundsp(i) = sqrt(max(e11(i),e22(i),g12(i),g23(i))/rho0(i))
210 evol(i) =one - exp(pla(i))
211 viscmax(i) = zero
212 uvar(i,1)= pla(i)
213 ENDDO
214C-------------------
215C CRITERE
216C-------------------
217 DO i=1,nel
218 jj(i) = 1
219 evl(i)=zero
220 ENDDO
221C
222 DO j=1, nfunc
223 DO i=1,nel
224C---first direction--
225 ipos1(i) = nint(uvar(i,1+j))
226 iad1(i) = npf(ifunc(i,j)) / 2 + 1
227 ilen1(i) = npf(ifunc(i,j)+1) / 2 - iad1(i) - ipos1(i)
228C--- second direction
229C IPOS2(I) = NINT(UVAR(I,3))
230C IAD2(I) = 0.5*NPF(IFUNC(I,2)) + 1
231C ILEN2(I) = 0.5*NPF(IFUNC(I,3)+1) - IAD2(I) - IPOS2(I)
232C--- third direction
233C IPOS3(I) = NINT(UVAR(I,4))
234C IAD3(I) = 0.5*NPF(IFUNC(I,3)) + 1
235C ILEN3(I) = 0.5*NPF(IFUNC(I,4)+1) - IAD3(I) - IPOS3(I)
236C....fourth direction
237C IPOS4(I) = NINT(UVAR(I,5))
238C IAD4(I) = 0.5*NPF(IFUNC(I,4)) + 1
239C ILE4(I) = 0.5*NPF(IFUNC(I,5)+1) - IAD4(I) - IPOS4(I)
240C... 45 direction curve
241C IPOS5(I) = NINT(UVAR(I,6))
242C IAD5(I) = 0.5*NPF(IFUNC(I,5)) + 1
243C ILE5(I) = 0.5*NPF(IFUNC(I,6)+1) - IAD5(I) - IPOS5(I)
244 ENDDO
245C extrapollation dans chaque direction
246C---first direction--
247 CALL vinter(tf,iad1,ipos1,ilen1,nel,evol,dydx1,y1)
248 CALL vinter(tf,iad1,ipos1,ilen1,nel,evl,dydx1,y2)
249C
250 IF(j==1) THEN
251 DO i=1, nel
252 IF(iflag(i)/=1)THEN
253 IF(evol(i)>zero)THEN
254 s1c(i)= y1(i)
255 s1t(i)= y2(i)
256 ELSE
257 s1c(i)= y2(i)
258 s1t(i)= y1(i)
259 ENDIF
260 ELSE
261 s1c(i)= y1(i)
262 s1t(i)= y1(i)
263 ENDIF
264 uvar(i,1+j)= ipos1(i)
265 ENDDO
266C second direction
267 ELSEIF(j==2)THEN
268C CALL vINTER(TF,IAD2,IPOS2,ILEN2,NEL,PLA,DYDX1,Y1)
269C CALL VINTER(TF,IAD2,IPOS2,ILEN2,NEL,EVL,DYDX1,Y2)
270C
271 DO i=1, nel
272 IF(iflag(i)/=1)THEN
273 IF(evol(i)>zero)THEN
274 s2c(i)= y1(i)
275 s2t(i)= y2(i)
276 ELSE
277 s2c(i)= y2(i)
278 s2t(i)= y1(i)
279 ENDIF
280 ELSE
281 s2c(i)= y1(i)
282 s2t(i)= y1(i)
283 ENDIF
284 uvar(i,1+j)= ipos1(i)
285 ENDDO
286CC third direction
287 ELSEIF(j==3)THEN
288C CALL VINTER(TF,IAD3,IPOS3,ILEN3,NEL,PLA,DYDX1,Y1)
289C CALL VINTER(TF,IAD3,IPOS3,ILEN3,NEL,EVL,DYDX1,Y2)
290C
291 DO i=1, nel
292 IF(iflag(i)/=1)THEN
293 IF(evol(i)>zero)THEN
294 s3c(i)= y1(i)
295 s3t(i)= y2(i)
296 ELSE
297 s3c(i)= y2(i)
298 s3t(i)= y1(i)
299 ENDIF
300 ELSE
301 s3c(i)= y1(i)
302 s3t(i)= y1(i)
303 ENDIF
304 uvar(i,1+j)= ipos1(i)
305 ENDDO
306CCC fourth direction
307 ELSEIF(j==4)THEN
308C CALL VINTER(TF,IAD4,IPOS4,ILEN4,NEL,PLA,DYDX1,Y1)
309C CALL VINTER(TF,IAD4,IPOS4,ILEN4,NEL,EVL,DYDX1,Y2)
310C
311 DO i=1, nel
312 IF(iflag(i)/=1)THEN
313 IF(evol(i)>zero)THEN
314 s4c(i)= y1(i)
315 s4t(i)= y2(i)
316 ELSE
317 s4c(i)= y2(i)
318 s4t(i)= y1(i)
319 ENDIF
320 ELSE
321 s4c(i)= y1(i)
322 s4t(i)= y1(i)
323 ENDIF
324 uvar(i,1+j)= ipos1(i)
325 ENDDO
326CCCC 45 direction
327 ELSEIF(j==5)THEN
328C CALL VINTER(TF,IAD5,IPOS5,ILEN5,NEL,PLA,DYDX1,Y1)
329C CALL VINTER(TF,IAD5,IPOS5,ILEN5,NEL,EVL,DYDX1,Y2)
330C
331 DO i=1, nel
332 IF(iflag(i)/=1)THEN
333 IF(evol(i)>zero)THEN
334 s45c(i)= y1(i)
335 s45t(i)= y2(i)
336 ELSE
337 s45c(i)= y2(i)
338 s45t(i)= y1(i)
339 ENDIF
340 ELSE
341 s45c(i)= y1(i)
342 s45t(i)= y1(i)
343 ENDIF
344 uvar(i,1+j)= ipos1(i)
345 ENDDO
346 ENDIF
347 ENDDO
348C--YEILD SURFACE COEFFICIENTS
349 DO i=1,nel
350 f1(i) = -one/s1c(i) + one/s1t(i)
351 f2(i) = -one/s2c(i) + one/s2t(i)
352 f11(i)= one/(s1c(i)*s1t(i))
353 f22(i)= one/(s2c(i)*s2t(i))
354 f44(i)= one/(s3c(i)*s3t(i))
355 f55(i)= one/(s4c(i)*s4t(i))
356 f12(i)=(-half)*sqrt(f11(i)*f22(i))
357 f23(i)=(-half)*f22(i)
358C... CORRECTION
359 IF(ifunc(i,5)>0)THEN
360 f12(i)= two/(s45c(i)*s45c(i)) -
361 . (half)*(f11(i)+ f22(i) + f44(i)) +
362 . (f1(i) + f2(i))/s45c(i)
363 ENDIF
364C... Check the yield condition according to Tsay Wu.
365 f(i) = f1(i)*signxx(i) + f2(i)* signyy(i) + f2(i)*signzz(i)+
366 . f11(i)*signxx(i)**2 + f22(i)*signyy(i)**2 + f22(i)*signzz(i)**2
367 .+f44(i)*signxy(i)**2 + f55(i)*signyz(i)**2 + f44(i)*signzx(i)**2
368 .+two*f12(i)*signxx(i)*signyy(i) + two*f23(i)*signyy(i)*signzz(i)
369 .+two*f12(i)*signzz(i)*signxx(i)
370C
371!! SEQ_OUTPUT(I) = F(I)
372C
373 IF(f(i)<=one) THEN
374 scale=one
375 ELSE !IF(F(I)>ONE)THEN
376 bb= f1(i)*signxx(i) + f2(i)*signyy(i)+ f2(i)*signzz(i)
377 cc = -one
378 aa =
379 . f11(i)*signxx(i)**2 + f22(i)*signyy(i)**2 + f22(i)*signzz(i)**2
380 .+f44(i)*signxy(i)**2 + f55(i)*signyz(i)**2 + f44(i)*signzx(i)**2
381 .+two*f12(i)*signxx(i)*signyy(i) + two*f23(i)*signyy(i)*signzz(i)
382 .+two*f12(i)*signzz(i)*signxx(i)
383C
384 dd= bb**2 - four*aa*cc
385 IF(dd<zero) THEN
386 CALL ancmsg(msgid=136,anmode=aninfo)
387 CALL arret(2)
388 ENDIF
389C...
390 ss1 = (-bb + sqrt(dd))/(two*aa)
391 ss2 = (-bb - sqrt(dd))/(two*aa)
392C
393 IF(ss1<=zero) scale = ss2
394 IF(ss2<=zero) scale = ss1
395C
396 IF(ss1>0.AND.ss2>0) THEN
397 WRITE(*,*)' TWO POSITIVE ROOTS IN STRANDFOAM'
398 IF(ss1<ss2) scale = ss1
399 IF(ss2<ss1) scale = ss2
400 ENDIF
401 ENDIF
402CCC STRESS SCALING
403C radial in stress space
404 signxx(i) = signxx(i) * scale
405 signyy(i) = signyy(i) * scale
406 signzz(i) = signzz(i) * scale
407 signxy(i) = signxy(i) * scale
408 signyz(i) = signyz(i) * scale
409 signzx(i) = signzx(i) * scale
410
411 ENDDO
412C-----
413 RETURN
my_real
#define my_real
Definition cppsort.cpp:32
max
#define max(a, b)
Definition macros.h:21
message_mod
Definition message_mod.F:1249
ancmsg
subroutine ancmsg(msgid, msgtype, anmode, i1, i2, i3, i4, i5, i6, i7, i8, i9, i10, i11, i12, i13, i14, i15, i16, i17, i18, i19, i20, r1, r2, r3, r4, r5, r6, r7, r8, r9, c1, c2, c3, c4, c5, c6, c7, c8, c9, prmode)
Definition message.F:889
arret
subroutine arret(nn)
Definition arret.F:87
vinter
subroutine vinter(tf, iad, ipos, ilen, nel, x, dydx, y)
Definition vinter.F:72