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