53
54
55
56 USE matparam_def_mod , ONLY : matparam_struct_
58
59
60
61#include "implicit_f.inc"
62#include "comlock.inc"
63
64
65
66#include "mvsiz_p.inc"
67
68
69
70#include "com01_c.inc"
71#include "param_c.inc"
72#include "com08_c.inc"
73#include "scr17_c.inc"
74#include "units_c.inc"
75 COMMON /tablesizf/ stf,snpc
76 INTEGER STF,SNPC
77
78
79
80 TYPE(MATPARAM_STRUCT_), INTENT(IN) :: MAT_PARAM
81 INTEGER, INTENT(IN) :: JCVT
82 INTEGER, INTENT(IN) :: JSPH
83 INTEGER MAT(MVSIZ),NGL(MVSIZ),NEL,EOSTYP
84 INTEGER,INTENT(IN) :: NVAREOS
86 . pm(npropm,*), off(*), sig(nel,6), eint(*), pla(*), sigf(*),
87 . rx(*), ry(*), rz(*), sx(*), sy(*), sz(*),
88 . epsf(*), dam(nel,5), epe(nel,3), epc(nel,3), a(mvsiz,6), vol(*),
89 . vnew(mvsiz), dvol(mvsiz), ssp(mvsiz),
90 . d1(mvsiz), d2(mvsiz), d3(mvsiz),d4(mvsiz),d5(mvsiz),d6(mvsiz),
91 . sold1(mvsiz), sold2(mvsiz), sold3(mvsiz), sold4(mvsiz),
92 . sold5(mvsiz), sold6(mvsiz),sigy(*),defp(*),seq_output(*),
93 . tsaiwu(mvsiz),vareos(nvareos*nel)
95 . rho0(nel), amu(nel) , amu2(nel), espe(nel), df
96 . psh(nel) , pnew(nel) , dpdm(nel), dpde(nel), rho(nel),
97 . temp(nel), bufmat(*),muold(nel)
98 INTEGER,INTENT(IN)::NPC(SNPC)
100 INTEGER,INTENT(IN) :: NVARTMP_EOS
101 INTEGER,INTENT(INOUt) :: VARTMP_EOS(NEL,NVARTMP_EOS)
102
103
104
105 INTEGER (MVSIZ), KD2(MVSIZ), KD3(MVSIZ),
106 . KD4(MVSIZ),ICC(MVSIZ),
107 . I, IDAM, KDX, MX, NINDX, INDEX(MVSIZ), J, ICC_1
110 . ax(mvsiz), ay(mvsiz), az(mvsiz), bx(mvsiz), by(mvsiz),
111 . bz(mvsiz), cx(mvsiz), cy(mvsiz),cz(mvsiz),
112 . wvec(mvsiz),s1(mvsiz), s2(mvsiz), s3(mvsiz),
113 . t1(mvsiz), t2(mvsiz), t3(mvsiz),t4(mvsiz),t5(mvsiz),t6(mvsiz),
114 . e1(mvsiz), e2(mvsiz), e3(mvsiz),e4(mvsiz),e5(mvsiz),e6(mvsiz),
115 . hard(mvsiz), sigmy(mvsiz),
alpha(mvsiz), efib(mvsiz),
116 . epsft(mvsiz), epsfc(mvsiz), sigeq(mvsiz), p(mvsiz),
117 . d11(mvsiz), d12(mvsiz), d13(mvsiz), d22(mvsiz),
118 . d23(mvsiz), d33(mvsiz), g12(mvsiz), g23(mvsiz), g31(mvsiz),
119 . c11(mvsiz), c22(mvsiz), c33(mvsiz), c12(mvsiz), c23(mvsiz),
120 . c13(mvsiz), f11(mvsiz), f22(mvsiz), f44(mvsiz), f55(mvsiz),
121 . f12(mvsiz), f23(mvsiz), f1(mvsiz), f2(mvsiz), f4(mvsiz),
122 . f5(mvsiz), delta(mvsiz), so1(mvsiz),
123 . so2(mvsiz), so3(mvsiz), so4(mvsiz), so5(mvsiz), so6(mvsiz),
124 . ds1(mvsiz), ds2(mvsiz), ds3(mvsiz), ds4(mvsiz), ds5(mvsiz),
125 . ds6(mvsiz), dp1(mvsiz), dp2(mvsiz), dp3(mvsiz), dp4(mvsiz),
126 . dp5(mvsiz), dp6(mvsiz), lamda(mvsiz), coef(mvsiz), plas(mvsiz),
127 . cn(mvsiz), cb(mvsiz), cnn(mvsiz), sigt1(mvsiz), sigt2(mvsiz),
128 . sigt3(mvsiz),cc(mvsiz),epdr(mvsiz), f3(mvsiz), f33(mvsiz),
129 . f6(mvsiz),f66(mvsiz),f13(mvsiz)
131 . deve1(mvsiz), deve2(mvsiz), deve3(mvsiz), dav(mvsiz), pold(mvsiz),
132 . ener(mvsiz)
134 . fib, ca, sigmx, epsp, dt5 ,
135 . sigym,ca_1, cb_1, cn_1, cc_1,
136 . epdr_1,alpha_1,f1_1,f2_1,f3_1,
137 . f4_1,f5_1,f6_1,f11_1,f22_1,
138 . f33_1,f44_1,f55_1,f66_1,f12_1,
139 . f23_1,f13_1,c11_1,c22_1,c33_1,
140 . c12_1,c23_1,c13_1,d11_1,d12_1,
141 . d13_1,d22_1,d23_1,d33_1,g12_1,
142 . g23_1,g31_1,wplaref,dwpla
143
144
145 IF(ncycle==0 .AND. irun==1) THEN
146 DO i=1,nel
147 kd1(i) = 0
148 kd2(i) = 0
149 kd3(i) = 0
150 kd4(i) = 0
151 ENDDO
152 ELSE
153 DO i=1,nel
154 idam=int(dam(i,5))-10000
155 kd1(i) = idam/1000
156 kdx = idam - kd1(i)*1000
157 kd2(i) = kdx/100
158 kdx = kdx - kd2(i)*100
159 kd3(i) = kdx/10
160 kd4(i) = kdx - kd3(i)*10
161
162 ENDDO
163 ENDIF
164
165
166
167 fib=zero
168 mx =mat(1)
169 alpha_1 =pm(39,mx)
170 DO i=1,nel
172 fib=fib+alpha_1
173 ENDDO
174
175
176
177
179 1 pm, a, rx, ry,
180 2 rz, sx, sy, sz,
181 3 ax, ay, az, bx,
182 4 by, bz, cx, cy,
183 5 cz, nel, jcvt, jsph)
184 CALL m14gtf(sig,ax ,ay ,az ,bx ,by,
185 2 bz ,cx ,cy ,cz ,d1 ,d2,
186 3 d3 ,d4 ,d5 ,d6 ,t1 ,t2,
187 4 t3 ,t4 ,t5 ,t6 ,e1 ,e2,
188 5 e3 ,e4 ,e5 ,e6 ,nel)
189
190
191
192
193 nindx=0
194 mx = mat(1)
195 ca_1 = pm(25,mx)
196 cb_1 = pm(26,mx)
197 cn_1 = pm(27,mx)
198 cc_1 = pm(50,mx)
199 epdr_1 = pm(51,mx)
200 icc_1 = nint(pm(52,mx))
201 DO i=1,nel
202 ca = ca_1
203 cb(i) = cb_1
204 cn(i) = cn_1
205 cc(i) = cc_1
206 epdr(i) = epdr_1
207 icc(i) = icc_1
208 epsp =
max(abs(e1(i)),abs(e2(i)),abs(e3(i)),
209 . half*abs(e4(i)),half*abs(e5(i)),half*abs(e6(i)))
210 IF(epsp>epdr(i).AND.cc(i)/=zero) THEN
211 epsp = one + cc(i) * log(epsp/epdr(i))
212 ELSE
213 epsp = one
214 ENDIF
215 IF(icc(i)==1)THEN
216 sigmx =pm(28,mx)*epsp
217 ELSEIF(icc(i)==2)THEN
218 sigmx =pm(28,mx)
219 ELSEIF(icc(i)==3)THEN
220 sigmx =pm(28,mx)*epsp
221 ELSEIF(icc(i)==4)THEN
222 sigmx =pm(28,mx)
223 ENDIF
224 cb(i) =cb(i)*epsp
225 ca =ca*epsp
226 sigmy(i)=
min(sigmx,ca+cb(i)*pla(i)**cn(i))
227 IF(sigmy(i)==sigmx .AND. off(i)==one) THEN
228 off(i)=zep99
229 kd4(i)=2
230
231 nindx=nindx+1
232 index(nindx)=i
233 ENDIF
234 sigt1(i)= pm(81,mx)
235 sigt2(i)= pm(82,mx)
236 sigt3(i)= pm(83,mx)
237 ssp(i) = pm(49,mx)
238 delta(i)= pm(79,mx)
239 ENDDO
240
241 IF(nindx/=0)THEN
242 DO j=1,nindx
243 i=index(j)
244#include "lockon.inc"
245 WRITE(iout,1000) ngl(i)
246#include "lockoff.inc"
247 END DO
248 END IF
249
250 DO i=1,nel
251 e1(i)=e1(i)*dt1
252 e2(i)=e2(i)*dt1
253 e3(i)=e3(i)*dt1
254 e4(i)=e4(i)*dt1
255 e5(i)=e5(i)*dt1
256 e6(i)=e6(i)*dt1
257 ENDDO
258
259 DO i=1,nel
260 epe(i,1)=epe(i,1)+e1(i)
261 epe(i,2)=epe(i,2)+e2(i)
262 epe(i,3)=epe(i,3)+e3(i)
263 ENDDO
264
265 mx =mat(1)
266 f1_1 =pm(53,mx)
267 f2_1 =pm(54,mx)
268 f3_1 =pm(55,mx)
269 f4_1 =pm(56,mx)
270 f5_1 =pm(57,mx)
271 f6_1 =pm(58,mx)
272 f11_1 =pm(59,mx)
273 f22_1 =pm(60,mx)
274 f33_1 =pm(61,mx)
275 f44_1 =pm(62,mx)
276 f55_1 =pm(63,mx)
277 f66_1 =pm(64,mx)
278 f12_1 =pm(65,mx)
279 f23_1 =pm(66,mx)
280 f13_1 =pm(67,mx)
281
282 c11_1 =pm(73,mx)
283 c22_1 =pm(74,mx)
284 c33_1 =pm(75,mx)
285 c12_1 =pm(76,mx)
286 c23_1 =pm(77,mx)
287 c13_1 =pm(78,mx)
288
289 d11_1 =pm(40,mx)
290 d12_1 =pm(41,mx)
291 d13_1 =pm(42,mx)
292 d22_1 =pm(43,mx)
293 d23_1 =pm(44,mx)
294 d33_1 =pm(45,mx)
295 g12_1 =pm(46,mx)
296 g23_1 =pm(47,mx)
297 g31_1 =pm(48,mx)
298
299 DO i=1,nel
300
301 f1(i) = f1_1
302 f2(i) = f2_1
303 f3(i) = f3_1
304 f4(i) = f4_1
305 f5(i) = f5_1
306 f6(i) = f6_1
307 f11(i) = f11_1
308 f22(i) = f22_1
309 f33(i) = f33_1
310 f44(i) = f44_1
311 f55(i) = f55_1
312 f66(i) = f66_1
313 f12(i) = f12_1
314 f23(i) = f23_1
315 f13(i) = f13_1
316
317 c11(i) = c11_1
318 c22(i) = c22_1
319 c33(i) = c33_1
320 c12(i) = c12_1
321 c23(i) = c23_1
322 c13(i) = c13_1
323
324 d11(i) = d11_1
325 d12(i) = d12_1
326 d13(i) = d13_1
327 d22(i) = d22_1
328 d23(i) = d23_1
329 d33(i) = d33_1
330 g12(i) = g12_1
331 g23(i) = g23_1
332 g31(i) = g31_1
333 ENDDO
334
335
336
337 DO i=1,nel
338 so1(i)=t1(i)
339 so2(i)=t2(i)
340 so3(i)=t3(i)
341 so4(i)=t4(i)
342 so5(i)=t5(i)
343 so6(i)=t6(i)
344 ENDDO
345
346 IF (eostyp == 18) THEN
347 DO i=1,nel
348 t1(i)=t1(i)+d11(i)*e1(i)+d12(i)*e2(i)+d13(i)*e3(i)
349 t2(i)=t2(i)+d12(i)*e1(i)+d22(i)*e2(i)+d23(i)*e3(i)
350 t3(i)=t3(i)+d13(i)*e1(i)+d23(i)*e2(i)+d33(i)*e3(i)
351 t4(i)=t4(i)+g12(i)*e4(i)
352 t5(i)=t5(i)+g23(i)*e5(i)
353 t6(i)=t6(i)+g31(i)*e6(i)
354 ENDDO
355 ELSE
356 DO i=1,nel
357 dav(i) = third*(e1(i)+e2(i)+e3(i))
358 pold(i) = third*(t1(i)+t2(i)+t3(i))
359 ener(i) = eint(i)
360 ENDDO
361 DO i=1,nel
362 deve1(i) = e1(i) - dav(i)
363 deve2(i) = e2(i) - dav(i)
364 deve3(i) = e3(i) - dav(i)
365 ENDDO
366 DO i=1,nel
367 t1(i)=t1(i)+d11(i)*deve1(i)+d12(i)*deve2(i)+d13(i)*deve3(i)-pold(i)
368 t2(i)=t2(i)+d12(i)*deve1(i)+d22(i)*deve2(i)+d23(i)*deve3(i)-pold(i)
369 t3(i)=t3(i)+d13(i)*deve1(i)+d23(i)*deve2(i)+d33(i)*deve3(i)-pold(i)
370 ENDDO
371
372 muold(1:nel) =( rho(1:nel)/rho0)-one
373 sig(1:nel,1:6) = zero
374 CALL eosmain(0 ,nel ,eostyp ,pm ,off , ener ,
375 2 rho ,rho0 ,amu ,amu2 ,espe ,
376 3 dvol ,df ,vnew ,mat ,psh ,
377 4 pnew ,dpdm ,dpde ,temp ,
378 5 bufmat ,sig ,muold ,12 ,
379 6 npc ,tf ,vareos ,nvareos, mat_param, bid, nvartmp_eos,vartmp_eos)
380 CALL eosmain(1 ,nel ,eostyp ,pm ,off , ener ,
381 2 rho ,rho0 ,amu ,amu2 ,espe ,
382 3 dvol ,df ,vnew ,mat ,psh ,
383 4 pnew ,dpdm ,dpde ,temp ,
384 5 bufmat ,sig ,muold ,12 ,
385 6 npc ,tf ,vareos ,nvareos, mat_param, bid, nvartmp_eos,vartmp_eos)
386 DO i=1,nel
387 t1(i)=t1(i)-pnew(i)
388 t2(i)=t2(i)-pnew(i)
389 t3(i)=t3(i)-pnew(i)
390 ENDDO
391 DO i=1,nel
392 ssp(i) = sqrt(abs(dpdm(i))/rho0(i))
393 ENDDO
394
395 ENDIF
396
397
398
399 IF (fib>0) THEN
400 mx = mat(1)
401 DO i=1,nel
402 efib(i) =pm(38,mx)
403 epsft(i)=pm(84,mx)
404 epsfc(i)=pm(85,mx)
405 epsf(i) = epsf(i)+ e1(i)
406 sigf(i) = efib(i)*epsf(i)
407 ENDDO
408
409
410
411 ENDIF
412
413
414
415 DO i=1,nel
416 IF(off(i) < em01) off(i)=zero
417 IF(off(i) < one) off(i)=off(i)*four_over_5
418 ENDDO
419
420
421
422 DO i=1,nel
423
424
425
426 wvec(i)=(one-dam(i,1))*sigt1(i)
427 IF(t1(i) > wvec(i)) THEN
428
429 IF(epc(i,1) == zero) THEN
430 epc(i,1)=
max(epe(i,1),zero)
431 ELSE
432 epc(i,1)=
max(epc(i,1)+e1(i),zero)
433 ENDIF
434
435 t1(i)=wvec(i)
436 t2(i)=t2(i)-d12(i)*e1(i)*dam(i,1)
437 t3(i)=t3(i)-d13(i)*e1(i)*dam(i,1)
438 IF(kd1(i)==0) THEN
439 kd1(i)=1
440 ENDIF
441 dam(i,1)=
min((dam(i,1)+delta(i)),one)
442 IF(dam(i,1)>=one .AND. kd1(i)/=2) THEN
443 kd1(i)=2
444 ENDIF
445 ENDIF
446
447 IF(e1(i) < zero .AND. dam(i,1) > zero)
448 . epc(i,1)=
max(epc(i,1)+e1(i),zero)
449
450
451
452 wvec(i)=(one-dam(i,2))*sigt2(i)
453 IF(t2(i) > wvec(i)) THEN
454
455 IF(epc(i,2)==zero) THEN
456 epc(i,2)=
max(epe(i,2),zero)
457 ELSE
458 epc(i,2)=
max(epc(i,2)+e2(i),zero)
459 ENDIF
460
461 t1(i)=t1(i)-d12(i)*e2(i)*dam(i,2)
462 t2(i)=wvec(i)
463 t3(i)=t3(i)-d23(i)*e2(i)*dam(i,2)
464 IF(kd2(i) == 0) THEN
465 kd2(i)=1
466 ENDIF
467 dam(i,2)=
min((dam(i,2)+delta(i)),one)
468 IF(dam(i,2) >= one .AND. kd2(i) /= 2) THEN
469 kd2(i)=2
470 ENDIF
471 ENDIF
472
473 IF(e2(i) < zero .AND. dam(i,2) > zero)
474 . epc(i,2)=
max(epc(i,2)+e2(i),zero)
475
476
477
478 wvec(i)=(one-dam(i,3))*sigt3(i)
479 IF(t3(i) > wvec(i)) THEN
480
481 IF(epc(i,3)==zero) THEN
482 epc(i,3)=
max(epe(i,3),zero)
483 ELSE
484 epc(i,3)=
max(epc(i,3)+e3(i),zero)
485 ENDIF
486
487 t1(i)=t1(i)-d13(i)*e3(i)*dam(i,3)
488 t2(i)=t2(i)-d23(i)*e3(i)*dam(i,3)
489 t3(i)=wvec(i)
490 IF(kd3(i)==0) THEN
491 kd3(i)=1
492 ENDIF
493 dam(i,3)=
min((dam(i,3)+delta(i)),one)
494 IF(dam(i,3) >= one .AND. kd3(i) /= 2) THEN
495 kd3(i)=2
496 ENDIF
497 ENDIF
498
499 IF(e3(i)<zero .AND. dam(i,3)>zero)
500 . epc(i,3)=
max(epc(i,3)+e3(i),zero)
501
502 ENDDO
503
504
505
506 DO i=1,nel
507 IF(t1(i) < zero.AND. epc(i,1) > zero) THEN
508 t1(i)=zero
509 t2(i)=t2(i)-d12(i)*e1(i)*dam(i,1)
510 t3(i)=t3(i)-d13(i)*e1(i)*dam(i,1)
511 ENDIF
512 IF(t2(i) < zero.AND. epc(i,2) > zero) THEN
513 t1(i)=t1(i)-d12(i)*e2(i)*dam(i,2)
514 t2(i)=zero
515 t3(i)=t3(i)-d23(i)*e2(i)*dam(i,2)
516 ENDIF
517 IF(t3(i) < zero.AND.epc(i,3) > zero) THEN
518 t3(i)=zero
519 t1(i)=t1(i)-d13(i)*e3(i)*dam(i,3)
520 t2(i)=t2(i)-d23(i)*e3(i)*dam(i,3)
521 ENDIF
522 ENDDO
523
524
525
526 DO i=1,nel
527 wvec(i)= f1(i)*t1(i) + f2(i)*t2(i) + f3(i)*t3(i)
528 . + f11(i)*t1(i)*t1(i) + f22(i)*t2(i)*t2(i) + f33(i)*t3(i)*t3(i)
529 . + f55(i)*t5(i)*t5(i) + f44(i)*t4(i)*t4(i) + f66(i)*t6(i)*t6(i)
530 . + two*f12(i)*t1(i)*t2(i) + two*f13(i)*t1(i)*t3(i)
531 . + two*f23(i)*t2(i)*t3(i)
532 dam(i,4)=wvec(i)
533 tsaiwu(i)=
max(
min(wvec(i)/sigmy(i),one),tsaiwu(i))
534
535
536 ENDDO
537
538 DO i=1,nel
539 coef(i)= zero
540 IF (wvec(i) > sigmy(i) .and .off(i) == one) THEN
541 coef(i)=one
542 IF(kd4(i)==0) THEN
543 kd4(i)=1
544 ENDIF
545 ENDIF
546 ENDDO
547
548 DO i=1,nel
549 dp1(i)=f1(i) + two*f11(i)*so1(i) + two*f12(i)*so2(i)
550 . + two*f13(i)*so3(i)
551 dp2(i)=f2(i) + two*f22(i)*so2(i) + two*f12(i)*so1(i)
552 . + two*f23(i)*so3(i)
553 dp3(i)=f3(i) + two*f33(i)*so3(i) + two*f13(i)*so1(i)
554 . + two*f23(i)*so2(i)
555 dp4(i)=two*f44(i)*so4(i)
556 dp5(i)=two*f55(i)*so5(i)
557 dp6(i)=two*f66(i)*so6(i)
558 ENDDO
559
560 3101 FORMAT(/' 205 WVEC SIGMY T1-T6 ',2e11.4/6e11.4)
561 3102 FORMAT(' F1 F2 F11 F22 F12 F23 F44 F55',3e11.4/5e11.4)
562
563 DO i=1,nel
564 ds1(i)=t1(i)-so1(i)
565 ds2(i)=t2(i)-so2(i)
566 ds3(i)=t3(i)-so3(i)
567 ds4(i)=t4(i)-so4(i)
568 ds5(i)=t5(i)-so5(i)
569 ds6(i)=t6(i)-so6(i)
570 ENDDO
571
572 DO i=1,nel
573 lamda(i)=(dp1(i)*ds1(i)+dp2(i)*ds2(i)+dp3(i)*ds3(i)
574 . +dp4(i)*ds4(i)+dp5(i)*ds5(i)+dp6(i)*ds6(i))*coef(i)
575 ENDDO
576
577 3103 FORMAT(' 207 LAMDA DS1-DS6 ',e11.4/6e11.4)
578
579 DO i=1,nel
580 cnn(i)=cn(i)-one
581 plas(i)= one
582 IF(pla(i) > zero) plas(i)=pla(i)**cnn(i)
583 ENDDO
584
585 DO i=1,nel
586 IF(lamda(i)/= zero) THEN
587 lamda(i)=lamda(i)*coef(i)/
588 . (dp1(i)*(d11(i)*dp1(i)+d12(i)*dp2(i)+d13(i)*dp3(i))+
589 . dp2(i)*(d12(i)*dp1(i)+d22(i)*dp2(i)+d23(i)*dp3(i))+
590 . dp3(i)*(d13(i)*dp1(i)+d23(i)*dp2(i)+d33(i)*dp3(i))+
591 . two*dp4(i)*g12(i)*dp4(i)+
592 . two*dp5(i)*g23(i)*dp5(i)+
593 . two*dp6(i)*g31(i)*dp6(i)+
594 . (so1(i)*dp1(i)+so2(i)*dp2(i)+so3(i)*dp3(i)+
595 . two*so4(i)*dp4(i)+two*so5(i)*dp5(i)+two*so6(i)*dp6(i))
596 . *cn(i)*cb(i)*plas(i) )
597
598 ENDIF
599 ENDDO
600
601 3104 FORMAT(' 208 LAMDA ',e11.4)
602
603 DO i=1,nel
604 dp1(i)=lamda(i)*dp1(i)
605 dp2(i)=lamda(i)*dp2(i)
606 dp3(i)=lamda(i)*dp3(i)
607 dp4(i)=lamda(i)*dp4(i)
608 dp5(i)=lamda(i)*dp5(i)
609 dp6(i)=lamda(i)*dp6(i)
610 ENDDO
611 3105 FORMAT(' 209 PLA DP1-DP6',e11.4/6e11.4)
612
613 DO i=1,nel
614 epe(i,1)=epe(i,1)-dp1(i)
615 epe(i,2)=epe(i,2)-dp2(i)
616 epe(i,3)=epe(i,3)-dp3(i)
617 ENDDO
618
619 DO i=1,nel
620 t1(i)=t1(i)-d11(i)*dp1(i)-d12(i)*dp2(i)-d13(i)*dp3(i)
621 t2(i)=t2(i)-d12(i)*dp1(i)-d22(i)*dp2(i)-d23(i)*dp3(i)
622 t3(i)=t3(i)-d13(i)*dp1(i)-d23(i)*dp2(i)-d33(i)*dp3(i)
623 t4(i)=t4(i)-g12(i)*dp4(i)*two
624 t5(i)=t5(i)-g23(i)*dp5(i)*two
625 t6(i)=t6(i)-g31(i)*dp6(i)*two
626 ENDDO
627
628 mx = mat(1)
629 wplaref = pm(68 ,mx)
630 DO i=1,nel
631 dwpla = half*
632 . (dp1(i)*(t1(i)+so1(i))+
633 . dp2(i)*(t2(i)+so2(i))+
634 . dp3(i)*(t3(i)+so3(i))+
635 . two*dp4(i)*(t4(i)+so4(i))+
636 . two*dp5(i)*(t5(i)+so5(i))+
637 . two*dp6(i)*(t6(i)+so6(i)))
638 dwpla =
max(dwpla ,zero) / wplaref
639 pla(i) = pla(i) + dwpla
640 pla(i)=
max(pla(i),zero)
641 ENDDO
642
643
644
645
646
647
648 CALL m14ftg(sig,ax ,ay ,az ,bx ,by ,
649 2 bz ,cx ,cy ,cz ,t1 ,t2 ,
650 3 t3 ,t4 ,t5 ,t6 ,nel)
651
652 DO i=1,nel
653 sig(i,1)=sig(i,1)*off(i)
654 sig(i,2)=sig(i,2)*off(i)
655 sig(i,3)=sig(i,3)*off(i)
656 sig(i,4)=sig(i,4)*off(i)
657 sig(i,5)=sig(i,5)*off(i)
658 sig(i,6)=sig(i,6)*off(i)
659 ENDDO
660
661
662
663 dt5=half*dt1
664 DO i=1,nel
665 eint(i)=eint(i)+dt5*vnew(i)*
666 . ( d1(i)*(sold1(i)+sig(i,1))
667 . + d2(i)*(sold2(i)+sig(i,2))
668 . + d3(i)*(sold3(i)+sig(i,3))
669 . + d4(i)*(sold4(i)+sig(i,4))
670 . + d5(i)*(sold5(i)+sig(i,5))
671 . + d6(i)*(sold6(i)+sig(i,6)))
672 eint(i)=eint(i)/vol(i)
673 dam(i,5)=kd1(i)*1000 + kd2(i)*100 + kd3(i)*10 + kd4(i) + 10000
674 ENDDO
675
676 DO i=1,nel
677 sigym =
max(em20,half*(one/f11(i)+one/f22(i)))
678 sigy(i)=sigmy(i)*sqrt(sigym)
679 defp(i)=pla(i)
680 ENDDO
681
682 1000 FORMAT(1x,'RUPTURE OF SOLID ELEMENT NUMBER ',i10)
683 6001 FORMAT(i6,' FAIL(14) E(',i5,') M(',a3,') T(',6e10.3,')')
684 RETURN
subroutine m14ama(pm, a, rx, ry, rz, sx, sy, sz, ax, ay, az, bx, by, bz, cx, cy, cz, nel, jcvt, jsph)
subroutine m14ftg(sig, ax, ay, az, bx, by, bz, cx, cy, cz, t1, t2, t3, t4, t5, t6, nel)
subroutine m14gtf(sig, ax, ay, az, bx, by, bz, cx, cy, cz, d1, d2, d3, d4, d5, d6, t1, t2, t3, t4, t5, t6, e1, e2, e3, e4, e5, e6, nel)
subroutine eosmain(iflag, nel, eostyp, pm, off, eint, rho, rho0, mu, mu2, espe, dvol, df, vnew, mat, psh, pnew, dpdm, dpde, theta, bufmat, sig, mu_bak, mlw, npf, tf, vareos, nvareos, mat_param, bfrac, nvartmp, vartmp)
1.15.0