OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
sigeps43g.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!|| sigeps43g ../engine/source/materials/mat/mat043/sigeps43g.F
25!||--- called by ------------------------------------------------------
26!|| mulawglc ../engine/source/materials/mat_share/mulawglc.F
27!||--- calls -----------------------------------------------------
28!|| finter ../engine/source/tools/curve/finter.F
29!|| vinter ../engine/source/tools/curve/vinter.F
30!||====================================================================
31 SUBROUTINE sigeps43g(
32 1 NEL ,NUVAR ,NPF ,
33 2 TF ,TIME ,UPARAM ,RHO0 ,
34 3 AREA ,EINT ,THK0 ,
35 4 EPSPXX ,EPSPYY ,EPSPXY ,EPSPYZ ,EPSPZX ,
36 5 DEPSXX ,DEPSYY ,DEPSXY ,DEPSYZ ,DEPSZX ,
37 5 DEPBXX ,DEPBYY ,DEPBXY ,
38 6 EPSXX ,EPSYY ,EPSXY ,EPSYZ ,EPSZX ,
39 7 SIGOXX ,SIGOYY ,SIGOXY ,SIGOYZ ,SIGOZX ,
40 7 MOMOXX ,MOMOYY ,MOMOXY ,
41 8 SIGNXX ,SIGNYY ,SIGNXY ,SIGNYZ ,SIGNZX ,
42 8 MOMNXX ,MOMNYY ,MOMNXY ,
43 9 SIGVXX ,SIGVYY ,SIGVXY ,SIGVYZ ,SIGVZX ,
44 A SOUNDSP,VISCMAX,THK ,PLA ,UVAR ,
45 B OFF ,NGL ,IPM ,MAT ,ETSE ,
46 C GS ,SIGY ,SHF ,SEQ_OUTPUT,EPSP )
47C-----------------------------------------------
48C I M P L I C I T T Y P E S
49C-----------------------------------------------
50#include "implicit_f.inc"
51C-----------------------------------------------
52C G L O B A L P A R A M E T E R S
53C-----------------------------------------------
54#include "mvsiz_p.inc"
55#include "param_c.inc"
56#include "com01_c.inc"
57C-----------------------------------------------
58C I N P U T A R G U M E N T S
59C-----------------------------------------------
60 INTEGER NEL, NUVAR, NGL(NEL), MAT(NEL),IPM(NPROPMI,*)
61 my_real
62 . TIME,UPARAM(*),
63 . AREA(NEL),RHO0(NEL),EINT(NEL,2),
64 . THK0(NEL),PLA(NEL),
65 . EPSPXX(NEL),EPSPYY(NEL),
66 . EPSPXY(NEL),EPSPYZ(NEL),EPSPZX(NEL),
67 . DEPSXX(NEL),DEPSYY(NEL),DEPSXY(NEL),
68 . DEPBXX(NEL),DEPBYY(NEL),DEPBXY(NEL),
69 . DEPSYZ(NEL),DEPSZX(NEL),
70 . EPSXX(NEL) ,EPSYY(NEL) ,
71 . EPSXY(NEL) ,EPSYZ(NEL) ,EPSZX(NEL) ,
72 . SIGOXX(NEL),SIGOYY(NEL),SIGOXY(NEL),
73 . MOMOXX(NEL),MOMOYY(NEL),MOMOXY(NEL),
74 . SIGOYZ(NEL),SIGOZX(NEL),GS(NEL),SHF(NEL),SEQ_OUTPUT(NEL)
75C-----------------------------------------------
76C O U T P U T A R G U M E N T S
77C-----------------------------------------------
78 my_real
79 . signxx(nel),signyy(nel),signxy(nel),
80 . momnxx(nel),momnyy(nel),momnxy(nel),
81 . signyz(nel),signzx(nel),epsp(nel),
82 . sigvxx(nel),sigvyy(nel),sigy(nel),
83 . sigvxy(nel),sigvyz(nel),sigvzx(nel),
84 . soundsp(nel),viscmax(nel),etse(nel)
85C-----------------------------------------------
86C I N P U T O U T P U T A R G U M E N T S
87C-----------------------------------------------
88 my_real uvar(nel,nuvar), off(nel),thk(nel)
89C-----------------------------------------------
90C VARIABLES FOR FUNCTION INTERPOLATION
91C-----------------------------------------------
92 INTEGER NPF(*)
93 my_real finter ,tf(*)
94C-----------------------------------------------
95C L O C A L V A R I A B L E S
96C-----------------------------------------------
97 INTEGER I,J,K,J1,J2,N,NINDX,NMAX,IADBUF,NFUNC,NS,
98 . NRATE,IAD1(MVSIZ),IPOS1(MVSIZ),ILEN1(MVSIZ),
99 . IAD2(MVSIZ),IPOS2(MVSIZ),ILEN2(MVSIZ),OPTE,
100 . JJ(MVSIZ),INDEX(MVSIZ),IFUNC(100), IFUNCE,MX,ISRATE
101 my_real
102 . E(MVSIZ),NU,A,B,C,FAC,DEZZ,S1,S2,S12,S3,
103 . dpla,epst,a1(mvsiz),a2(mvsiz),g(mvsiz),g3(mvsiz),
104 . dydx1(mvsiz),dydx2(mvsiz),rate(mvsiz,2),svm(mvsiz),
105 . yld(mvsiz),y1(mvsiz),y2(mvsiz),dr,f1,f2,
106 . yfac(mvsiz,2),nnu1,nu1(mvsiz),
107 . nu2(mvsiz),nu3(mvsiz),nu4(mvsiz),nu5(mvsiz),dpla_i,
108 . fail(mvsiz),epsmax,epsr1,epsr2,
109 . err,f,df,yld_i,dpla_j(mvsiz),
110 . aa1,aa2,aa3,fn(3),fm(3),fnm(3),sn1,sn2,sm1,sm2,
111 . pn(3),p_m(3),dfn(3),dfm(3),dfnm(3),djac(3),jac(3),c1,
112 . xp(3),pnm1(3),pnm2(3),jq(mvsiz),jq2(mvsiz),s(mvsiz),
113 . gm(mvsiz),cm(mvsiz),qtier(mvsiz),cnm(mvsiz),
114 . num1(mvsiz),num2(mvsiz),an(mvsiz),bn(mvsiz),
115 . am1(mvsiz),am2(mvsiz),gama(mvsiz),gama2(mvsiz),
116 . h(mvsiz),einf,ce,dydxe(mvsiz),
117 . escale(mvsiz)
118 my_real
119 . a_1,a_2,a_3,amn_xy(mvsiz),jac_i(3),jac_2(3),
120 . a01,a02,a03,a12,
121 . anxx(mvsiz),anyy(mvsiz),anxy(mvsiz),an_xy(mvsiz),
122 . amxx(mvsiz),amyy(mvsiz),amxy(mvsiz),am_xy(mvsiz),
123 . anmxx(mvsiz),anmyy(mvsiz),anmxy(mvsiz),anm_xy(mvsiz),
124 . sn11(mvsiz),sn22(mvsiz),sm11(mvsiz),sm22(mvsiz),
125 . b_1(mvsiz),b_2(mvsiz),b_3(mvsiz),q12(mvsiz),q21(mvsiz),
126 . sfn,sfm,sfnm,dsfn,dsfm,dsfnm,tol
127C
128 DATA nmax/3/,ns/10/
129 tol =em6
130C-----------------------------------------------
131C USER VARIABLES INITIALIZATION
132C-----------------------------------------------
133 mx = mat(1)
134 nfunc = ipm(10,mx)
135 DO j=1,nfunc
136 ifunc(j) = ipm(10+j,mx)
137 ENDDO
138 iadbuf = ipm(7,mx)-1
139 israte = ipm(3,mx)
140 nu = uparam(iadbuf+6)
141 a01 = uparam(iadbuf+7)
142 a02 = uparam(iadbuf+8)
143 a03 = uparam(iadbuf+9)
144 a12 = uparam(iadbuf+10)
145 nrate = nint(uparam(iadbuf+1))
146 epsmax=uparam(iadbuf+ns+2*nrate+1)
147C
148 IF(epsmax==zero)THEN
149 IF(tf(npf(ifunc(1)+1)-1)==zero)THEN
150 epsmax=tf(npf(ifunc(1)+1)-2)
151 ELSE
152 epsmax= ep30
153 ENDIF
154 ENDIF
155 nnu1 = nu / (one - nu)
156 epsr1 =uparam(iadbuf+ns+2*nrate+2)
157 epsr2 =uparam(iadbuf+ns+2*nrate+3)
158c---------------------
159 opte = uparam(iadbuf+ns+2*nrate + 10)
160 einf = uparam(iadbuf+ns+2*nrate+ 11)
161 ce = uparam(iadbuf+ns+2*nrate+ 12)
162c--------------------
163 DO i=1,nel
164 e(i) = uparam(iadbuf+2)
165 a1(i) = uparam(iadbuf+3)
166 a2(i) = uparam(iadbuf+4)
167 g(i) = uparam(iadbuf+5)
168 g3(i) = three*g(i)
169C 011 C1=THK(I)*DOUZ
170 c1=thk0(i)*one_over_12
171 am1(i) = a1(i)*c1
172 am2(i) = a2(i)*c1
173 gm(i) = g(i)*c1
174 ENDDO
175C
176 IF (opte == 1)THEN
177 ifunce = uparam(iadbuf+ns+2*nrate+ 9)
178 DO i=1,nel
179 IF(pla(i) > zero)THEN
180 escale(i) = finter(ifunc(ifunce),pla(i),npf,tf,dydxe(i))
181 e(i) = escale(i)* e(i)
182 g(i) = half*e(i)/(one+nu)
183 gs(i) = g(i)*shf(i)
184 g3(i) = three*g(i)
185 a1(i) = e(i)/(one - nu*nu)
186 a2(i) = nu*a1(i)
187 am1(i) = a1(i)*c1
188 am2(i) = a2(i)*c1
189 gm(i) = g(i)*c1
190 ENDIF
191 ENDDO
192 ELSEIF ( ce /= zero) THEN
193 DO i=1,nel
194 IF(pla(i) > zero)THEN
195 e(i) = e(i)-(e(i)-einf)*(one-exp(-ce*pla(i)))
196 g(i) = half*e(i)/(one+nu)
197 gs(i) = g(i)*shf(i)
198 g3(i) = three*g(i)
199 a1(i) = e(i)/(one - nu*nu)
200 a2(i) = nu*a1(i)
201 am1(i) = a1(i)*c1
202 am2(i) = a2(i)*c1
203 gm(i) = g(i)*c1
204 ENDIF
205 ENDDO
206 ENDIF
207C
208 IF (isigi==0) THEN
209 IF(time==0.0)THEN
210 DO i=1,nel
211 uvar(i,1)=zero
212 uvar(i,2)=zero
213 DO j=1,nrate
214 uvar(i,j+2)=zero
215 ENDDO
216 ENDDO
217 ENDIF
218 ENDIF
219C
220 DO i=1,nel
221 signxx(i)=sigoxx(i)+a1(i)*depsxx(i)+a2(i)*depsyy(i)
222 signyy(i)=sigoyy(i)+a2(i)*depsxx(i)+a1(i)*depsyy(i)
223 signxy(i)=sigoxy(i)+g(i) *depsxy(i)
224 momnxx(i)=momoxx(i)+am1(i)*depbxx(i)+am2(i)*depbyy(i)
225 momnyy(i)=momoyy(i)+am2(i)*depbxx(i)+am1(i)*depbyy(i)
226 momnxy(i)=momoxy(i)+gm(i) *depbxy(i)
227 signyz(i)=sigoyz(i)+gs(i) *depsyz(i)
228 signzx(i)=sigozx(i)+gs(i) *depszx(i)
229C
230 soundsp(i) = sqrt(a1(i)/rho0(i))
231 viscmax(i) = zero
232 etse(i) = one
233C-------------------
234C STRAIN RATE
235C-------------------
236 IF (israte == 0) THEN
237 epsp(i) = half*( abs(epspxx(i)+epspyy(i))
238 . + sqrt( (epspxx(i)-epspyy(i))*(epspxx(i)-epspyy(i))
239 . + epspxy(i)*epspxy(i) ) )
240 ENDIF
241C-------------------
242C STRAIN
243C-------------------
244 epst = half*( epsxx(i)+epsyy(i)
245 . + sqrt( (epsxx(i)-epsyy(i))*(epsxx(i)-epsyy(i))
246 . + epsxy(i)*epsxy(i) ) )
247 fail(i) = max(zero,min(one,(epsr2-epst)/(epsr2-epsr1)))
248C
249 ENDDO
250C-------------------
251C HARDENING LAW
252C-------------------
253 DO i=1,nel
254 jj(i) = 1
255 ENDDO
256 DO i=1,nel
257 iadbuf = ipm(7,mat(i))-1
258 DO j=2,nrate-1
259 IF(epsp(i)>=uparam(iadbuf+ns+j)) jj(i) = j
260 ENDDO
261 rate(i,1)=uparam(iadbuf+ns+jj(i))
262 rate(i,2)=uparam(iadbuf+ns+jj(i)+1)
263 yfac(i,1)=uparam(iadbuf+ns+nrate+jj(i))
264 yfac(i,2)=uparam(iadbuf+ns+nrate+jj(i)+1)
265 ENDDO
266 DO i=1,nel
267 j1 = jj(i)
268 j2 = j1+1
269 ipos1(i) = nint(uvar(i,j1))
270 iad1(i) = npf(ifunc(j1)) / 2 + 1
271 ilen1(i) = npf(ifunc(j1)+1) / 2 - iad1(i) - ipos1(i)
272 ipos2(i) = nint(uvar(i,j2))
273 iad2(i) = npf(ifunc(j2)) / 2 + 1
274 ilen2(i) = npf(ifunc(j2)+1) / 2 - iad2(i) - ipos2(i)
275 ENDDO
276C
277 CALL vinter(tf,iad1,ipos1,ilen1,nel,pla,dydx1,y1)
278 CALL vinter(tf,iad2,ipos2,ilen2,nel,pla,dydx2,y2)
279C
280 DO i=1,nel
281 j1 = jj(i)
282 j2 = j1+1
283 y1(i)=y1(i)*yfac(i,1)
284 y2(i)=y2(i)*yfac(i,2)
285 fac = (epsp(i) - rate(i,1))/(rate(i,2) - rate(i,1))
286 yld(i) = fail(i)*(y1(i) + fac*(y2(i)-y1(i)))
287 yld(i) = max(yld(i),em20)
288 sigy(i) = yld(i)
289 h(i) = fail(i)*(dydx1(i) + fac*(dydx2(i)-dydx1(i)))
290 uvar(i,j1) = ipos1(i)
291 uvar(i,j2) = ipos2(i)
292 ENDDO
293C-------------------------
294C HILL CRITERION
295C-------------------------
296 DO i=1,nel
297 c1=pla(i)*e(i)
298 gama(i)=three_half*(c1+yld(i))/(three_half*c1+yld(i))
299 gama2(i)=gama(i)*gama(i)
300 cm(i)= sixteen*gama2(i)
301 cnm(i)= sqr16_3*gama(i)
302 qtier(i)=four_over_3*gama2(i)
303 h(i) = max(zero,h(i))
304C
305 s1=a01*(signxx(i)*signxx(i)+cm(i)*momnxx(i)*momnxx(i))
306 s2=a02*(signyy(i)*signyy(i)+cm(i)*momnyy(i)*momnyy(i))
307 s3=a03*(signyy(i)*signxx(i)+cm(i)*momnxx(i)*momnyy(i))
308 anxy(i)=a12*signxy(i)*signxy(i)
309 amxy(i)=a12*momnxy(i)*momnxy(i)*cm(i)
310 f1=s1+s2-s3+anxy(i)+amxy(i)
311 s1=a01*(signxx(i)*momnxx(i))
312 s2=a02*(signyy(i)*momnyy(i))
313 s3=a03*(signxx(i)*momnyy(i)+signyy(i)*momnxx(i))*half
314 anmxy(i)=a12*signxy(i)*momnxy(i)
315 f2=cnm(i)*(s1+s2-s3+anmxy(i))
316 anmxy(i)=anmxy(i)*cnm(i)
317 svm(i)=sqrt(f1+abs(f2))
318!! SEQ_OUTPUT(I) = SVM(I)
319 dezz =-(depsxx(i)+depsyy(i))*nnu1
320 thk(i) = thk(i)+thk(i)*dezz*off(i)
321 ENDDO
322C-------------------------
323C GATHER PLASTIC FLOW
324C-------------------------
325 nindx=0
326 DO i=1,nel
327 IF(svm(i)>yld(i).AND.off(i)==one) THEN
328 nindx=nindx+1
329 index(nindx)=i
330 ENDIF
331 ENDDO
332 IF(nindx==0) RETURN
333C---------------------------
334C DEP EN CONTRAINTE PLANE
335C---------------------------
336 DO j=1,nindx
337 i=index(j)
338 nu2(i) = 1.-nu*nu
339 nu3(i) = nu*half
340 nu4(i) = half -nu3(i)
341 nu5(i) = one - nnu1
342 dpla_j(i)=(svm(i)-yld(i))/(g3(i)*qtier(i)+h(i))
343 etse(i)= h(i)/(h(i)+e(i))
344 s1=a01*nu*two - a03
345 s2=a02*nu*two - a03
346 s12=a03-nu*(a01+a02)
347 s3=sqrt(nu2(i)*(a01-a02)*(a01-a02)+s12*s12)
348 IF (abs(s1)<em20) THEN
349 q12(i)=zero
350 ELSE
351 q12(i)=-(a01-a02+s3)/s1
352 ENDIF
353 IF (abs(s2)<em20) THEN
354 q21(i)=zero
355 ELSE
356 q21(i)=(a01-a02+s3)/s2
357 ENDIF
358 jq(i)=one/(one - q12(i)*q21(i))
359 jq2(i)=jq(i)*jq(i)
360 a=a01*q12(i)
361 b=a02*q21(i)
362 a_1=(a01+a03*q21(i)+b*q21(i))*jq2(i)
363 a_2=(a02+a03*q12(i)+a*q12(i))*jq2(i)
364 a_3=(a+b)*jq2(i)*2.0+a03*(jq2(i)*2.0-jq(i))
365 sn11(i)=signxx(i)+signyy(i)*q12(i)
366 sn22(i)=q21(i)*signxx(i)+signyy(i)
367 sm11(i)=momnxx(i)+momnyy(i)*q12(i)
368 sm22(i)=q21(i)*momnxx(i)+momnyy(i)
369 anxx(i)=a_1*sn11(i)*sn11(i)
370 anyy(i)=a_2*sn22(i)*sn22(i)
371 an_xy(i)=a_3*sn11(i)*sn22(i)
372 amxx(i)=a_1*sm11(i)*sm11(i)*cm(i)
373 amyy(i)=a_2*sm22(i)*sm22(i)*cm(i)
374 am_xy(i)=a_3*sm11(i)*sm22(i)*cm(i)
375 anmxx(i)=a_1*sn11(i)*sm11(i)*cnm(i)
376 anmyy(i)=a_2*sn22(i)*sm22(i)*cnm(i)
377 anm_xy(i)=a_3*sn11(i)*sm22(i)*cnm(i)*half
378 amn_xy(i)=a_3*sm11(i)*sn22(i)*cnm(i)*half
379 a=a03*nu3(i)
380 b=s3*jq(i)
381 b_1(i)=a02-a-b
382 b_2(i)=a01-a+b
383 b_3(i)=a12*nu4(i)
384 ENDDO
385C-------------------------------
386C TIENT COMPTE DU COUPLAGE
387C-------------------------------
388 DO n=1,nmax
389 DO j=1,nindx
390 i=index(j)
391 dpla_i=dpla_j(i)
392 yld_i =yld(i)+h(i)*dpla_i
393 dr =a1(i)*dpla_i/yld_i
394 xp(1) =b_1(i)*dr
395 xp(2) =b_2(i)*dr
396 xp(3) =b_3(i)*dr
397 c1=one + qtier(i)
398 b=twop444*gama2(i)
399 DO k=1,3
400 djac(k)=c1+b*xp(k)
401 jac(k)=one+(djac(k)+c1)*xp(k)*half
402 jac_i(k)=one/jac(k)
403 jac_2(k)= jac_i(k)*jac_i(k)
404 a=xp(k)*zep444*gama2(i)
405 dfn(k)=fivep5+sixteenp5*a
406 fn(k)=one+(dfn(k)+fivep5)*a*half
407 a=onep8333*xp(k)
408 dfm(k)=onep8333+a
409 fm(k)=one +(dfm(k)*xp(k)+a)*half
410 dfnm(k)=-b*xp(k)
411 fnm(k)=one+dfnm(k)*xp(k)*half
412 ENDDO
413 a=jac_2(1)*fn(1)
414 b=jac_2(2)*fn(2)
415 c=jac_2(3)*fn(3)
416 sfn=a*(anxx(i)-an_xy(i))+b*(anyy(i)-an_xy(i))+c*anxy(i)
417 a=jac_2(1)*fm(1)
418 b=jac_2(2)*fm(2)
419 c=jac_2(3)*fm(3)
420 sfm=a*(amxx(i)-am_xy(i))+b*(amyy(i)-am_xy(i))+c*amxy(i)
421 a=jac_2(1)*fnm(1)
422 b=jac_2(2)*fnm(2)
423 c=jac_2(3)*fnm(3)
424 sfnm=a*(anmxx(i)-anm_xy(i))+b*(anmyy(i)-amn_xy(i))
425 . +c*anmxy(i)
426 c1=abs(sfnm)/max(sfn,sfm)
427 IF (c1<tol) THEN
428 s(i)=zero
429 ELSEIF (sfnm<zero) THEN
430 s(i)=-one
431 ELSE
432 s(i)=one
433 ENDIF
434 f =sfn+sfm+s(i)*sfnm-yld_i*yld_i
435 c1=zep444*gama2(i)
436 s12=jac_2(1)*b_1(i)
437 s1=s12*c1
438 aa1=two*jac_i(1)*djac(1)*s12
439 s12=jac_2(2)*b_2(i)
440 s2=s12*c1
441 aa2=two*jac_i(2)*djac(2)*s12
442 s12=jac_2(3)*b_3(i)
443 s3=s12*c1
444 aa3=two*jac_i(3)*djac(3)*s12
445 a=s1*dfn(1)-aa1*fn(1)
446 b=s2*dfn(2)-aa2*fn(2)
447 c=s3*dfn(3)-aa3*fn(3)
448 dsfn=a*(anxx(i)-an_xy(i))+b*(anyy(i)-an_xy(i))+c*anxy(i)
449 a=s1*dfm(1)-aa1*fm(1)
450 b=s2*dfm(2)-aa2*fm(2)
451 c=s3*dfm(3)-aa3*fm(3)
452 dsfm=a*(amxx(i)-am_xy(i))+b*(amyy(i)-am_xy(i))+c*amxy(i)
453 a=s1*dfnm(1)-aa1*fnm(1)
454 b=s2*dfnm(2)-aa2*fnm(2)
455 c=s3*dfnm(3)-aa3*fnm(3)
456 dsfnm=a*(anmxx(i)-anm_xy(i))+b*(anmyy(i)-amn_xy(i))
457 . +c*anmxy(i)
458 df =(dsfn+dsfm+s(i)*dsfnm)*
459 . (a1(i)-dr*h(i))/yld_i-two*h(i)*yld_i
460 dpla_j(i)=max(zero,dpla_i-f/df)
461 ENDDO
462 ENDDO
463C------------------------------------------
464C CONTRAINTES PLASTIQUEMENT ADMISSIBLES
465C------------------------------------------
466 DO j=1,nindx
467 i=index(j)
468 pla(i) = pla(i) + dpla_j(i)
469 yld_i =yld(i)+h(i)*dpla_j(i)
470 dr =a1(i)*dpla_j(i)/yld_i
471 xp(1) =b_1(i)*dr
472 xp(2) =b_2(i)*dr
473 xp(3) =b_3(i)*dr
474 c1=one+qtier(i)
475 b=twop444*gama2(i)
476 DO k=1,3
477 a=b*xp(k)
478 jac(k)=1.+c1*xp(k)+a*xp(k)
479 jac_i(k)=one/jac(k)
480 jac_2(k)= jac_i(k)*jac_i(k)
481 fnm(k)=one-a*xp(k)*half
482 a=xp(k)*jac_i(k)
483 pn(k)=jac_i(k)+qtier(i)*a
484 p_m(k)=jac_i(k)+a
485 pnm1(k)=-sqr4_3*gama(i)*a
486 pnm2(k)=pnm1(k)*one_over_12
487 ENDDO
488 a=jac_2(1)*fnm(1)
489 b=jac_2(2)*fnm(2)
490 c=jac_2(3)*fnm(3)
491 sfnm=a*(anmxx(i)-anm_xy(i))+b*(anmyy(i)-amn_xy(i))
492 . +c*anmxy(i)
493 sn1=sn11(i)*pn(1)+sm11(i)*pnm1(1)*s(i)
494 sn2=sn22(i)*pn(2)+sm22(i)*pnm1(2)*s(i)
495 s3=signxy(i)*pn(3)+momnxy(i)*pnm1(3)*s(i)
496 sm1=sm11(i)*p_m(1)+sn11(i)*pnm2(1)*s(i)
497 sm2=sm22(i)*p_m(2)+sn22(i)*pnm2(2)*s(i)
498 momnxy(i)=momnxy(i)*p_m(3)+signxy(i)*pnm2(3)*s(i)
499 signxx(i)=jq(i)*(sn1-sn2*q12(i))
500 signyy(i)=jq(i)*(sn2-sn1*q21(i))
501 signxy(i)=s3
502 momnxx(i)=jq(i)*(sm1-sm2*q12(i))
503 momnyy(i)=jq(i)*(sm2-sm1*q21(i))
504 s1=a01*signxx(i)+a02*signyy(i)
505 . -a03*(signxx(i)+signyy(i))*half
506 dezz = - nu5(i)*dpla_j(i)*s1/yld_i
507 thk(i) = thk(i)+thk(i)* dezz*off(i)
508 ENDDO
509C
510 DO i=1,nel
511 IF(pla(i)>epsmax.AND.off(i)==one)off(i)=four_over_5
512 ENDDO
513C
514 RETURN
515 END
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
subroutine sigeps43g(nel, nuvar, npf, tf, time, uparam, rho0, area, eint, thk0, epspxx, epspyy, epspxy, epspyz, epspzx, depsxx, depsyy, depsxy, depsyz, depszx, depbxx, depbyy, depbxy, epsxx, epsyy, epsxy, epsyz, epszx, sigoxx, sigoyy, sigoxy, sigoyz, sigozx, momoxx, momoyy, momoxy, signxx, signyy, signxy, signyz, signzx, momnxx, momnyy, momnxy, sigvxx, sigvyy, sigvxy, sigvyz, sigvzx, soundsp, viscmax, thk, pla, uvar, off, ngl, ipm, mat, etse, gs, sigy, shf, seq_output, epsp)
Definition sigeps43g.F:47
subroutine vinter(tf, iad, ipos, ilen, nel, x, dydx, y)
Definition vinter.F:72