OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
r4def3.F
Go to the documentation of this file.
1Copyright> OpenRadioss
2Copyright> Copyright (C) 1986-2025 Altair Engineering Inc.
3Copyright>
4Copyright> This program is free software: you can redistribute it and/or modify
5Copyright> it under the terms of the GNU Affero General Public License as published by
6Copyright> the Free Software Foundation, either version 3 of the License, or
7Copyright> (at your option) any later version.
8Copyright>
9Copyright> This program is distributed in the hope that it will be useful,
10Copyright> but WITHOUT ANY WARRANTY; without even the implied warranty of
11Copyright> MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12Copyright> GNU Affero General Public License for more details.
13Copyright>
14Copyright> You should have received a copy of the GNU Affero General Public License
15Copyright> along with this program. If not, see <https://www.gnu.org/licenses/>.
16Copyright>
17Copyright>
18Copyright> Commercial Alternative: Altair Radioss Software
19Copyright>
20Copyright> As an alternative to this open-source version, Altair also offers Altair Radioss
21Copyright> software under a commercial license. Contact Altair to discuss further if the
22Copyright> commercial version may interest you: https://www.altair.com/radioss/.
23!||====================================================================
24!|| r4def3 ../engine/source/elements/spring/r4def3.f
25!||--- called by ------------------------------------------------------
26!|| rforc3 ../engine/source/elements/spring/rforc3.F
27!||--- calls -----------------------------------------------------
28!|| redef3 ../engine/source/elements/spring/redef3.F90
29!|| repla3 ../engine/source/elements/spring/repla3.f
30!||--- uses -----------------------------------------------------
31!|| python_funct_mod ../common_source/modules/python_mod.F90
32!|| redef3_mod ../engine/source/elements/spring/redef3.F90
33!||====================================================================
34 SUBROUTINE r4def3(PYTHON,
35 1 SKEW, GEO, FX, FY,
36 2 FZ, E, DX, DY,
37 3 DZ, NPF, TF, OFF,
38 4 DPX, DPY, DPZ, DPX2,
39 5 DPY2, DPZ2, FXEP, FYEP,
40 6 FZEP, X0, Y0, Z0,
41 7 XMOM, YMOM, ZMOM, RX,
42 8 RY, RZ, RPX, RPY,
43 9 RPZ, XMEP, YMEP, ZMEP,
44 A RPX2, RPY2, RPZ2, ANIM,
45 B POSX, POSY, POSZ, POSXX,
46 C POSYY, POSZZ, FR_WAVE, E6,
47 D NEL, EXX2, EYX2, EZX2,
48 E EXY2, EYY2, EZY2, EXZ2,
49 F EYZ2, EZZ2, AL2DP, IGEO,
50 G CRIT_NEW,X0_ERR, ALDP, YIELDX,
51 H YIELDY, YIELDZ, YIELDX2, YIELDY2,
52 I YIELDZ2, NGL, MGN, EXX,
53 J EYX, EZX, EXY, EYY,
54 K EZY, EXZ, EYZ, EZZ,
55 L XCR, RX1, RY1, RZ1,
56 M RX2, RY2, RZ2, XIN,
57 N AK, XM, XKM, XCM,
58 O XKR, VX1, VX2, VY1,
59 P VY2, VZ1, VZ2, NUVAR,
60 Q UVAR, DX0, DY0, DZ0,
61 R RX0, RY0, RZ0, FX0,
62 S FY0, FZ0, XMOM0, YMOM0,
63 T ZMOM0, NFT, STF, SANIN,
64 U IRESP, SNPC)
65 USE python_funct_mod
66 USE redef3_mod
67C-----------------------------------------------
68C I m p l i c i t T y p e s
69C-----------------------------------------------
70#include "implicit_f.inc"
71#include "comlock.inc"
72C-----------------------------------------------
73C G l o b a l P a r a m e t e r s
74C-----------------------------------------------
75#include "mvsiz_p.inc"
76C-----------------------------------------------
77C C o m m o n B l o c k s
78C-----------------------------------------------
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"
87C-----------------------------------------------
88C D u m m y A r g u m e n t s
89C-----------------------------------------------
90 type(python_) :: PYTHON
91 INTEGER, INTENT(IN) :: NFT
92 INTEGER, INTENT(IN) :: STF !< Size of TF
93 INTEGER, INTENT(IN) :: SANIN !< Size of ANIM
94 INTEGER, INTENT(IN) :: IRESP !< Single Precision flag
95 INTEGER, INTENT(IN) :: SNPC !< Size of NPF
96 INTEGER NPF(SNPC),IGEO(NPROPGI,*),NEL,NGL(*),MGN(*),NUVAR
97C REAL
98 my_real
99 . SKEW(LSKEW,*), GEO(NPROPG,*), FX(*), FY(*), FZ(*), E(*), DX(*),
100 . DY(*), DZ(*), TF(STF), OFF(*), DPX(*), DPY(*), DPZ(*), FXEP(*),
101 . FYEP(*), FZEP(*), X0(*), Y0(*), Z0(*), XMOM(*), YMOM(*),
102 . ZMOM(*), RX(*), RY(*), RZ(*), RPX(*), RPY(*), RPZ(*), XMEP(*),
103 . YMEP(*), ZMEP(*), DPX2(*), DPY2(*), DPZ2(*),RPX2(*), RPY2(*),
104 . RPZ2(*),ANIM(SANIN),POSX(*),POSY(*),POSZ(*),POSXX(*),
105 . POSYY(*),POSZZ(*),FR_WAVE(*),E6(NEL,6),
106 . EXX2(MVSIZ), EYX2(MVSIZ), EZX2(MVSIZ),
107 . EXY2(MVSIZ), EYY2(MVSIZ), EZY2(MVSIZ),
108 . EXZ2(MVSIZ), EYZ2(MVSIZ), EZZ2(MVSIZ),
109 . CRIT_NEW(*), X0_ERR(MVSIZ),YIELDX(*),YIELDY(*),
110 . YIELDZ(*),YIELDX2(*),YIELDY2(*),YIELDZ2(*),
111 . EXX(MVSIZ), EYX(MVSIZ), EZX(MVSIZ), EXY(MVSIZ),
112 . EYY(MVSIZ), EZY(MVSIZ), EXZ(MVSIZ), EYZ(MVSIZ),
113 . EZZ(MVSIZ), XCR(MVSIZ), RX1(MVSIZ), RX2(MVSIZ),
114 . RY1(MVSIZ), RY2(MVSIZ), RZ1(MVSIZ), RZ2(MVSIZ),
115 . XIN(MVSIZ),AK(MVSIZ),XM(MVSIZ),XKM(MVSIZ),XCM(MVSIZ),
116 . XKR(MVSIZ),VX1(MVSIZ),VX2(MVSIZ),VY1(MVSIZ),VY2(MVSIZ),
117 . VZ1(MVSIZ),VZ2(MVSIZ),UVAR(NUVAR,*),DX0(*),DY0(*),DZ0(*),
118 . RX0(*),RY0(*),RZ0(*),FX0(*),FY0(*),FZ0(*),XMOM0(*),YMOM0(*),ZMOM0(*)
119 DOUBLE PRECISION ALDP(MVSIZ),AL2DP(MVSIZ)
120 TARGET :: uvar
121C-----------------------------------------------
122C L o c a l V a r i a b l e s
123C-----------------------------------------------
124 INTEGER INDX(MVSIZ),
125 . iecrou(mvsiz), ifunc(mvsiz), ifv(mvsiz), ifunc2(mvsiz),
126 . i, ileng, j, kk, ifail(mvsiz),ifail2(mvsiz),
127 . nindx,israte, ifunc3(mvsiz)
128C REAL
129 my_real
130 . xk(mvsiz), yk(mvsiz),
131 . zk(mvsiz),xc(mvsiz), yc(mvsiz), zc(mvsiz),xh(mvsiz),
132 . xhr(mvsiz),dxold(mvsiz), dyold(mvsiz), dzold(mvsiz),
133 . b(mvsiz), d(mvsiz), epla(mvsiz),
134 . dv(mvsiz),vrt(mvsiz),vrr(mvsiz),
135 . dmn(mvsiz),dmx(mvsiz),xl0(mvsiz),crit(mvsiz),
136 . xn(mvsiz),ff(mvsiz),lscale(mvsiz),ee(mvsiz),gf3(mvsiz),
137 . hx(mvsiz), hy(mvsiz), hz(mvsiz)
138 my_real
139 . at,dt05,xka,yka,zka,cc,cn,xa,xb,dlim,vfail,
140 . x21, y21, z21, epxy, epxz,
141 . vx21, vy21, vz21, ryavp, rzavp,eyzp,exzp,
142 . ryav, rzav,den, c, cp, exyp,
143 . x21phi, y21phi, z21phi, vx21phi, vy21phi, vz21phi,
144 . ryav1, rzav1, ryav1p, rzav1p,asrate,not_used,not_used2(2)
145 DOUBLE PRECISION X0DP(MVSIZ)
146 my_real ,DIMENSION(:), POINTER :: XX_OLD
147 TARGET :: NOT_USED2
148C-----------------------------------------------
149C
150 not_used = zero
151 not_used2 = zero
152C
153 DO i=1,nel
154 epla(i)=zero
155 xm(i)=geo(1,mgn(i))
156 xk(i)=geo(3,mgn(i))
157 xc(i)=geo(4,mgn(i))
158 yk(i)=geo(10,mgn(i))
159 yc(i)=geo(11,mgn(i))
160 zk(i)=geo(15,mgn(i))
161 zc(i)=geo(16,mgn(i))
162 xka=xk(i)*geo(41,mgn(i))
163 yka=yk(i)*geo(45,mgn(i))
164 zka=zk(i)*geo(49,mgn(i))
165 xkm(i)= max(xka,yka,zka)
166 hx(i) = geo(141,mgn(i))
167 hy(i) = geo(142,mgn(i))
168 hz(i) = geo(143,mgn(i))
169 xh(i)= max(hx(i),hy(i),hz(i))
170 xcm(i)= max(xc(i),yc(i),zc(i))
171 xcm(i)= xcm(i)+xh(i)
172
173 xkr(i)= max(yka,zka) * aldp(i) * aldp(i)
174 xcr(i)= (max(yc(i),zc(i)) + max(hy(i),hz(i))) * aldp(i) * aldp(i)
175 vrt(i) = geo(101,mgn(i))
176 vrr(i) = geo(102,mgn(i))
177 ifail(i) = nint(geo(79,mgn(i)))
178 ifail2(i)= nint(geo(95,mgn(i)))
179 ENDDO
180C
181 IF (inispri /= 0 .and. tt == zero) THEN
182 DO i=1,nel
183 xl0(i)= x0(i)
184! if not initialized length
185 IF (xl0(i) == zero) xl0(i) = aldp(i)
186 ENDDO
187 ENDIF
188C
189 IF (tt == zero) THEN
190 DO i=1,nel
191 x0(i)= aldp(i) ! cast double to My_real
192 ENDDO
193 ENDIF
194C
195 IF (scodver >= 101) THEN
196 IF (tt == zero) THEN
197 DO i=1,nel
198 x0_err(i)= aldp(i)-x0(i) ! difference between double and my_real
199 ENDDO
200 ENDIF
201 ENDIF
202C
203 IF ( inispri /= 0 .and. tt == zero) THEN
204 DO i=1,nel
205 x0(i)= xl0(i)
206 ENDDO
207 ENDIF
208C
209 DO i=1,nel
210 x0dp(i)= x0(i) ! cast double to My_real
211 ENDDO
212C
213 IF (scodver >= 101) THEN
214 DO i=1,nel
215 x0dp(i)= x0(i) + x0_err(i) ! difference between double and my_real
216 ENDDO
217 ENDIF
218C---------------------
219C TRANSLATIONS
220C---------------------
221 DO i=1,nel
222 dxold(i)=dx(i)
223 dyold(i)=dy(i)
224 dzold(i)=dz(i)
225 ENDDO
226!
227 IF (inispri /= 0 .and. tt == zero) THEN
228 DO i=1,nel
229 dxold(i)=dx0(i)
230 dyold(i)=dy0(i)
231 dzold(i)=dz0(i)
232 ENDDO
233 ENDIF
234C
235 dt05=half*dt1
236 IF (ismdisp > 0) THEN
237 DO i=1,nel
238 vx21 = vx2(i)-vx1(i)
239 vy21 = vy2(i)-vy1(i)
240 vz21 = vz2(i)-vz1(i)
241 dx(i) = dxold(i)+(vx21*exx(i)+vy21*eyx(i)+vz21*ezx(i))*dt1
242 dy(i) = dyold(i)+(vx21*exy(i)+vy21*eyy(i)+vz21*ezy(i))*dt1
243 dz(i) = dzold(i)+(vx21*exz(i)+vy21*eyz(i)+vz21*ezz(i))*dt1
244C
245 x21 = (rx2(i)+rx1(i))
246 y21 = (ry2(i)+ry1(i))
247 z21 = (rz2(i)+rz1(i))
248C
249 ryav1 = (x21*exy2(i)+y21*eyy2(i)+z21*ezy2(i))
250 rzav1 = (x21*exz2(i)+y21*eyz2(i)+z21*ezz2(i))
251C
252 ryav = dt05 * ryav1
253 rzav = dt05 * rzav1
254C
255 dy(i) = dy(i) - rzav * al2dp(i)
256 dz(i) = dz(i) + ryav * al2dp(i)
257C
258 crit(i) = zero
259 ENDDO
260 ELSE
261 DO i=1,nel
262 vx21 = vx2(i)-vx1(i)
263 vy21 = vy2(i)-vy1(i)
264 vz21 = vz2(i)-vz1(i)
265C
266 epxy = (vx21*exy2(i)+vy21*eyy2(i)+vz21*ezy2(i))*dt05
267 epxz = (vx21*exz2(i)+vy21*eyz2(i)+vz21*ezz2(i))*dt05
268C
269 x21 = (rx2(i)+rx1(i))
270 y21 = (ry2(i)+ry1(i))
271 z21 = (rz2(i)+rz1(i))
272C
273 ryav1 = (x21*exy2(i)+y21*eyy2(i)+z21*ezy2(i))
274 rzav1 = (x21*exz2(i)+y21*eyz2(i)+z21*ezz2(i))
275C
276 at=epxz/max(al2dp(i),em30)
277 at=atan(at)
278 ryav = dt05 * (ryav1) + two * at
279 at=epxy/max(al2dp(i),em30)
280 at=atan(at)
281 rzav = dt05 * (rzav1) - two * at
282C
283 dx(i) = aldp(i) - x0dp(i)
284 dy(i) = dyold(i) - rzav * al2dp(i)
285 dz(i) = dzold(i) + ryav * al2dp(i)
286C
287 crit(i) = zero
288 ENDDO
289 ENDIF !(ISMDISP > 0) THEN
290C
291 DO i=1,nel
292 ileng=nint(geo(93,mgn(i)))
293 IF (ileng /= 0) THEN
294 xl0(i)=x0dp(i)
295 ELSE
296 xl0(i)=one
297 ENDIF
298 ENDDO
299C-------------------------------
300 nindx = 0
301 DO i=1,nel
302 ifunc(i) = igeo(101,mgn(i))
303 ifv(i) = igeo(102,mgn(i))
304 ifunc2(i)= igeo(103,mgn(i))
305 ifunc3(i)= igeo(119,mgn(i))
306 iecrou(i)= nint(geo(7,mgn(i)))
307 ak(i) = geo(41,mgn(i))
308 b(i) = geo(42,mgn(i))
309 d(i) = geo(43,mgn(i))
310 ee(i) = geo(40 ,mgn(i))
311 gf3(i) = geo(132 ,mgn(i))
312 ff(i) = geo(44,mgn(i))
313 lscale(i)= geo(39 ,mgn(i))
314 dmn(i) = geo(65,mgn(i))
315 dmx(i) = geo(66,mgn(i))
316 ENDDO
317 IF (nuvar > 0) THEN
318 xx_old => uvar(1,1:nel)
319 ELSE
320 xx_old => not_used2
321 ENDIF
322 CALL redef3(python,
323 1 fx, xk, dx, fxep,
324 2 dxold, dpx, tf, npf,
325 3 xc, off, e6(1,1), dpx2,
326 4 anim, anim_fe(11),posx,
327 5 xl0, dmn, dmx, dv,
328 6 ff, lscale, ee, gf3,
329 7 ifunc3, yieldx, aldp, ak,
330 8 b, d, iecrou, ifunc,
331 9 ifv, ifunc2, epla, xx_old,
332 a nel, nft, stf, sanin ,
333 b dt1, iresp, impl_s, idyna,
334 c snpc)
335C
336 dlim = zero
337 DO i=1,nel
338 cc = geo(103,mgn(i))
339 cn = geo(109,mgn(i))
340 xa = geo(115,mgn(i))
341 xb = geo(121,mgn(i))
342 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i)/= zero) THEN
343 IF (ifail2(i) == 0) THEN
344 xa = one
345 xb = two
346 IF (dx(i) > zero) THEN
347 dlim = dx(i) / (xl0(i)*dmx(i))
348 ELSE
349 dlim = dx(i) / (xl0(i)*dmn(i))
350 ENDIF
351 ELSE
352 vfail = cc * (abs(dv(i)/vrt(i)))**cn
353 IF (ifail2(i) == 1) THEN
354 IF (dx(i) > zero) THEN
355 dlim = dx(i) / (xl0(i)*dmx(i) + vfail)
356 ELSE
357 dlim = dx(i) / (xl0(i)*dmn(i) - vfail)
358 ENDIF
359 ELSEIF (ifail2(i) == 2) THEN
360 IF (fx(i) > zero) THEN
361 dlim = fx(i) / (dmx(i) + vfail)
362 ELSE
363 dlim = fx(i) / (dmn(i) - vfail)
364 ENDIF
365 ELSEIF (ifail2(i) == 3) THEN
366 dlim = max(zero,e6(i,1)) / (dmx(i) + vfail)
367 ENDIF
368 ENDIF
369 IF (ifail(i) == 0) THEN
370! Uniaxial failure
371 crit(i) = max(crit(i),xa*dlim)
372 ELSE
373! Multiaxial failure
374 crit(i)= crit(i) + xa*dlim**xb
375 ENDIF
376 ENDIF
377 ENDDO
378C
379 DO i=1,nel
380 ifunc(i) = igeo(104,mgn(i))
381 ifv(i) = igeo(105,mgn(i))
382 ifunc2(i)= igeo(106,mgn(i))
383 ifunc3(i)= igeo(120,mgn(i))
384 iecrou(i)=nint(geo(14,mgn(i)))
385 ak(i) =geo(45,mgn(i))
386 b(i) =geo(46,mgn(i))
387 d(i) =geo(47,mgn(i))
388 ee(i) =geo(180,mgn(i))
389 gf3(i) =geo(133,mgn(i))
390 ff(i) =geo(48,mgn(i))
391 lscale(i)= geo(174,mgn(i))
392 dmn(i) =geo(67,mgn(i))
393 dmx(i) =geo(68,mgn(i))
394 ENDDO
395C
396 kk = 1 + numelr * anim_fe(11)
397 IF (nuvar > 0) xx_old => uvar(2,1:nel)
398 CALL redef3(python,
399 1 fy, yk, dy, fyep,
400 2 dyold, dpy, tf, npf,
401 3 yc, off, e6(1,2), dpy2,
402 4 anim(kk), anim_fe(12),posy,
403 5 xl0, dmn, dmx, dv,
404 6 ff, lscale, ee, gf3,
405 7 ifunc3, yieldy, aldp, ak,
406 8 b, d, iecrou, ifunc,
407 9 ifv, ifunc2, epla, xx_old,
408 a nel, nft, stf, sanin,
409 b dt1, iresp, impl_s, idyna,
410 c snpc)
411
412 DO i=1,nel
413 cc = geo(104,mgn(i))
414 cn = geo(110,mgn(i))
415 xa = geo(116,mgn(i))
416 xb = geo(122,mgn(i))
417 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
418 IF (ifail2(i) == 0) THEN
419 xa = one
420 xb = two
421 IF (dy(i) > zero) THEN
422 dlim = dy(i) / (xl0(i)*dmx(i))
423 ELSE
424 dlim = dy(i) / (xl0(i)*dmn(i))
425 ENDIF
426 ELSE
427 vfail = cc * (abs(dv(i)/vrt(i)))**cn
428 IF (ifail2(i) == 1) THEN
429 IF (dy(i) > zero) THEN
430 dlim = dy(i) / (xl0(i)*dmx(i) + vfail)
431 ELSE
432 dlim = dy(i) / (xl0(i)*dmn(i) - vfail)
433 ENDIF
434 ELSEIF (ifail2(i) == 2) THEN
435 IF (fy(i) > zero) THEN
436 dlim = fy(i) / (dmx(i) + vfail)
437 ELSE
438 dlim = fy(i) / (dmn(i) - vfail)
439 ENDIF
440 ELSEIF (ifail2(i) == 3) THEN
441 dlim = max(zero,e6(i,2)) / (dmx(i) + vfail)
442 ENDIF
443 ENDIF
444 IF (ifail(i) == 0) THEN
445! Uniaxial failure
446 crit(i) = max(crit(i),xa*dlim)
447 ELSE
448! Multiaxial failure
449 crit(i)= crit(i) + xa*dlim**xb
450 ENDIF
451 ENDIF
452 ENDDO
453C
454 DO i=1,nel
455 ifunc(i) = igeo(107,mgn(i))
456 ifv(i) = igeo(108,mgn(i))
457 ifunc2(i)= igeo(109,mgn(i))
458 ifunc3(i)= igeo(121,mgn(i))
459 iecrou(i)=nint(geo(18,mgn(i)))
460 ak(i) =geo(49,mgn(i))
461 b(i) =geo(50,mgn(i))
462 d(i) =geo(51,mgn(i))
463 ee(i) =geo(181,mgn(i))
464 gf3(i) =geo(134,mgn(i))
465 ff(i) =geo(52,mgn(i))
466 lscale(i)=geo(175,mgn(i))
467 dmn(i) =geo(69,mgn(i))
468 dmx(i) =geo(77,mgn(i))
469 ENDDO
470C
471 kk = 1 + numelr * (anim_fe(11)+anim_fe(12))
472 IF (nuvar > 0) xx_old => uvar(3,1:nel)
473 CALL redef3(python,
474 1 fz, zk, dz, fzep,
475 2 dzold, dpz, tf, npf,
476 3 zc, off, e6(1,3), dpz2,
477 4 anim(kk), anim_fe(13),posz,
478 5 xl0, dmn, dmx, dv,
479 6 ff, lscale, ee, gf3,
480 7 ifunc3, yieldz, aldp, ak,
481 8 b, d, iecrou, ifunc,
482 9 ifv, ifunc2, epla, xx_old,
483 a nel, nft, stf, sanin ,
484 b dt1, iresp, impl_s, idyna,
485 c snpc)
486C
487 DO i=1,nel
488 cc = geo(105,mgn(i))
489 cn = geo(111,mgn(i))
490 xa = geo(117,mgn(i))
491 xb = geo(123,mgn(i))
492 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
493 IF (ifail2(i) == 0) THEN
494 xa = one
495 xb = two
496 IF (dz(i) > zero)THEN
497 dlim = dz(i) / (xl0(i)*dmx(i))
498 ELSE
499 dlim = dz(i) / (xl0(i)*dmn(i))
500 ENDIF
501 ELSE
502 vfail = cc * (abs(dv(i)/vrt(i)))**cn
503 IF (ifail2(i) == 1) THEN
504 IF (dz(i) > zero) THEN
505 dlim = dz(i) / (xl0(i)*dmx(i) + vfail)
506 ELSE
507 dlim = dz(i) / (xl0(i)*dmn(i) - vfail)
508 ENDIF
509 ELSEIF (ifail2(i) == 2) THEN
510 IF (fz(i) > zero) THEN
511 dlim = fz(i) / (dmx(i) + vfail)
512 ELSE
513 dlim = fz(i) / (dmn(i) - vfail)
514 ENDIF
515 ELSEIF (ifail2(i) == 3) THEN
516 dlim = max(zero,e6(i,3)) / (dmx(i) + vfail)
517 ENDIF
518 ENDIF
519 IF (ifail(i) == 0) THEN
520! Uniaxial failure
521 crit(i) = max(crit(i),xa*dlim)
522 ELSE
523! Multiaxial failure
524 crit(i)= crit(i) + xa *dlim**xb
525 ENDIF
526 ENDIF
527 ENDDO
528C---------------------
529C ROTATIONS
530C---------------------
531 DO i=1,nel
532 xin(i)=geo(9,mgn(i))
533 xk(i)=geo(19,mgn(i))
534 xc(i)=geo(20,mgn(i))
535 yk(i)=geo(23,mgn(i))
536 yc(i)=geo(24,mgn(i))
537 zk(i)=geo(27,mgn(i))
538 zc(i)=geo(28,mgn(i))
539 hx(i) = geo(144,mgn(i))
540 hy(i) = geo(145,mgn(i))
541 hz(i) = geo(146,mgn(i))
542
543 xhr(i)= max(hx(i),hy(i),hz(i))
544
545 xkr(i)= max(xk(i)*geo(53,mgn(i)),
546 . yk(i)*geo(57,mgn(i)),
547 . zk(i)*geo(61,mgn(i)))+xkr(i)
548 xcr(i)= max(xc(i),yc(i),zc(i)) + xhr(i) +xcr(i)+xh(i)
549 ENDDO
550C
551 DO i=1,nel
552 dxold(i)=rx(i)
553 dyold(i)=ry(i)
554 dzold(i)=rz(i)
555 ENDDO
556!
557 IF ( inispri /= 0 .AND. tt == zero) THEN
558 DO i=1,nel
559 dxold(i)=rx0(i)
560 dyold(i)=ry0(i)
561 dzold(i)=rz0(i)
562 ENDDO
563 ENDIF
564C
565 DO i=1,nel
566 x21 = (rx2(i)-rx1(i))*dt1
567 y21 = (ry2(i)-ry1(i))*dt1
568 z21 = (rz2(i)-rz1(i))*dt1
569 rx(i) = dxold(i) + x21*exx2(i)+y21*eyx2(i)+z21*ezx2(i)
570 ry(i) = dyold(i) + x21*exy2(i)+y21*eyy2(i)+z21*ezy2(i)
571 rz(i) = dzold(i) + x21*exz2(i)+y21*eyz2(i)+z21*ezz2(i)
572 ENDDO
573C-------------------------------
574 DO i=1,nel
575 ifunc(i) = igeo(110,mgn(i))
576 ifv(i) = igeo(111,mgn(i))
577 ifunc2(i)= igeo(112,mgn(i))
578 ifunc3(i)= igeo(122,mgn(i))
579 iecrou(i)=nint(geo(22,mgn(i)))
580 ak(i) =geo(53,mgn(i))
581 b(i) =geo(54,mgn(i))
582 d(i) =geo(55,mgn(i))
583 ee(i) =geo(182,mgn(i))
584 gf3(i) =geo(135,mgn(i))
585 ff(i)= geo(56,mgn(i))
586 lscale(i)= geo(176,mgn(i))
587 dmn(i) =geo(71,mgn(i))
588 dmx(i) =geo(72,mgn(i))
589 ENDDO
590 IF (nuvar > 0) xx_old => uvar(4,1:nel)
591 CALL redef3(python,
592 1 xmom, xk, rx, xmep,
593 2 dxold, rpx, tf, npf,
594 3 xc, off, e6(1,4), rpx2,
595 4 anim, 0, posxx,
596 5 xl0, dmn, dmx, dv,
597 6 ff, lscale, ee, gf3,
598 7 ifunc3, yieldx2, aldp, ak,
599 8 b, d, iecrou, ifunc,
600 9 ifv, ifunc2, epla, xx_old,
601 a nel, nft, stf, sanin,
602 b dt1, iresp, impl_s, idyna,
603 c snpc)
604C
605 DO i=1,nel
606 cc = geo(106,mgn(i))
607 cn = geo(112,mgn(i))
608 xa = geo(118,mgn(i))
609 xb = geo(124,mgn(i))
610 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
611 IF (ifail2(i) == 0) THEN
612 xa = one
613 xb = two
614 IF (rx(i) > zero) THEN
615 dlim = rx(i) / (xl0(i)*dmx(i))
616 ELSE
617 dlim = rx(i) / (xl0(i)*dmn(i))
618 ENDIF
619 ELSE
620 vfail = cc * (abs(dv(i)/vrr(i)))**cn
621 IF (ifail2(i) == 1) THEN
622 IF (rx(i) > zero) THEN
623 dlim = rx(i) / (xl0(i)*dmx(i) + vfail)
624 ELSE
625 dlim = rx(i) / (xl0(i)*dmn(i) - vfail)
626 ENDIF
627 ELSEIF (ifail2(i) == 2) THEN
628 IF(xmom(i)>0.)THEN
629 dlim = xmom(i)/(dmx(i) + vfail)
630 ELSE
631 dlim = xmom(i)/(dmn(i) - vfail)
632 ENDIF
633 ELSEIF (ifail2(i) == 3) THEN
634 dlim = max(zero,e6(i,4)) / (dmx(i) + vfail)
635 ENDIF
636 ENDIF
637 IF (ifail(i) == 0) THEN
638! Uniaxial failure
639 crit(i) = max(crit(i),xa*dlim)
640 ELSE
641! Multiaxial failure
642 crit(i)= crit(i) + xa *dlim**xb
643 ENDIF
644 ENDIF
645 ENDDO
646C
647 DO i=1,nel
648 ifunc(i) = igeo(113,mgn(i))
649 ifv(i) = igeo(114,mgn(i))
650 ifunc2(i)= igeo(115,mgn(i))
651 ifunc3(i)= igeo(123,mgn(i))
652 iecrou(i)=nint(geo(26,mgn(i)))
653 ak(i) =geo(57,mgn(i))
654 b(i) =geo(58,mgn(i))
655 d(i) =geo(59,mgn(i))
656 ee(i)= geo(183,mgn(i))
657 gf3(i)= geo(136,mgn(i))
658 ff(i)= geo(60,mgn(i))
659 lscale(i)= geo(177,mgn(i))
660 dmn(i) =geo(73,mgn(i))
661 dmx(i) =geo(74,mgn(i))
662 ENDDO
663 IF (nuvar > 0) xx_old => uvar(5,1:nel)
664 CALL redef3(python,
665 1 ymom, yk, ry, ymep,
666 2 dyold, rpy, tf, npf,
667 3 yc, off, e6(1,5), rpy2,
668 4 anim, 0, posyy,
669 5 xl0, dmn, dmx, dv,
670 6 ff, lscale, ee, gf3,
671 7 ifunc3, yieldy2, aldp, ak,
672 8 b, d, iecrou, ifunc,
673 9 ifv, ifunc2, epla, xx_old,
674 a nel, nft, stf, sanin,
675 b dt1, iresp, impl_s, idyna,
676 c snpc)
677C
678 DO i=1,nel
679 cc = geo(107,mgn(i))
680 cn = geo(113,mgn(i))
681 xa = geo(119,mgn(i))
682 xb = geo(125,mgn(i))
683 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
684 IF (ifail2(i) == 0) THEN
685 xa = one
686 xb = two
687 IF (ry(i) > zero) THEN
688 dlim = ry(i) / (xl0(i)*dmx(i))
689 ELSE
690 dlim = ry(i) / (xl0(i)*dmn(i))
691 ENDIF
692 ELSE
693 vfail = cc * (abs(dv(i)/vrr(i)))**cn
694 IF (ifail2(i) == 1) THEN
695 IF (ry(i) > zero) THEN
696 dlim = ry(i) / (xl0(i)*dmx(i) + vfail)
697 ELSE
698 dlim = ry(i) / (xl0(i)*dmn(i) - vfail)
699 ENDIF
700 ELSEIF (ifail2(i) == 2) THEN
701 IF (ymom(i) > zero)THEN
702 dlim = ymom(i)/(dmx(i) + vfail)
703 ELSE
704 dlim = ymom(i)/(dmn(i) - vfail)
705 ENDIF
706 ELSEIF (ifail2(i) == 3) THEN
707 dlim = max(zero,e6(i,5)) / (dmx(i) + vfail)
708 ENDIF
709 ENDIF
710 IF (ifail(i) == 0) THEN
711! Uniaxial failure
712 crit(i) = max(crit(i),xa*dlim)
713 ELSE
714! Multiaxial failure
715 crit(i)= crit(i) + xa *dlim**xb
716 ENDIF
717 ENDIF
718 ENDDO
719C
720 DO i=1,nel
721 ifunc(i) = igeo(116,mgn(i))
722 ifv(i) = igeo(117,mgn(i))
723 ifunc2(i)= igeo(118,mgn(i))
724 ifunc3(i)= igeo(124,mgn(i))
725 iecrou(i)=nint(geo(30,mgn(i)))
726 ak(i) =geo(61,mgn(i))
727 b(i) =geo(62,mgn(i))
728 d(i) =geo(63,mgn(i))
729 ee(i) =geo(184,mgn(i))
730 gf3(i) =geo(137,mgn(i))
731 ff(i) =geo(64,mgn(i))
732 lscale(i)= geo(178,mgn(i))
733 dmn(i) =geo(75,mgn(i))
734 dmx(i) =geo(76,mgn(i))
735 ENDDO
736 IF (nuvar > 0) xx_old => uvar(6,1:nel)
737 CALL redef3(python,
738 1 zmom, zk, rz, zmep,
739 2 dzold, rpz, tf, npf,
740 3 zc, off, e6(1,6), rpz2,
741 4 anim, 0, poszz,
742 5 xl0, dmn, dmx, dv,
743 6 ff, lscale, ee, gf3,
744 7 ifunc3, yieldz2, aldp, ak,
745 8 b, d, iecrou, ifunc,
746 9 ifv, ifunc2, epla, xx_old,
747 a nel, nft, stf, sanin,
748 b dt1, iresp, impl_s, idyna,
749 c snpc)
750C
751 DO i=1,nel
752 cc = geo(108,mgn(i))
753 cn = geo(114,mgn(i))
754 xa = geo(120,mgn(i))
755 xb = geo(126,mgn(i))
756 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
757 IF (ifail2(i) == 0) THEN
758 xa = one
759 xb = two
760 IF (rz(i) > zero) THEN
761 dlim = rz(i) / (xl0(i)*dmx(i))
762 ELSE
763 dlim = rz(i) / (xl0(i)*dmn(i))
764 ENDIF
765 ELSE
766 vfail = cc * (abs(dv(i)/vrr(i)))**cn
767 IF (ifail2(i) == 1) THEN
768 IF (rz(i) > zero)THEN
769 dlim = rz(i) / (xl0(i)*dmx(i) + vfail)
770 ELSE
771 dlim = rz(i) / (xl0(i)*dmn(i) - vfail)
772 ENDIF
773 ELSEIF (ifail2(i) == 2) THEN
774 IF (zmom(i) > zero)THEN
775 dlim = zmom(i)/(dmx(i) + vfail)
776 ELSE
777 dlim = zmom(i)/(dmn(i) - vfail)
778 ENDIF
779 ELSEIF (ifail2(i) == 3) THEN
780 dlim = max(zero,e6(i,6)) / (dmx(i) + vfail)
781 ENDIF
782 ENDIF
783 IF (ifail(i) == 0) THEN
784! Uniaxial failure
785 crit(i) = max(crit(i),xa*dlim)
786 ELSE
787! Multiaxial failure
788 crit(i)= crit(i) + xa *dlim**xb
789 ENDIF
790 ENDIF
791 ENDDO
792C
793 DO i=1,nel
794 e(i) = e6(i,1)+e6(i,2)+e6(i,3)+e6(i,4)+e6(i,5)+e6(i,6)
795 ENDDO
796C-------------------------------
797C COUPLED FAILURE
798C-------------------------------
799 DO i=1,nel
800 israte = int(geo(127,mgn(i)))
801C---- smoothing factor alpha = 2PI*fc*dt/(2PI*fc*dt+1) ---
802 asrate = (2*pi*geo(128,mgn(i))*dt1)/(one+2*pi*geo(128,mgn(i))*dt1)
803 IF (israte /= 0) THEN
804 IF (crit_new(i) < one) THEN
805 crit(i) = min(crit(i),one+em3)
806 crit(i) = asrate*crit(i) + (one - asrate)*crit_new(i)
807 crit_new(i) = min(crit(i),one)
808 ELSE
809 crit_new(i) = one
810 ENDIF
811 ELSE
812 IF (crit_new(i) < one) THEN
813 crit_new(i) = min(crit(i),one)
814 ELSE
815 crit_new(i) = one
816 ENDIF
817 ENDIF
818 IF (off(i) == one) THEN
819 IF (crit(i) >= one) THEN
820 off(i)=zero
821 nindx = nindx + 1
822 indx(nindx) = i
823 idel7nok = 1
824 ENDIF
825 ENDIF
826 ENDDO
827C
828 DO j=1,nindx
829 i = indx(j)
830#include "lockon.inc"
831 WRITE(iout, 1000) ngl(i)
832 WRITE(istdo,1100) ngl(i),tt
833#include "lockoff.inc"
834 ENDDO
835C-------------------------------
836C COUPLED PLASTICITY
837C-------------------------------
838 CALL repla3(
839 1 xk, rpx, tf, npf,
840 2 iecrou, ifunc, ifv, epla,
841 3 nel)
842 CALL repla3(
843 1 yk, rpy, tf, npf,
844 2 iecrou, ifunc, ifv, epla,
845 3 nel)
846 CALL repla3(
847 1 zk, rpz, tf, npf,
848 2 iecrou, ifunc, ifv, epla,
849 3 nel)
850C
851 DO i=1,nel
852 xk(i)=geo(3,mgn(i))
853 yk(i)=geo(10,mgn(i))
854 zk(i)=geo(15,mgn(i))
855 ENDDO
856C
857 CALL repla3(
858 1 xk, dpx, tf, npf,
859 2 iecrou, ifunc, ifv, epla,
860 3 nel)
861 CALL repla3(
862 1 yk, dpy, tf, npf,
863 2 iecrou, ifunc, ifv, epla,
864 3 nel)
865 CALL repla3(
866 1 zk, dpz, tf, npf,
867 2 iecrou, ifunc, ifv, epla,
868 3 nel)
869 DO i=1,nel
870 xm(i)=xm(i)*xl0(i)
871 xkm(i)=xkm(i)/xl0(i)
872 xcm(i)=xcm(i)/xl0(i)
873 xin(i)=xin(i)*xl0(i)
874 xkr(i)=xkr(i)/xl0(i)
875 xcr(i)=xcr(i)/xl0(i)
876 ENDDO
877C---
878 1000 FORMAT(1x,'-- rupture of spring element number ',I10)
879 1100 FORMAT(1X,'-- rupture of spring element :',I10,' at time :',G11.4)
880C---
881 RETURN
882 END
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
r4def3
subroutine r4def3(python, skew, geo, fx, fy, fz, e, dx, dy, dz, npf, tf, off, dpx, dpy, dpz, dpx2, dpy2, dpz2, fxep, fyep, fzep, x0, y0, z0, xmom, ymom, zmom, rx, ry, rz, rpx, rpy, rpz, xmep, ymep, zmep, rpx2, rpy2, rpz2, anim, posx, posy, posz, posxx, posyy, poszz, fr_wave, e6, nel, exx2, eyx2, ezx2, exy2, eyy2, ezy2, exz2, eyz2, ezz2, al2dp, igeo, crit_new, x0_err, aldp, yieldx, yieldy, yieldz, yieldx2, yieldy2, yieldz2, ngl, mgn, exx, eyx, ezx, exy, eyy, ezy, exz, eyz, ezz, xcr, rx1, ry1, rz1, rx2, ry2, rz2, xin, ak, xm, xkm, xcm, xkr, vx1, vx2, vy1, vy2, vz1, vz2, nuvar, uvar, dx0, dy0, dz0, rx0, ry0, rz0, fx0, fy0, fz0, xmom0, ymom0, zmom0, nft, stf, sanin, iresp, snpc)
Definition r4def3.F:65
repla3
subroutine repla3(xk, dpx, tf, npf, iecrou, ifunc, ifv, epla, nel)
Definition repla3.F:39