OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
sigeps90.F File Reference
#include "implicit_f.inc"
#include "scr05_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine sigeps90 (nel, nuparam, nuvar, nfunc, ifunc, 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, sig0xx, sig0yy, sig0zz, sig0xy, sig0yz, sig0zx, signxx, signyy, signzz, signxy, signyz, signzx, soundsp, viscmax, uvar, off, ngl, ismstr)

Function/Subroutine Documentation

◆ sigeps90()

subroutine sigeps90 ( integer, intent(in) nel,
integer, intent(in) nuparam,
integer, intent(in) nuvar,
integer, intent(in) nfunc,
integer, dimension(nfunc), intent(in) ifunc,
integer, dimension(*), intent(in) npf,
dimension(*), intent(in) tf,
intent(in) time,
intent(in) timestep,
dimension(*), intent(in) uparam,
dimension(nel), intent(in) rho0,
dimension(nel), intent(in) rho,
dimension(nel), intent(in) volume,
dimension(nel), intent(in) eint,
dimension(nel), intent(in) epspxx,
dimension(nel), intent(in) epspyy,
dimension(nel), intent(in) epspzz,
dimension(nel), intent(in) epspxy,
dimension(nel), intent(in) epspyz,
dimension(nel), intent(in) epspzx,
dimension(nel), intent(in) depsxx,
dimension(nel), intent(in) depsyy,
dimension(nel), intent(in) depszz,
dimension(nel), intent(in) depsxy,
dimension(nel), intent(in) depsyz,
dimension(nel), intent(in) depszx,
dimension(nel), intent(in) epsxx,
dimension(nel), intent(in) epsyy,
dimension(nel), intent(in) epszz,
dimension(nel), intent(in) epsxy,
dimension(nel), intent(in) epsyz,
dimension(nel), intent(in) epszx,
dimension(nel), intent(in) sig0xx,
dimension(nel), intent(in) sig0yy,
dimension(nel), intent(in) sig0zz,
dimension(nel), intent(in) sig0xy,
dimension(nel), intent(in) sig0yz,
dimension(nel), intent(in) sig0zx,
dimension(nel), intent(out) signxx,
dimension(nel), intent(out) signyy,
dimension(nel), intent(out) signzz,
dimension(nel), intent(out) signxy,
dimension(nel), intent(out) signyz,
dimension(nel), intent(out) signzx,
dimension(nel), intent(out) soundsp,
dimension(nel), intent(out) viscmax,
dimension(nel,nuvar), intent(inout) uvar,
dimension(nel), intent(inout) off,
integer, dimension(nel), intent(in) ngl,
integer, intent(in) ismstr )

