OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
sigeps71.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!|| sigeps71 ../engine/source/materials/mat/mat071/sigeps71.F
25!||--- called by ------------------------------------------------------
26!|| mulaw ../engine/source/materials/mat_share/mulaw.F90
27!||--- calls -----------------------------------------------------
28!|| valpvec_v ../engine/source/materials/mat/mat033/sigeps33.F
29!|| valpvecdp_v ../engine/source/materials/mat/mat033/sigeps33.F
30!||====================================================================
31 SUBROUTINE sigeps71(
32 1 NEL ,NUPARAM,NUVAR ,TIME ,TIMESTEP,UPARAM ,
33 2 RHO0 ,RHO ,VOLUME ,EINT ,
34 3 EPSPXX ,EPSPYY ,EPSPZZ ,EPSPXY ,EPSPYZ ,EPSPZX ,
35 4 DEPSXX ,DEPSYY ,DEPSZZ ,DEPSXY ,DEPSYZ ,DEPSZX ,
36 5 EPSXX ,EPSYY ,EPSZZ ,EPSXY ,EPSYZ ,EPSZX ,
37 6 SIGOXX ,SIGOYY ,SIGOZZ ,SIGOXY ,SIGOYZ ,SIGOZX ,
38 7 SIGNXX ,SIGNYY ,SIGNZZ ,SIGNXY ,SIGNYZ ,SIGNZX ,
39 8 SIGVXX ,SIGVYY ,SIGVZZ ,SIGVXY ,SIGVYZ ,SIGVZX ,
40 9 SOUNDSP,VISCMAX,UVAR ,OFF ,NGL ,IPM ,
41 A MAT ,JTHE ,TEMP ,ISMSTR ,ETSE)
42C-----------------------------------------------
43c Law for SMA (Shape memory alloy - NiTinol)
44c based on Auricchio 1997
45C-----------------------------------------------
46C I m p l i c i t T y p e s
47C-----------------------------------------------
48#include "implicit_f.inc"
49#include "comlock.inc"
50C-----------------------------------------------
51C G l o b a l P a r a m e t e r s
52C-----------------------------------------------
53#include "mvsiz_p.inc"
54C---------+---------+---+---+--------------------------------------------
55C VAR | SIZE |TYP| RW| DEFINITION
56C---------+---------+---+---+--------------------------------------------
57C NEL | 1 | I | R | SIZE OF THE ELEMENT GROUP NEL
58C NUPARAM | 1 | I | R | SIZE OF THE USER PARAMETER ARRAY
59C NUVAR | 1 | I | R | NUMBER OF USER ELEMENT VARIABLES
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 "param_c.inc"
94#include "scr05_c.inc"
95C---------+---------+---+---+--------------------------------------------
96C
97 INTEGER ,INTENT(IN) :: NEL,NUPARAM,NUVAR
98 INTEGER ,INTENT(IN) :: JTHE
99 INTEGER NGL(NEL),MAT(NEL),IPM(NPROPMI,*)
100 my_real
101 . TIME,TIMESTEP,UPARAM(NUPARAM),
102 . RHO(NEL),RHO0(NEL),VOLUME(NEL),EINT(NEL),
103 . EPSPXX(NEL),EPSPYY(NEL),EPSPZZ(NEL),
104 . EPSPXY(NEL),EPSPYZ(NEL),EPSPZX(NEL),
105 . DEPSXX(NEL),DEPSYY(NEL),DEPSZZ(NEL),
106 . DEPSXY(NEL),DEPSYZ(NEL),DEPSZX(NEL),
107 . epsxx(nel) ,epsyy(nel) ,epszz(nel) ,
108 . epsxy(nel) ,epsyz(nel) ,epszx(nel) ,
109 . sigoxx(nel),sigoyy(nel),sigozz(nel),
110 . sigoxy(nel),sigoyz(nel),sigozx(nel)
111 my_real, DIMENSION(NEL) ,INTENT(INOUT) :: temp
112C-----------------------------------------------
113C O U T P U T A r g u m e n t s
114C-----------------------------------------------
115 my_real
116 . signxx(nel),signyy(nel),signzz(nel),
117 . signxy(nel),signyz(nel),signzx(nel),
118 . sigvxx(nel),sigvyy(nel),sigvzz(nel),
119 . sigvxy(nel),sigvyz(nel),sigvzx(nel),
120 . soundsp(nel),etse(nel),viscmax(nel)
121C-----------------------------------------------
122C I N P U T O U T P U T A r g u m e n t s
123C-----------------------------------------------
124 my_real, DIMENSION(NEL) ,INTENT(INOUT) :: off
125 my_real, DIMENSION(NEL,NUVAR) ,INTENT(INOUT) ::uvar
126C-----------------------------------------------
127C L o c a l V a r i a b l e s
128C-----------------------------------------------
129 INTEGER I,J,J1,J2,I1,I2,KK,IADBUF,EFLAG,ISMSTR,IFUNC(100)
130 my_real
131 . E,EMART,NU,G,G2,WAVE,SQDT,A,B,C,FCT,FCTP, DFTR,UNMXN,DB,
132 . ALPHA,EPSL,AA1,FM,DFMSA,DFMAS,UXX,NN,BETA,GM,KM,
133 . CB,CC,CAAS,CBAS,POLD,DGT,DKT,CP,TINI,INVE, n1,n2,n3,
134 . k,p,sxx,cas,csa,tsas,tfas, tssa,tfsa,rv_pui,dfs,
135 . syy,szz,sxy,syz,szx,fass,fsas,fasf,fsaf,rsas,rfas,
136 . sv,fs,fs0,yld_ass,yld_asf,yld_sas,yld_saf,rssa,rfsa,
137 . dfm, fss,dsxx,dsyy,dsxy,dsyz,dszx,dszz,var,h,et,
138 . pm,delta,x1,x2,test,test2,ftest,gnew,knew,betan,
139 . nx2,ny2 ,nz2,nxy2,nyz2,nzx2,ne,dnx,dny,dnz,dnxy,dnyz,dnzx,
140 . nxx(nel),nyy(nel),nzz(nel),nxy(nel),nyz(nel),nzx(nel),
141 . e1(nel),e2(nel),e3(nel),e4(nel),e5(nel),e6(nel),trde(nel),
142 . de1(nel),de2(nel),de3(nel),gt(nel),kt(nel),
143 . ee1(nel),ee2(nel),ee3(nel),pp(nel),nne(nel),det(nel),
144 . sigxx(nel), sigyy(nel), sigzz(nel)
145 my_real
146 . ev(mvsiz,3)
147 my_real
148 . av(mvsiz,6),evv(mvsiz,3),dirprv(mvsiz,3,3)
149
150C-----------------------------------------------
151C USER VARIABLES INITIALIZATION
152C-----------------------------------------------
153
154 e = uparam(1)
155 nu = uparam(2)
156 g = uparam(3)
157 k = uparam(4)
158 aa1 = uparam(5)
159 yld_ass = uparam(6)
160 yld_asf = uparam(7)
161 yld_sas = uparam(8)
162 yld_saf = uparam(9)
163 alpha = uparam(10)
164 epsl = uparam(11)/(sqrt(two_third)+alpha)
165 emart = uparam(14)
166 eflag = int(uparam(15))
167 gm = uparam(16)
168 km = uparam(17)
169 g2 = two*g
170 beta = epsl*(g2+nine*k*alpha*alpha)
171 sqdt = sqrt(two/three)
172 cas = uparam(18)
173 csa = uparam(19)
174 tsas = uparam(20)
175 tfas = uparam(21)
176 tssa = uparam(22)
177 tfsa = uparam(23)
178 cp = uparam(24)
179 tini = uparam(25)
180c
181C
182 DO i=1,nel
183 av(i,1)=epsxx(i)
184 av(i,2)=epsyy(i)
185 av(i,3)=epszz(i)
186 av(i,4)=epsxy(i) * half
187 av(i,5)=epsyz(i) * half
188 av(i,6)=epszx(i) * half
189 ENDDO
190C Eigenvalues needed to be calculated in double precision
191C for a simple precision executing
192 IF (iresp==1) THEN
193 CALL valpvecdp_v(av,evv,dirprv,nel)
194 ELSE
195 CALL valpvec_v(av,evv,dirprv,nel)
196 ENDIF
197
198 IF(ismstr== 0.OR. ismstr==2.OR. ismstr==4) THEN
199 DO i=1,nel
200C ---- (STRAIN IS LOGARITHMIC)
201 ev(i,1)=exp(evv(i,1))
202 ev(i,2)=exp(evv(i,2))
203 ev(i,3)=exp(evv(i,3))
204 ENDDO
205 ELSEIF(ismstr==10.OR.ismstr==12) THEN
206 DO i =1,nel
207 ev(i,1)=sqrt(evv(i,1)+ one)! = lambda
208 ev(i,2)=sqrt(evv(i,2)+ one)
209 ev(i,3)=sqrt(evv(i,3)+ one)
210 ENDDO
211 ELSE
212C ---- STRAIN IS ENGINEERING)
213 DO i=1,nel
214 ev(i,1)=evv(i,1)+ one
215 ev(i,2)=evv(i,2)+ one
216 ev(i,3)=evv(i,3)+ one
217 ENDDO
218 ENDIF
219 DO i =1,nel
220 !DET A (A = GRADIENT OF DEFORMATION)
221 det(i) =ev(i,1)*ev(i,2)*ev(i,3)
222 IF(det(i)/=zero) THEN
223 trde(i) = log(det(i)) !TRACE OF DEFORMATION
224 rv_pui = exp((-third)*trde(i))
225 ELSE
226 rv_pui = zero
227 trde(i)= zero
228 ENDIF
229 ee1(i) = log(ev(i,1)*rv_pui) !DEVIATOR OF DEFORMATION
230 ee2(i) = log(ev(i,2)*rv_pui)
231 ee3(i) = log(ev(i,3)*rv_pui)
232 ENDDO
233C-----------------------------------------------
234C Temperature update in adiabatic conditions
235C--------------------
236 IF (jthe == 0) THEN
237 DO i=1,nel
238 temp(i) = tini + eint(i) / rho0(i)/cp/max(em15,volume(i))
239 ENDDO
240 ENDIF
241C-----------------------------------------------
242
243 IF (eflag > zero)THEN
244C=======================================================================
245 rsas = yld_ass *(sqdt+alpha)-cas*tsas
246 rfas = yld_asf *(sqdt+alpha)-cas*tfas
247 rssa = yld_sas *(sqdt+alpha)-csa*tssa
248 rfsa = yld_saf *(sqdt+alpha)-csa*tfsa
249 DO i = 1,nel
250C
251 fm = uvar(i,1) ! fraction of Martensite
252 gt(i) = g + fm * (gm - g) !G_n
253 kt(i) = k + fm * (km - k) !K_n
254 !pressure
255 p = kt(i) * (trde(i) - three*alpha*epsl*fm)
256 ! n= e/ norm(e)
257 ne = sqrt( ee1(i)**2 + ee2(i)**2 + ee3(i)**2)
258 nxx(i) =ee1(i)/max(ne,em20)
259 nyy(i) =ee2(i)/max(ne,em20)
260 nzz(i) =ee3(i)/max(ne,em20)
261! Estimation dev(sigma_n+1)
262 sxx= two*gt(i)*(ee1(i) -epsl*fm*nxx(i))
263 syy= two*gt(i)*(ee2(i) -epsl*fm*nyy(i))
264 szz= two*gt(i)*(ee3(i) -epsl*fm*nzz(i))
265
266 sv = sqrt( sxx*sxx + syy*syy + szz*szz )
267c sound velocity
268 soundsp(i) = sqrt(aa1/rho0(i))
269
270 viscmax(i) = zero
271C-------------------
272C transformation
273C-------------------
274 dfmsa = zero
275 dfmas = zero
276c loading function
277 fs = sv + three*alpha*p - cas*temp(i)
278C----------------------
279C Check Austenite -----> martensite
280 fass = fs - rsas
281 fasf = fs - rfas
282 fs0 = uvar(i,2)
283 beta = epsl*(two*gt(i)+nine*kt(i)*alpha*alpha)
284 IF((fs - fs0) > zero .AND. fass > zero.AND. fasf < zero .AND. fm < one )THEN
285 ! DFMAS = MIN(ONE, -(FS-FS0)*(ONE-FM)/(FASF-BETA*(ONE-FM) ) ) ! (should be positive)
286 db = (two * (gm-g) +nine*alpha*alpha*(km-k)) *epsl
287 unmxn = one - fm
288 dftr = two*ne*(gm-g) + three*alpha*trde(i)*(km-k)
289 dfmas = min(one, -(fs-fs0)*(one-fm)/(fasf-beta*(one-fm) ) )
290 a = unmxn *db
291 b = (rfas-fs+unmxn*(beta-dftr))
292 c = unmxn*(fs0 - fs)
293 DO kk = 1,3
294 fct = dfmas*dfmas *a+ dfmas* b +c
295 fctp = two*dfmas *a+ b
296 dfmas = dfmas - fct / fctp
297 ENDDO
298 dfmas = min(one,dfmas ) ! (should be positive)
299 ENDIF
300C----------------------
301C Check marteniste -----> austenite
302c Unloading function
303 fs = sv + three*alpha*p - csa*temp(i)
304 fsas = fs -rssa
305 fsaf = fs -rfsa
306 fs0 = uvar(i,3)
307
308 IF((fs - fs0) < zero .AND. fsas < zero.AND. fsaf > zero )THEN
309 !DFMSA = MAX(-ONE , FM * (FS - FS0)/ (FSAF+BETA*FM) )
310 db = (two * (gm-g) +nine*alpha*alpha*(km-k))*epsl
311 dftr = two * (gm-g)*ne+ three*alpha*(km-k)*trde(i)
312 dfmsa = zero
313 a = fm *db
314 b = -(rfsa-fs+fm*(dftr-beta))
315 c = -fm*(fs - fs0)
316 DO kk = 1,3
317 fct = dfmsa*dfmsa *a+ dfmsa* b +c
318 fctp = two*dfmsa *a+ b
319 dfmsa = dfmsa - fct / fctp
320 ENDDO
321 dfmsa = max(-one , dfmsa )
322 ENDIF
323C--------------------------------
324C new martensite fraction
325 dfm = dfmas + dfmsa
326 IF(dfm < zero .AND. fm == zero) dfm = zero
327C----------------------
328 !UPDATE
329C----------------------
330
331 dgt = dfm * (gm - g)
332 dkt = dfm * (km - k)
333
334 sxx = sxx -two*gt(i) * epsl*nxx(i)*dfm
335 . + two*dgt* (ee1(i)-epsl*nxx(i)*dfm)
336 syy = syy -two*gt(i)* epsl*dfm*nyy(i)
337 . + two*dgt* (ee2(i)-epsl*nyy(i)*dfm)
338 szz = szz -two*gt(i)* epsl*dfm*nzz(i)
339 . + two*dgt* (ee3(i)-epsl*nzz(i)*dfm)
340
341 p = p - kt(i)*epsl*three*alpha*dfm
342 . + dkt *(trde(i) -epsl*three*alpha*dfm)
343 !KIRCHHOFF STRESS
344 sigxx(i)= sxx + p
345 sigyy(i)= syy + p
346 sigzz(i)= szz + p
347 sv = sqrt( sxx*sxx + syy*syy + szz*szz )
348 fs = sv + three*alpha*p
349 !CAUCHY STRESS
350 IF(det(i)/=zero)THEN
351 inve = one/det(i)
352 ELSE
353 inve = zero
354 ENDIF
355 sigxx(i)= sigxx(i) *inve
356 sigyy(i)= sigyy(i) *inve
357 sigzz(i)= sigzz(i) *inve
358C TRANSFORM PRINCIPAL CAUCHY STRESSES TO GLOBAL DIRECTIONS
359 signxx(i) = dirprv(i,1,1)*dirprv(i,1,1)*sigxx(i)
360 . + dirprv(i,1,2)*dirprv(i,1,2)*sigyy(i)
361 . + dirprv(i,1,3)*dirprv(i,1,3)*sigzz(i)
362 signyy(i) = dirprv(i,2,2)*dirprv(i,2,2)*sigyy(i)
363 . + dirprv(i,2,3)*dirprv(i,2,3)*sigzz(i)
364 . + dirprv(i,2,1)*dirprv(i,2,1)*sigxx(i)
365 signzz(i) = dirprv(i,3,3)*dirprv(i,3,3)*sigzz(i)
366 . + dirprv(i,3,1)*dirprv(i,3,1)*sigxx(i)
367 . + dirprv(i,3,2)*dirprv(i,3,2)*sigyy(i)
368 signxy(i) = dirprv(i,1,1)*dirprv(i,2,1)*sigxx(i)
369 . + dirprv(i,1,2)*dirprv(i,2,2)*sigyy(i)
370 . + dirprv(i,1,3)*dirprv(i,2,3)*sigzz(i)
371 signyz(i) = dirprv(i,2,2)*dirprv(i,3,2)*sigyy(i)
372 . + dirprv(i,2,3)*dirprv(i,3,3)*sigzz(i)
373 . + dirprv(i,2,1)*dirprv(i,3,1)*sigxx(i)
374 signzx(i) = dirprv(i,3,3)*dirprv(i,1,3)*sigzz(i)
375 . + dirprv(i,3,1)*dirprv(i,1,1)*sigxx(i)
376 . + dirprv(i,3,2)*dirprv(i,1,2)*sigyy(i)
377c
378 uvar(i,1) = uvar(i,1) + dfm
379 uvar(i,1) = max(zero, uvar(i,1))
380 uvar(i,1) = min(one, uvar(i,1))
381 uvar(i,2) = fs- cas*temp(i)
382 uvar(i,3) = fs- csa*temp(i)
383 IF (dfmas /= zero) dfs = abs(uvar(i,2) - fs0)
384 IF (dfmsa /= zero) dfs = abs(uvar(i,3) - fs0)
385 IF (dfs /= zero .AND. epsl /= zero .AND. dfm/= zero) THEN
386 h = dfs/epsl/dfm
387 etse(i) = h *(one+nu) /( (e + uvar(i,1)*(emart-e)) + h)
388 ELSE
389 etse(i) = one
390 ENDIF
391 uvar(i,10) = epsxx(i)
392 ENDDO
393C=======================================================================
394 ELSE !EFLAG = 0
395 rsas = yld_ass *(sqdt+alpha)-cas*tsas
396 rfas = yld_asf *(sqdt+alpha)-cas*tfas
397 rssa = yld_sas *(sqdt+alpha)-csa*tssa
398 rfsa = yld_saf *(sqdt+alpha)-csa*tfsa
399 DO i = 1,nel
400C
401 fm = uvar(i,1) ! fraction of Martensite
402 !pressure
403 p = k * (trde(i) - three*alpha*epsl*fm)
404 ! n= e/ norm(e)
405 ne = sqrt( ee1(i)**2 + ee2(i)**2 + ee3(i)**2)
406 nxx(i) =ee1(i)/max(ne,em20)
407 nyy(i) =ee2(i)/max(ne,em20)
408 nzz(i) =ee3(i)/max(ne,em20)
409! Estimation dev(sigma_n+1)
410 sxx= g2*(ee1(i) -epsl*fm*nxx(i))
411 syy= g2*(ee2(i) -epsl*fm*nyy(i))
412 szz= g2*(ee3(i) -epsl*fm*nzz(i))
413
414 sv = sqrt( sxx*sxx + syy*syy + szz*szz )
415c sound velocity
416 soundsp(i) = sqrt(aa1/rho0(i))
417 viscmax(i) = zero
418C-------------------
419C transformation
420C-------------------
421 dfmsa = zero
422 dfmas = zero
423
424c loading function
425 fs = sv + three*alpha*p - cas*temp(i) ! loading function A--> M
426 !FSS = FS - YLD_ASS
427C----------------------
428C Check Austenite -----> martensite
429 fass = fs -rsas
430 fasf = fs -rfas
431 fs0 = uvar(i,2)
432 IF((fs - fs0) > zero .AND. fass > zero.AND. fasf < zero.AND. fm < one ) then!
433 dfmas = min(one, -(fs-fs0)*(one-fm)/(fasf-beta*(one-fm) ) ) ! (should be positive)
434 ENDIF
435C----------------------
436C Check marteniste -----> austenite
437 fs = sv + three*alpha*p - csa*temp(i) ! Unloading function M--> A
438 fsas = fs - rssa
439 fsaf = fs - rfsa
440 fs0 = uvar(i,3)
441 IF((fs - fs0) < zero .AND. fsas < zero .AND. fsaf > zero .AND. fm > zero ) THEN
442 dfmsa = max(-one , fm * (fs - fs0)/ (fsaf+beta*fm) ) !compute marteniste fraction
443 ENDIF
444C----------------------
445C new martensite fraction
446 dfm = dfmas + dfmsa
447 IF(dfm < zero .AND. fm == zero) dfm = zero
448C----------------------
449 !UPDATE
450 sxx = sxx -g2* epsl*dfm*nxx(i)
451 syy = syy -g2* epsl*dfm*nyy(i)
452 szz = szz -g2* epsl*dfm*nzz(i)
453 p = p - k*epsl*three*alpha*dfm
454 !KIRCHHOFF STRESS
455 sigxx(i)= sxx + p
456 sigyy(i)= syy + p
457 sigzz(i)= szz + p
458 sv = sqrt( sxx*sxx + syy*syy + szz*szz )
459 fs = sv + three*alpha*p
460 !CAUCHY STRESS
461 IF(det(i)/=zero)THEN
462 inve = one/det(i)
463 ELSE
464 inve = zero
465 ENDIF
466 sigxx(i)= sigxx(i) *inve
467 sigyy(i)= sigyy(i) *inve
468 sigzz(i)= sigzz(i) *inve
469C TRANSFORM PRINCIPAL CAUCHY STRESSES TO GLOBAL DIRECTIONS
470 signxx(i) = dirprv(i,1,1)*dirprv(i,1,1)*sigxx(i)
471 . + dirprv(i,1,2)*dirprv(i,1,2)*sigyy(i)
472 . + dirprv(i,1,3)*dirprv(i,1,3)*sigzz(i)
473 signyy(i) = dirprv(i,2,2)*dirprv(i,2,2)*sigyy(i)
474 . + dirprv(i,2,3)*dirprv(i,2,3)*sigzz(i)
475 . + dirprv(i,2,1)*dirprv(i,2,1)*sigxx(i)
476 signzz(i) = dirprv(i,3,3)*dirprv(i,3,3)*sigzz(i)
477 . + dirprv(i,3,1)*dirprv(i,3,1)*sigxx(i)
478 . + dirprv(i,3,2)*dirprv(i,3,2)*sigyy(i)
479 signxy(i) = dirprv(i,1,1)*dirprv(i,2,1)*sigxx(i)
480 . + dirprv(i,1,2)*dirprv(i,2,2)*sigyy(i)
481 . + dirprv(i,1,3)*dirprv(i,2,3)*sigzz(i)
482 signyz(i) = dirprv(i,2,2)*dirprv(i,3,2)*sigyy(i)
483 . + dirprv(i,2,3)*dirprv(i,3,3)*sigzz(i)
484 . + dirprv(i,2,1)*dirprv(i,3,1)*sigxx(i)
485 signzx(i) = dirprv(i,3,3)*dirprv(i,1,3)*sigzz(i)
486 . + dirprv(i,3,1)*dirprv(i,1,1)*sigxx(i)
487 . + dirprv(i,3,2)*dirprv(i,1,2)*sigyy(i)
488 uvar(i,1) = uvar(i,1) + dfm
489 uvar(i,1) = max(zero, uvar(i,1))
490 uvar(i,1) = min(one, uvar(i,1))
491 uvar(i,2) = fs- cas*temp(i)
492 uvar(i,3) = fs- csa*temp(i)
493 dfs = zero
494 IF (dfmas /= zero) dfs = abs(uvar(i,2) - fs0)
495 IF (dfmsa /= zero) dfs = abs(uvar(i,3) - fs0)
496 IF (dfs /= zero .AND. epsl /= zero.AND.dfm/= zero) THEN
497 h = dfs/epsl/dfm
498 etse(i) = h * (one+nu)/(e + h)
499 ELSE
500 etse(i) = one
501 ENDIF
502 uvar(i,4) = epsl*dfm ! transformation strain
503 uvar(i,7) = uvar(i,7)+epsl*dfm
504 uvar(i,10) = epsxx(i)
505 ENDDO
506C
507 ENDIF ! EFLAG
508 RETURN
509 END
510C
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
valpvecdp_v
subroutine valpvecdp_v(sig, val, vec, nel)
Definition sigeps33.F:2231
valpvec_v
subroutine valpvec_v(sig, val, vec, nel)
Definition sigeps33.F:1902
sigeps71
subroutine sigeps71(nel, nuparam, nuvar, 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, ipm, mat, jthe, temp, ismstr, etse)
Definition sigeps71.F:42