40
41
42
43#include "implicit_f.inc"
44
45
46
47#include "mvsiz_p.inc"
48
49
50
51#include "impl1_c.inc"
52
53
54
55 INTEGER JFT, JLT, IPLA,IRTY,ISRATE,NEL
56 INTEGER,INTENT(IN) :: VP
57
59 . ezz(*), off(*), pla(*),temp(*),z3,z4,
60 . etse(*),gs(*),epsp(*)
62 . a1, a2, g, ymax0,
63 . ca0, cb0, cn, young, yld(mvsiz),
64 . cc, epdr(mvsiz), nu, epchk(mvsiz),tstar(mvsiz),
65 . dpla(mvsiz),fisokin,gama_imp(*),signor(mvsiz,5),
66 . hardm(*),depsxx(mvsiz),depsyy(mvsiz),depsxy(mvsiz),depsyz(mvsiz),
67 . depszx(mvsiz),signxx(nel),signyy(nel),signxy(nel),
68 . signyz(nel),signzx(nel),sigbakxx(nel),sigbakyy(nel),sigbakxy(nel),
69 . sigoxx(nel),sigoyy(nel),
70 . sigoxy(nel),sigoyz(nel),sigozx(nel)
71
72
73
74 INTEGER ICC
75 INTEGER I, J, N, NINDX, INDEX(MVSIZ), NMAX,IKFLG
76
78 . a(mvsiz), b(mvsiz) , dpla_i(mvsiz), dpla_j(mvsiz), dr(mvsiz),
79 . h(mvsiz), nu1(mvsiz), nu2(mvsiz) , p(mvsiz) , q(mvsiz),
80 . svm(mvsiz),ca(mvsiz), cb(mvsiz) ,
ymax(mvsiz) ,
81 . svm2(mvsiz),yld2(mvsiz),hi(mvsiz) ,hk(mvsiz) ,logep(mvsiz),
82 . nu11(mvsiz),nu21(mvsiz),anu1(mvsiz),bnu2(mvsiz),h2(mvsiz),
83 . err,f,df,pla_i,p2,q2,r,s1,s2,s3,yld_i,nnu1,nnu2,
84 . nu3,nu4,nu5,nu6,sigz,pp,aa,bb,c,s11,s22,s12,s1s2,s122,umr,
85 . vm2,qq,small,mt,tm,beta,aaa ,hkin,
alpha,plap1
86
87 DATA nmax/3/
88
89 nindx = 0
90 small = em7
91 DO i=jft,jlt
92 etse(i) = one
93
94 ca(i) = ca0
95 cb(i) = cb0
97 h(i) = zero
98
99 signxx(i)=sigoxx(i)
100 signyy(i)=sigoyy(i)
101 signxy(i)=sigoxy(i)
102 signyz(i)=sigoyz(i)
103 signzx(i)=sigozx(i)
104 ENDDO
105
106
107
108 ikflg=0
109 IF (fisokin > 0) THEN
110 DO i=jft,jlt
111 signxx(i)=signxx(i)-sigbakxx(i)
112 signyy(i)=signyy(i)-sigbakyy(i)
113 signxy(i)=signxy(i)-sigbakxy(i)
114 ENDDO
115 ikflg=jlt
116 END IF
117
118
119
120 DO i=jft,jlt
121 signxx(i)=signxx(i)+a1*depsxx(i)+a2*depsyy(i)
122 signyy(i)=signyy(i)+a2*depsxx(i)+a1*depsyy(i)
123 signxy(i)=signxy(i)+g*depsxy(i)
124 signyz(i)=signyz(i)+gs(i)*depsyz(i)
125 signzx(i)=signzx(i)+gs(i)*depszx(i)
126 ENDDO
127 DO i=jft,jlt
128 logep(i) = zero
129 ENDDO
130
131
132
133
134 IF (cc /= zero) THEN
135 IF (irty == 0) THEN
136 DO i=jft,jlt
137 IF (israte == 0.AND.vp == 2) epsp(i) =
138 .
max( abs(depsxx(i)), abs(depsyy(i)), half*abs(depsxy(i)
139 epsp(i) =
max(epsp(i),epdr(i))
140 logep(i) = log(epsp(i)/epdr(i))
141 ENDDO
142#include "vectorize.inc"
143 DO i=jft,jlt
145 IF (tstar(i) == zero) THEN
146 q(i) = (one + cc * logep(i))
147 ELSE
148 q(i) = (one + cc * logep(i))*(one-exp(mt*log(tstar(i))))
149 ENDIF
150 q(i) =
max(q(i),em20)
151 ca(i) = ca(i) * q(i)
152 cb(i) = cb(i) * q(i)
153 IF (icc == 1)
ymax(i) =
ymax(i) * q(i)
154 ENDDO
155 ELSEIF (irty == 1) THEN
156 DO i=jft,jlt
157 IF (israte == 0.AND.vp == 2) epsp(i) =
max( abs(depsxx(i)),
158 . abs(depsyy(i)), half*abs(depsxy(i)))
159 epsp(i) =
max(epsp(i),em20)
160 logep(i) = log(epsp(i)/epdr(i))
161 ENDDO
162 DO i=jft,jlt
163 q(i) = logep(i)
164 q(i) = cc*exp((-z3+z4 * q(i))*temp(i))
165 IF (icc == 1)
ymax(i)=
ymax(i) + q(i)
166 ca(i) = ca(i) + q(i)
167 ENDDO
168 ENDIF
169 ELSE
170 IF (irty == 0) THEN
172 DO i=jft,jlt
173 IF (tstar(i) /= zero) THEN
174 q(i) = one - exp(mt*log(tstar(i)))
175 q(i) =
max(q(i),em20)
176 ca(i) = ca(i) * q(i)
177 cb(i) = cb(i) * q(i)
178 END IF
179 ENDDO
180 ENDIF
181 ENDIF
182
183
184
185 DO i=jft,jlt
186 dpla(i) = zero
187 IF (pla(i) == zero) THEN
188 yld(i)= ca(i)
189 ELSE
190 beta = cb(i)*(one-fisokin)
191 yld(i)= ca(i)+beta*exp(cn*log(pla(i)))
192 ENDIF
194 ENDDO
195
196
197
198 IF (ipla == 0) THEN
199
200
201
202 DO i=jft,jlt
203 svm(i)=sqrt(signxx(i)*signxx(i)
204 . +signyy(i)*signyy(i)
205 . -signxx(i)*signyy(i)
206 . +three*signxy(i)*signxy(i))
207 ENDDO
208
209#include "vectorize.inc"
210 DO i=jft,jlt
211 r =
min(one,yld(i)/(svm(i)+em15))
212 IF (r < one) THEN
213 signxx(i)=signxx(i)*r
214 signyy(i)=signyy(i)*r
215 signxy(i)=signxy(i)*r
216 dpla(i) = off(i) *
max(zero,(svm(i)-yld(i))/young)
217 s1=half*(signxx(i)+signyy(i))
218 ezz(i) = dpla(i) * s1 /yld(i)
219 pla(i) = pla(i) + dpla(i)
220 epchk(i) =
max(pla(i),epchk(i))
221 IF (yld(i) >=
ymax(i))
THEN
222 h(i)=zero
223 ELSE
224 h(i)=cn*cb(i)*exp((cn-one)*log(pla(i)+small))
225 ENDIF
226 etse(i)= h(i)/(h(i)+young)
227 ENDIF
228 ENDDO
229
230 ELSEIF (ipla == 1) THEN
231
232
233
234 DO i=jft,jlt
235 s1=signxx(i)+signyy(i)
236 s2=signxx(i)-signyy(i)
237 s3=signxy(i)
238 a(i)=fourth*s1*s1
239 b(i)=three_over_4*s2*s2+three*s3*s3
240 svm(i)=sqrt(a(i)+b(i))
241 ENDDO
242
243
244
245 nindx=0
246 DO i=jft,jlt
247 IF (svm(i) > yld(i) .AND. off(i) == one) THEN
248 nindx=nindx+1
249 index(nindx)=i
250 ENDIF
251 ENDDO
252 IF (nindx == 0) THEN
253 IF (ikflg > 0) THEN
254 DO i=jft,jlt
255 signxx(i)=signxx(i) + sigbakxx(i)
256 signyy(i)=signyy(i) + sigbakyy(i)
257 signxy(i)=signxy(i) + sigbakxy(i)
258 ENDDO
259 END IF
260 IF (impl_s > 0.AND.ikt > 0) THEN
261 DO i = jft,jlt
262 gama_imp(i) = zero
263 END DO
264 END IF
265 RETURN
266 END IF
267
268
269
270#include "vectorize.inc"
271 DO j=1,nindx
272 i=index(j)
273 nu1(i)=one/(one-nu)
274 nu2(i)=one/(one+nu)
275 IF (yld(i) >=
ymax(i))
THEN
276 h(i)=zero
277 ELSE
278 h(i)=cn*cb(i)*exp((cn-one)*log(pla(i)+small))
279 ENDIF
280 dpla_j(i)=(svm(i)-yld(i))/(three*g+h(i))
281 etse(i)= h(i)/(h(i)+young)
282 anu1(i) = a(i)*nu1(i)
283 bnu2(i) = three*b(i)*nu2(i)
284 h2(i)= two*h(i)
285 ENDDO
286
287 IF (ikflg == 0) THEN
288 DO n=1,nmax
289#include "vectorize.inc"
290 DO j=1,nindx
291 i=index(j)
292 dpla_i(i)=dpla_j(i)
293 pla_i =pla(i)+dpla_i(i)
294 dpla(i) = dpla_j(i)
295 IF (pla_i == zero) THEN
297 ELSE
298 yld_i =
min(
ymax(i),ca(i)+cb(i)*exp(cn*log(pla_i)))
299 ENDIF
300 dr(i) =half*young*dpla_i(i)/yld_i
301 p(i) =one/(one+dr(i)*nu1(i))
302 q(i) =one/(one+three*dr(i)*nu2(i))
303 p2 =p(i)*p(i)
304 q2 =q(i)*q(i)
305 f =a(i)*p2+b(i)*q2-yld_i*yld_i
306 df =-(anu1(i)*p2*p(i)+bnu2(i)*q2*q(i))
307 . *(young-dr(i)*h2(i))/yld_i
308 . -h2(i)*yld_i
309 IF (dpla_i(i) > zero) THEN
311 ELSE
312 dpla_j(i)=zero
313 ENDIF
314 ENDDO
315 ENDDO
316
317 ELSE
318#include "vectorize.inc"
319 DO j=1,nindx
320 i=index(j)
321 beta = h(i)*fisokin
322 hi(i) = h(i)-beta
323 hk(i) = two_third*beta
324 aaa = three*hk(i)/young
325 nu11(i) = nu1(i) + aaa
326 nu21(i) = three*nu2(i) + aaa
327 anu1(i) = a(i)*nu11(i)
328 bnu2(i) = b(i)*nu21(i)
329 h2(i)= two*hi(i)
330 ENDDO
331
332 DO n=1,nmax
333#include "vectorize.inc"
334 DO j=1,nindx
335 i=index(j)
336 dpla_i(i)=dpla_j(i)
337 pla_i =pla(i)+dpla_i(i)
338 dpla(i) = dpla_j(i)
339 beta = one-fisokin
340 IF (pla_i == zero) THEN
342 ELSE
343 yld_i =
min(
ymax(i),ca(i)+beta*cb(i)*exp(cn*log(pla_i)))
344 ENDIF
345 dr(i) =half*young*dpla_i(i)/yld_i
346 p(i) =one/(one+dr(i)*nu11(i))
347 q(i) =one/(one+dr(i)*nu21(i))
348 p2 =p(i)*p(i)
349 q2 =q(i)*q(i)
350 f =a(i)*p2+b(i)*q2-yld_i*yld_i
351 df =-(anu1(i)*p2*p(i)+bnu2(i)*q2*q(i))
352 . *(young-dr(i)*h2(i))/yld_i
353 . -h2(i)*yld_i
354 IF (dpla_i(i) > zero) THEN
355 dpla_j(i)=
max(zero,dpla_i(i)-f/df)
356 ELSE
357 dpla_j(i)=zero
358 ENDIF
359 ENDDO
360 ENDDO
361 ENDIF
362
363
364
365#include "vectorize.inc"
366 DO j=1,nindx
367 i=index(j)
368 pla(i) = pla(i) + dpla_i(i)
369 epchk(i) =
max(pla(i),epchk(i))
370 s1=(signxx(i)+signyy(i))*p(i)
371 s2=(signxx(i)-signyy(i))*q(i)
372 signxx(i)=half*(s1+s2)
373 signyy(i)=half*(s1-s2)
374 signxy(i)=signxy(i)*q(i)
375 ezz(i) = dr(i)*s1/young
376 ENDDO
377
378 ELSEIF (ipla == 2) THEN
379
380
381
382 DO i=jft,jlt
383 svm2(i)= signxx(i)*signxx(i)
384 . +signyy(i)*signyy(i)
385 . -signxx(i)*signyy(i)
386 . +three*signxy(i)*signxy(i)
387 svm(i)=sqrt(svm2(i))
388 END DO
389
390
391
392 nindx=0
393 DO i=jft,jlt
394 yld2(i)=yld(i)*yld(i)
395 IF (svm2(i) > yld2(i) .AND. off(i) == one) THEN
396 nindx=nindx+1
397 index(nindx)=i
398 END IF
399 END DO
400
401 IF (nindx /= 0) THEN
402
403
404
405#include "vectorize.inc"
406 DO j=1,nindx
407 i=index(j)
408 IF (yld(i) >=
ymax(i))
THEN
409 h(i)=zero
410 ELSE
411 h(i)=cn*cb(i)*exp((cn-one)*log(pla(i)+small))
412 ENDIF
413 etse(i)= h(i)/(h(i)+young)
414
415 aa=(svm2(i)-yld2(i))
416 . /(five*svm2(i)+three*(-signxx(i)*signyy(i)+signxy(i)*signxy(i)))
417 s1=(one-two*aa)*signxx(i)+ aa*signyy(i)
418 s2=aa*signxx(i)+(one-two*aa)*signyy(i)
419 s3=(one-three*aa)*signxy(i)
420 signxx(i)=s1
421 signyy(i)=s2
422 signxy(i)=s3
423 dpla(i) = off(i)*(svm(i)-yld(i))/(three*g+h(i))
424 pla(i) = pla(i) + dpla(i)
425
426 yld(i) =yld(i)+h(i)*dpla(i)
427 END DO
428
429#include "vectorize.inc"
430 DO j=1,nindx
431 i=index(j)
432 svm(i)=sqrt( signxx(i)*signxx(i)
433 . +signyy(i)*signyy(i)
434 . -signxx(i)*signyy(i)
435 . +three*signxy(i)*signxy(i))
436 r =
min(one,yld(i)/
max(em20,svm(i)))
437 signxx(i)=signxx(i)*r
438 signyy(i)=signyy(i)*r
439 signxy(i)=signxy(i)*r
440 ezz(i) = dpla(i) * half*(signxx(i)+signyy(i)) / yld(i)
441 END DO
442 END IF
443 ENDIF
444
445
446
447 IF (impl_s > 0) THEN
448 IF (ikt > 0) THEN
449 DO i = jft,jlt
450 IF (dpla(i) > zero) THEN
451
452 pla_i =pla(i)
453 beta = one-fisokin
454 yld(i) =
min(
ymax(i),ca(i)+beta*cb(i)*exp(cn*log(pla_i)))
455 gama_imp(i)= three_half*dpla(i)/yld(i)
456
457 signor(i,4)=fisokin*h(i)
458 signor(i,5)=(one-fisokin)*h(i)
459
460 signor(i,1)=third*(two*signxx(i)-signyy(i))
461 signor(i,2)=third*(two*signyy(i)-signxx(i))
462 signor(i,3)=two*signxy(i)
463 ELSE
464 gama_imp(i) = zero
465 END IF
466 END DO
467 END IF
468 END IF
469
470
471
472 IF (ikflg > 0) THEN
473 IF (ipla == 1 )THEN
474#include "vectorize.inc"
475 DO j=1,nindx
476 i=index(j)
477 pla_i =pla(i)
478 beta = one-fisokin
479
480 IF (pla_i == zero) THEN
481 yld(i) =ca(i)
482 ELSE
483 yld(i) =
min(
ymax(i),ca(i)+beta*cb(i)*exp(cn*log(pla_i)))
484 ENDIF
485 END DO
486 END IF
487 DO i=jft,jlt
488 hkin = fisokin*h(i)
490 sigbakxx(i) = sigbakxx(i) +
alpha*signxx(i)
491 sigbakyy(i) = sigbakyy(i) +
alpha*signyy(i)
492 sigbakxy(i) = sigbakxy(i) +
alpha*signxy(i)
493
494 signxx(i)=signxx(i) + sigbakxx(i)
495 signyy(i)=signyy(i) + sigbakyy(i)
496 signxy(i)=signxy(i) + sigbakxy(i)
497 ENDDO
498
499 END IF
500
501
502
503
504 DO i=jft,jlt
505 hardm(i) = h(i)
506 ENDDO
507
508 RETURN
subroutine ymax(idn, fac, npc, pld, stiffmin, stiffmax, stiffini, stiffavg)