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

Go to the source code of this file.

Functions/Subroutines

subroutine asso_plas76 (nel, nuparam, nuvar, nfunc, ifunc, nvartmp, npf, tf, time, timestep, uparam, vartmp, rho0, pla, dpla, et, numtabl, table, depsxx, depsyy, depszz, depsxy, depsyz, depszx, sigoxx, sigoyy, sigozz, sigoxy, sigoyz, sigozx, signxx, signyy, signzz, signxy, signyz, signzx, soundsp, uvar, off, epsd, yld)

Function/Subroutine Documentation

◆ asso_plas76()

subroutine asso_plas76 ( integer nel,
integer nuparam,
integer nuvar,
integer nfunc,
integer, dimension(nfunc) ifunc,
integer, intent(in) nvartmp,
integer, dimension(snpc) npf,
tf,
intent(in) time,
intent(in) timestep,
dimension(nuparam), intent(in) uparam,
integer, dimension(nel,nvartmp), intent(inout) vartmp,
intent(in) rho0,
intent(inout) pla,
intent(out) dpla,
intent(out) et,
integer numtabl,
type(table_4d_), dimension(numtabl), intent(in) table,
intent(in) depsxx,
intent(in) depsyy,
intent(in) depszz,
intent(in) depsxy,
intent(in) depsyz,
intent(in) depszx,
intent(in) sigoxx,
intent(in) sigoyy,
intent(in) sigozz,
intent(in) sigoxy,
intent(in) sigoyz,
intent(in) sigozx,
intent(out) signxx,
intent(out) signyy,
intent(out) signzz,
intent(out) signxy,
intent(out) signyz,
intent(out) signzx,
intent(out) soundsp,
intent(inout) uvar,
intent(inout) off,
intent(inout) epsd,
intent(inout) yld )

Definition at line 34 of file asso_plas76.F.

