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

Go to the source code of this file.

Functions/Subroutines

subroutine mat104_nodam_nice (nel, ngl, nuparam, nuvar, volume, fheat, time, timestep, uparam, uvar, jthe, off, rho0, rho, pla, dpla, epsd, soundsp, depsxx, depsyy, depszz, depsxy, depsyz, depszx, sigoxx, sigoyy, sigozz, sigoxy, sigoyz, sigozx, signxx, signyy, signzz, signxy, signyz, signzx, sigy, et, temp, seq, inloc, dplanl, jlag)

Function/Subroutine Documentation

◆ mat104_nodam_nice()

subroutine mat104_nodam_nice ( integer nel,
integer, dimension(nel), intent(in) ngl,
integer nuparam,
integer nuvar,
intent(in) volume,
intent(inout) fheat,
time,
timestep,
intent(in) uparam,
intent(inout) uvar,
integer jthe,
intent(inout) off,
intent(in) rho0,
intent(in) rho,
intent(inout) pla,
intent(inout) dpla,
intent(inout) epsd,
intent(out) soundsp,
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) sigy,
intent(out) et,
intent(inout) temp,
intent(inout) seq,
integer inloc,
intent(in) dplanl,
integer, intent(in) jlag )

Definition at line 28 of file mat104_nodam_nice.F.

