29 SUBROUTINE cdleni(PM ,GEO ,STIFN ,STIFR,IXC,
30 . PX1 ,PX2 ,PY1 ,PY2 ,THK,
31 . IGEO,DT ,SH4TREE,ALDT,UPARAM ,
32 . IPM ,NLAY,PM_STACK,ISUBSTACK,STRC,
34 . X2L ,X3L ,X4L ,Y2L ,Y3L ,Y4L ,
35 . IGEO_STACK ,GROUP_PARAM)
43#include "implicit_f.inc"
52#include "remesh_c.inc"
53#include "vect01_c.inc"
58 INTEGER IXC(NIXC,*), IGEO(NPROPGI,*), SH4TREE(KSH4TREE,*),
59 . (NPROPMI,*),NLAY,ISUBSTACK,IGEO_STACK(4*NPT_STACK+2,*)
61 . PM(NPROPM,*), GEO(NPROPG,*),STIFN(*),STIFR(*),UPARAM(*),
62 . PX1(*),PX2(*),PY1(*),PY2(*),THK(*),DT(*),ALDT(*),PM_STACK(20,*),
63 .
area(mvsiz), strc(*),
64 . x2l(mvsiz),x3l(mvsiz),x4l(mvsiz),y2l(mvsiz),y3l(mvsiz),y4l(mvsiz)
65 TYPE (GROUP_PARAM_) :: GROUP_PARAM
69 INTEGER I,N, IMT, IPMAT, IGTYP,IPPID,IADB,
70 . I1,I3,IPTHK,IPPOS,I2,MATLY,IGMAT,IPGMAT,,NIP,MLAWLY
72 . ssp(mvsiz), al1(mvsiz),al(mvsiz), almin(mvsiz),
73 . al2(mvsiz), al3(mvsiz), al4(mvsiz), al5(mvsiz), al6(mvsiz)
75 . viscmx,a11,a11r,a12,b1,b2,c1,vv,sti,stir,viscdef,dtdyn,rho,
76 . young,nu,gmax,thkly,posly,fac,z0
77 my_real,
DIMENSION(MVSIZ) :: zoffset
79 igtyp = nint(geo(12,iprop))
80 igmat = igeo(98,iprop)
84 zoffset(lft:llt) = zero
94 zoffset(i) = z0 - half*thk(i)
96 ELSEIF (ipos== 3 .OR. ipos == 4)
THEN
103 zoffset(lft:llt) = zero
106 IF ((igtyp == 11 .AND. igmat < 0) .OR. igtyp == 16)
THEN
111 imt = igeo(ipmat+n,iprop)
112 ssp(i)=
max(ssp(i),pm(27,imt))
115 ELSEIF (mtn == 42)
THEN
118 imt = igeo(ipmat+n,iprop)
122 a11 = gmax*(one + nu)/(one - nu**2)
123 ssp(i)=
max(ssp(i), sqrt(a11/rho))
126 ELSEIF (mtn == 69)
THEN
129 imt = igeo(ipmat+n,iprop)
132 gmax = uparam(iadb+1)*uparam(iadb+6)
133 . + uparam(iadb+2)*uparam(iadb+7)
134 . + uparam(iadb+3)*uparam(iadb+8)
135 . + uparam(iadb+4)*uparam(iadb+9)
136 . + uparam(iadb+5)*uparam(iadb+10)
138 a11 = gmax*(one + nu)/(one - nu**2)
139 ssp(i)=
max(ssp(i), sqrt(a11/rho))
142 ELSEIF (mtn == 65)
THEN
145 imt = igeo(ipmat+n,iprop)
148 ssp(i)=
max(ssp(i), sqrt(young/rho))
154 imt = igeo(ipmat+n,iprop)
158 ssp(i)=
max(ssp(i), sqrt(young/(one-nu*nu)/rho))
162 ELSEIF(igtyp == 11 .AND. igmat > 0)
THEN
164 ssp(i) = geo(ipgmat +9 ,iprop)
166 ELSEIF(igtyp == 52 .OR.
167 . ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0))
THEN
169 ssp(i) = pm_stack(9 ,isubstack)
171 ELSEIF(igtyp == 17 .AND. igmat < 0)
THEN
178 imt = igeo_stack(ipmat + n,isubstack)
179 ssp(i)=
max(ssp(i),pm(27,imt))
182 ELSEIF (mtn == 42)
THEN
185 imt = igeo_stack(ipmat + n,isubstack)
189 a11 = gmax*(one + nu)/(one - nu**2)
190 ssp(i)=
max(ssp(i), sqrt(a11/rho))
193 ELSEIF (mtn == 69)
THEN
196 imt = igeo_stack(ipmat + n,isubstack)
199 gmax = uparam(iadb+1)*uparam(iadb+6)
200 . + uparam(iadb+2)*uparam(iadb+7)
201 . + uparam(iadb+3)*uparam(iadb+8)
202 . + uparam(iadb+4)*uparam(iadb+9)
203 . + uparam(iadb+5)*uparam(iadb+10)
205 a11 = gmax*(one + nu)/(one - nu**2)
206 ssp(i)=
max(ssp(i), sqrt(a11/rho))
209 ELSEIF (mtn == 69)
THEN
212 imt = igeo_stack(ipmat + n,isubstack)
215 gmax = uparam(iadb+1)*uparam(iadb+6)
216 . + uparam(iadb+2)*uparam(iadb+7)
217 . + uparam(iadb+3)*uparam(iadb+8)
218 . + uparam(iadb+4)*uparam(iadb+9)
219 . + uparam(iadb+5)*uparam(iadb+10)
221 a11 = gmax*(one + nu)/(one - nu**2)
222 ssp(i)=
max(ssp(i), sqrt(a11/rho))
228 imt = igeo_stack(ipmat + n,isubstack)
231 ssp(i)=
max(ssp(i), sqrt(young/rho))
237 imt = igeo_stack(ipmat + n,isubstack)
241 ssp(i)=
max(ssp(i), sqrt(young/(one-nu*nu)/rho))
245 ELSEIF(igtyp == 51 .AND. igmat < 0)
THEN
250 imt = igeo_stack(ipmat + n,isubstack)
251 mlawly = nint(pm(19,imt))
252 IF (mlawly <= 28)
THEN
253 ssp(i)=
max(ssp(i),pm(27,imt))
254 ELSEIF (mlawly == 42)
THEN
258 a11 = gmax*(one + nu)/(one - nu**2)
259 ssp(i)=
max(ssp(i), sqrt(a11/rho))
260 ELSEIF (mlawly == 69)
THEN
263 gmax = uparam(iadb+1)*uparam(iadb+6)
264 . + uparam(iadb+2)*uparam(iadb+7)
265 . + uparam(iadb+3)*uparam(iadb+8)
266 . + uparam(iadb+4)*uparam(iadb+9)
267 . + uparam(iadb+5)*uparam(iadb+10)
269 a11 = gmax*(one + nu)/(one - nu**2)
270 ssp(i)=
max(ssp(i), sqrt(a11/rho))
271 ELSEIF (mlawly == 65)
THEN
274 ssp(i)=
max(ssp(i), sqrt(young/rho))
279 ssp(i)=
max(ssp(i), sqrt(young/(one-nu*nu)/rho))
288 ELSEIF (mtn == 42)
THEN
293 a11 = gmax*(one + nu)/(one - nu**2)
294 ssp(i)=
max(ssp(i), sqrt(a11/rho))
296 ELSEIF (mtn == 69)
THEN
300 gmax = uparam(iadb+1)*uparam(iadb+6)
301 . + uparam(iadb+2)*uparam(iadb+7)
302 . + uparam(iadb+3)*uparam(iadb+8)
303 . + uparam(iadb+4)*uparam(iadb+9)
304 . + uparam(iadb+5)*uparam(iadb+10)
306 a11 = gmax*(one + nu)/(one - nu**2)
307 ssp(i)=
max(ssp(i), sqrt(a11/rho))
309 ELSEIF (mtn == 65)
THEN
313 ssp(i)=sqrt(young/rho)
320 ssp(i)=sqrt(young/(one-nu*nu)/rho)
325 al1(i)= x2l(i) * x2l(i) + y2l(i) * y2l(i)
326 al2(i)=(x3l(i)-x2l(i))*(x3l(i)-x2l(i))+(y3l(i)-y2l(i))*(y3l(i)-y2l(i))
327 al3(i)=(x4l(i)-x3l(i))*(x4l(i)-x3l(i))+(y4l(i)-y3l(i))*(y4l(i)-y3l(i))
328 al4(i)= x4l(i) * x4l(i) + y4l(i) * y4l(i)
329 al5(i)=(x4l(i)-x2l(i))*(x4l(i)-x2l(i))+(y4l(i)-y2l(i))*(y4l(i)-y2l(i))
330 al6(i)= x3l(i) * x3l(i) + y3l(i) * y3l(i)
334 al(i)=
min(al1(i),al2(i),al3(i),al4(i),al5(i),al6(i))
335 IF(al3(i) == zero) al(i)=
min(al1(i),al2(i),al4(i))
341 ELSEIF(mtn == 25.OR.mtn == 27)
THEN
347 viscmx = group_param%VISC_DM
348 IF (viscmx == zero) viscmx = viscdef
349 IF (mtn == 1 .OR.mtn == 2.OR.mtn == 3.OR.
350 . mtn == 22.OR.mtn == 23) viscmx=zero
351 viscmx = sqrt(one + viscmx*viscmx) - viscmx
353 dtdyn =
area(i)/sqrt(
max(al5(i),al6(i)))
354 aldt(i) =
max(dtdyn,almin(i))
355 dt(i) = aldt(i)*viscmx/ssp(i)
362 IF (igtyp == 11 .AND. igmat > 0)
THEN
364 a11 =geo(ipgmat + 5,iprop)
365 a11r =geo(ipgmat + 7,iprop)
366 b1 = px1(i)*px1(i)+py1(i)*py1(i)
367 b2 = px2(i)*px2(i)+py2(i)*py2(i)
369 fac =
max(b1,b2) / (
area(i) * vv)
370 sti = fac * thk(i) * a11
371 stir = fac*a11r * thk(i)*(thk(i)**2 +
area(i))*one_over_12
372 stifn(ixc(2,i))=stifn(ixc(2,i))+sti
373 stifn(ixc(3,i))=stifn(ixc(3,i))+sti
374 stifn(ixc(4,i))=stifn(ixc(4,i))+sti
375 stifn(ixc(5,i))=stifn(ixc(5,i))+sti
376 stifr(ixc(2,i))=stifr(ixc(2,i))+stir
377 stifr(ixc(3,i))=stifr(ixc(3,i))+stir
378 stifr(ixc(4,i))=stifr(ixc(4,i))+stir
379 stifr(ixc(5,i))=stifr(ixc(5,i))+stir
382 ELSEIF(igtyp == 52 .OR.
383 . ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0 ))
THEN
385 a11 = pm_stack(5 ,isubstack)
386 a11r = pm_stack(7 ,isubstack)
387 b1 = px1(i)*px1(i)+py1(i)*py1(i)
388 b2 = px2(i)*px2(i)+py2(i)*py2(i)
390 fac =
max(b1,b2) / (
area(i) * vv)
391 sti = fac * thk(i) * a11
392 stir = fac*a11r * thk(i)*((thk(i)**2 +
area(i))*one_over_12 +
393 . zoffset(i)*zoffset(i) )
394 stifn(ixc(2,i))=stifn(ixc(2,i))+sti
395 stifn(ixc(3,i))=stifn(ixc(3,i))+sti
397 stifn(ixc(5,i))=stifn(ixc(5,i))+sti
398 stifr(ixc(2,i))=stifr(ixc(2,i))+stir
399 stifr(ixc(3,i))=stifr(ixc(3,i))+stir
401 stifr(ixc(5,i))=stifr(ixc(5,i))+stir
407 b1 = px1(i)*px1(i)+py1(i)*py1(i)
408 b2 = px2(i)*px2(i)+py2(i)*py2(i)
411 . * thk(i) * a11 / (
area(i) * vv)
412 stir = sti * (thk(i)*thk(i) +
area(i)) / 12.
413 stifn(ixc(2,i))=stifn(ixc(2,i))+sti
414 stifn(ixc(3,i))=stifn(ixc(3,i))+sti
415 stifn(ixc(4,i))=stifn(ixc(4,i))+sti
416 stifn(ixc(5,i))=stifn(ixc(5,i))+sti
417 stifr(ixc(2,i))=stifr(ixc(2,i))+stir
418 stifr(ixc(3,i))=stifr(ixc(3,i))+stir
419 stifr(ixc(4,i))=stifr(ixc(4,i))+stir
420 stifr(ixc(5,i))=stifr(ixc(5,i))+stir
425 IF(igtyp == 11 .AND. igmat > 0)
THEN
428 IF(sh4tree(3,n) >= 0)
THEN
429 a11 =geo(ipgmat + 5,iprop)
430 a11r =geo(ipgmat + 7,iprop)
431 b1 = px1(i)*px1(i)+py1(i)*py1(i)
432 b2 = px2(i)*px2(i)+py2(i)*py2(i)
434 fac =
max(b1,b2) / (
area(i) * vv)
435 sti = fac * thk(i) * a11
436 stir = fac * a11r * thk(i)*(thk(i)**2 +
area(i))*one_over_12
437 stifn(ixc(2,i))=stifn(ixc(2,i))+sti
438 stifn(ixc(3,i))=stifn(ixc(3,i))+sti
439 stifn(ixc(4,i))=stifn(ixc(4,i))+sti
440 stifn(ixc(5,i))=stifn(ixc(5,i))+sti
441 stifr(ixc(2,i))=stifr(ixc(2,i))+stir
442 stifr(ixc(3,i))=stifr(ixc(3,i))+stir
443 stifr(ixc(4,i))=stifr(ixc(4,i))+stir
444 stifr(ixc(5,i))=stifr(ixc(5,i))+stir
448 ELSEIF(igtyp == 52 .OR.
449 . ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0 ))
THEN
452 IF(sh4tree(3,n) >= 0)
THEN
453 a11 = pm_stack(5 ,isubstack)
454 a11r = pm_stack(7 ,isubstack)
455 b1 = px1(i)*px1(i)+py1(i)*py1(i)
456 b2 = px2(i)*px2(i)+py2(i)*py2(i)
458 fac =
max(b1,b2) / (
area(i) * vv)
459 sti = fac * thk(i) * a11
460 stir = fac * a11r * thk(i)*((thk(i)**2 +
area(i))*one_over_12 +
461 . zoffset(i)*zoffset(i) )
462 stifn(ixc(2,i))=stifn(ixc(2,i))+sti
463 stifn(ixc(3,i))=stifn(ixc(3,i))+sti
464 stifn(ixc(4,i))=stifn(ixc(4,i))+sti
465 stifn(ixc(5,i))=stifn(ixc(5,i))+sti
466 stifr(ixc(2,i))=stifr(ixc(2,i))+stir
467 stifr(ixc(3,i))=stifr(ixc(3,i))+stir
468 stifr(ixc(4,i))=stifr(ixc(4,i))+stir
469 stifr(ixc(5,i))=stifr(ixc(5,i))+stir
477 IF(sh4tree(3,n) >= 0)
THEN
479 b1 = px1(i)*px1(i)+py1(i)*py1(i)
480 b2 = px2(i)*px2(i)+py2(i)*py2(i)
483 . * thk(i) * a11 / (
area(i) * vv)
484 stir = sti * (thk(i)*thk(i) +
area(i)) / 12.
485 stifn(ixc(2,i))=stifn(ixc(2,i))+sti
486 stifn(ixc(3,i))=stifn(ixc(3,i))+sti
487 stifn(ixc(4,i))=stifn(ixc(4,i))+sti
488 stifn(ixc(5,i))=stifn(ixc(5,i))+sti
489 stifr(ixc(2,i))=stifr(ixc(2,i))+stir
490 stifr(ixc(3,i))=stifr(ixc(3,i))+stir
491 stifr(ixc(4,i))=stifr(ixc(4,i))+stir
492 stifr(ixc(5,i))=stifr(ixc(5,i))+stir
501 IF(geo(5,iprop)/=zero)geo(5,iprop)=
min(geo(5,iprop),dt(i))