42C-----------------------------------------------
43C M o d u l e s
44C-----------------------------------------------
45 USE table4d_mod
46 USE table_mat_vinterp_mod
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 "tabsiz_c.inc"
55C=======================================================================
56C I N P U T A r g u m e n t s
57C-----------------------------------------------
58 INTEGER NEL,NUPARAM,NUVAR,NFUNC,NUMTABL
59 INTEGER ,INTENT(IN) :: NVARTMP
60 my_real,INTENT(IN) :: time,timestep
61 my_real,DIMENSION(NUPARAM),INTENT(IN) :: uparam(nuparam)
62 my_real,DIMENSION(NEL),INTENT(IN) ::
63 . rho0,depsxx,depsyy,depszz,
64 . depsxy,depsyz,depszx,sigoxx,sigoyy,sigozz,
65 . sigoxy,sigoyz,sigozx
66C-----------------------------------------------
67C O U T P U T A r g u m e n t s
68C-----------------------------------------------
69 INTEGER ,DIMENSION(NEL,NVARTMP) ,INTENT(INOUT) :: VARTMP
70 my_real,DIMENSION(NEL),INTENT(OUT) ::
71 . signxx,signyy,signzz,signxy,signyz,signzx,
72 . soundsp,dpla,et
73C-----------------------------------------------
74C I N P U T O U T P U T A r g u m e n t s
75C-----------------------------------------------
76 my_real,DIMENSION(NEL),INTENT(INOUT) ::
77 . off,yld,pla,epsd
78 my_real,DIMENSION(NEL,NUVAR),INTENT(INOUT) ::
79 . uvar
80 TYPE(TABLE_4D_),DIMENSION(NUMTABL),INTENT(IN) :: TABLE
81C-----------------------------------------------
82C VARIABLES FOR FUNCTION INTERPOLATION
83C-----------------------------------------------
84 INTEGER :: NPF(SNPC), IFUNC(NFUNC)
85 my_real :: tf(stf)
86C-----------------------------------------------
87C L o c a l V a r i a b l e s
88C-----------------------------------------------
89 INTEGER :: I,J,II,ICONV,NINDX,ICAS,ITER
90 INTEGER ,PARAMETER :: FUNC_TRAC = 1
91 INTEGER ,PARAMETER :: FUNC_COMP = 2
92 INTEGER ,PARAMETER :: FUNC_SHEAR = 3
93 INTEGER ,DIMENSION(NEL) :: INDX,IAD,ILEN
94 my_real :: g,aa1,aa2,c1,nupc,xfac
95 my_real ,DIMENSION(NEL) :: plat,plac,plas,epspt,epspc,epsps,sigt,
96 . sigs,dsigs_dp,sigc,dsigc_dp,sdxx,sdyy,sdzz,sdxy,sdyz,sdzx,dsigt_dp,
97 . nup,p,svm,a0,a1,a2,phi,dplac,dplat,dplas,df
98 my_real
99 . normxx,normyy,normzz,normxy,normyz,normzx,norm,aa,cb,
100 . normxx_n,normyy_n,normzz_n,normxy_n,normyz_n,normzx_n,
101 . dsigt_dlam,dsigc_dlam,dsigs_dlam,da1_dsigs,da1_dsigt,
102 . da1_dsigc,da2_dsigs,da2_dsigt,da2_dsigc,da0_dlam,da1_dlam,
103 . da2_dlam,dfdsig2,dphi_dlam,dpxx,dpyy,dpzz,dpxy,dpyz,dpzx,
104 . dpdt_c,dpdt_s,dpdt_t,dpdt,dlam,cc,epdt_min,
105 . epdt_max,epdc_max,epdc_min,epds_min,epds_max,asrate
106 my_real ,DIMENSION(NEL,2) :: xvec
107 my_real ,DIMENSION(NEL,1) :: xvec1
108 my_real ,DIMENSION(NUMTABL) :: tfac
109 my_real ,DIMENSION(NFUNC) :: yfac
110 LOGICAL :: CONV(NEL)
111 my_real, PARAMETER :: sfac = 1.05d0 ! Security factor of ICONV
112C
113 !=======================================================================
114 ! - INITIALISATION OF COMPUTATION ON TIME STEP
115 !=======================================================================
116 ! Recovering model parameters
117 ! -> Elastic parameters
118 g = uparam(4) ! Shear modulus
119 aa1 = uparam(6) ! First component of elastic matrix
120 aa2 = uparam(7) ! Second component of elastic matrix
121 c1 = uparam(8) ! Bulk modulus
122 ! -> Plastic parameters
123 nupc = uparam(9) ! Plastic Poisson ratio
124 ! -> Flags
125 iconv = nint(uparam(15)) ! Flag to ensure convexity
126 asrate = min(one,uparam(16)*timestep)
127 icas = nint(uparam(17)) ! Flag for tabulated case
128 !icas ifunt | ifunc | ifuncs
129 ! -1 1 | 1 | 1
130 ! 0 1 | 0 | 0
131 ! 1 1 | 1 | 0
132 ! 2 1 | 0 | 1
133 xfac = uparam(18) ! Strain-rate scale factor
134 epdt_min = uparam(19)
135 epdt_max = uparam(20)
136 epdc_min = uparam(21)
137 epdc_max = uparam(22)
138 epds_min = uparam(23)
139 epds_max = uparam(24)
140 tfac(1) = uparam(25)
141 tfac(2) = uparam(26)
142 tfac(3) = uparam(27)
143 yfac(1) = uparam(28)
144 yfac(2) = uparam(29)
145c
146 ! Recovering internal variables
147 DO i=1,nel
148 IF(off(i) < em01) off(i) = zero
149 IF(off(i) < one) off(i) = off(i)*four_over_5
150 plat(i) = uvar(i,1) ! Uniaxial tension plastic strain
151 plac(i) = uvar(i,2) ! Uniaxial compression plastic strain
152 plas(i) = uvar(i,3) ! Shear plastic strain
153 epspt(i) = uvar(i,4) ! Uniaxial tension plastic strain-rate
154 epspc(i) = uvar(i,5) ! Uniaxial compression plastic strain-rate
155 epsps(i) = uvar(i,6) ! Shear plastic strain-rate
156 epspt(i) = min(epdt_max, max(epspt(i),epdt_min))
157 epspc(i) = min(epdc_max, max(epspc(i),epdc_min))
158 epsps(i) = min(epds_max, max(epsps(i),epds_min))
159 dpla(i) = zero
160 dplac(i) = zero
161 dplat(i) = zero
162 dplas(i) = zero
163 ENDDO
164c Initialize plastic Poisson ratio
165 IF (time == zero) THEN
166 nup(1:nel) = nupc
167 uvar(1:nel,7) = nupc
168 ELSE
169 nup(1:nel) = uvar(1:nel,7)
170 END IF
171C
172 ! Computation of yield stresses
173 xvec(1:nel,1) = plat(1:nel)
174 xvec(1:nel,2) = epspt(1:nel)*xfac
175 CALL table_mat_vinterp(table(func_trac),nel,nel,vartmp(1,1),xvec,sigt,dsigt_dp)
176 sigt = sigt*tfac(1)
177 dsigt_dp = dsigt_dp*tfac(1)
178 IF (table(func_comp)%NOTABLE > 0) THEN
179 xvec(1:nel,1) = plac(1:nel)
180 xvec(1:nel,2) = epspc(1:nel)*xfac
181 CALL table_mat_vinterp(table(func_comp),nel,nel,vartmp(1,3),xvec,sigc,dsigc_dp)
182 sigc = sigc*tfac(2)
183 dsigc_dp = dsigc_dp*tfac(2)
184 ENDIF
185 IF (table(func_shear)%NOTABLE > 0) THEN
186 xvec(1:nel,1) = plas(1:nel)
187 xvec(1:nel,2) = epsps(1:nel)*xfac
188 CALL table_mat_vinterp(table(func_shear),nel,nel,vartmp(1,5),xvec,sigs,dsigs_dp)
189 sigs = sigs*tfac(3)
190 dsigs_dp = dsigs_dp*tfac(3)
191 ENDIF
192 ! Select case for tabulated yield stresses
193 IF (icas == 0) THEN
194 sigc(1:nel) = sigt(1:nel)
195 sigs(1:nel) = sigt(1:nel)/sqr3
196 ELSEIF (icas == 1) THEN
197 DO i = 1,nel
198 sigs(i) = sqrt(sigc(i)*sigt(i)/three)
199 ENDDO
200 ENDIF
201C
202 !========================================================================
203 ! - COMPUTATION OF TRIAL VALUES
204 !========================================================================
205 DO i=1,nel
206C
207 ! Computation of the trial stress tensor
208 signxx(i) = sigoxx(i) + (aa1*depsxx(i)
209 . + aa2*(depsyy(i) + depszz(i)))
210 signyy(i) = sigoyy(i) + (aa1*depsyy(i)
211 . + aa2*(depsxx(i) + depszz(i)))
212 signzz(i) = sigozz(i) + (aa1*depszz(i)
213 . + aa2*(depsxx(i) + depsyy(i)))
214 signxy(i) = sigoxy(i) + g*depsxy(i)
215 signyz(i) = sigoyz(i) + g*depsyz(i)
216 signzx(i) = sigozx(i) + g*depszx(i)
217C
218 ! Computation of the pressure of the trial stress tensor
219 p(i) = -third*(signxx(i) + signyy(i) + signzz(i))
220 ! Computation of the Von Mises equivalent stress of the trial stress tensor
221 sdxx(i) = signxx(i) + p(i)
222 sdyy(i) = signyy(i) + p(i)
223 sdzz(i) = signzz(i) + p(i)
224 sdxy(i) = signxy(i)
225 sdyz(i) = signyz(i)
226 sdzx(i) = signzx(i)
227 svm(i) = half*(sdxx(i)**2 + sdyy(i)**2 + sdzz(i)**2)
228 . + (sdxy(i)**2 + sdzx(i)**2 + sdyz(i)**2)
229 svm(i) = sqrt(three*svm(i))
230 ENDDO
231c
232 !========================================================================
233 ! - COMPUTATION OF YIELD FONCTION
234 !========================================================================
235 ! Compute yield criterion parameters
236 IF (iconv == 1) THEN
237 ! Ensured convexity
238 DO i = 1,nel
239 conv(i) = .false.
240 IF (sigs(i) < sfac*sqrt(sigc(i)*sigt(i)/three)) THEN
241 sigs(i) = sfac*sqrt(sigc(i)*sigt(i)/three)
242 conv(i) = .true.
243 ENDIF
244 ENDDO
245 ENDIF
246 ! Compute yield criterion parameters A0,A1,A2
247 DO i = 1,nel
248 aa = one/sigc(i)/sigt(i)
249 a0(i) = three*(sigs(i)**2)
250 a1(i) = nine*(sigs(i)**2)*(sigc(i) - sigt(i))*aa
251 a2(i) = nine*(sigc(i)*sigt(i) - three*(sigs(i)**2))*aa
252 ENDDO
253c
254 ! -> Checking plastic behavior for all elements
255 nindx = 0
256 DO i=1,nel
257 phi(i) = (svm(i)**2) - a0(i) - a1(i)*p(i) - a2(i)*p(i)*p(i)
258 IF (phi(i) >= zero .AND. off(i) == one) THEN
259 nindx = nindx + 1
260 indx(nindx) = i
261 ENDIF
262 ENDDO
263c
264 !====================================================================
265 ! - PLASTIC RETURN MAPPING WITH CUTTING PLANE METHOD
266 !====================================================================
267 IF (nindx > 0) THEN
268c
269 ! Loop over the iterations
270 DO iter = 1,3
271c
272 ! Loop over yielding elements
273 DO ii = 1,nindx
274 i = indx(ii)
275c
276 ! Note : in this part, the purpose is to compute for each iteration
277 ! a plastic multiplier allowing to update internal variables to satisfy
278 ! the consistency condition using the cutting plane algorithm
279 ! Its expression at each iteration is : DLAMBDA = - PHI/DPHI_DLAMBDA
280 ! -> PHI : current value of yield function (known)
281 ! -> DPHI_DLAMBDA : derivative of PHI with respect to DLAMBDA by taking
282 ! into account of internal variables kinetic :
283 ! plasticity ...
284c
285 ! 1 - Computation of DPHI_DSIG the normal to the yield surface
286 !-------------------------------------------------------------
287C
288 ! Computation of normal to yield surface
289 normxx = three*sdxx(i) + third*(a1(i) + two*a2(i)*p(i))
290 normyy = three*sdyy(i) + third*(a1(i) + two*a2(i)*p(i))
291 normzz = three*sdzz(i) + third*(a1(i) + two*a2(i)*p(i))
292 normxy = three*sdxy(i)
293 normyz = three*sdyz(i)
294 normzx = three*sdzx(i)
295C
296 ! Computation of the normalized normal to the yield surface
297 cb = a1(i) + two*a2(i)*p(i)
298 norm = max(em20,sqrt(six*(svm(i)**2) + third*cb*cb))
299 normxx_n = normxx/norm
300 normyy_n = normyy/norm
301 normzz_n = normzz/norm
302 normxy_n = normxy/norm
303 normyz_n = normyz/norm
304 normzx_n = normzx/norm
305C
306 ! 2 - Computation of DPHI_DLAMBDA
307 !---------------------------------------------------------
308c
309 ! a) Derivative with respect stress increments tensor DSIG
310 ! --------------------------------------------------------
311 dfdsig2 = normxx*(aa1*normxx_n + aa2*(normyy_n + normzz_n)) +
312 . normyy*(aa1*normyy_n + aa2*(normxx_n + normzz_n)) +
313 . normzz*(aa1*normzz_n + aa2*(normxx_n + normyy_n)) +
314 . two*normxy*two*g*normxy_n +
315 . two*normyz*two*g*normyz_n +
316 . two*normzx*two*g*normzx_n
317c
318 ! b) Derivative of yield criterion parameters
319 ! --------------------------------------------
320C
321 ! Derivative of yield surfaces with respect to yield criterion parameter A0,A1,A2
322 dsigt_dlam = dsigt_dp(i)*(svm(i)/norm)*(three/(one + nup(i)))
323 IF (table(func_comp)%NOTABLE > 0) dsigc_dlam = dsigc_dp(i)*(svm(i)/norm)*(three/(one + nup(i)))
324 IF (table(func_shear)%NOTABLE > 0) dsigs_dlam = dsigs_dp(i)*sqr3*(svm(i)/norm)
325 IF (icas == 0) THEN
326 dsigc_dlam = dsigt_dlam
327 dsigs_dlam = (one/sqr3)*dsigt_dlam
328 ELSEIF (icas == 1) THEN
329 dsigs_dlam = (one/sqr3)*(one/(two*sqrt(sigt(i)*sigc(i))))*
330 . (dsigc_dlam*sigt(i) + sigc(i)*dsigt_dlam)
331 ENDIF
332 IF (iconv == 1) THEN
333 IF (conv(i)) THEN
334 dsigs_dlam = sfac*(one/sqr3)*(one/(two*sqrt(sigt(i)*sigc(i))))*
335 . (dsigc_dlam*sigt(i) + sigc(i)*dsigt_dlam)
336 ENDIF
337 ENDIF
338C
339 ! -> A1 derivatives
340 cc = sigs(i)/sigc(i)/sigt(i)
341 ! -> With respect to SIGS
342 da1_dsigs = eighteen*(sigc(i) - sigt(i))*cc
343 ! -> With respect to SIGC
344 da1_dsigc = nine*(sigs(i)/sigc(i))**2
345 ! -> With respect to SIGT
346 da1_dsigt = -nine*(sigs(i)/sigt(i))**2
347C
348 ! -> A2 derivatives
349 ! -> With respect to SIGS
350 da2_dsigs = -cinquante4*cc
351 ! -> With respect to SIGC
352 da2_dsigc = twenty7*cc*sigs(i)/sigc(i)
353 ! -> With respect to SIGT
354 da2_dsigt = twenty7*cc*sigs(i)/sigt(i)
355c
356 ! -> A parameters derivatives with respect to lambda
357 ! -> A0 with respect to lambda
358 da0_dlam = six*sigs(i)*dsigs_dlam
359 ! -> A1 with respect to lambda
360 da1_dlam = da1_dsigs*dsigs_dlam + da1_dsigt*dsigt_dlam + da1_dsigc*dsigc_dlam
361 ! -> A2 with respect to lambda
362 da2_dlam = da2_dsigs*dsigs_dlam + da2_dsigt*dsigt_dlam + da2_dsigc*dsigc_dlam
363c
364c
365 ! 3 - Computation of plastic multiplier and variables update
366 !----------------------------------------------------------
367c
368 ! Derivative of PHI with respect to DLAM
369 dphi_dlam = - dfdsig2 - da0_dlam - p(i)*da1_dlam - (p(i)**2)*da2_dlam
370 dphi_dlam = sign(max(abs(dphi_dlam),em20),dphi_dlam)
371 dlam = -phi(i)/dphi_dlam
372c
373 ! Plastic strains tensor update
374 dpxx = dlam*normxx_n
375 dpyy = dlam*normyy_n
376 dpzz = dlam*normzz_n
377 dpxy = dlam*normxy_n
378 dpyz = dlam*normyz_n
379 dpzx = dlam*normzx_n
380c
381 ! Elasto-plastic stresses update
382 signxx(i) = signxx(i) - (aa1*dpxx + aa2*(dpyy + dpzz))
383 signyy(i) = signyy(i) - (aa1*dpyy + aa2*(dpxx + dpzz))
384 signzz(i) = signzz(i) - (aa1*dpzz + aa2*(dpxx + dpyy))
385 signxy(i) = signxy(i) - two*g*dpxy
386 signyz(i) = signyz(i) - two*g*dpyz
387 signzx(i) = signzx(i) - two*g*dpzx
388c
389 ! Cumulated plastic strains update
390 plat(i) = plat(i) + dlam*(svm(i)/norm)*three/(one + nup(i))
391 plac(i) = plat(i)
392 plas(i) = plas(i) + sqr3*(svm(i)/norm)*dlam
393 pla(i) = pla(i) + two*dlam*(svm(i)/norm)
394 ! Plastic strain increments update
395 dplat(i) = dplat(i) + dlam*(svm(i)/norm)*three/(one + nup(i))
396 dplac(i) = dplat(i)
397 dplas(i) = dplas(i) + sqr3*(svm(i)/norm)*dlam
398 dpla(i) = dpla(i) + two*dlam*(svm(i)/norm)
399C
400 ! Update pressure
401 p(i) = -third*(signxx(i) + signyy(i) + signzz(i))
402 ! Update Von Mises stress
403 sdxx(i) = signxx(i) + p(i)
404 sdyy(i) = signyy(i) + p(i)
405 sdzz(i) = signzz(i) + p(i)
406 sdxy(i) = signxy(i)
407 sdyz(i) = signyz(i)
408 sdzx(i) = signzx(i)
409 svm(i) = half*(sdxx(i)**2 + sdyy(i)**2 + sdzz(i)**2)
410 . + (sdxy(i)**2 + sdzx(i)**2 + sdyz(i)**2)
411 svm(i) = sqrt(three*svm(i))
412 ENDDO
413C
414 ! Update yield stresses
415 xvec(1:nel,1) = plat(1:nel)
416 xvec(1:nel,2) = epspt(1:nel)*xfac
417 CALL table_mat_vinterp(table(func_trac),nel,nel,vartmp(1,1),xvec,sigt,dsigt_dp)
418 sigt = sigt*tfac(1)
419 dsigt_dp = dsigt_dp*tfac(1)
420 IF (table(func_comp)%NOTABLE > 0) THEN
421 xvec(1:nel,1) = plac(1:nel)
422 xvec(1:nel,2) = epspc(1:nel)*xfac
423 CALL table_mat_vinterp(table(func_comp),nel,nel,vartmp(1,3),xvec,sigc,dsigc_dp)
424 sigc(1:nel) = sigc*tfac(2)
425 dsigc_dp(1:nel) = dsigc_dp*tfac(2)
426 ENDIF
427 IF (table(func_shear)%NOTABLE > 0) THEN
428 xvec(1:nel,1) = plas(1:nel)
429 xvec(1:nel,2) = epsps(1:nel)*xfac
430 CALL table_mat_vinterp(table(func_shear),nel,nel,vartmp(1,5),xvec,sigs,dsigs_dp)
431 sigs(1:nel) = sigs*tfac(3)
432 dsigs_dp(1:nel) = dsigs_dp*tfac(3)
433 ENDIF
434 ! Select case for tabulated yield stresses
435 IF (icas == 0) THEN
436 DO ii = 1,nindx
437 i = indx(ii)
438 sigc(i) = sigt(i)
439 sigs(i) = sigt(i)/sqr3
440 ENDDO
441 ELSEIF (icas == 1) THEN
442 DO ii = 1,nindx
443 i = indx(ii)
444 sigs(i) = sqrt(sigc(i)*sigt(i)/three)
445 ENDDO
446 ENDIF
447C
448 ! Update yield function value
449 IF (iconv == 1) THEN
450 DO ii = 1,nindx
451 i = indx(ii)
452 conv(i) = .false.
453 IF (sigs(i) < sfac*sqrt(sigc(i)*sigt(i)/three)) THEN
454 sigs(i) = sfac*sqrt(sigc(i)*sigt(i)/three)
455 conv(i) = .true.
456 ENDIF
457 ENDDO
458 ENDIF
459 DO ii = 1,nindx
460 i = indx(ii)
461 aa = one/sigc(i)/sigt(i)
462 a0(i) = three*(sigs(i)**2)
463 a1(i) = nine*(sigs(i)**2)*(sigc(i) - sigt(i))*aa
464 a2(i) = nine*(sigc(i)*sigt(i) - three*(sigs(i)**2))*aa
465 phi(i) = svm(i)**2 - a0(i) - a1(i)*p(i) - a2(i)*p(i)*p(i)
466 ENDDO
467 ENDDO
468 ENDIF
469c-----------------------------------------------
470c Update plastic Poisson coefficient
471c-----------------------------------------------
472 IF (ifunc(1) > 0) THEN
473 iad(1:nel) = npf(ifunc(1)) / 2 + 1
474 ilen(1:nel) = npf(ifunc(1)+1) / 2 - iad(1:nel) - vartmp(1:nel,8)
475!
476 CALL vinter(tf,iad,vartmp(1:nel,8),ilen,nel,pla,df,nup)
477!
478 uvar(1:nel,7) = yfac(1) * nup(1:nel)
479 uvar(1:nel,7) = max(zero, min(nup(1:nel), half))
480 END IF
481c
482 !====================================================================
483 ! - STORING NEW VALUES
484 !====================================================================
485 DO i=1,nel
486 uvar(i,1) = plat(i)
487 uvar(i,2) = plac(i)
488 uvar(i,3) = plas(i)
489 dpdt_t = dplat(i)/max(timestep,em20)
490 uvar(i,4) = asrate*dpdt_t + (one-asrate)*epspt(i)
491 dpdt_c = dplac(i)/max(timestep,em20)
492 uvar(i,5) = asrate*dpdt_c + (one-asrate)*epspc(i)
493 dpdt_s = dplas(i)/max(timestep,em20)
494 uvar(i,6) = asrate*dpdt_s + (one-asrate)*epsps(i)
495 dpdt = dpla(i)/max(timestep,em20)
496 epsd(i) = asrate*dpdt + (one-asrate)*epsd(i)
497 ! Computation of soundspeed
498 soundsp(i) = sqrt((c1+four*g/three)/rho0(i))
499 ! Computation of the yield stress
500 yld(i) = a0(i) + a1(i)*p(i) + a2(i)*p(i)*p(i)
501 ! Computation of the hourglass coefficient
502 et(i) = half
503 ENDDO
504c
my_real
#define my_real
Definition cppsort.cpp:32
norm
norm(diag(diag(diag(inv(mat))) -id.SOL), 2) % destroy mumps instance id.JOB
func_comp
subroutine func_comp(table, ntable, nptmax, npt_trac, npt_shear, npt_comp, x_comp, y_comp, nup)
Definition law76_upd.F:344
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
table4d_mod
Definition table4d_mod.F:52
table_mat_vinterp_mod
Definition table_mat_vinterp.F:32
table_mat_vinterp_mod::table_mat_vinterp
subroutine table_mat_vinterp(table, dimx, nel, ipos, xx, yy, dydx)
Definition table_mat_vinterp.F:45
vinter
subroutine vinter(tf, iad, ipos, ilen, nel, x, dydx, y)
Definition vinter.F:72