OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
m2law.F File Reference
#include "implicit_f.inc"
#include "comlock.inc"
#include "mvsiz_p.inc"
#include "com08_c.inc"
#include "param_c.inc"
#include "impl1_c.inc"
#include "units_c.inc"
#include "lockon.inc"
#include "lockoff.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine m2law (pm, off, sig, eint, rho, qold, epxe, epsd, vol, stifn, dt2t, neltst, ityptst, offg, geo, pid, amu, vol_avg, mumax, mat, ngl, ssp, dvol, aire, vnew, vd2, deltax, vis, d1, d2, d3, d4, d5, d6, pnew, psh, qnew, ssp_eq, sold1, sold2, sold3, sold4, sold5, sold6, sigy, defp, dpla, epsp, tstar, etse, mssa, dmels, tempel, sigbak, al_imp, signor, conde, dtel, g_dt, nel, ipm, rhoref, rhosp, ipg, dmg, ity, jtur, jthe, jsph, ismstr, jsms, plap, npg, ieos, dpdm, fheat, glob_therm, jlag)

Function/Subroutine Documentation

◆ m2law()

subroutine m2law ( pm,
off,
sig,
eint,
rho,
qold,
epxe,
epsd,
vol,
stifn,
dt2t,
integer neltst,
integer ityptst,
offg,
geo,
integer, dimension(nel) pid,
amu,
vol_avg,
mumax,
integer, dimension(nel) mat,
integer, dimension(nel) ngl,
ssp,
dvol,
aire,
vnew,
vd2,
deltax,
vis,
d1,
d2,
d3,
d4,
d5,
d6,
pnew,
psh,
qnew,
ssp_eq,
sold1,
sold2,
sold3,
sold4,
sold5,
sold6,
sigy,
defp,
dpla,
epsp,
tstar,
etse,
mssa,
dmels,
tempel,
sigbak,
al_imp,
signor,
conde,
dtel,
integer g_dt,
integer nel,
integer, dimension(npropmi,*) ipm,
rhoref,
rhosp,
integer ipg,
dmg,
integer, intent(in) ity,
integer, intent(in) jtur,
integer, intent(in) jthe,
integer, intent(in) jsph,
integer, intent(in) ismstr,
integer, intent(in) jsms,
intent(inout) plap,
integer, intent(in) npg,
integer, intent(in) ieos,
intent(in) dpdm,
intent(inout) fheat,
type (glob_therm_), intent(inout) glob_therm,
integer, intent(in) jlag )

Definition at line 35 of file m2law.F.

