OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
sigeps70.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!|| sigeps70 ../engine/source/materials/mat/mat070/sigeps70.F
25!||--- called by ------------------------------------------------------
26!|| mulaw ../engine/source/materials/mat_share/mulaw.f90
27!||--- calls -----------------------------------------------------
28!|| finter ../engine/source/tools/curve/finter.F
29!|| table_mat_vinterp ../engine/source/materials/tools/table_mat_vinterp.F
30!|| vinter2 ../engine/source/tools/curve/vinter.F
31!||--- uses -----------------------------------------------------
32!|| matparam_def_mod ../common_source/modules/mat_elem/matparam_def_mod.F90
33!|| table_mat_vinterp_mod ../engine/source/materials/tools/table_mat_vinterp.F
34!||====================================================================
35 SUBROUTINE sigeps70(
36 1 NEL, NUPARAM, NUVAR, MFUNC,
37 2 KFUNC, NPF, TF, TIME,
38 3 TIMESTEP,UPARAM, RHO0, RHO,
39 4 OFFG, RHOREF, NVARTMP, VARTMP,
40 5 RHOSP, DEPSXX,
41 7 DEPSYY, DEPSZZ, DEPSXY, DEPSYZ,
42 8 DEPSZX, EPSXX, EPSYY, EPSZZ,
43 9 EPSXY, EPSYZ, EPSZX, SIG0XX,
44 A SIG0YY, SIG0ZZ, SIG0XY, SIG0YZ,
45 B SIG0ZX, SIGNXX, SIGNYY, SIGNZZ,
46 C SIGNXY, SIGNYZ, SIGNZX,
47 E SOUNDSP, VISCMAX, UVAR,
48 F OFF, NGL, IPM,
49 G MAT, EPSP, ET, ISMSTR,
50 H IHET, JSMS, MATPARAM)
51C-----------------------------------------------
52C M o d u l e s
53C-----------------------------------------------
54 USE matparam_def_mod
56C-----------------------------------------------
57C I m p l i c i t T y p e s
58C-----------------------------------------------
59#include "implicit_f.inc"
60#include "comlock.inc"
61C-----------------------------------------------
62C G l o b a l P a r a m e t e r s
63C-----------------------------------------------
64#include "mvsiz_p.inc"
65C---------+---------+---+---+--------------------------------------------
66C VAR | SIZE |TYP| RW| DEFINITION
67C---------+---------+---+---+--------------------------------------------
68C NEL | 1 | I | R | SIZE OF THE ELEMENT GROUP NEL
69C NUPARAM | 1 | I | R | SIZE OF THE USER PARAMETER ARRAY
70C NUVAR | 1 | I | R | NUMBER OF USER ELEMENT VARIABLES
71C---------+---------+---+---+--------------------------------------------
72C MFUNC | 1 | I | R | NUMBER FUNCTION USED FOR THIS USER LAW not used
73C KFUNC | NFUNC | I | R | FUNCTION INDEX not used
74C NPF | * | I | R | FUNCTION ARRAY
75C TF | * | F | R | FUNCTION ARRAY
76C---------+---------+---+---+--------------------------------------------
77C TIME | 1 | F | R | CURRENT TIME
78C TIMESTEP| 1 | F | R | CURRENT TIME STEP
79C UPARAM | NUPARAM | F | R | USER MATERIAL PARAMETER ARRAY
80C RHO0 | NEL | F | R | INITIAL DENSITY
81C RHO | NEL | F | R | DENSITY
82C ... | | | |
83C DEPSXX | NEL | F | R | STRAIN INCREMENT XX
84C DEPSYY | NEL | F | R | STRAIN INCREMENT YY
85C ... | | | |
86C EPSXX | NEL | F | R | STRAIN XX
87C EPSYY | NEL | F | R | STRAIN YY
88C ... | | | |
89C SIG0XX | NEL | F | R | OLD ELASTO PLASTIC STRESS XX
90C SIG0YY | NEL | F | R | OLD ELASTO PLASTIC STRESS YY
91C ... | | | |
92C---------+---------+---+---+--------------------------------------------
93C SIGNXX | NEL | F | W | NEW ELASTO PLASTIC STRESS XX
94C SIGNYY | NEL | F | W | NEW ELASTO PLASTIC STRESS YY
95C ... | | | |
96C SOUNDSP | NEL | F | W | SOUND SPEED (NEEDED FOR TIME STEP)
97C VISCMAX | NEL | F | W | MAXIMUM DAMPING MODULUS(NEEDED FOR TIME STEP)
98C---------+---------+---+---+--------------------------------------------
99C UVAR |NEL*NUVAR| F |R/W| USER ELEMENT VARIABLE ARRAY
100C OFF | NEL | F |R/W| DELETED ELEMENT FLAG (=1. ON, =0. OFF)
101C---------+---------+---+---+--------------------------------------------
102#include "param_c.inc"
103#include "sms_c.inc"
104#include "impl1_c.inc"
105C-----------------------------------------------
106 INTEGER ,INTENT(IN) :: NVARTMP
107 INTEGER ,INTENT(INOUT) :: VARTMP(NEL,NVARTMP)
108 INTEGER, INTENT(IN) :: ISMSTR, IHET, JSMS
109 INTEGER NEL, NUPARAM, NUVAR,
110 . NGL(NEL),MAT(NEL),IPM(NPROPMI,*)
111 my_real
112 . TIME,TIMESTEP,UPARAM(*),
113 . RHO(NEL),RHO0(NEL),
114 . DEPSXX(NEL),DEPSYY(NEL),DEPSZZ(NEL),
115 . DEPSXY(NEL),DEPSYZ(NEL),DEPSZX(NEL),
116 . EPSXX(NEL) ,EPSYY(NEL) ,EPSZZ(NEL) ,
117 . EPSXY(NEL) ,EPSYZ(NEL) ,EPSZX(NEL) ,
118 . SIG0XX(NEL),SIG0YY(NEL),SIG0ZZ(NEL),
119 . SIG0XY(NEL),SIG0YZ(NEL),SIG0ZX(NEL),
120 . EPSP(NEL)
121 TYPE(matparam_struct_) ,INTENT(IN) :: MATPARAM
122C-----------------------------------------------
123C O U T P U T A r g u m e n t s
124C-----------------------------------------------
125 my_real
126 . signxx(nel),signyy(nel),signzz(nel),
127 . signxy(nel),signyz(nel),signzx(nel),
128 . soundsp(nel),viscmax(nel),et(nel),
129 . rhosp(mvsiz)
130C-----------------------------------------------
131C I N P U T O U T P U T A r g u m e n t s
132C-----------------------------------------------
133 my_real
134 . uvar(nel,nuvar), off(nel), pla(nel), offg(nel), rhoref(mvsiz)
135C-----------------------------------------------
136C VARIABLES FOR FUNCTION INTERPOLATION
137C-----------------------------------------------
138 INTEGER NPF(*), MFUNC, KFUNC(MFUNC)
139 my_real
140 . tf(*),finter
141 EXTERNAL finter
142C EXTERNAL FINTER
143C Y = FINTER(IFUNC(J),X,NPF,TF,DYDX)
144C Y : y = f(x)
145C X : x
146C DYDX : f'(x) = dy/dx
147C IFUNC(J): FUNCTION INDEX
148C J : FIRST(J=1), SECOND(J=2) .. FUNCTION USED FOR THIS LAW
149C NPF,TF : FUNCTION PARAMETER
150C-----------------------------------------------
151C L o c a l V a r i a b l e s
152C-----------------------------------------------
153 INTEGER :: I,II,J,J1,J2,IADBUF,NFUNC,
154 . nindex_epssm,nindex_epst,nindex_mu,nload,nunload,mx,
155 . iload0,idamage,itens,numtabl
156 INTEGER, DIMENSION(100) :: IFUNC
157 INTEGER, DIMENSION(NEL) :: JJ,IE_CST,ILOAD
158 INTEGER, DIMENSION(NEL) :: INDEX_EPSSM,INDEX_EPST,INDEX_MU
159 INTEGER, DIMENSION(NEL) :: IPOS1,IADP1,ILENP1
160 INTEGER, DIMENSION(NEL,2) :: IPOS
161 my_real :: YLD_EMAX,R,FAC,YP1,YP2,EPSP0,EPSPMAX,E,DELTA,SVM2,SVM,
162 . AA,E0,NU,EMAX,EPSSMAX,DSIG,P,EXPO,HYS,DE,YFACC,MU,ALPHA
163 my_real, DIMENSION(NEL) :: ALPHA_1,AA1,AA2,EPST,YLDMIN,YLDMAX,YLDELAS,
164 . G,SIG0,EPS0,EPSS,TAB_LOC1,TAB_LOC2,YLD,DFC
165 my_real, DIMENSION(NEL,2) :: XVEC
166C=======================================================================
167 MX = mat(1)
168 nfunc = ipm(10,mx)
169 DO j=1,nfunc
170 ifunc(j)=ipm(10+j,mx)
171 ENDDO
172C
173 iadbuf = ipm(7,mx)-1
174 e0 = uparam(iadbuf + 2)
175 aa = uparam(iadbuf + 3)
176 epssmax = uparam(iadbuf + 4)
177 g(:) = uparam(iadbuf + 5)
178 nu = uparam(iadbuf + 6)
179 epsp0 = uparam(iadbuf + 9) ! static strain rate
180 nload = nint(uparam(iadbuf + 7))
181 nunload = nint(uparam(iadbuf + 8))
182 idamage = nint(uparam(iadbuf + 2*nfunc + 9))
183 itens = nint(uparam(iadbuf + 2*nfunc + 13))
184 yfacc = uparam(iadbuf + 2*nfunc + 8)
185 expo = uparam(iadbuf + 2*nfunc + 10)
186 hys = uparam(iadbuf + 2*nfunc + 11)
187 emax = uparam(iadbuf + 2*nfunc + 12)
188 yld_emax= uparam(iadbuf + 2*nfunc + 15)
189 numtabl = matparam%NTABLE
190C-----------------------------------------------
191C epst_spherique
192C-------------------
193 DO i=1,nel
194 epst(i) = epsxx(i)**2 + epsyy(i)**2 + epszz(i)**2
195 . + half*(epsxy(i)**2+epsyz(i)**2+epszx(i)**2)
196 epst(i) = sqrt(epst(i))
197 eps0(i) = uvar(i,1)
198 sig0(i) = uvar(i,2)
199 ENDDO
200C-------------------
201c Criterion
202C-------------------
203 nindex_epst = 0
204 nindex_epssm = 0
205c
206c static yield at epssmax
207 DO i=1,nel
208 ipos(i,1) = vartmp(i,1)
209 ipos(i,2) = 1
210 xvec(i,1) = epst(i)
211 xvec(i,2) = epsp0
212 ENDDO
213c
214 CALL table_mat_vinterp(matparam%TABLE(1),nel,nel,ipos,xvec,yld,dfc)
215 vartmp(1:nel,1)=ipos(1:nel,1)
216c
217 DO i=1,nel
218 IF (epst(i) >= epssmax) THEN
219 nindex_epssm = nindex_epssm + 1
220 index_epssm(nindex_epssm) = i
221 yldelas(i) = yld_emax + emax*(epst(i) - epssmax)
222 ELSE
223 nindex_epst = nindex_epst + 1
224 index_epst(nindex_epst) = i
225 yldelas(i) = yld(i)
226 ENDIF
227 ENDDO
228c=======================================================================
229c Loading case
230c=======================================================================
231 IF (nload > 1) THEN ! strain rate dependent
232 epspmax = matparam%TABLE(1)%X(2)%VALUES(nload) ! max strain rate
233#include "vectorize.inc"
234 DO i=1,nel
235 ipos(i,1) = vartmp(i,2)
236 ipos(i,2) = vartmp(i,5)
237 xvec(i,1) = min(epst(i), epssmax)
238 xvec(i,2) = epsp(i)
239c XVEC(I,2) = MIN(EPSP(I), EPSPMAX)
240c XVEC(I,2) = MAX(EPSP(I), EPSP0)
241 ENDDO
242c
243 CALL table_mat_vinterp(matparam%TABLE(1),nel,nel,ipos,xvec,yld,dfc)
244c
245 ELSE IF (nload == 1) THEN ! static only
246c
247 ipos(:,1) = vartmp(:,2)
248 ipos(:,2) = vartmp(:,5)
249 DO i=1,nel
250 xvec(i,1) = min(epst(i), epssmax)
251 ENDDO
252c
253 CALL table_mat_vinterp(matparam%TABLE(1),nel,nel,ipos,xvec,yld,dfc)
254c
255 END IF
256c
257 vartmp(1:nel,2) = ipos(1:nel,1)
258 vartmp(1:nel,5) = ipos(1:nel,2)
259
260 yldmax(:) = yld(:)
261c-----------------------------------------------------------------------
262c case EPST(I) >= EPSSMAX
263c-----------------------------------------------------------------------
264 IF (nindex_epssm > 0) THEN
265#include "vectorize.inc"
266 DO ii=1,nindex_epssm
267 i = index_epssm(ii)
268 yldmax(i) = yldmax(i) + emax*(epst(i) - epssmax)
269 ENDDO
270 END IF
271c=======================================================================
272c Unloading case
273c=======================================================================
274 IF (nunload > 0) THEN
275 IF (nunload > 1) THEN
276 epspmax = matparam%TABLE(1)%X(2)%VALUES(nunload) ! max strain rate
277#include "vectorize.inc"
278 DO i=1,nel
279 ipos(i,1) = vartmp(i,3)
280 ipos(i,2) = vartmp(i,6)
281 xvec(i,1) = min(epst(i), epssmax)
282 xvec(i,2) = epsp(i) ! MAX(MIN(EPSP(I), EPSPMAX), EPSP0)
283 ENDDO
284c
285 CALL table_mat_vinterp(matparam%TABLE(2),nel,nel,ipos,xvec,yld,dfc)
286c
287 ELSE IF (nunload == 1) THEN
288 ipos(:,1) = vartmp(:,3)
289 ipos(:,2) = vartmp(:,6)
290 DO i=1,nel
291 xvec(i,1) = epst(i) ! MAX(MIN(EPST(I), EPSSMAX),EPSP0)
292 ENDDO
293c
294 CALL table_mat_vinterp(matparam%TABLE(2),nel,nel,ipos,xvec,yld,dfc)
295c
296 END IF
297
298 vartmp(1:nel,3)=ipos(1:nel,1)
299 vartmp(1:nel,6)=ipos(1:nel,2)
300c
301 yldmin(:) = yld(:)
302c-----------------------------------------------------------------------
303c case EPST(I) >= EPSSMAX
304c-----------------------------------------------------------------------
305 IF (nindex_epssm > 0) THEN
306#include "vectorize.inc"
307 DO ii=1,nindex_epssm
308 i = index_epssm(ii)
309 yldmin(i) = yldmin(i) + emax*(epst(i) - epssmax)
310 ENDDO
311 END IF
312c
313 END IF ! UNLOAD
314c=======================================================================
315c end loading / unloading
316c=======================================================================
317 IF (itens > 0) THEN
318 DO i=1,nel
319 alpha = one/max(em20,uvar(i,10))
320 sig0xx(i) = alpha*sig0xx(i)
321 sig0yy(i) = alpha*sig0yy(i)
322 sig0zz(i) = alpha*sig0zz(i)
323 sig0xy(i) = alpha*sig0xy(i)
324 sig0zx(i) = alpha*sig0zx(i)
325 sig0yz(i) = alpha*sig0yz(i)
326 ENDDO
327 ENDIF
328C =====================================================
329 DO i = 1,nel
330 iload0 = int(uvar(i,5))
331 ie_cst(i)= 0
332 delta = epst(i) - uvar(i,4)
333 e = uvar(i,3)
334C
335 IF (delta >= zero ) THEN
336 yld(i) = yldmax(i)
337 iload(i) = 1
338 ELSEIF (delta < zero) THEN
339 yld(i) = yldmin(i)
340 iload(i) = -1
341 IF (idamage /= 0 ) THEN
342 yld(i) = yldmax(i)
343 ENDIF
344 ENDIF
345C
346 epss(i) = epst(i) - yld(i)/ e
347 epss(i) = max(zero, epss(i))
348 de = aa*(epss(i) - uvar(i,1))
349 IF (iload(i) == 1) THEN
350 e = e + max(de, zero)
351 IF(iload0 == -1) e= uvar(i,3)
352 uvar(i,1) = max(uvar(i,1), epss(i))
353 ELSE
354 e = e + min(de ,zero)
355 IF(iload0 == 1) e= uvar(i,3)
356 uvar(i,1) = min(epss(i),uvar(i,1))
357 ENDIF
358 e = min(e, emax)
359 e = max(e, e0)
360 uvar(i,3) = e
361C ----
362 aa1(i) = e*(one-nu)/(one + nu)/(one - two*nu)
363 aa2(i) = aa1(i)*nu/(one - nu)
364 g(i) = half*e/(one + nu)
365C ----
366 signxx(i)= aa1(i)*depsxx(i) + aa2(i)*(depsyy(i) + depszz(i))
367 signyy(i)= aa1(i)*depsyy(i) + aa2(i)*(depsxx(i) + depszz(i))
368 signzz(i)= aa1(i)*depszz(i) + aa2(i)*(depsxx(i) + depsyy(i))
369 signxy(i)= g(i) *depsxy(i)
370 signyz(i)= g(i) *depsyz(i)
371 signzx(i)= g(i) *depszx(i)
372 dsig = (signxx(i)**2 + signyy(i)**2 + signzz(i)**2)
373 . + two*(signxy(i)**2 + signzx(i)**2 + signyz(i)**2)
374 dsig =sqrt(dsig)
375C estimate stress
376 signxx(i)=sig0xx(i) + aa1(i)*depsxx(i)
377 . + aa2(i)*(depsyy(i) + depszz(i))
378 signyy(i)=sig0yy(i) + aa1(i)*depsyy(i)
379 . + aa2(i)*(depsxx(i) + depszz(i))
380 signzz(i)=sig0zz(i) + aa1(i)*depszz(i)
381 . + aa2(i)*(depsxx(i) + depsyy(i))
382 signxy(i)=sig0xy(i) + g(i) *depsxy(i)
383 signyz(i)=sig0yz(i) + g(i) *depsyz(i)
384 signzx(i)=sig0zx(i) + g(i) *depszx(i)
385
386 svm2 = (signxx(i)**2 + signyy(i)**2 + signzz(i)**2)
387 . + two*(signxy(i)**2 + signzx(i)**2 + signyz(i)**2)
388 svm =sqrt(svm2)
389C
390 IF (idamage == 0 ) THEN
391 IF(delta >= zero ) THEN
392 yld(i) = min(yldmax(i),svm)
393 ELSEIF(delta < zero) THEN
394 yld(i) = max(yldmin(i),svm)
395 yld(i) = min(yld(i),uvar(i,6))
396 ie_cst(i) = 1
397 IF (dsig > yldmin(i) .AND. dsig > svm)THEN
398 yld(i) = yldmin(i)
399 ie_cst(i) = 0
400 ENDIF
401 ENDIF
402 ELSE
403 yld(i) = yldmax(i)
404 uvar(i,8) = uvar(i,8) + half*(yldelas(i) + uvar(i,9))*delta
405 uvar(i,8) = max(zero, uvar(i,8))
406 uvar(i,2) = max(uvar(i,2) , uvar(i,8))
407 ENDIF
408 uvar(i,9) = yldelas(i)
409 ENDDO
410C
411C =====================================================
412C sound speed
413C =====================================================
414 IF (idtmins/=2.OR.jsms==0)THEN ! compatibilite ascendante
415 DO i = 1,nel
416 !-----------------------
417 ! P = C1 mu, mu = rho/rho0-1
418 ! d(rho)/d(mu) = rho0
419 !-----------------------
420 soundsp(i) = sqrt(aa1(i)/rho0(i))
421 rhosp(i) = rho0(i)
422 viscmax(i) = zero
423 ENDDO
424 ELSEIF (ismstr==1 .OR. ismstr==11)THEN ! Small Strain
425 DO i = 1,nel
426 !-----------------------
427 ! P = C1 mu, mu = rho/rho0-1
428 ! d(rho)/d(mu) = rho0
429 !-----------------------
430 soundsp(i) = sqrt(aa1(i)/rhoref(i))
431 rhosp(i) = rhoref(i)
432 viscmax(i) = zero
433 ENDDO
434 ELSE
435 DO i = 1,nel
436 IF (abs(offg(i)) <= one)THEN ! Large Strain
437 !-----------------------
438 ! P = C1 mu, mu = ln(rho/rho0), rho= rho0 exp(mu)
439 ! d(rho)/d(mu) = rho
440 !-----------------------
441 soundsp(i) = sqrt(aa1(i)/rho(i))
442 rhosp(i) = rho(i)
443 ELSE ! Small Strain
444 !-----------------------
445 ! P = C1 mu, mu = rho/rho0-1
446 ! d(rho)/d(mu) = rho0
447 !-----------------------
448 soundsp(i) = sqrt(aa1(i)/rhoref(i))
449 rhosp(i) = rhoref(i)
450 END IF
451 viscmax(i) = zero
452 ENDDO
453 END IF
454C
455C-------------------
456C projection spherique
457C-------------------
458 IF (idamage == 0 ) THEN
459 DO i=1,nel
460 signxx(i)= aa1(i)*epsxx(i) + aa2(i)*(epsyy(i) + epszz(i))
461 signyy(i)= aa1(i)*epsyy(i) + aa2(i)*(epsxx(i) + epszz(i))
462 signzz(i)= aa1(i)*epszz(i) + aa2(i)*(epsxx(i) + epsyy(i))
463 signxy(i)= g(i) *epsxy(i)
464 signyz(i)= g(i) *epsyz(i)
465 signzx(i)= g(i) *epszx(i)
466C
467 svm2 = (signxx(i)**2 + signyy(i)**2 + signzz(i)**2)
468 . + two*(signxy(i)**2 + signzx(i)**2 + signyz(i)**2)
469 svm =sqrt(svm2)
470 r = yld(i)/max(em20,svm)
471 et(i) = r
472 signxx(i)=signxx(i)*r
473 signyy(i)=signyy(i)*r
474 signzz(i)=signzz(i)*r
475 signxy(i)=signxy(i)*r
476 signyz(i)=signyz(i)*r
477 signzx(i)=signzx(i)*r
478c E = E0(I) + AA(I)*EPS0(I)
479c EPSS(I) = EPST(I) - YLD(I)/E
480c EPSS(I) = MAX(EPSS(I), ZERO)
481c UVAR(I,1) = EPSS(I)
482c E = E0(I) + AA(I)*EPSS(I)
483c UVAR(I,3) = E
484C
485 IF (ie_cst(i) == 1) THEN
486 iload0 = int(uvar(i,5))
487 IF(iload0 /= iload(i))THEN
488 iload(i) = iload0
489 uvar(i,1) = eps0(i)
490 ENDIF
491 ENDIF
492 uvar(i,4) = epst(i)
493 uvar(i,2) = svm*r
494 uvar(i,5) = iload(i)
495 uvar(i,6) = yld(i)
496 uvar(i,7) = epsp(i)
497 ENDDO
498c
499 ELSEIF (idamage == 1) THEN
500c
501 DO i=1,nel
502 r = one
503 signxx(i)= aa1(i)*epsxx(i) + aa2(i)*(epsyy(i) + epszz(i))
504 signyy(i)= aa1(i)*epsyy(i) + aa2(i)*(epsxx(i) + epszz(i))
505 signzz(i)= aa1(i)*epszz(i) + aa2(i)*(epsxx(i) + epsyy(i))
506 signxy(i)= g(i) *epsxy(i)
507 signyz(i)= g(i) *epsyz(i)
508 signzx(i)= g(i) *epszx(i)
509C
510 svm2 = (signxx(i)**2 + signyy(i)**2 + signzz(i)**2)
511 . + two*(signxy(i)**2 + signzx(i)**2 + signyz(i)**2)
512 svm =sqrt(svm2)
513 r = (yld(i)/max(em20,svm))
514 et(i) = r
515 signxx(i)=signxx(i)*r
516 signyy(i)=signyy(i)*r
517 signzz(i)=signzz(i)*r
518 signxy(i)=signxy(i)*r
519 signyz(i)=signyz(i)*r
520 signzx(i)=signzx(i)*r
521C
522 IF (iload(i) == -1) THEN
523 r = (yldmin(i)/max(em20,yldelas(i)))
524 p = third*(signxx(i) + signyy(i) + signzz(i))
525 signxx(i)=(signxx(i) - p)*r + p
526 signyy(i)=(signyy(i) - p)*r + p
527 signzz(i)=(signzz(i) - p)*r + p
528 signxy(i)=signxy(i)*r
529 signyz(i)=signyz(i)*r
530 signzx(i)=signzx(i)*r
531 ENDIF
532C
533 svm2 = (signxx(i)**2 + signyy(i)**2 + signzz(i)**2)
534 . + two*(signxy(i)**2 + signzx(i)**2 + signyz(i)**2)
535 svm =sqrt(svm2)
536 uvar(i,4) = epst(i)
537 uvar(i,5) = iload(i)
538 uvar(i,6) = yld(i)
539 uvar(i,7) = epsp(i)
540 ENDDO
541c
542 ELSEIF(idamage == 2 ) THEN
543c
544 DO i=1,nel
545 signxx(i)= aa1(i)*epsxx(i) + aa2(i)*(epsyy(i) + epszz(i))
546 signyy(i)= aa1(i)*epsyy(i) + aa2(i)*(epsxx(i) + epszz(i))
547 signzz(i)= aa1(i)*epszz(i) + aa2(i)*(epsxx(i) + epsyy(i))
548 signxy(i)= g(i) *epsxy(i)
549 signyz(i)= g(i) *epsyz(i)
550 signzx(i)= g(i) *epszx(i)
551C
552 svm2 = (signxx(i)**2 + signyy(i)**2 + signzz(i)**2)
553 . + two*(signxy(i)**2 + signzx(i)**2 + signyz(i)**2)
554 svm =sqrt(svm2)
555 r = yld(i)/max(em20,svm)
556 et(i) = r
557 signxx(i)=signxx(i)*r
558 signyy(i)=signyy(i)*r
559 signzz(i)=signzz(i)*r
560 signxy(i)=signxy(i)*r
561 signyz(i)=signyz(i)*r
562 signzx(i)=signzx(i)*r
563C
564 IF (iload(i) == -1) THEN
565 r = yldmin(i)/max(em20,yldelas(i))
566 signxx(i)=signxx(i)*r
567 signyy(i)=signyy(i)*r
568 signzz(i)=signzz(i)*r
569 signxy(i)=signxy(i)*r
570 signyz(i)=signyz(i)*r
571 signzx(i)=signzx(i)*r
572 et(i) = r*et(i)
573 ENDIF
574C
575 svm2 = (signxx(i)**2 + signyy(i)**2 + signzz(i)**2)
576 . + two*(signxy(i)**2 + signzx(i)**2 + signyz(i)**2)
577 svm =sqrt(svm2)
578 uvar(i,4) = epst(i)
579 uvar(i,5) = iload(i)
580 uvar(i,6) = yld(i)
581 uvar(i,7) = epsp(i)
582 ENDDO
583c
584 ELSEIF (idamage == 3) THEN
585c
586 DO i=1,nel
587 signxx(i)= aa1(i)*epsxx(i) + aa2(i)*(epsyy(i) + epszz(i))
588 signyy(i)= aa1(i)*epsyy(i) + aa2(i)*(epsxx(i) + epszz(i))
589 signzz(i)= aa1(i)*epszz(i) + aa2(i)*(epsxx(i) + epsyy(i))
590 signxy(i)= g(i) *epsxy(i)
591 signyz(i)= g(i) *epsyz(i)
592 signzx(i)= g(i) *epszx(i)
593C
594 svm2 = (signxx(i)**2 + signyy(i)**2 + signzz(i)**2)
595 . + two*(signxy(i)**2 + signzx(i)**2 + signyz(i)**2)
596 svm =sqrt(svm2)
597 r = (yld(i)/max(em20,svm))
598 signxx(i)=signxx(i)*r
599 signyy(i)=signyy(i)*r
600 signzz(i)=signzz(i)*r
601 signxy(i)=signxy(i)*r
602 signyz(i)=signyz(i)*r
603 signzx(i)=signzx(i)*r
604 et(i) = r
605C
606 IF (iload(i)==-1 .AND.uvar(i,2)/=zero)THEN
607 r = one - (uvar(i,8)/uvar(i,2))**expo
608 r = one - (one - hys)*r
609 p = third*(signxx(i) + signyy(i) + signzz(i))
610 signxx(i)=(signxx(i) - p)*r + p
611 signyy(i)=(signyy(i) - p)*r + p
612 signzz(i)=(signzz(i) - p)*r + p
613 signxy(i)=signxy(i)*r
614 signyz(i)=signyz(i)*r
615 signzx(i)=signzx(i)*r
616 et(i) = r*et(i)
617 ENDIF
618C
619 uvar(i,4) = epst(i)
620 uvar(i,5) = iload(i)
621 uvar(i,6) = yld(i)
622 uvar(i,7) = epsp(i)
623 ENDDO
624c
625 ELSEIF (idamage == 4) THEN
626c
627 DO i=1,nel
628 signxx(i)= aa1(i)*epsxx(i) + aa2(i)*(epsyy(i) + epszz(i))
629 signyy(i)= aa1(i)*epsyy(i) + aa2(i)*(epsxx(i) + epszz(i))
630 signzz(i)= aa1(i)*epszz(i) + aa2(i)*(epsxx(i) + epsyy(i))
631 signxy(i)= g(i) *epsxy(i)
632 signyz(i)= g(i) *epsyz(i)
633 signzx(i)= g(i) *epszx(i)
634C
635 svm2 = (signxx(i)**2 + signyy(i)**2 + signzz(i)**2)
636 . + two*(signxy(i)**2 + signzx(i)**2 + signyz(i)**2)
637 svm =sqrt(svm2)
638 r = (yld(i)/max(em20,svm))
639 et(i) = r
640 signxx(i)=signxx(i)*r
641 signyy(i)=signyy(i)*r
642 signzz(i)=signzz(i)*r
643 signxy(i)=signxy(i)*r
644 signyz(i)=signyz(i)*r
645 signzx(i)=signzx(i)*r
646C
647 IF (iload(i)==-1 .AND.uvar(i,2)/= zero)THEN
648 r = one - (uvar(i,8)/uvar(i,2))**expo
649 r = one - (one - hys)*r
650 signxx(i)=signxx(i)*r
651 signyy(i)=signyy(i)*r
652 signzz(i)=signzz(i)*r
653 signxy(i)=signxy(i)*r
654 signyz(i)=signyz(i)*r
655 signzx(i)=signzx(i)*r
656 ENDIF
657 uvar(i,4) = epst(i)
658 uvar(i,5) = iload(i)
659 uvar(i,6) = yld(i)
660 uvar(i,7) = epsp(i)
661 ENDDO
662c
663 ENDIF ! case of IDAMAGE
664C-------------------------------------------------------
665 IF (itens > 0 ) THEN
666 alpha_1(1:nel) = one
667 nindex_mu = 0
668C
669 DO i=1,nel
670 mu = 1 - rho(i)/rho0(i)
671 IF (mu > 0) THEN
672 nindex_mu = nindex_mu + 1
673 index_mu(nindex_mu) = i
674 ipos1(nindex_mu) = vartmp(i,4)
675 iadp1(nindex_mu) = npf(ifunc(nfunc)) / 2 + 1
676 ilenp1(nindex_mu) = npf(ifunc(nfunc)+1) / 2 - iadp1(nindex_mu) - ipos1(nindex_mu)
677 tab_loc1(nindex_mu) = mu
678 ENDIF
679 ENDDO
680C
681 IF ( nindex_mu > 0 ) THEN
682
683 CALL vinter2(tf,iadp1,ipos1,ilenp1,nindex_mu,tab_loc1,dfc,yld)
684
685 DO ii=1,nindex_mu
686 i = index_mu(ii)
687 alpha_1(i) = yfacc*yld(ii)
688 alpha_1(i) = max(zero, alpha_1(i))
689 vartmp(i,4) = ipos1(ii)
690 ENDDO
691 ENDIF
692C
693 DO i=1,nel
694 signxx(i) = alpha_1(i)*signxx(i)
695 signyy(i) = alpha_1(i)*signyy(i)
696 signzz(i) = alpha_1(i)*signzz(i)
697 signxy(i) = alpha_1(i)*signxy(i)
698 signzx(i) = alpha_1(i)*signzx(i)
699 signyz(i) = alpha_1(i)*signyz(i)
700 uvar(i,10) = alpha_1(i)
701 et(i) = alpha_1(i)*et(i)
702 ENDDO
703 ENDIF ! ITENS > 0
704C-----------------------------------------------------
705 IF (impl_s > 0 .OR. ihet > 1) THEN
706 DO i=1,nel
707 et(i) = min(one , et(i))
708 ENDDO
709 ENDIF
710C------------------------------------
711 RETURN
712 END SUBROUTINE
subroutine sigeps70(nel, nuparam, nuvar, mfunc, kfunc, npf, tf, time, timestep, uparam, rho0, rho, offg, rhoref, nvartmp, vartmp, rhosp, 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, ipm, mat, epsp, et, ismstr, ihet, jsms, matparam)
Definition sigeps70.F:51
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
subroutine mulaw(lft, llt, nft, mtn, jcvt, pm, off, sig, eint, rho, vol, strain, gama, uvar, bufmat, tf, npf, imat, ngl, nuvar, nvartmp, vartmp, geo, pid, epsd, wxx, wyy, wzz, jsph, ssp, voln, vis, d1, d2, d3, d4, d5, d6, dvol, sold1, sold2, sold3, sold4, sold5, sold6, rx, ry, rz, sx, sy, sz, tx, ty, tz, ismstr, mfxx, mfxy, mfxz, mfyx, mfyy, mfyz, mfzx, mfzy, mfzz, ipm, isorth, nel, matparam)
Definition mulaw.F:54
subroutine table_mat_vinterp(table, dimx, nel, ipos, xx, yy, dydx)
subroutine vinter2(tf, iad, ipos, ilen, nel0, x, dydx, y)
Definition vinter.F:143