Definition at line 32 of file sigeps90.F.

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-----------------------------------------------
49C---------+---------+---+---+--------------------------------------------
50C VAR | SIZE |TYP| RW| DEFINITION
51C---------+---------+---+---+--------------------------------------------
52C NEL | 1 | I | R | SIZE OF THE ELEMENT GROUP NEL
53C NUPARAM | 1 | I | R | SIZE OF THE USER PARAMETER ARRAY
54C NUVAR | 1 | I | R | NUMBER OF USER ELEMENT VARIABLES
55C---------+---------+---+---+--------------------------------------------
56C MFUNC | 1 | I | R | NUMBER FUNCTION USED FOR THIS USER LAW not used
57C KFUNC | NFUNC | I | R | FUNCTION INDEX not used
58C NPF | * | I | R | FUNCTION ARRAY
59C TF | * | F | R | FUNCTION ARRAY
60C---------+---------+---+---+--------------------------------------------
61C TIME | 1 | F | R | CURRENT TIME
62C TIMESTEP| 1 | F | R | CURRENT TIME STEP
63C UPARAM | NUPARAM | F | R | USER MATERIAL PARAMETER ARRAY
64C RHO0 | NEL | F | R | INITIAL DENSITY
65C RHO | NEL | F | R | DENSITY
66C VOLUME | NEL | F | R | VOLUME
67C EINT | NEL | F | R | TOTAL INTERNAL ENERGY
68C EPSPXX | NEL | F | R | STRAIN RATE XX
69C EPSPYY | NEL | F | R | STRAIN RATE YY
70C ... | | | |
71C DEPSXX | NEL | F | R | STRAIN INCREMENT XX
72C DEPSYY | NEL | F | R | STRAIN INCREMENT YY
73C ... | | | |
74C EPSXX | NEL | F | R | STRAIN XX
75C EPSYY | NEL | F | R | STRAIN YY
76C ... | | | |
77C SIG0XX | NEL | F | R | OLD ELASTO PLASTIC STRESS XX
78C SIG0YY | NEL | F | R | OLD ELASTO PLASTIC STRESS YY
79C ... | | | |
80C---------+---------+---+---+--------------------------------------------
81C SIGNXX | NEL | F | W | NEW ELASTO PLASTIC STRESS XX
82C SIGNYY | NEL | F | W | NEW ELASTO PLASTIC STRESS YY
83C ... | | | |
84C SIGVXX | NEL | F | W | VISCOUS STRESS XX
85C SIGVYY | NEL | F | W | VISCOUS STRESS YY
86C ... | | | |
87C SOUNDSP | NEL | F | W | SOUND SPEED (NEEDED FOR TIME STEP)
88C VISCMAX | NEL | F | W | MAXIMUM DAMPING MODULUS(NEEDED FOR TIME STEP)
89C---------+---------+---+---+--------------------------------------------
90C UVAR |NEL*NUVAR| F |R/W| USER ELEMENT VARIABLE ARRAY
91C OFF | NEL | F |R/W| DELETED ELEMENT FLAG (=1. ON, =0. OFF)
92C---------+---------+---+---+--------------------------------------------
93#include "scr05_c.inc"
94C
95 INTEGER, INTENT(IN) :: NEL, NUPARAM, NUVAR, NGL(NEL),ISMSTR
96
97 my_real, INTENT(IN) ::
98 . time,timestep,uparam(*),
99 . rho(nel),rho0(nel),volume(nel),eint(nel),
100 . epspxx(nel),epspyy(nel),epspzz(nel),
101 . epspxy(nel),epspyz(nel),epspzx(nel),
102 . depsxx(nel),depsyy(nel),depszz(nel),
103 . depsxy(nel),depsyz(nel),depszx(nel),
104 . epsxx(nel) ,epsyy(nel) ,epszz(nel) ,
105 . epsxy(nel) ,epsyz(nel) ,epszx(nel) ,
106 . sig0xx(nel),sig0yy(nel),sig0zz(nel),
107 . sig0xy(nel),sig0yz(nel),sig0zx(nel)
108C-----------------------------------------------
109C O U T P U T A r g u m e n t s
110C-----------------------------------------------
111 my_real,INTENT(OUT)::
112 . signxx(nel),signyy(nel),signzz(nel),
113 . signxy(nel),signyz(nel),signzx(nel),
114 . soundsp(nel),viscmax(nel)
115C-----------------------------------------------
116C I N P U T O U T P U T A r g u m e n t s
117C-----------------------------------------------
118 my_real,INTENT(INOUT) ::
119 . uvar(nel,nuvar), off(nel)
120C-----------------------------------------------
121C VARIABLES FOR FUNCTION INTERPOLATION
122C-----------------------------------------------
123 INTEGER, INTENT(IN) :: NPF(*), NFUNC, IFUNC(NFUNC)
124 my_real, INTENT(IN) :: tf(*)
125 my_real, EXTERNAL :: finter
126C EXTERNAL FINTER
127C Y = FINTER(IFUNC(J),X,NPF,TF,DYDX)
128C Y : y = f(x)
129C X : x
130C DYDX : f'(x) = dy/dx
131C IFUNC(J): FUNCTION INDEX
132C J : FIRST(J=1), SECOND(J=2) .. FUNCTION USED FOR THIS LAW
133C NPF,TF : FUNCTION PARAMETER
134C-----------------------------------------------
135C L o c a l V a r i a b l e s
136C-----------------------------------------------
137 INTEGER
138 . II,I,J,K,I1,J1,J2,IFLAG,ILOAD(NEL,3),ILOADE(NEL),
139 . IDAM,INDX_L(NEL),INDX_UNL(NEL),NE_L,NE_UNL
140 my_real
141 . e0,aa,epsmax,g,nu,shape,hys,emax,
142 . yfac,yfacj1,yfacj2,ratej1,ratej2, epse,ep1,
143 . ep2,ep3,ep4,ep5,ep6,ert11,ert12,ert13,ert21,
144 . ert22,ert23,ert31,ert32,ert33,sj1,sj2,fac,t1,t2,t3,
145 . dam,epe,emin,e_min(nel),e1,e2,tmp,quasi_eint
146 my_real
147 . av(6,nel), evv(3,nel),ev(nel,3),strain(nel,3),
148 . strainrate(nel,3),s(nel,3),sqstat(nel,3),df(3),
149 . dirprv(3,3,nel),epsp(3),ecurent(nel),e(nel),deint,
150 . rateeps
151C-----------------------------------------------
152C USER VARIABLES INITIALIZATION
153C-----------------------------------------------
154 e0 = uparam(1)
155 g = uparam(4)
156 nu = uparam(5)
157 shape = uparam(6)
158 hys = uparam(7)
159 iflag = uparam(9)
160 idam = uparam(10)
161C
162C-----------------------------------------------
163C
164 DO i=1,nel
165 av(1,i) = epsxx(i)
166 av(2,i) = epsyy(i)
167 av(3,i) = epszz(i)
168 av(4,i) = half*epsxy(i)
169 av(5,i) = half*epsyz(i)
170 av(6,i) = half*epszx(i)
171 ENDDO
172CEigenvalues needed to be calculated in double precision
173C for a simple precision executing*
174 IF (iresp==1) THEN
175 CALL valpvecdp(av,evv,dirprv,nel)
176 ELSE
177 CALL valpvec(av,evv,dirprv,nel)
178 ENDIF
179C-ISMSTR=0-NO SMALL STRAIN OPTION:STRAINS ARE LOGARITHMIC, STRESS IS CAUCHY
180C-ISMSTR=1-SMALL STRAIN OPTION:STRAINS ARE ENGINEERING, STRESS IS CAUCHY
181C-ISMSTR=2-SMALL STRAIN OPTION:STRAINS ARE ENGINEERING, STRESS IS BIOT
182C-ISMSTR=3-NO SMALL STRAIN OPTION:STRESS IS BIOT
183 IF(ismstr==0.OR.ismstr==2.OR.ismstr==4) THEN
184 DO i=1,nel
185C ---- (STRAIN IS LOGARITHMIC)
186 ev(i,1)=exp(evv(1,i))
187 ev(i,2)=exp(evv(2,i))
188 ev(i,3)=exp(evv(3,i))
189 ENDDO
190 ELSEIF(ismstr==10.OR.ismstr==12) THEN
191 DO i =1,nel
192 ev(i,1)=sqrt(evv(1,i) + one )
193 ev(i,2)=sqrt(evv(2,i) + one )
194 ev(i,3)=sqrt(evv(3,i) + one )
195 ENDDO
196 ELSE
197C ---- STRAIN IS ENGINEERING)
198 DO i=1,nel
199 ev(i,1)=evv(1,i) + one
200 ev(i,2)=evv(2,i) + one
201 ev(i,3)=evv(3,i) + one
202 ENDDO
203 ENDIF
204
205C engineering strain and strain rate
206 DO i=1,nel
207C engineering strain e = lambda-1 , according the input curve
208C e=1-lambda (e > 0 compression and e < 0 traction)
209 strain(i,1) = one - ev(i,1)
210 strain(i,2) = one - ev(i,2)
211 strain(i,3) = one - ev(i,3)
212C
213 ep1 = epspxx(i)
214 ep2 = epspyy(i)
215 ep3 = epspzz(i)
216 ep4 = half*epspxy(i)
217 ep5 = half*epspyz(i)
218 ep6 = half*epspzx(i)
219C phi_trans*L*phi_t
220 ert11 =dirprv(1,1,i)*ep1 + dirprv(2,1,i)*ep4 + dirprv(3,1,i)*ep6
221 ert12 =dirprv(1,2,i)*ep1 + dirprv(2,2,i)*ep4 + dirprv(3,2,i)*ep6
222 ert13 =dirprv(1,3,i)*ep1 + dirprv(2,3,i)*ep4 + dirprv(3,3,i)*ep6
223
224 ert21 =dirprv(1,1,i)*ep4 + dirprv(2,1,i)*ep2 + dirprv(3,1,i)*ep5
225 ert22 =dirprv(1,2,i)*ep4 + dirprv(2,2,i)*ep2 + dirprv(3,2,i)*ep5
226 ert23 =dirprv(1,3,i)*ep4 + dirprv(2,3,i)*ep2 + dirprv(3,3,i)*ep5
227
228 ert31 =dirprv(1,1,i)*ep6 + dirprv(2,1,i)*ep5 + dirprv(3,1,i)*ep3
229 ert32 =dirprv(1,2,i)*ep6 + dirprv(2,2,i)*ep5 + dirprv(3,2,i)*ep3
230 ert33 =dirprv(1,3,i)*ep6 + dirprv(2,3,i)*ep5 + dirprv(3,3,i)*ep3
231C
232 epsp(1) = dirprv(1,1,i)*ert11 + dirprv(2,1,i)*ert21
233 . + dirprv(3,1,i)*ert31
234 epsp(2) = dirprv(1,2,i)*ert12 + dirprv(2,2,i)*ert22
235 . + dirprv(3,2,i)*ert32
236 epsp(3) = dirprv(1,3,i)*ert13 + dirprv(2,3,i)*ert23
237 . + dirprv(3,3,i)*ert33
238 strainrate(i,1) = abs(epsp(1))
239 strainrate(i,2) = abs(epsp(2))
240 strainrate(i,3) = abs(epsp(3))
241C computing energy increase
242!! YFAC = UPARAM(NFUNC + 11)
243!! QUASI_EINT= ZERO
244!! DO K=1,3
245!! EPSE = STRAIN(I,K)
246!! SQSTAT(I,K) = YFAC*FINTER(IFUNC(1),EPSE,NPF,TF,DF(K))
247C compute current energy
248!! QUASI_EINT= QUASI_EINT + HALF*STRAIN(I,K)*SQSTAT(I,K)
249!! ENDDO
250!! DEINT = QUASI_EINT - UVAR(I,9)
251!! UVAR(I,9) = QUASI_EINT
252C -----------
253C check loading and unloading.
254 iloade(i) = 1
255 DO k=1,3
256 epe = epsp(k)*strain(i,k)
257 iload(i,k) = 1
258 IF(epe > em10 )iload(i,k) = -1
259 IF(iload(i,k) == -1 ) iloade(i) = -1
260 ENDDO
261C
262 DO k=1,3
263 IF(iload(i,k) == -1) strainrate(i,k) = zero
264 ENDDO
265
266 uvar(i,3) = strainrate(i,1)
267 uvar(i,4) = strainrate(i,2)
268 uvar(i,5) = strainrate(i,3)
269C
270 rateeps = max(strainrate(i,1),strainrate(i,2),strainrate(i,3))
271
272!! EPSD(I) = RATEEPS
273 ENDDO
274C sous groupe
275 indx_l(1:nel) = 0
276 indx_unl(1:nel) = 0
277 ne_l = 0
278 ne_unl = 0
279C
280 DO i=1,nel
281 IF(iloade(i) == 1) THEN
282 ne_l = ne_l + 1
283 indx_l(ne_l) = i
284 uvar(i,8) = one
285 ELSEIF(iloade(i) == -1 ) THEN
286 ne_unl = ne_unl + 1
287 indx_unl(ne_unl) = i
288 uvar(i,8) = -one
289 ENDIF
290 ENDDO
291C case with unloading stress-strain curve only the quasistatic curve is used
292 IF(iflag == 1) THEN
293 DO i=1,nel
294 yfac = uparam(nfunc + 11)
295 emax = zero
296 DO k=1,3
297 epse = strain(i,k)
298 s(i,k) = yfac*finter(ifunc(1),epse,npf,tf,df(k))
299 emax = max(emax, yfac*df(k))
300 ENDDO
301 e(i) = max(e0,emax)
302 ENDDO
303 ENDIF
304C unloading with damage based on the energy.
305 IF(iflag == 2) THEN
306 DO i=1,nel
307 yfac = uparam(nfunc + 11)
308 ecurent(i)= zero
309 emax = zero
310 emin = ep20
311 DO k=1,3
312 epse = strain(i,k)
313 sqstat(i,k) = yfac*finter(ifunc(1),epse,npf,tf,df(k))
314 s(i,k) = sqstat(i,k)
315 emax = max(emax, yfac*df(k))
316 emin = min(emin, yfac*df(k))
317C computecurent energy
318 ecurent(i)= ecurent(i) + half*strain(i,k)*sqstat(i,k)
319!! ECURENT(I)= ECURENT(I) + HALF*EV(I,K)*ABS(SQSTAT(I,K))
320 ENDDO
321 uvar(i,2) = max(uvar(i,2),ecurent(i))
322 e(i) = max(e0,emax)
323 e_min(i) =emin
324 ENDDO
325C flag is a hidden flag, only idam=0 is activated, Idam > o not tested.
326 IF(idam == 0) THEN
327 DO ii=1,ne_unl
328 i = indx_unl(ii)
329 uvar(i,1) = ecurent(i)
330 IF(uvar(i,2) > zero) THEN
331 dam = one - (ecurent(i)/uvar(i,2))**shape
332 dam = one - (one - hys)*dam
333C global
334 DO k=1,3
335 s(i,k)= dam*s(i,k)
336 ENDDO
337 ENDIF
338 ENDDO ! NE_UNL
339 ELSE ! not tested Idam/=0
340C damage by direction to be tested for
341 DO ii=1,ne_unl
342 i = indx_unl(ii)
343 uvar(i,1) = ecurent(i)
344 IF(uvar(i,2) > zero) THEN
345 dam = one - (ecurent(i)/uvar(i,2))**shape
346 dam = one - (one - hys)*dam
347 DO k=1,3
348 IF(iload(i,k) < 0)s(i,k) = dam*s(i,k)
349 ENDDO
350 ENDIF
351 ENDDO ! nel
352 ENDIF ! IDAM
353 ENDIF ! iflag=2
354C
355C =====================================================
356 DO i = 1,nel
357C S > 0 for compression - curve definition
358C S < 0 for traction
359 t1 = -s(i,1)/ev(i,2)/ev(i,3)
360 t2 = -s(i,2)/ev(i,1)/ev(i,3)
361 t3 = -s(i,3)/ev(i,1)/ev(i,2)
362C
363C cauchy to glabale
364C
365 signxx(i) = dirprv(1,1,i)*dirprv(1,1,i)*t1
366 . + dirprv(1,2,i)*dirprv(1,2,i)*t2
367 . + dirprv(1,3,i)*dirprv(1,3,i)*t3
368
369 signyy(i) = dirprv(2,2,i)*dirprv(2,2,i)*t2
370 . + dirprv(2,3,i)*dirprv(2,3,i)*t3
371 . + dirprv(2,1,i)*dirprv(2,1,i)*t1
372
373 signzz(i) = dirprv(3,3,i)*dirprv(3,3,i)*t3
374 . + dirprv(3,1,i)*dirprv(3,1,i)*t1
375 . + dirprv(3,2,i)*dirprv(3,2,i)*t2
376
377 signxy(i) = dirprv(1,1,i)*dirprv(2,1,i)*t1
378 . + dirprv(1,2,i)*dirprv(2,2,i)*t2
379 . + dirprv(1,3,i)*dirprv(2,3,i)*t3
380
381 signyz(i) = dirprv(2,2,i)*dirprv(3,2,i)*t2
382 . + dirprv(2,3,i)*dirprv(3,3,i)*t3
383 . + dirprv(2,1,i)*dirprv(3,1,i)*t1
384
385 signzx(i) = dirprv(3,3,i)*dirprv(1,3,i)*t3
386 . + dirprv(3,1,i)*dirprv(1,1,i)*t1
387 . + dirprv(3,2,i)*dirprv(1,2,i)*t2
388C==================================================
389 aa = e(i)*(one-nu)/(one + nu)/(one - two*nu)
390 soundsp(i) = sqrt(aa/rho0(i))
391 viscmax(i) = zero
392C
393 ENDDO
394CE0
395C------------------------------------
396 RETURN
my_real
#define my_real
Definition cppsort.cpp:32
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
valpvecdp
subroutine valpvecdp(sig, val, vec, nel)
Definition matrix.F:533
valpvec
subroutine valpvec(sig, val, vec, nel)
Definition matrix.F:215