37 !=======================================================================
38 ! Implicit types
39 !=======================================================================
40#include "implicit_f.inc"
41 !=======================================================================
42 ! Common
43 !=======================================================================
44 !=======================================================================
45 ! dummy arguments
46 !=======================================================================
47 INTEGER NEL,NUPARAM,NUVAR,JTHE,INLOC
48 INTEGER ,INTENT(IN) :: JLAG
49 INTEGER ,DIMENSION(NEL), INTENT(IN) :: NGL
50 my_real
51 . time,timestep
52 my_real,DIMENSION(NUPARAM), INTENT(IN) ::
53 . uparam
54 my_real, DIMENSION(NEL), INTENT(IN) :: volume
55 my_real ,DIMENSION(NEL), INTENT(INOUT) :: fheat
56 my_real,DIMENSION(NEL), INTENT(IN) ::
57 . rho0,rho,dplanl,
58 . depsxx,depsyy,depszz,depsxy,depsyz,depszx,
59 . sigoxx,sigoyy,sigozz,sigoxy,sigoyz,sigozx
60c
61 my_real ,DIMENSION(NEL), INTENT(OUT) ::
62 . soundsp,sigy,et,
63 . signxx,signyy,signzz,signxy,signyz,signzx
64c
65 my_real ,DIMENSION(NEL), INTENT(INOUT) ::
66 . pla,dpla,epsd,off,temp,seq
67 my_real ,DIMENSION(NEL,NUVAR), INTENT(INOUT) ::
68 . uvar
69c
70 !=======================================================================
71 ! Local Variables
72 !=======================================================================
73 INTEGER I,K,II,NSP,NINDX,INDEX(NEL),NICE
74c
75 my_real
76 . young,bulk,lam,g,g2,nu,cdr,kdr,hard,yld0,qvoce,bvoce,jcc,
77 . epsp0,mtemp,tini,tref,eta,cp,dpis,dpad,asrate,afiltr
78 my_real
79 . h,ldav,trdep,trdfds,sigvm,sigdr2,yld2i,omega,
80 . dtemp,fcosh,fsinh,dpdt,dlam,ddep
81 my_real
82 . dpxx,dpyy,dpzz,dpxy,dpyz,dpzx,dsdrdj2,dsdrdj3,
83 . dj3dsxx,dj3dsyy,dj3dszz,dj3dsxy,dj3dsyz,dj3dszx,
84 . dj2dsxx,dj2dsyy,dj2dszz,dj2dsxy,dj2dsyz,dj2dszx,
85 . dfdsxx,dfdsyy,dfdszz,dfdsxy,dfdsyz,dfdszx,
86 . normxx,normyy,normzz,normxy,normyz,normzx,
87 . denom,dphi,dphi_dtrsig,sig_dfdsig,dfdsig2,
88 . dphi_dsig,dphi_dyld,dphi_dfdr,dpdt_nl,
89 . dyld_dpla,dyld_dtemp,dtemp_dlam,normsig
90c
91 my_real, DIMENSION(NEL) ::
92 . dsigxx,dsigyy,dsigzz,dsigxy,dsigyz,dsigzx,trsig,
93 . sxx,syy,szz,sxy,syz,szx,sigm,j2,j3,sigdr,yld,weitemp,
94 . hardp,fhard,frate,ftherm,fdr,phi0,phi,dpla_dlam
95 !=======================================================================
96 ! DRUCKER - VOCE - JOHNSON-COOK MATERIAL LAW
97 !=======================================================================
98 !UVAR(1) PHI YIELD FUNCTION VALUE
99 !=======================================================================
100c
101 !=======================================================================
102 ! - INITIALISATION OF COMPUTATION ON TIME STEP
103 !=======================================================================
104 ! Recovering model parameters
105 ! Elastic parameters
106 young = uparam(1) ! Young modulus
107 bulk = uparam(2) ! Bulk modulus
108 g = uparam(3) ! Shear modulus
109 g2 = uparam(4) ! 2*Shear modulus
110 lam = uparam(5) ! Lambda Hooke parameter
111 nu = uparam(6) ! Poisson ration
112 ! Plastic criterion and hardening parameters [Drucker, 1948]
113 nice = nint(uparam(11)) ! Choice of the Nice method
114 cdr = uparam(12) ! Drucker coefficient
115 kdr = uparam(13) ! Drucker 1/K coefficient
116 tini = uparam(14) ! Initial temperature
117 hard = uparam(15) ! Linear hardening
118 yld0 = uparam(16) ! Initial yield stress
119 qvoce = uparam(17) ! 1st Voce parameter
120 bvoce = uparam(18) ! 2nd Voce parameter
121 ! Strain-rate dependence parameters
122 jcc = uparam(20) ! Johnson-Cook strain rate coefficient
123 epsp0 = uparam(21) ! Johnson-Cook reference strain rate
124 ! Thermal softening and self-heating parameters
125 mtemp = uparam(22) ! Thermal softening slope
126 tref = uparam(23) ! Reference temperature
127 eta = uparam(24) ! Taylor-Quinney coefficient
128 cp = uparam(25) ! Thermal mass capacity
129 dpis = uparam(26) ! Isothermal plastic strain rate
130 dpad = uparam(27) ! Adiabatic plastic strain rate
131 ! Plastic strain-rate filtering parameters
132 asrate = uparam(28) ! Plastic strain rate filtering frequency
133 afiltr = asrate*timestep/(asrate*timestep + one)
134c
135 ! Recovering internal variables
136 DO i=1,nel
137 IF (off(i) < em01) off(i) = zero
138 IF (off(i) < one) off(i) = off(i)*four_over_5
139 ! User inputs
140 phi0(i) = uvar(i,1) ! Previous yield function value
141 ! Standard inputs
142 dpla(i) = zero ! Initialization of the plastic strain increment
143 et(i) = one ! Initialization of hourglass coefficient
144 hardp(i) = zero ! Initialization of hardening modulus
145 ENDDO
146c
147 ! Initialization of temperature and self-heating weight factor
148 IF (time == zero .and. jthe == 0) THEN
149 temp(1:nel) = tini
150 ENDIF
151 IF (cp > zero .OR. jthe /= 0) THEN
152 IF (inloc == 0 .OR. jthe /= 0) THEN
153 DO i=1,nel
154 IF (epsd(i) < dpis) THEN
155 weitemp(i) = zero
156 ELSEIF (epsd(i) > dpad) THEN
157 weitemp(i) = one
158 ELSE
159 weitemp(i) = ((epsd(i)-dpis)**2 )
160 . * (three*dpad - two*epsd(i) - dpis)
161 . / ((dpad-dpis)**3)
162 ENDIF
163 ENDDO
164 ENDIF
165 ENDIF
166c
167 ! Computation of the initial yield stress
168 DO i = 1,nel
169 ! a) - Hardening law
170 fhard(i) = yld0 + hard*pla(i) + qvoce*(one-exp(-bvoce*pla(i)))
171 ! b) - Correction factor for strain-rate dependence (Johnson-Cook)
172 frate(i) = one
173 IF (epsd(i) > epsp0) frate(i) = frate(i) + jcc*log(epsd(i)/epsp0)
174 ! c) - Correction factor for thermal softening
175 ftherm(i) = one - mtemp * (temp(i) - tref)
176 ! d) - Computation of the yield function and value check
177 yld(i) = fhard(i)*frate(i)*ftherm(i)
178 ! e) - Checking values
179 yld(i) = max(em10, yld(i))
180 ENDDO
181c
182 !========================================================================
183 ! - COMPUTATION OF TRIAL VALUES
184 !========================================================================
185 DO i=1,nel
186 ! Computation of the trial stress tensor
187 ldav = (depsxx(i) + depsyy(i) + depszz(i)) * lam
188 signxx(i) = sigoxx(i) + depsxx(i)*g2 + ldav
189 signyy(i) = sigoyy(i) + depsyy(i)*g2 + ldav
190 signzz(i) = sigozz(i) + depszz(i)*g2 + ldav
191 signxy(i) = sigoxy(i) + depsxy(i)*g
192 signyz(i) = sigoyz(i) + depsyz(i)*g
193 signzx(i) = sigozx(i) + depszx(i)*g
194 ! Computation of the trace of the trial stress tensor
195 trsig(i) = signxx(i) + signyy(i) + signzz(i)
196 sigm(i) = -trsig(i) * third
197 ! Computation of the deviatoric trial stress tensor
198 sxx(i) = signxx(i) + sigm(i)
199 syy(i) = signyy(i) + sigm(i)
200 szz(i) = signzz(i) + sigm(i)
201 sxy(i) = signxy(i)
202 syz(i) = signyz(i)
203 szx(i) = signzx(i)
204 ! Second deviatoric invariant
205 j2(i) = half*(sxx(i)**2 + syy(i)**2 + szz(i)**2 )
206 . + sxy(i)**2 + syz(i)**2 + szx(i)**2
207 ! Third deviatoric invariant
208 j3(i) = sxx(i) * syy(i) * szz(i)
209 . + sxy(i) * syz(i) * szx(i) * two
210 . - sxx(i) * syz(i)**2
211 . - syy(i) * szx(i)**2
212 . - szz(i) * sxy(i)**2
213 ! Drucker equivalent stress
214 fdr(i) = j2(i)**3 - cdr*(j3(i)**2)
215 ! Checking equivalent stress values
216 IF (fdr(i) > zero) THEN
217 sigdr(i) = kdr * exp((one/six)*log(fdr(i))) ! FDR(I)**(1/6)
218 ELSE
219 sigdr(i) = zero
220 ENDIF
221 ENDDO
222c
223 !========================================================================
224 ! - COMPUTATION OF YIELD FONCTION
225 !========================================================================
226 phi(1:nel) = (sigdr(1:nel) / yld(1:nel))**2 - one
227c
228 ! Checking plastic behavior for all elements
229 nindx = 0
230 DO i=1,nel
231 IF ((phi(i) >= zero).AND.(off(i) == one)) THEN
232 nindx=nindx+1
233 index(nindx)=i
234 ENDIF
235 ENDDO
236c
237 !====================================================================
238 ! - PLASTIC CORRECTION WITH NICE - EXPLICIT METHOD
239 !====================================================================
240#include "vectorize.inc"
241 DO ii = 1, nindx
242c
243 ! Number of the element with plastic behaviour
244 i = index(ii)
245c
246 ! Computation of the trial stress increment
247 ldav = (depsxx(i) + depsyy(i) + depszz(i)) * lam
248 dsigxx(i) = depsxx(i)*g2 + ldav
249 dsigyy(i) = depsyy(i)*g2 + ldav
250 dsigzz(i) = depszz(i)*g2 + ldav
251 dsigxy(i) = depsxy(i)*g
252 dsigyz(i) = depsyz(i)*g
253 dsigzx(i) = depszx(i)*g
254 dlam = zero
255c
256 ! Computation of Drucker equivalent stress of the previous stress tensor
257 trsig(i) = sigoxx(i) + sigoyy(i) + sigozz(i)
258 sigm(i) = -trsig(i) * third
259 sxx(i) = sigoxx(i) + sigm(i)
260 syy(i) = sigoyy(i) + sigm(i)
261 szz(i) = sigozz(i) + sigm(i)
262 sxy(i) = sigoxy(i)
263 syz(i) = sigoyz(i)
264 szx(i) = sigozx(i)
265 j2(i) = half*(sxx(i)**2 + syy(i)**2 + szz(i)**2 )
266 . + sxy(i)**2 + syz(i)**2 + szx(i)**2
267 j3(i) = sxx(i)*syy(i)*szz(i) + sxy(i)*syz(i)*szx(i) * two
268 . - sxx(i)*syz(i)**2 - syy(i)*szx(i)**2 - szz(i)*sxy(i)**2
269 fdr(i) = j2(i)**3 - cdr*(j3(i)**2)
270 sigdr(i) = kdr * exp((one/six)*log(fdr(i)))
271c
272 ! Note : in this part, the purpose is to compute for each iteration
273 ! a plastic multiplier allowing to update internal variables to satisfy
274 ! the consistency condition.
275 ! Its expression is : DLAMBDA = - (PHI_OLD + DPHI) / DPHI_DLAMBDA
276 ! -> PHI_OLD : old value of yield function (known)
277 ! -> DPHI : yield function prediction (to compute)
278 ! -> DPHI_DLAMBDA : derivative of PHI with respect to DLAMBDA by taking
279 ! into account of internal variables kinetic :
280 ! plasticity, temperature and damage (to compute)
281c
282 ! 1 - Computation of DPHI_DSIG the normal to the yield surface
283 !-------------------------------------------------------------
284c
285 ! Derivative with respect to the equivalent stress
286 yld2i = one/(yld(i)**2)
287 dphi_dsig = two*sigdr(i)*yld2i
288c
289 ! Computation of the Eulerian norm of the stress tensor
290 normsig = sqrt(sigoxx(i)*sigoxx(i)
291 . + sigoyy(i)*sigoyy(i)
292 . + sigozz(i)*sigozz(i)
293 . + two*sigoxy(i)*sigoxy(i)
294 . + two*sigoyz(i)*sigoyz(i)
295 . + two*sigozx(i)*sigozx(i))
296 normsig = max(normsig,one)
297c
298 ! Derivative with respect to Fdr
299 fdr(i) = (j2(i)/(normsig**2))**3 - cdr*((j3(i)/(normsig**3))**2)
300 dphi_dfdr = dphi_dsig*kdr*(one/six)*exp(-(five/six)*log(fdr(i)))
301 dsdrdj2 = dphi_dfdr*three*(j2(i)/(normsig**2))**2
302 dsdrdj3 = -dphi_dfdr*two*cdr*(j3(i)/(normsig**3))
303 ! dJ3/dSig
304 dj3dsxx = two_third*(syy(i)*szz(i)-syz(i)**2)/(normsig**2)
305 . - third*(sxx(i)*szz(i)-szx(i)**2)/(normsig**2)
306 . - third*(sxx(i)*syy(i)-sxy(i)**2)/(normsig**2)
307 dj3dsyy = - third*(syy(i)*szz(i)-syz(i)**2)/(normsig**2)
308 . + two_third*(sxx(i)*szz(i)-szx(i)**2)/(normsig**2)
309 . - third*(sxx(i)*syy(i)-sxy(i)**2)/(normsig**2)
310 dj3dszz = - third*(syy(i)*szz(i)-syz(i)**2)/(normsig**2)
311 . - third*(sxx(i)*szz(i)-szx(i)**2)/(normsig**2)
312 . + two_third*(sxx(i)*syy(i)-sxy(i)**2)/(normsig**2)
313 dj3dsxy = two*(sxx(i)*sxy(i) + sxy(i)*syy(i) + szx(i)*syz(i))/(normsig**2)
314 dj3dsyz = two*(sxy(i)*szx(i) + syy(i)*syz(i) + syz(i)*szz(i))/(normsig**2)
315 dj3dszx = two*(sxx(i)*szx(i) + sxy(i)*syz(i) + szx(i)*szz(i))/(normsig**2)
316 ! dPhi/dSig
317 normxx = dsdrdj2*sxx(i)/normsig + dsdrdj3*dj3dsxx
318 normyy = dsdrdj2*syy(i)/normsig + dsdrdj3*dj3dsyy
319 normzz = dsdrdj2*szz(i)/normsig + dsdrdj3*dj3dszz
320 normxy = two*dsdrdj2*sxy(i)/normsig + dsdrdj3*dj3dsxy
321 normyz = two*dsdrdj2*syz(i)/normsig + dsdrdj3*dj3dsyz
322 normzx = two*dsdrdj2*szx(i)/normsig + dsdrdj3*dj3dszx
323c
324 ! Restoring previous value of the yield function
325 phi(i) = phi0(i)
326c
327 ! Computation of yield surface trial increment DPHI
328 dphi = normxx * dsigxx(i)
329 . + normyy * dsigyy(i)
330 . + normzz * dsigzz(i)
331 . + normxy * dsigxy(i)
332 . + normyz * dsigyz(i)
333 . + normzx * dsigzx(i)
334c
335 ! 2 - Computation of DPHI_DLAMBDA
336 !---------------------------------------------------------
337c
338 ! a) Derivative with respect stress increments tensor DSIG
339 ! --------------------------------------------------------
340 dfdsig2 = normxx * normxx * g2
341 . + normyy * normyy * g2
342 . + normzz * normzz * g2
343 . + normxy * normxy * g
344 . + normyz * normyz * g
345 . + normzx * normzx * g
346c
347 ! b) Derivatives with respect to plastic strain P
348 ! ------------------------------------------------
349c
350 ! i) Derivative of the yield stress with respect to plastic strain dYLD / dPLA
351 ! ----------------------------------------------------------------------------
352 hardp(i) = hard + qvoce*bvoce*exp(-bvoce*pla(i))
353 dyld_dpla = hardp(i)*frate(i)*ftherm(i)
354c
355 ! ii) Derivative of dPLA with respect to DLAM
356 ! -------------------------------------------
357 sig_dfdsig = signxx(i) * normxx
358 . + signyy(i) * normyy
359 . + signzz(i) * normzz
360 . + signxy(i) * normxy
361 . + signyz(i) * normyz
362 . + signzx(i) * normzx
363 dpla_dlam(i) = sig_dfdsig / yld(i)
364c
365 ! c) Derivatives with respect to the temperature TEMP
366 ! ---------------------------------------------------
367 IF (jthe == 0 .AND. cp > zero .AND. inloc == 0) THEN
368 ! i) Derivative of the yield stress with respect to temperature dYLD / dTEMP
369 ! ---------------------------------------------------------------------------
370 dyld_dtemp = -fhard(i)*frate(i)*mtemp
371 ! ii) Derivative of the temperature TEMP with respect to DLAM
372 ! -----------------------------------------------------------
373 dtemp_dlam = weitemp(i)*(eta/(rho0(i)*cp))*sig_dfdsig
374 ELSE
375 dyld_dtemp = zero
376 dtemp_dlam = zero
377 ENDIF
378c
379 ! d) Derivative with respect to the yield stress
380 ! ----------------------------------------------
381 sigdr2 = sigdr(i)**2
382 dphi_dyld = -two*sigdr2/(yld(i)**3)
383c
384 ! e) Derivative of PHI with respect to DLAM ( = -DENOM)
385 denom = dfdsig2 - (dphi_dyld * dyld_dpla * dpla_dlam(i))
386 IF (jthe == 0 .AND. cp > zero .AND. inloc == 0) THEN
387 denom = denom - (dphi_dyld * dyld_dtemp * dtemp_dlam)
388 ENDIF
389 denom = sign(max(abs(denom),em20) ,denom)
390c
391 ! 3 - Computation of plastic multiplier and variables update
392 !----------------------------------------------------------
393c
394 ! Computation of the plastic multiplier
395 dlam = (dphi + phi(i)) / denom
396 dlam = max(dlam, zero)
397c
398 ! Plastic strains tensor update
399 dpxx = dlam * normxx
400 dpyy = dlam * normyy
401 dpzz = dlam * normzz
402 dpxy = dlam * normxy
403 dpyz = dlam * normyz
404 dpzx = dlam * normzx
405c
406 ! Elasto-plastic stresses update
407 signxx(i) = sigoxx(i) + dsigxx(i) - dpxx*g2
408 signyy(i) = sigoyy(i) + dsigyy(i) - dpyy*g2
409 signzz(i) = sigozz(i) + dsigzz(i) - dpzz*g2
410 signxy(i) = sigoxy(i) + dsigxy(i) - dpxy*g
411 signyz(i) = sigoyz(i) + dsigyz(i) - dpyz*g
412 signzx(i) = sigozx(i) + dsigzx(i) - dpzx*g
413 trsig(i) = signxx(i) + signyy(i) + signzz(i)
414 sigm(i) = -trsig(i) * third
415 sxx(i) = signxx(i) + sigm(i)
416 syy(i) = signyy(i) + sigm(i)
417 szz(i) = signzz(i) + sigm(i)
418 sxy(i) = signxy(i)
419 syz(i) = signyz(i)
420 szx(i) = signzx(i)
421c
422 ! Cumulated plastic strain and strain rate update
423 ddep = (dlam/yld(i))*sig_dfdsig
424 dpla(i) = dpla(i) + ddep
425 pla(i) = pla(i) + ddep
426c
427 ! Drucker equivalent stress update
428 j2(i) = half*(sxx(i)**2 + syy(i)**2 + szz(i)**2 )
429 . + sxy(i)**2 + syz(i)**2 + szx(i)**2
430 j3(i) = sxx(i) * syy(i) * szz(i) + two * sxy(i) * syz(i) * szx(i)
431 . - sxx(i) * syz(i)**2 - syy(i) * szx(i)**2 - szz(i) * sxy(i)**2
432 fdr(i) = j2(i)**3 - cdr*(j3(i)**2)
433 sigdr(i) = kdr * exp((one/six)*log(fdr(i)))
434c
435 ! Temperature update
436 IF (jthe == 0 .AND. cp > zero .AND. inloc == 0) THEN
437 dtemp = weitemp(i)*yld(i)*ddep*eta/(rho0(i)*cp)
438 temp(i) = temp(i) + dtemp
439 ENDIF
440 ftherm(i) = one - mtemp* (temp(i) - tref)
441c
442 ! Hardening law update
443 fhard(i) = yld0 + hard * pla(i) + qvoce*(one-exp(-bvoce*pla(i)))
444c
445 ! Yield stress update
446 yld(i) = fhard(i) * frate(i) * ftherm(i)
447 yld(i) = max(yld(i), em10)
448c
449 ! Yield function value update
450 sigdr2 = sigdr(i)**2
451 yld2i = one / yld(i)**2
452 phi(i) = sigdr2 * yld2i - one
453 ENDDO
454 !===================================================================
455 ! - END OF PLASTIC CORRECTION WITH NICE - EXPLICIT METHOD
456 !===================================================================
457c
458 ! Storing new values
459 DO i=1,nel
460 ! USR Outputs & coefficient for hourglass
461 IF (dpla(i) > zero) THEN
462 uvar(i,1) = phi(i) ! Yield function value
463 et(i) = hardp(i)*frate(i) / (hardp(i)*frate(i) + young)
464 ELSE
465 uvar(i,1) = zero
466 et(i) = one
467 ENDIF
468 seq(i) = sigdr(i) ! SIGEQ
469 dpdt = dpla(i) / max(em20,timestep)
470 ! Plastic strain-rate (filtered)
471 epsd(i) = afiltr * dpdt + (one - afiltr) * epsd(i)
472 ! Computation of the sound speed
473 soundsp(i) = sqrt((bulk + four_over_3*g)/rho0(i))
474 ! storing the yield stress
475 sigy(i) = yld(i)
476 ! Non-local temperature
477 IF ((inloc > 0).AND.(off(i) == one)) THEN
478 IF (cp > zero .AND. jthe == 0) THEN
479 ! Computation of the weight factor
480 dpdt_nl = max(dplanl(i),zero)/max(timestep,em20)
481 IF (dpdt_nl < dpis) THEN
482 weitemp(i) = zero
483 ELSEIF (dpdt_nl > dpad) THEN
484 weitemp(i) = one
485 ELSE
486 weitemp(i) = ((dpdt_nl-dpis)**2 )
487 . * (three*dpad - two*dpdt_nl - dpis)
488 . / ((dpad-dpis)**3)
489 ENDIF
490 dtemp = weitemp(i)*dplanl(i)*yld(i)*eta/(rho0(i)*cp)
491 temp(i) = temp(i) + dtemp
492 ENDIF
493 ENDIF
494 ENDDO
495c
496C-----------------------------------------------
497C plastic work dissipated as heat - in case of /heat/mat
498C-----------------------------------------------
499 IF (jthe < 0 .AND. jlag /= 0) THEN
500 DO i=1,nel
501 fheat(i) = fheat(i) + eta*weitemp(i)*yld(i)*dpla(i)*volume(i)
502 ENDDO
503 END IF
504!-----------
505 RETURN
my_real
#define my_real
Definition cppsort.cpp:32
the
end diagonal values have been computed in the(sparse) matrix id.SOL
max
#define max(a, b)
Definition macros.h:21