56C-----------------------------------------------
57C M o d u l e s
58C-----------------------------------------------
59 use glob_therm_mod
60C-----------------------------------------------
61C I m p l i c i t T y p e s
62C-----------------------------------------------
63#include "implicit_f.inc"
64#include "comlock.inc"
65C-----------------------------------------------
66C G l o b a l P a r a m e t e r s
67C-----------------------------------------------
68#include "mvsiz_p.inc"
69C-----------------------------------------------
70C C o m m o n B l o c k s
71C-----------------------------------------------
72#include "com08_c.inc"
73#include "param_c.inc"
74#include "impl1_c.inc"
75#include "units_c.inc"
76C-----------------------------------------------
77C D u m m y A r g u m e n t s
78C-----------------------------------------------
79 INTEGER, INTENT(IN) :: NPG
80 INTEGER, INTENT(IN) :: ISMSTR
81 INTEGER, INTENT(IN) :: JSMS
82 INTEGER, INTENT(IN) :: ITY
83 INTEGER, INTENT(IN) :: JTUR
84 INTEGER, INTENT(IN) :: JTHE
85 INTEGER, INTENT(IN) :: JSPH
86 INTEGER, INTENT(IN) :: IEOS
87 INTEGER, INTENT(IN) :: JLAG
88C
89 INTEGER NELTST,ITYPTST,PID(NEL),G_DT,NEL,IPG
90 INTEGER MAT(NEL),NGL(NEL), IPM(NPROPMI,*)
91 my_real, DIMENSION(NEL) ,INTENT(INOUT) :: plap
92 my_real :: dt2t
93 my_real pm(npropm,*), off(nel), sig(nel,6), eint(nel), rho(nel), qold(nel),
94 . epxe(nel), epsd(nel), vol(nel), stifn(nel), offg(nel),geo(npropg,*),
95 . mumax(nel),amu(nel),vol_avg(nel)
96 my_real vnew(nel), vd2(nel), deltax(nel), ssp(nel), aire(nel), vis(nel),
97 . psh(nel), pnew(nel),qnew(nel) ,ssp_eq(nel), dvol(nel),
98 . d1(nel), d2(nel), d3(nel), d4(nel), d5(nel), d6(nel),
99 . sold1(mvsiz), sold2(mvsiz), sold3(mvsiz),
100 . sold4(mvsiz), sold5(mvsiz), sold6(mvsiz),
101 . tstar(mvsiz), dpla(mvsiz), epsp(mvsiz), dmg(nel)
102 my_real, DIMENSION(MVSIZ) ,INTENT(IN) :: dpdm
103 my_real, DIMENSION(NEL) ,INTENT(INOUT) :: fheat ! Heat due to plastic work for Heat Equation with lagrangian framework
104 my_real sigy(nel) ,defp(nel),etse(nel), mssa(nel), dmels(nel), tempel(nel),
105 . sigbak(nel,6),al_imp(nel),signor(mvsiz,6),conde(nel),dtel(nel),
106 . rhoref(nel) ,rhosp(nel)
107 type (glob_therm_) ,intent(inout) :: glob_therm
108C-----------------------------------------------
109C L o c a l V a r i a b l e s
110C-----------------------------------------------
111 INTEGER :: I, J, ICC,IRTY,ISRATE,VP,IDEV,MX,IBID,NINDX
112 INTEGER :: INDX(NEL)
113 my_real rho0, plap1,pold,pshift,bulk,p, epmx,cc,cn,epdr, cmx,
114 . z3,z4,rhocp,rhocpi,t_ref,t_melt,asrate,
115 . e1, e2, e3, e4, e5,e6, einc, g0,g1,g2,g3,g3h,
116 . scale, bid1, bid2, bid3, dta,mt,
117 . dsxx,dsyy,dszz,dsxy,dsyz,dszx,alpha,hkin,
118 . facq0,fisokin,beta,vm,vm_1,norm_1,ca0,cb,sigm0
119 my_real epd(mvsiz),g(mvsiz), ak(mvsiz),pc(mvsiz), qh(mvsiz),
120 . aj2(mvsiz), dav(mvsiz),ca(mvsiz),sigmx(mvsiz),
121 . sigexx(mvsiz),sigeyy(mvsiz),sigezz(mvsiz),
122 . sigexy(mvsiz),sigeyz(mvsiz),sigezx(mvsiz)
123 my_real, DIMENSION(MVSIZ) :: bidmvsiz
124C-----------------------------------------------
125C B o d y
126C-----------------------------------------------
127 facq0 = one
128C
129 mx = mat(1)
130 rho0 = pm( 1,mx)
131 bulk = pm(32,mx)
132 cn = pm(40,mx)
133 cc = pm(43,mx)
134 icc = nint(pm(49,mx))
135 epmx = pm(41,mx)
136 epdr = pm(44,mx)
137 irty = nint(pm(50,mx))
138 z3 = pm(51,mx)
139 cmx = pm(45,mx)
140 israte = ipm(3,mx)
141 asrate = min(one,pm(9,mx)*dt1)
142 fisokin = pm(55,mx)
143 g0 = pm(22,mx)
144 ca0 = pm(38,mx)
145 cb = pm(39,mx)
146 sigm0 = pm(42,mx)
147 vp = ipm(255,mx)
148 pshift = pm(88,mx)
149!
150 IF (irty == 1) THEN ! Zerilli
151 rhocpi = pm(53,mx)
152 t_ref = pm(54,mx)
153 z4 = pm(52,mx)
154 ELSE ! Johnson-Cook
155 rhocp = pm(69,mx)
156 IF (rhocp > zero) rhocpi = one / rhocp
157 t_ref = pm(79,mx)
158 t_melt = pm(80,mx)
159 END IF
160!
161 DO i=1,nel
162 g(i) = g0*off(i)
163 ca(i) = ca0
164 sigmx(i)= sigm0
165 etse(i) = one
166 ENDDO
167C
168C------------------------------------------
169C ECROUISSAGE CINE
170C------------------------------------------
171 IF (fisokin > zero) THEN
172 DO i=1,nel
173 sig(i,1)=sig(i,1)-sigbak(i,1)
174 sig(i,2)=sig(i,2)-sigbak(i,2)
175 sig(i,3)=sig(i,3)-sigbak(i,3)
176 sig(i,4)=sig(i,4)-sigbak(i,4)
177 sig(i,5)=sig(i,5)-sigbak(i,5)
178 sig(i,6)=sig(i,6)-sigbak(i,6)
179 ENDDO
180 ENDIF
181C
182 DO i=1,nel
183 p =-third*(sig(i,1)+sig(i,2)+sig(i,3))
184 dav(i)=-third*(d1(i)+d2(i)+d3(i))
185 g1=dt1*g(i)
186 g2=two*g1
187 ssp(i)=sqrt((onep333*g(i)+bulk)/rho0)
188c -- deviatoric elastic stress predictor
189 sig(i,1)=sig(i,1)+p+g2*(d1(i)+dav(i))
190 sig(i,2)=sig(i,2)+p+g2*(d2(i)+dav(i))
191 sig(i,3)=sig(i,3)+p+g2*(d3(i)+dav(i))
192 sig(i,4)=sig(i,4)+g1*d4(i)
193 sig(i,5)=sig(i,5)+g1*d5(i)
194 sig(i,6)=sig(i,6)+g1*d6(i)
195C
196c -- von Mises stress at elastic predictor
197 aj2(i)=half*(sig(i,1)**2+sig(i,2)**2+sig(i,3)**2)
198 . + sig(i,4)**2+sig(i,5)**2+sig(i,6)**2
199 aj2(i)=sqrt(three*aj2(i))
200 ENDDO
201C
202c -- storing elastic predictors for stiffness computation
203 IF (impl_s>0 .OR. fisokin>zero) THEN
204 DO i=1,nel
205 sigexx(i) = sig(i,1)
206 sigeyy(i) = sig(i,2)
207 sigezz(i) = sig(i,3)
208 sigexy(i) = sig(i,4)
209 sigeyz(i) = sig(i,5)
210 sigezx(i) = sig(i,6)
211 ENDDO
212 ENDIF
213C-------------
214C STRAIN RATE (JOHNSON-COOK, ZERILLI-ARMSTRONG)
215C-------------
216 idev = vp - 2
217 CALL mstrain_rate(nel ,israte ,asrate ,epsd ,idev ,
218 . d1, d2, d3, d4, d5, d6)
219
220 epsp(1:nel) = epsd(1:nel)
221 epd(1:nel) = one
222c
223 IF (cc /= zero) THEN
224 IF(vp == 1)THEN
225 DO i=1,nel
226 epd(i)= max(plap(i),epdr)
227 epd(i)= log(epd(i)/epdr)
228 ENDDO
229 ELSE
230 DO i=1,nel
231 epd(i)= max(epsd(i),em15)
232 epd(i)= log(epd(i)/epdr)
233 ENDDO
234 ENDIF ! VP
235!
236 IF (irty == 0) THEN ! J-C
237 DO i=1,nel
238 mt=max(em15,z3)
239 epd(i)= max(zero,epd(i))
240 epd(i)= (one +cc * epd(i))*(one -tstar(i)**mt)
241 IF(icc==1) sigmx(i)= sigm0*epd(i)
242 ENDDO
243 ELSEIF (irty==1) THEN
244 DO i=1,nel
245 epd(i)= cc*exp((-z3+z4 * epd(i))*tempel(i))
246 IF (icc==1) sigmx(i) = sigm0 + epd(i)
247 ca(i) = ca0 + epd(i)
248 epd(i)=one
249 ENDDO
250 END IF
251 ELSEIF (irty == 0) THEN
252 mt = max(em15,z3)
253 epd(1:nel) = one - tstar(1:nel)**mt
254 ENDIF
255C-------------
256C CRITERE
257C-------------
258 IF (fisokin == zero) THEN ! isotropic hardening
259 IF(cn==one) THEN
260 ak(1:nel)= ca(1:nel)+cb*epxe(1:nel)
261 qh(1:nel)= cb*epd(1:nel)
262 ELSE
263 DO i=1,nel
264 IF(epxe(i)>zero) THEN
265 ak(i)=ca(i)+cb*epxe(i)**cn
266 IF(cn>one) THEN
267 qh(i) = (cb*cn*epxe(i)**(cn - one))*epd(i)
268 ELSE
269 qh(i) = (cb*cn/epxe(i)**(one -cn))*epd(i)
270 ENDIF
271 ELSE
272 ak(i)=ca(i)
273 qh(i)=zero
274 ENDIF
275 ENDDO
276 ENDIF
277 DO i=1,nel
278 ak(i)=ak(i)*epd(i)
279 IF(sigmx(i)<ak(i))THEN
280 ak(i)=sigmx(i)
281 qh(i)=zero
282 ENDIF
283 sigy(i) = ak(i)
284 IF(epxe(i)>epmx)THEN
285 ak(i)=zero
286 qh(i)=zero
287 ENDIF
288 ENDDO
289!
290 ELSE ! kinematic hardening
291!
292 DO i=1,nel
293 beta = one-fisokin
294 ! SIGY is used for hourglass stress compute--
295 IF(cn==one) THEN
296 sigy(i) = ca(i)+cb*epxe(i)
297 ak(i)= ca(i)+beta*cb*epxe(i)
298 qh(i)= cb*epd(i)
299 ELSE
300 IF(epxe(i)>zero) THEN
301 sigy(i)=ca(i)+cb*epxe(i)**cn
302 ak(i)=ca(i)+beta*cb*epxe(i)**cn
303 IF(cn>one) THEN
304 qh(i)= (cb*cn*epxe(i)**(cn - one))*epd(i)
305 ELSE
306 qh(i)= (cb*cn/epxe(i)**(one -cn))*epd(i)
307 ENDIF
308 ELSE
309 ak(i)=ca(i)
310 sigy(i)=ca(i)
311 qh(i)=zero
312 ENDIF
313 ENDIF
314 ak(i)=ak(i)*epd(i)
315 sigy(i)=sigy(i)*epd(i)
316 IF(sigmx(i)<ak(i))THEN
317 ak(i)=sigmx(i)
318 qh(i)=zero
319 ENDIF
320 sigy(i)=min(sigy(i),sigmx(i))
321 IF(epxe(i)>epmx)THEN
322 ak(i)=zero
323 qh(i)=zero
324 ENDIF
325 END DO ! I=1,NEL
326 END IF ! (FISOKIN == 0 ) THEN
327!---------------------------------------------------------
328 IF (impl_s==0) THEN
329 DO i=1,nel
330 scale= min(one,ak(i)/ max(aj2(i),em15))
331C plastic strain increment.
332 dpla(i)=(one-scale)*aj2(i)/max(three*g(i)+qh(i),em15)
333C actual yield stress.
334 ak(i)=ak(i)+(one - fisokin)*dpla(i)*qh(i)
335 scale= min(one,ak(i)/ max(aj2(i),em15))
336 sig(i,1)=scale*sig(i,1)
337 sig(i,2)=scale*sig(i,2)
338 sig(i,3)=scale*sig(i,3)
339 sig(i,4)=scale*sig(i,4)
340 sig(i,5)=scale*sig(i,5)
341 sig(i,6)=scale*sig(i,6)
342 epxe(i)=epxe(i)+dpla(i)
343 ENDDO
344!
345 ELSE ! implicit : nonlinear hardening requires iterations in radial return
346 CALL m2iter_imp(sig, epxe, aj2, g,
347 1 ca, cb, cn, epd,
348 2 sigmx, epmx, dpla, ak,
349 3 qh, sigy, fisokin, nel)
350!
351 END IF ! IMPL_S
352c-----------------------------------------
353 IF (fisokin > zero) THEN
354 DO i=1,nel
355 dsxx = sigexx(i) - sig(i,1)
356 dsyy = sigeyy(i) - sig(i,2)
357 dszz = sigezz(i) - sig(i,3)
358 dsxy = sigexy(i) - sig(i,4)
359 dsyz = sigeyz(i) - sig(i,5)
360 dszx = sigezx(i) - sig(i,6)
361C
362 hkin = two_third*fisokin*qh(i)
363 alpha = hkin/max(two*g(i)+hkin,em15)
364C ..updates back stresses
365 sigbak(i,1) = sigbak(i,1) + alpha*dsxx
366 sigbak(i,2) = sigbak(i,2) + alpha*dsyy
367 sigbak(i,3) = sigbak(i,3) + alpha*dszz
368 sigbak(i,4) = sigbak(i,4) + alpha*dsxy
369 sigbak(i,5) = sigbak(i,5) + alpha*dsyz
370 sigbak(i,6) = sigbak(i,6) + alpha*dszx
371C ..gets stresses from shifted stresses and back stresses
372 sig(i,1)=sig(i,1) + sigbak(i,1)
373 sig(i,2)=sig(i,2) + sigbak(i,2)
374 sig(i,3)=sig(i,3) + sigbak(i,3)
375 sig(i,4)=sig(i,4) + sigbak(i,4)
376 sig(i,5)=sig(i,5) + sigbak(i,5)
377 sig(i,6)=sig(i,6) + sigbak(i,6)
378 ENDDO
379 END IF !(FISOKIN > 0 ) THEN
380c-----------------------------------------
381!
382 bidmvsiz(1:mvsiz) = zero
383 IF (jsph==0) THEN
384 CALL mqviscb(
385 1 pm, off, rho, bidmvsiz,
386 2 bidmvsiz,ssp, bidmvsiz,stifn,
387 3 dt2t, neltst, ityptst, aire,
388 4 offg, geo, pid, vnew,
389 5 vd2, deltax, vis, d1,
390 6 d2, d3, pnew, psh,
391 7 mat, ngl, qnew, ssp_eq,
392 8 vol, mssa, dmels, ibid,
393 9 facq0, conde, dtel, g_dt,
394 a ipm, rhoref, rhosp, nel,
395 b ity, ismstr, jtur, jthe,
396 c jsms, npg , glob_therm)
397 ELSE
398 CALL mdtsph(
399 1 pm, off, rho, bidmvsiz,
400 2 bidmvsiz,bidmvsiz,stifn, dt2t,
401 3 neltst, ityptst, offg, geo,
402 4 pid, mumax, ssp, vnew,
403 5 vd2, deltax, vis, d1,
404 6 d2, d3, pnew, psh,
405 7 mat, ngl, qnew, ssp_eq,
406 8 g_dt, dtel, nel, ity,
407 9 jtur, jthe)
408 ENDIF
409C
410 dta = half*dt1
411C
412 nindx = 0
413 indx(1:nel) = 0
414 DO i=1,nel
415 IF ((epxe(i) > epmx).AND.(dmg(i)==zero)) THEN
416 nindx = nindx + 1
417 indx(nindx) = i
418 dmg(i) = one
419 ENDIF
420 ENDDO
421
422 IF (ieos == 0) THEN ! add pressure to the deviatoric stress
423 DO i=1,nel
424 pnew(i) = bulk*amu(i)
425 sig(i,1) =(sig(i,1) - pnew(i))*off(i)
426 sig(i,2) =(sig(i,2) - pnew(i))*off(i)
427 sig(i,3) =(sig(i,3) - pnew(i))*off(i)
428 sig(i,4) = sig(i,4) * off(i)
429 sig(i,5) = sig(i,5) * off(i)
430 sig(i,6) = sig(i,6) * off(i)
431 ENDDO
432!
433 DO i=1,nel
434 e1 = d1(i)*(sold1(i)+sig(i,1))
435 e2 = d2(i)*(sold2(i)+sig(i,2))
436 e3 = d3(i)*(sold3(i)+sig(i,3))
437 e4 = d4(i)*(sold4(i)+sig(i,4))
438 e5 = d5(i)*(sold5(i)+sig(i,5))
439 e6 = d6(i)*(sold6(i)+sig(i,6))
440 einc = vol_avg(i)*(e1+e2+e3+e4+e5+e6)*dta - half*dvol(i)*(qold(i)+qnew(i))
441 eint(i) = (eint(i)+einc*off(i)) / max(em15,vol(i))
442 ENDDO
443
444 ELSE
445 ! if EOS is used, material law calculates only deviatoric stress tensor
446 ! sound speed depends on pressure derivative over volume change
447 ! calculated in EOS
448 DO i = 1, nel
449 ssp(i) = sqrt((dpdm(i) + four*g(i)/three)/rho0)
450!
451 pold = (sold1(i)+sold2(i)+sold3(i)) * third
452 sold1(i) = sold1(i) - pold
453 sold2(i) = sold2(i) - pold
454 sold3(i) = sold3(i) - pold
455 e1 = d1(i) * (sold1(i)+sig(i,1))
456 e2 = d2(i) * (sold2(i)+sig(i,2))
457 e3 = d3(i) * (sold3(i)+sig(i,3))
458 e4 = d4(i) * (sold4(i)+sig(i,4))
459 e5 = d5(i) * (sold5(i)+sig(i,5))
460 e6 = d6(i) * (sold6(i)+sig(i,6))
461 einc = vol_avg(i) * (e1+e2+e3+e4+e5+e6) * dta
462 eint(i) = eint(i) + (einc+half*dvol(i)*(pold-pshift-qold(i)-qnew(i)))*off(i)
463 ENDDO
464 END IF
465C
466 DO i=1,nel
467 qold(i) = qnew(i)
468 defp(i) = epxe(i)
469 sigy(i) = max(sigy(i),ak(i))
470 ENDDO
471 DO i=1,nel
472 IF (dpla(i) > 0) etse(i)= half*qh(i)*off(i)/max(g(i),em15)
473 ENDDO
474!----------------------------------------------------------------
475! implicit
476!----------------------------------------------------------------
477 IF (impl_s>0) THEN
478 IF (ikt==0) RETURN
479 IF (fisokin == zero) THEN
480 DO i=1,nel
481 IF (dpla(i)>0)THEN
482 IF(cn==one) THEN
483 qh(i)= cb*epd(i)
484 ELSEIF(cn>one) THEN
485 qh(i)= (cb*cn*epxe(i)**(cn - one))*epd(i)
486 ELSE
487 qh(i)= (cb*cn/epxe(i)**(one -cn))*epd(i)
488 ENDIF
489 ENDIF
490 ENDDO
491 ENDIF
492C
493 DO i = 1,nel
494 IF (dpla(i)>zero) THEN
495
496c ..... Von Mises stress at the elastic predictor (point B)
497c ..... SIGEXX, etc. are deviatoric stresses
498 vm =half*(sigexx(i)**2+sigeyy(i)**2+sigezz(i)**2)
499 1 +sigexy(i)**2+sigeyz(i)**2+sigezx(i)**2
500 vm_1 =one/sqrt(three*vm)
501 g3 = three*g(i)
502 g3h = max(g3+qh(i),em15)
503 scale = max(zero,one-g3h*dpla(i)*vm_1)
504c ..... NORM_1 normalizes deviatoric stresses, includes consistent
505c ..... stiffness matrix parameter beta, von Mises at B, and two_pmi
506 norm_1=g3*vm_1*sqrt(scale/g3h)
507c ..... Deviatoric stresses "normalized"
508 signor(i,1)=sigexx(i)*norm_1
509 signor(i,2)=sigeyy(i)*norm_1
510 signor(i,3)=sigezz(i)*norm_1
511 signor(i,4)=sigexy(i)*norm_1
512 signor(i,5)=sigeyz(i)*norm_1
513 signor(i,6)=sigezx(i)*norm_1
514
515c ..... Parameter alpha of consistent matrix
516 al_imp(i)= one - g3*dpla(i)*vm_1
517 ELSE
518 al_imp(i)=one
519 ENDIF
520 ENDDO
521 ENDIF
522!---------------------------------------------------------
523 ! end implicit
524!---------------------------------------------------------
525! update and filter plastic strain rate for VP=1
526!---------------------------------------------------------
527!
528 IF (vp== 1) THEN
529 DO i=1,nel
530 plap1 = dpla(i)/max(em20,dt1)
531 plap(i) = asrate * plap1 + (one - asrate) * plap(i)
532 ENDDO
533 ENDIF
534C----------------------------------------------
535C TEMPERATURE (Heating due to plastic work)
536C----------------------------------------------
537 IF (jthe < 0 .AND. jlag /= 0) THEN
538 DO i=1,nel
539 fheat(i) = fheat(i) + sigy(i)*dpla(i)*vol(i)
540 ENDDO
541 ELSEIF(rhocp > zero)THEN
542 DO i=1,nel
543 tempel(i) = tempel(i) + sigy(i)*dpla(i) / rhocp
544 ! temperature and internal energy must be incremented consistantly
545 ! internal energy is incremented later in parent subroutine (mmain)
546 ! with total energy deformation which already includes plastic work
547 ! Edef = 0.5 * VOL * sum ( sig.eps_dot , i=1..6)
548 ! so internal energy and temperature remain consistant
549 ENDDO
550 END IF
551!---------------------------------------------------------
552 ! Printout element deletion
553 IF (nindx > 0) THEN
554 DO j=1,nindx
555#include "lockon.inc"
556 WRITE(iout, 1000) ngl(indx(j)),ipg
557 WRITE(istdo,1100) ngl(indx(j)),ipg,tt
558#include "lockoff.inc"
559 ENDDO
560 ENDIF
561!---------------------------------------------------------------
562 1000 FORMAT(1x,'EXCEEDED EPS_MAX ON SOLID ELEMENT NUMBER ',i10,
563 . ': DEVIATORIC STRESS SET TO 0 ON INTEGRATION POINT ',i5 )
564 1100 FORMAT(1x,'EXCEEDED EPS_MAX ON SOLID ELEMENT NUMBER ',i10,
565 . ': DEVIATORIC STRESS SET TO 0 ON INTEGRATION POINT ',i5 ,
566 . ' AT TIME :',g11.4)
567!-----------
568 RETURN
#define my_real
Definition cppsort.cpp:32
subroutine dtel(ssp, pm, geo, pid, mat, rho0, vis, deltax, aire, vol, dtx)
Definition dtel.F:46
#define alpha
Definition eval.h:35
subroutine m2iter_imp(sig, epxe, aj2, g, ca, cb, cn, epd, sigmx, epmx, dpla1, ak, qh, sigy, fisokin, nel)
Definition m2iter_imp.F:33
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
subroutine mdtsph(pm, off, rho, rk, t, re, sti, dt2t, neltst, ityptst, offg, geo, pid, mumax, ssp, vol, vd2, deltax, vis, d1, d2, d3, pnew, psh, mat, ngl, qvis, ssp_eq, g_dt, dtsph, nel, ity, jtur, jthe)
Definition mdtsph.F:47
subroutine mqviscb(pm, off, rho, rk, temp, ssp, re, sti, dt2t, neltst, ityptst, aire, offg, geo, pid, vol, vd2, deltax, vis, d1, d2, d3, pnew, psh, mat, ngl, qvis, ssp_eq, vol0, mssa, dmels, igeo, facq0, conde, dtel, g_dt, ipm, rhoref, rhosp, nel, ity, ismstr, jtur, jthe, jsms, npg, glob_therm)
Definition mqviscb.F:56
subroutine mstrain_rate(nel, israte, asrate, epsd, idev, ep1, ep2, ep3, ep4, ep5, ep6)