68
69 USE redef3_mod
70 USE python_funct_mod
71
72
73
74#include "implicit_f.inc"
75#include "comlock.inc"
76
77
78
79#include "mvsiz_p.inc"
80
81
82
83#include "param_c.inc"
84#include "com04_c.inc"
85#include "com08_c.inc"
86#include "scr14_c.inc"
87#include "scr17_c.inc"
88#include "units_c.inc"
89#include "com01_c.inc"
90#include "impl1_c.inc"
91
92
93
94 type(python_), intent(inout) :: PYTHON
95 INTEGER, INTENT(IN) :: NFT
96 INTEGER NPF(*),IGEO(NPROPGI,*),NEL,NGL(*),PID(*),MID(*),NUVAR,
97 . IPM(NPROPMI,*)
98 INTEGER, INTENT(IN) :: STF
99 INTEGER, INTENT(IN) :: SANIN
100 INTEGER, INTENT(IN) :: IRESP
101 INTEGER, INTENT(IN) :: SNPC
102 INTEGER, INTENT(IN) :: SZYIELD_COMP
103 INTEGER, INTENT(IN) :: SZXXOLD_COMP
104
106 . skew(lskew,*), geo(npropg,*), fx(*), fy(*), fz(*), e(*), dx(*),
107 . dy(*), dz(*), tf(*), off(*), dpx(*), dpy(*), dpz(*), fxep(*),
108 . fyep(*), fzep(*), x0(*), y0(*), z0(*), xmom(*), ymom(*),
109 . zmom(*), rx(*), ry(*), rz(*), rpx(*), rpy(*), rpz(*), xmep(*),
110 . ymep(*), zmep(*), dpx2(*), dpy2(*), dpz2(*),rpx2(*), rpy2(*),
111 . rpz2(*),anim(*),fr_wave(*),e6(nel,6),
112 . exx2(mvsiz), eyx2(mvsiz), ezx2(mvsiz),
113 . exy2(mvsiz), eyy2(mvsiz), ezy2(mvsiz),
114 . exz2(mvsiz), eyz2(mvsiz), ezz2(mvsiz),
115 . crit_new(*), x0_err(mvsiz),yieldx(*),yieldy(*),
116 . yieldz(*),yieldx2(*),yieldy2(*),yieldz2(*),
117 . exx(mvsiz), eyx(mvsiz), ezx(mvsiz), exy(mvsiz),
118 . eyy(mvsiz), ezy(mvsiz), exz(mvsiz), eyz(mvsiz),
119 . ezz(mvsiz), xcr(mvsiz), rx1(mvsiz), rx2(mvsiz),
120 . ry1(mvsiz), ry2(mvsiz), rz1(mvsiz), rz2(mvsiz),
121 . xin(mvsiz),ak(mvsiz),xm(mvsiz),xkm(mvsiz),xcm(mvsiz),
122 . xkr(mvsiz),vx1(mvsiz),vx2(mvsiz),vy1(mvsiz),vy2(mvsiz),
123 . vz1(mvsiz),vz2(mvsiz),uvar(nuvar,*),uparam(*),mass(*),
124 . dx0(*),dy0(*),dz0(*),rx0(*),ry0(*),rz0(*)
125 my_real,
DIMENSION(6,NEL),
INTENT(INOUT) :: posx,posy,posz,posxx,posyy,poszz
126 my_real,
INTENT(INOUT) :: yieldxc(szyield_comp),yieldyc(szyield_comp),
127 . yieldzc(szyield_comp),yieldrxc(szyield_comp),yieldryc(szyield_comp),
128 . yieldrzc(szyield_comp)
129 my_real,
INTENT(INOUT) :: dxoldc(szxxold_comp),dyoldc(szxxold_comp),
130 . dzoldc(szxxold_comp),drxoldc(szxxold_comp),dryoldc(szxxold_comp),
131 . drzoldc(szxxold_comp)
132 TARGET :: uvar
133 DOUBLE PRECISION ALDP(MVSIZ),AL2DP(MVSIZ)
134
135
136
137 INTEGER INDX(MVSIZ),
138 . IECROU(MVSIZ), IFUNC(MVSIZ), IFV(MVSIZ), IFUNC2(MVSIZ),
139 . I, ILENG, J, KK, IFAIL(MVSIZ),IFAIL2(MVSIZ),
140 . NINDX,ISRATE, IFUNC3(MVSIZ),I1,I2,I3,I4,I5,I6,I7,I8,
141 . I9,I10,I11,I12,I13,I14,I15,IF1,IF2,IF3,IF4,IADBUF,NUPAR,
142 . IF5,IFUNC4(MVSIZ)
143
145 . xk(mvsiz), yk(mvsiz),
146 . zk(mvsiz),xc(mvsiz), yc(mvsiz), zc(mvsiz),xh(mvsiz),
147 . xhr(mvsiz),dxold(mvsiz), dyold(mvsiz), dzold(mvsiz),
148 . b(mvsiz), d(mvsiz), epla(mvsiz),
149 . dv(mvsiz),vrt(mvsiz),vrr(mvsiz),
150 . dmn(mvsiz),dmx(mvsiz),xl0(mvsiz),crit(mvsiz),
151 . xn(mvsiz),ff(mvsiz),lscale(mvsiz),ee(mvsiz),gf3(mvsiz),
152 . hx(mvsiz), hy(mvsiz), hz(mvsiz)
154 . at,dt05,xka,yka,zka,cc,cn,xa,xb,dlim,vfail,
155 . x21, y21, z21, epxy, epxz,
156 . vx21, vy21, vz21, ryavp, rzavp,eyzp,exzp,
157 . ryav, rzav,den, c, cp, exyp,
158 . x21phi, y21phi, z21phi, vx21phi, vy21phi, vz21phi,
159 . ryav1, rzav1, ryav1p, rzav1p,asrate,not_used(2)
161 DOUBLE PRECISION X0DP(MVSIZ)
162 my_real ,
DIMENSION(:),
POINTER :: coord_old
163 LOGICAL :: ANY_PYTHON_FUNCTION
164 TARGET :: not_used
165
166
167 not_used = zero
168
169 nupar = 4
170 i1 = nupar
171 i2 = i1 + 6
172 i3 = i2 + 6
173 i4 = i3 + 6
174 i5 = i4 + 6
175 i6 = i5 + 6
176 i7 = i6 + 6
177 i8 = i7 + 6
178 i9 = i8 + 6
179 i10 = i9 + 6
180 i11 = i10 + 6
181 i12 = i11 + 6
182 i13 = i12 + 6
183 i14 = i13 + 6
184 nupar = nupar + 84
185 i15 = i14 + 6 + 28
186 DO i=1,nel
187 iadbuf= ipm(7,mid(i)) - 1
188 epla(i)=zero
189 xm(i)=mass(i)
190 xk(i)=uparam(iadbuf + i11 + 1)
191 xc(i)=uparam(iadbuf + i12 + 1)
192 yk(i)=uparam(iadbuf + i11 + 2)
193 yc(i)=uparam(iadbuf + i12 + 2)
194 zk(i)=uparam(iadbuf + i11 + 3)
195 zc(i)=uparam(iadbuf + i12 + 3)
196 xka=xk(i)*uparam(iadbuf + i1 + 1)
197 yka=yk(i)*uparam(iadbuf + i1 + 2)
198 zka=zk(i)*uparam(iadbuf + i1 + 3)
199 xkm(i)=
max(xka,yka,zka)
200 hx(i) = uparam(iadbuf + i14 + 1)
201 hy(i) = uparam(iadbuf + i14 + 2)
202 hz(i) = uparam(iadbuf + i14 + 3)
203
204 xh(i)=
max(hx(i),hy(i),hz(i))
205 xcm(i)=
max(xc(i),yc(i),zc(i))
206 xcm(i)= xcm(i)+xh(i)
207
208 xkr(i)=
max(yka,zka) * aldp(i) * aldp(i)
209 xcr(i)= (
max(yc(i),zc(i)) +
max(hy(i),hz(i))) * aldp(i) * aldp(i)
210 vrt(i) = uparam(iadbuf + nupar + 1)
211 vrr(i) = uparam(iadbuf + nupar + 2)
212 ifail(i) = nint(uparam(iadbuf + 1 ))
213 ifail2(i)= nint(uparam(iadbuf + 3 ))
214 ENDDO
215
216 IF (inispri /= 0 .and. tt == zero) THEN
217 DO i=1,nel
218 xl0(i)= x0(i)
219
220 IF (xl0(i) == zero) xl0(i) = aldp(i)
221 ENDDO
222 ENDIF
223
224 IF (tt == zero) THEN
225 DO i=1,nel
226 x0(i)= aldp(i)
227 ENDDO
228 ENDIF
229
230 IF (scodver >= 101) THEN
231 IF (tt == zero) THEN
232 DO i=1,nel
233 x0_err(i)= aldp(i)-x0(i)
234 ENDDO
235 ENDIF
236 ENDIF
237
238 IF ( inispri /= 0 .and. tt == zero) THEN
239 DO i=1,nel
240 x0(i)= xl0(i)
241 ENDDO
242 ENDIF
243
244 DO i=1,nel
245 x0dp(i)= x0(i)
246 ENDDO
247
248 IF (scodver >= 101) THEN
249 DO i=1,nel
250 x0dp(i)= x0(i) + x0_err(i)
251 ENDDO
252 ENDIF
253
254
255
256 DO i=1,nel
257 dxold(i)=dx(i)
258 dyold(i)=dy(i)
259 dzold(i)=dz(i)
260 ENDDO
261
262 IF (inispri /= 0 .and. tt == zero) THEN
263 DO i=1,nel
264 dxold(i)=dx0(i)
265 dyold(i)=dy0(i)
266 dzold(i)=dz0(i)
267 ENDDO
268 ENDIF
269
270 dt05=half*dt1
271 IF (ismdisp > 0) THEN
272 DO i=1,nel
273 vx21 = vx2(i)-vx1(i)
274 vy21 = vy2(i)-vy1(i)
275 vz21 = vz2(i)-vz1(i)
276 dx(i) = dxold(i)+(vx21*exx(i)+vy21*eyx(i)+vz21*ezx(i))*dt1
277 dy(i) = dyold(i)+(vx21*exy(i)+vy21*eyy(i)+vz21*ezy(i))*dt1
278 dz(i) = dzold(i)+(vx21*exz(i)+vy21*eyz(i)+vz21*ezz(i))*dt1
279
280 x21 = (rx2(i)+rx1(i))
281 y21 = (ry2(i)+ry1(i))
282 z21 = (rz2(i)+rz1(i))
283
284 ryav1 = (x21*exy2(i)+y21*eyy2(i)+z21*ezy2(i))
285 rzav1 = (x21*exz2(i)+y21*eyz2(i)+z21*ezz2(i))
286
287 ryav = dt05 * ryav1
288 rzav = dt05 * rzav1
289
290 dy(i) = dy(i) - rzav * al2dp(i)
291 dz(i) = dz(i) + ryav * al2dp(i)
292
293 crit(i) = zero
294 ENDDO
295 ELSE
296 DO i=1,nel
297 vx21 = vx2(i)-vx1(i)
298 vy21 = vy2(i)-vy1(i)
299 vz21 = vz2(i)-vz1(i)
300
301 epxy = (vx21*exy2(i)+vy21*eyy2(i)+vz21*ezy2(i))*dt05
302 epxz = (vx21*exz2(i)+vy21*eyz2(i)+vz21*ezz2(i))*dt05
303
304 x21 = (rx2(i)+rx1(i))
305 y21 = (ry2(i)+ry1(i))
306 z21 = (rz2(i)+rz1(i))
307
308 ryav1 = (x21*exy2(i)+y21*eyy2(i)+z21*ezy2(i))
309 rzav1 = (x21*exz2(i)+y21*eyz2(i)+z21*ezz2(i))
310
311 at=epxz/
max(al2dp(i),em30)
312 at=atan(at)
313 ryav = dt05 * (ryav1) + two * at
314 at=epxy/
max(al2dp(i),em30)
315 at=atan(at)
316 rzav = dt05 * (rzav1) - two * at
317
318 dx(i) = aldp(i) - x0dp(i)
319 dy(i) = dyold(i) - rzav * al2dp(i)
320 dz(i) = dzold(i) + ryav * al2dp(i)
321
322 crit(i) = zero
323 ENDDO
324 ENDIF
325
326 DO i=1,nel
327 iadbuf = ipm(7,mid(i)) - 1
328 ileng = nint(uparam(iadbuf + 2))
329 IF (ileng /= 0) THEN
330 xl0(i)=x0dp(i)
331 ELSE
332 xl0(i)=one
333 ENDIF
334 ENDDO
335
336 nindx = 0
337 if1 = 0
338 if2 = 6
339 if3 = 12
340 if4 = 18
341 if5 = 24
342 DO i=1,nel
343 iadbuf = ipm(7,mid(i)) - 1
344 ifunc(i) = ipm(10 + if1 + 1,mid(i))
345 ifv(i) = ipm(10 + if2 + 1,mid(i))
346 ifunc2(i)= ipm(10 + if3 + 1,mid(i))
347 ifunc3(i)= ipm(10 + if4 + 1,mid(i))
348 ifunc4(i)= ipm(10 + if5 + 1,mid(i))
349 iecrou(i)= nint(uparam(iadbuf + i13 + 1))
350 ak(i) = uparam(iadbuf + i1 + 1)
351 b(i) = uparam(iadbuf + i2 + 1)
352 d(i) = uparam(iadbuf + i3 + 1)
353 ee(i) = uparam(iadbuf + i4 + 1)
354 gf3(i) = uparam(iadbuf + i5 + 1)
355 ff(i) = uparam(iadbuf + i6 + 1)
356 lscale(i)= uparam(iadbuf + i7 + 1)
357 dmn(i) = uparam(iadbuf + i8 + 1)
358 dmx(i) = uparam(iadbuf + i9 + 1)
359 ENDDO
360 IF (nuvar >=1) THEN
361 coord_old => uvar(1,1:nel)
362 ELSE
363 coord_old => not_used
364 ENDIF
365 CALL redef3(python,
366 . fx, xk, dx, fxep,
367 . dxold, dpx, tf, npf,
368 . xc, off, e6(1,1), dpx2,
369 . anim, anim_fe(11),posx,
370 . xl0, dmn, dmx, dv,
371 . ff, lscale, ee, gf3,
372 . ifunc3, yieldx, aldp, ak,
373 . b, d, iecrou, ifunc,
374 . ifv, ifunc2, epla, coord_old,
375 . nel, nft, stf, sanin,
376 . dt1, iresp, impl_s, idyna,
377 . snpc, yieldc=yieldxc, xx_oldc=dxoldc,
378 . fx0=uparam(iadbuf+i15+1),ifunc4=ifunc4,pos6=posx(6,1))
379
380 DO i=1,nel
381 cc = uparam(iadbuf + nupar + 3)
382 cn = uparam(iadbuf + nupar + 9)
383 xa = uparam(iadbuf + nupar + 15)
384 xb = uparam(iadbuf + nupar + 21)
385 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i)/= zero) THEN
386 IF (ifail2(i) == 0) THEN
387 xa = one
388 xb = two
389 IF (dx(i) > zero) THEN
390 dlim = dx(i) / (xl0(i)*dmx(i))
391 ELSE
392 dlim = dx(i) / (xl0(i)*dmn(i))
393 ENDIF
394 ELSE
395 vfail = cc * (abs(dv(i)/vrt(i)))**cn
396 IF (ifail2(i) == 1) THEN
397 IF (dx(i) > zero) THEN
398 dlim = dx(i) / (xl0(i)*dmx(i) + vfail)
399 ELSE
400 dlim = dx(i) / (xl0(i)*dmn
401 ENDIF
402 ELSEIF (ifail2(i) == 2) THEN
403 IF (fx(i) > zero) THEN
404 dlim = fx(i) / (dmx(i) + vfail)
405 ELSE
406 dlim = fx(i) / (dmn(i) - vfail)
407 ENDIF
408 ELSEIF (ifail2(i) == 3) THEN
409 dlim =
max(zero,e6(i,1)) / (dmx(i) + vfail)
410 ENDIF
411 ENDIF
412 IF (ifail(i) == 0) THEN
413
414 crit(i) =
max(crit(i),xa*dlim)
415 ELSE
416
417 crit(i)= crit(i) + xa*dlim**xb
418 ENDIF
419 ENDIF
420 ENDDO
421
422 DO i=1,nel
423 iadbuf = ipm(7,mid(i)) - 1
424 ifunc(i) = ipm(10 + if1 + 2,mid(i))
425 ifv(i) = ipm(10 + if2 + 2,mid(i))
426 ifunc2(i)= ipm(10 + if3 + 2,mid(i))
427 ifunc3(i)= ipm(10 + if4 + 2,mid(i))
428 ifunc4(i)= ipm(10 + if5 + 2,mid(i))
429 iecrou(i)= nint(uparam(iadbuf + i13 + 2))
430 ak(i) = uparam(iadbuf + i1 + 2)
431 b(i) = uparam(iadbuf + i2 + 2)
432 d(i) = uparam(iadbuf + i3 + 2)
433 ee(i) = uparam(iadbuf + i4 + 2)
434 gf3(i) = uparam(iadbuf + i5 + 2)
435 ff(i) = uparam(iadbuf + i6 + 2)
436 lscale(i)= uparam(iadbuf + i7 + 2)
437 dmn(i) = uparam(iadbuf + i8 + 2)
438 dmx(i) = uparam(iadbuf + i9 + 2)
439 ENDDO
440
441 kk = 1 + numelr * anim_fe(11)
442 IF (nuvar >= 2) THEN
443 coord_old => uvar(2,1:nel)
444 ELSE
445 coord_old => not_used
446 ENDIF
447 CALL redef3(python,
448 . fy, yk, dy, fyep
449 . dyold, dpy, tf, npf,
450 . yc, off, e6(1,2), dpy2,
451 . anim(kk), anim_fe(12),posy,
452 . xl0, dmn, dmx, dv,
453 . ff, lscale, ee, gf3,
454 . ifunc3, yieldy, aldp, ak,
455 . b, d, iecrou, ifunc,
456 . ifv, ifunc2, epla, coord_old,
457 . nel, nft, stf, sanin,
458 . dt1, iresp, impl_s, idyna,
459 . snpc, yieldc=yieldyc, xx_oldc=dyoldc,
460 . fx0=uparam(iadbuf+i15+2),ifunc4=ifunc4, pos6=posy(6,1))
461
462 DO i=1,nel
463 iadbuf = ipm(7,mid(i)) - 1
464 cc = uparam(iadbuf + nupar + 4)
465 cn = uparam(iadbuf + nupar + 10)
466 xa = uparam(iadbuf + nupar + 16)
467 xb = uparam(iadbuf + nupar + 22)
468 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
469 IF (ifail2(i) == 0) THEN
470 xa = one
471 xb = two
472 IF (dy(i) > zero) THEN
473 dlim = dy(i) / (xl0(i)*dmx(i))
474 ELSE
475 dlim = dy(i) / (xl0(i)*dmn(i))
476 ENDIF
477 ELSE
478 vfail = cc * (abs(dv(i)/vrt(i)))**cn
479 IF (ifail2(i) == 1) THEN
480 IF (dy(i) > zero) THEN
481 dlim = dy(i) / (xl0(i)*dmx(i) + vfail)
482 ELSE
483 dlim = dy(i) / (xl0(i)*dmn(i) - vfail)
484 ENDIF
485 ELSEIF (ifail2(i) == 2) THEN
486 IF (fy(i) > zero) THEN
487 dlim = fy(i) / (dmx(i) + vfail)
488 ELSE
489 dlim = fy(i) / (dmn(i) - vfail)
490 ENDIF
491 ELSEIF (ifail2(i) == 3) THEN
492 dlim =
max(zero,e6(i,2)) / (dmx(i) + vfail)
493 ENDIF
494 ENDIF
495 IF (ifail(i) == 0) THEN
496
497 crit(i) =
max(crit(i),xa*dlim)
498 ELSE
499
500 crit(i)= crit(i) + xa*dlim**xb
501 ENDIF
502 ENDIF
503 ENDDO
504
505 DO i=1,nel
506 iadbuf = ipm(7,mid(i)) - 1
507 ifunc(i) = ipm(10 + if1 + 3,mid(i))
508 ifv(i) = ipm(10 + if2 + 3,mid(i))
509 ifunc2(i)= ipm(10 + if3 + 3,mid(i))
510 ifunc3(i)= ipm(10 + if4 + 3,mid(i))
511 ifunc4(i)= ipm(10 + if5 + 3,mid(i))
512 iecrou(i)= nint(uparam(iadbuf + i13 + 3))
513 ak(i) = uparam(iadbuf + i1 + 3)
514 b(i) = uparam(iadbuf + i2 + 3)
515 d(i) = uparam(iadbuf + i3 + 3)
516 ee(i) = uparam(iadbuf + i4 + 3)
517 gf3(i) = uparam(iadbuf + i5 + 3)
518 ff(i) = uparam(iadbuf + i6 + 3)
519 lscale(i)= uparam(iadbuf + i7 + 3)
520 dmn(i) = uparam(iadbuf + i8 + 3)
521 dmx(i) = uparam(iadbuf + i9 + 3)
522 ENDDO
523
524 kk = 1 + numelr * (anim_fe(11)+anim_fe(12))
525 IF (nuvar >= 3) THEN
526 coord_old => uvar(3,1:nel)
527 ELSE
528 coord_old => not_used
529 ENDIF
530 CALL redef3(python,
531 . fz, zk, dz, fzep,
532 . dzold, dpz, tf, npf,
533 . zc, off, e6(1,3), dpz2,
534 . anim(kk), anim_fe(13),posz,
535 . xl0, dmn, dmx, dv,
536 . ff, lscale, ee, gf3,
537 . ifunc3, yieldz, aldp, ak,
538 . b, d, iecrou, ifunc,
539 . ifv, ifunc2, epla, coord_old,
540 . nel, nft, stf, sanin,
541 . dt1, iresp, impl_s, idyna,
542 . snpc, yieldc=yieldzc, xx_oldc=dzoldc,
543 . fx0=uparam(iadbuf+i15+3),ifunc4=ifunc4, pos6=posz(6,1))
544
545 DO i=1,nel
546 iadbuf = ipm(7,mid(i)) - 1
547 cc = uparam(iadbuf + nupar + 5)
548 cn = uparam(iadbuf + nupar + 11)
549 xa = uparam(iadbuf + nupar + 17)
550 xb = uparam(iadbuf + nupar + 23)
551 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
552 IF (ifail2(i) == 0) THEN
553 xa = one
554 xb = two
555 IF (dz(i) > zero)THEN
556 dlim = dz(i) / (xl0(i)*dmx(i))
557 ELSE
558 dlim = dz(i) / (xl0(i)*dmn(i))
559 ENDIF
560 ELSE
561 vfail = cc * (abs(dv(i)/vrt(i)))**cn
562 IF (ifail2(i) == 1) THEN
563 IF (dz(i) > zero) THEN
564 dlim = dz(i) / (xl0(i)*dmx(i) + vfail)
565 ELSE
566 dlim = dz(i) / (xl0(i)*dmn(i) - vfail)
567 ENDIF
568 ELSEIF (ifail2(i) == 2) THEN
569 IF (fz(i) > zero) THEN
570 dlim = fz(i) / (dmx(i) + vfail)
571 ELSE
572 dlim = fz(i) / (dmn(i) - vfail)
573 ENDIF
574 ELSEIF (ifail2(i) == 3) THEN
575 dlim =
max(zero,e6(i,3)) / (dmx(i) + vfail)
576 ENDIF
577 ENDIF
578 IF (ifail(i) == 0) THEN
579
580 crit(i) =
max(crit(i),xa*dlim)
581 ELSE
582
583 crit(i)= crit(i) + xa*dlim**xb
584 ENDIF
585 ENDIF
586 ENDDO
587
588
589
590 DO i=1,nel
591 iadbuf= ipm(7,mid(i)) - 1
592 xin(i)= geo(2,pid(i))
593 xk(i) = uparam(iadbuf + i11 + 4)
594 xc(i) = uparam(iadbuf + i12 + 4)
595 yk(i) = uparam(iadbuf + i11 + 5)
596 yc(i) = uparam(iadbuf + i12 + 5)
597 zk(i) = uparam(iadbuf + i11 + 6)
598 zc(i) = uparam(iadbuf + i12 + 6)
599 hx(i) = uparam(iadbuf + i14 + 4)
600 hy(i) = uparam(iadbuf + i14 + 5)
601 hz(i) = uparam(iadbuf + i14 + 6)
602
603 xhr(i)=
max(hx(i),hy(i),hz(i))
604
605 xkr(i)=
max(xk(i)*uparam(iadbuf + i1 + 4),
606 . yk(i)*uparam(iadbuf + i1 + 5),
607 . zk(i)*uparam(iadbuf + i1 + 6))+xkr(i)
608 xcr(i)=
max(xc(i),yc(i),zc(i)) + xhr(i) +xcr(i)+xh(i)
609 ENDDO
610
611 DO i=1,nel
612 dxold(i)=rx(i)
613 dyold(i)=ry(i)
614 dzold(i)=rz(i)
615 ENDDO
616
617 IF ( inispri /= 0 .AND. tt == zero) THEN
618 DO i=1,nel
619 dxold(i)=rx0(i)
620 dyold(i)=ry0(i)
621 dzold(i)=rz0(i)
622 ENDDO
623 ENDIF
624
625 DO i=1,nel
626 x21 = (rx2(i)-rx1(i))*dt1
627 y21 = (ry2(i)-ry1(i))*dt1
628 z21 = (rz2(i)-rz1(i))*dt1
629 rx(i) = dxold(i) + x21*exx2(i)+y21*eyx2(i)+z21*ezx2(i)
630 ry(i) = dyold(i) + x21*exy2(i)+y21*eyy2(i)+z21*ezy2(i)
631 rz(i) = dzold(i) + x21*exz2(i)+y21*eyz2(i)+z21*ezz2(i)
632 ENDDO
633
634 DO i=1,nel
635 iadbuf = ipm(7,mid(i)) - 1
636 ifunc(i) = ipm(10 + if1 + 4,mid(i))
637 ifv(i) = ipm(10 + if2 + 4,mid(i))
638 ifunc2(i)= ipm(10 + if3 + 4,mid(i))
639 ifunc3(i)= ipm(10 + if4 + 4,mid(i))
640 ifunc4(i)= ipm(10 + if5 + 4,mid(i))
641 iecrou(i)= nint(uparam(iadbuf + i13 + 4))
642 ak(i) = uparam(iadbuf + i1 + 4)
643 b(i) = uparam(iadbuf + i2 + 4)
644 d(i) = uparam(iadbuf + i3 + 4)
645 ee(i) = uparam(iadbuf + i4 + 4)
646 gf3(i) = uparam(iadbuf + i5 + 4)
647 ff(i) = uparam(iadbuf + i6 + 4)
648 lscale(i)= uparam(iadbuf + i7 + 4)
649 dmn(i) = uparam(iadbuf + i8 + 4)
650 dmx(i) = uparam(iadbuf + i9 + 4)
651 ENDDO
652 IF (nuvar >= 4) THEN
653 coord_old => uvar(4,1:nel)
654 ELSE
655 coord_old => not_used
656 ENDIF
657 CALL redef3(python,
658 . xmom, xk, rx, xmep,
659 . dxold, rpx, tf, npf,
660 . xc, off, e6(1,4), rpx2,
661 . anim, 0, posxx,
662 . xl0, dmn, dmx, dv,
663 . ff, lscale, ee, gf3,
664 . ifunc3, yieldx2, aldp, ak,
665 . b, d, iecrou, ifunc,
666 . ifv, ifunc2, epla, coord_old,
667 . nel, nft, stf, sanin,
668 . dt1, iresp, impl_s, idyna,
669 . snpc, yieldc=yieldrxc, xx_oldc=drxoldc,
670 . fx0=uparam(iadbuf+i15+4),ifunc4=ifunc4, pos6=posxx(6,1))
671
672 DO i=1,nel
673 iadbuf= ipm(7,mid(i)) - 1
674 cc = uparam(iadbuf + nupar + 6)
675 cn = uparam(iadbuf + nupar + 12)
676 xa = uparam(iadbuf + nupar + 18)
677 xb = uparam(iadbuf + nupar + 24)
678 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
679 IF (ifail2(i) == 0) THEN
680 xa = one
681 xb = two
682 IF (rx(i) > zero) THEN
683 dlim = rx(i) / (xl0(i)*dmx(i))
684 ELSE
685 dlim = rx(i) / (xl0(i)*dmn(i))
686 ENDIF
687 ELSE
688 vfail = cc * (abs(dv(i)/vrr(i)))**cn
689 IF (ifail2(i) == 1) THEN
690 IF (rx(i) > zero) THEN
691 dlim = rx(i) / (xl0(i)*dmx(i) + vfail)
692 ELSE
693 dlim = rx(i) / (xl0(i)*dmn(i) - vfail)
694 ENDIF
695 ELSEIF (ifail2(i) == 2) THEN
696 IF(xmom(i)>0.)THEN
697 dlim = xmom(i)/(dmx(i) + vfail)
698 ELSE
699 dlim = xmom(i)/(dmn(i) - vfail)
700 ENDIF
701 ELSEIF (ifail2(i) == 3) THEN
702 dlim =
max(zero,e6(i,4)) / (dmx(i) + vfail)
703 ENDIF
704 ENDIF
705 IF (ifail(i) == 0) THEN
706
707 crit(i) =
max(crit(i),xa*dlim)
708 ELSE
709
710 crit(i)= crit(i) + xa*dlim**xb
711 ENDIF
712 ENDIF
713 ENDDO
714
715 DO i=1,nel
716 iadbuf = ipm(7,mid(i)) - 1
717 ifunc(i) = ipm(10 + if1 + 5,mid(i))
718 ifv(i) = ipm(10 + if2 + 5,mid(i))
719 ifunc2(i)= ipm(10 + if3 + 5,mid(i))
720 ifunc3(i)= ipm(10 + if4 + 5,mid(i))
721 ifunc4(i)= ipm(10 + if5 + 5,mid(i))
722 iecrou(i)= nint(uparam(iadbuf + i13 + 5))
723 ak(i) = uparam(iadbuf + i1 + 5)
724 b(i) = uparam(iadbuf + i2 + 5)
725 d(i) = uparam(iadbuf + i3 + 5)
726 ee(i) = uparam(iadbuf + i4 + 5)
727 gf3(i) = uparam(iadbuf + i5 + 5)
728 ff(i) = uparam(iadbuf + i6 + 5)
729 lscale(i)= uparam(iadbuf + i7 + 5)
730 dmn(i) = uparam(iadbuf + i8 + 5)
731 dmx(i) = uparam(iadbuf + i9 + 5)
732 ENDDO
733 IF (nuvar >= 5) THEN
734 coord_old => uvar(5,1:nel)
735 ELSE
736 coord_old => not_used
737 ENDIF
738 CALL redef3(python,
739 . ymom, yk, ry, ymep,
740 . dyold, rpy, tf, npf,
741 . yc, off, e6(1,5), rpy2,
742 . anim, 0, posyy,
743 . xl0, dmn, dmx, dv,
744 . ff, lscale, ee, gf3,
745 . ifunc3, yieldy2, aldp, ak,
746 . b, d, iecrou, ifunc,
747 . ifv, ifunc2, epla, coord_old,
748 . nel, nft, stf, sanin,
749 . dt1, iresp, impl_s, idyna,
750 . snpc, yieldc=yieldryc, xx_oldc=dryoldc,
751 . fx0=uparam(iadbuf+i15+5),ifunc4=ifunc4, pos6=posyy(6,1))
752
753 DO i=1,nel
754 iadbuf= ipm(7,mid(i)) - 1
755 cc = uparam(iadbuf + nupar + 7)
756 cn = uparam(iadbuf + nupar + 13)
757 xa = uparam(iadbuf + nupar + 19)
758 xb = uparam(iadbuf + nupar + 25)
759 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
760 IF (ifail2(i) == 0) THEN
761 xa = one
762 xb = two
763 IF (ry(i) > zero) THEN
764 dlim = ry(i) / (xl0(i)*dmx(i))
765 ELSE
766 dlim = ry(i) / (xl0(i)*dmn(i))
767 ENDIF
768 ELSE
769 vfail = cc * (abs(dv(i)/vrr(i)))**cn
770 IF (ifail2(i) == 1) THEN
771 IF (ry(i) > zero) THEN
772 dlim = ry(i) / (xl0(i)*dmx(i) + vfail)
773 ELSE
774 dlim = ry(i) / (xl0(i)*dmn(i) - vfail)
775 ENDIF
776 ELSEIF (ifail2(i) == 2) THEN
777 IF (ymom(i) > zero)THEN
778 dlim = ymom(i)/(dmx(i) + vfail)
779 ELSE
780 dlim = ymom(i)/(dmn(i) - vfail)
781 ENDIF
782 ELSEIF (ifail2(i) == 3) THEN
783 dlim =
max(zero,e6(i,5)) / (dmx(i) + vfail)
784 ENDIF
785 ENDIF
786 IF (ifail(i) == 0) THEN
787
788 crit(i) =
max(crit(i),xa*dlim)
789 ELSE
790
791 crit(i)= crit(i) + xa*dlim**xb
792 ENDIF
793 ENDIF
794 ENDDO
795
796 DO i=1,nel
797 iadbuf = ipm(7,mid(i)) - 1
798 ifunc(i) = ipm(10 + if1 + 6,mid(i))
799 ifv(i) = ipm(10 + if2 + 6,mid(i))
800 ifunc2(i)= ipm(10 + if3 + 6,mid(i))
801 ifunc3(i)= ipm(10 + if4 + 6,mid(i))
802 ifunc4(i)= ipm(10 + if5 + 6,mid(i))
803 iecrou(i)= nint(uparam(iadbuf + i13 + 6))
804 ak(i) = uparam(iadbuf + i1 + 6)
805 b(i) = uparam(iadbuf + i2 + 6)
806 d(i) = uparam(iadbuf + i3 + 6)
807 ee(i) = uparam(iadbuf + i4 + 6)
808 gf3(i) = uparam(iadbuf + i5 + 6)
809 ff(i) = uparam(iadbuf + i6 + 6)
810 lscale(i)= uparam(iadbuf + i7 + 6)
811 dmn(i) = uparam(iadbuf + i8 + 6)
812 dmx(i) = uparam(iadbuf + i9 + 6)
813 ENDDO
814 IF (nuvar >= 6) THEN
815 coord_old => uvar(6,1:nel)
816 ELSE
817 coord_old => not_used
818 ENDIF
819 CALL redef3(python,
820 . zmom, zk, rz, zmep,
821 . dzold, rpz, tf, npf,
822 . zc, off, e6(1,6), rpz2,
823 . anim, 0, poszz,
824 . xl0, dmn, dmx, dv,
825 . ff, lscale, ee, gf3,
826 . ifunc3, yieldz2, aldp, ak,
827 . b, d, iecrou, ifunc,
828 . ifv, ifunc2, epla, coord_old,
829 . nel, nft, stf, sanin,
830 . dt1, iresp, impl_s, idyna,
831 . snpc, yieldc=yieldrzc, xx_oldc=drzoldc,
832 . fx0=uparam(iadbuf+i15+6),ifunc4=ifunc4,pos6=poszz(6,1))
833
834 DO i=1,nel
835 iadbuf= ipm(7,mid(i)) - 1
836 cc = uparam(iadbuf + nupar + 8)
837 cn = uparam(iadbuf + nupar + 14)
838 xa = uparam(iadbuf + nupar + 20)
839 xb = uparam(iadbuf + nupar + 26)
840 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
841 IF (ifail2(i) == 0) THEN
842 xa = one
843 xb = two
844 IF (rz(i) > zero) THEN
845 dlim = rz(i) / (xl0(i)*dmx(i))
846 ELSE
847 dlim = rz(i) / (xl0(i)*dmn(i))
848 ENDIF
849 ELSE
850 vfail = cc * (abs(dv(i)/vrr(i)))**cn
851 IF (ifail2(i) == 1) THEN
852 IF (rz(i) > zero)THEN
853 dlim = rz(i) / (xl0(i)*dmx(i) + vfail)
854 ELSE
855 dlim = rz(i) / (xl0(i)*dmn(i) - vfail)
856 ENDIF
857 ELSEIF (ifail2(i) == 2) THEN
858 IF (zmom(i) > zero)THEN
859 dlim = zmom(i)/(dmx(i) + vfail)
860 ELSE
861 dlim = zmom(i)/(dmn(i) - vfail)
862 ENDIF
863 ELSEIF (ifail2(i) == 3) THEN
864 dlim =
max(zero,e6(i,6)) / (dmx(i) + vfail)
865 ENDIF
866 ENDIF
867 IF (ifail(i) == 0) THEN
868! uniaxial failure
869 crit(i) =
max(crit(i),xa*dlim)
870 ELSE
871
872 crit(i)= crit(i) + xa *dlim**xb
873 ENDIF
874 ENDIF
875 ENDDO
876
877 DO i=1,nel
878 e(i) = e6(i,1)+e6(i,2)+e6(i,3)+e6(i,4)+e6(i,5)+e6(i,6)
879 ENDDO
880
881
882
883 DO i=1,nel
884 iadbuf = ipm(7,mid(i)) - 1
885 israte = nint(uparam(iadbuf + nupar + 27))
886
887 asrate = (2*pi*uparam(iadbuf + nupar + 28)*dt1)/(one+2*pi*uparam(iadbuf + nupar
888 IF (israte /= 0) THEN
889 IF (crit_new(i) < one) THEN
890 crit(i) =
min(crit(i),one+em3)
891 crit(i) = asrate*crit(i) + (one - asrate)*crit_new(i)
892 crit_new(i) =
min(crit(i),one)
893 ELSE
894 crit_new(i) = one
895 ENDIF
896 ELSE
897 IF (crit_new(i) < one) THEN
898 crit_new(i) =
min(crit(i),one)
899 ELSE
900 crit_new(i) = one
901 ENDIF
902 ENDIF
903 IF (off(i) == one) THEN
904 IF (crit(i) >= one) THEN
905 off(i)=zero
906 nindx = nindx + 1
907 indx(nindx) = i
908 idel7nok = 1
909 ENDIF
910 ENDIF
911 ENDDO
912
913 DO j=1,nindx
914 i = indx(j)
915#include "lockon.inc"
916 WRITE(iout, 1000) ngl(i)
917 WRITE(istdo,1100) ngl(i),tt
918#include "lockoff.inc"
919 ENDDO
920
921
922
924 1 xk, rpx, tf, npf,
925 2 iecrou, ifunc, ifv, epla,
926 3 nel)
928 1 yk, rpy, tf, npf,
929 2 iecrou, ifunc, ifv, epla,
930 3 nel)
932 1 zk, rpz, tf, npf,
933 2 iecrou, ifunc, ifv, epla,
934 3 nel)
935
936 DO i=1,nel
937 iadbuf= ipm(7,mid(i)) - 1
938 xk(i)=uparam(iadbuf + i11 + 1)
939 yk(i)=uparam(iadbuf + i11 + 2)
940 zk(i)=uparam(iadbuf + i11 + 3)
941 ENDDO
942
944 1 xk, dpx, tf, npf,
945 2 iecrou, ifunc, ifv, epla,
946 3 nel)
948 1 yk, dpy, tf, npf,
949 2 iecrou, ifunc, ifv, epla,
950 3 nel)
952 1 zk, dpz, tf, npf,
953 2 iecrou, ifunc, ifv, epla,
954 3 nel)
955 DO i=1,nel
956 xkm(i)=xkm(i)/xl0(i)
957 xcm(i)=xcm(i)/xl0(i)
958 xin(i)=xin(i)*xl0(i)
959 xkr(i)=xkr(i)/xl0(i)
960 xcr(i)=xcr(i)/xl0(i)
961 ENDDO
962
963 1000 FORMAT(1x,'-- RUPTURE OF SPRING ELEMENT NUMBER ',i10)
964 1100 FORMAT(1x,'-- RUPTURE OF SPRING ELEMENT :',i10,' AT TIME :',g11.4)
965
966 RETURN
subroutine repla3(xk, dpx, tf, npf, iecrou, ifunc, ifv, epla, nel)