OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
r23l114def3.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!|| r23l114def3 ../engine/source/elements/spring/r23l114def3.F
25!||--- called by ------------------------------------------------------
26!|| r23law114 ../engine/source/elements/spring/r23law114.F
27!||--- calls -----------------------------------------------------
28!|| material_flow ../engine/source/tools/seatbelts/material_flow.f
29!|| redef3 ../engine/source/elements/spring/redef3.F90
30!|| redef_seatbelt ../engine/source/tools/seatbelts/redef_seatbelt.F90
31!|| repla3 ../engine/source/elements/spring/repla3.F
32!||--- uses -----------------------------------------------------
33!|| python_funct_mod ../common_source/modules/python_mod.F90
34!|| redef3_mod ../engine/source/elements/spring/redef3.F90
35!|| redef_seatbelt_mod ../engine/source/tools/seatbelts/redef_seatbelt.F90
36!|| seatbelt_mod ../common_source/modules/seatbelt_mod.F
37!|| sensor_mod ../common_source/modules/sensor_mod.F90
38!||====================================================================
39 SUBROUTINE r23l114def3( PYTHON,
40 1 SKEW, IPM, IGEO, MID,
41 2 PID, GEO, UPARAM, FX,
42 3 FY, FZ, E, DX,
43 4 DY, DZ, NPF, TF,
44 5 OFF, DPX, DPY, DPZ,
45 6 DPX2, DPY2, DPZ2, FXEP,
46 7 FYEP, FZEP, X0, Y0,
47 8 Z0, XMOM, YMOM, ZMOM,
48 9 RX, RY, RZ, RPX,
49 A RPY, RPZ, XMEP, YMEP,
50 B ZMEP, RPX2, RPY2, RPZ2,
51 C ANIM, POSX, POSY, POSZ,
52 D POSXX, POSYY, POSZZ, FR_WAVE,
53 E E6, NEL, EXX2, EYX2,
54 F EZX2, EXY2, EYY2, EZY2,
55 G EXZ2, EYZ2, EZZ2, AL2DP,
56 H NGL, CRIT_NEW, X0_ERR, ALDP,
57 I YIELDX, YIELDY, YIELDZ, YIELDX2,
58 J YIELDY2, YIELDZ2, EXX, EYX,
59 K EZX, EXY, EYY, EZY,
60 L EXZ, EYZ, EZZ, XCR,
61 M RX1, RY1, RZ1, RX2,
62 N RY2, RZ2, XIN, AK,
63 O XM, XKM, XCM, XKR,
64 P VX1, VX2, VY1, VY2,
65 Q VZ1, VZ2, NUVAR, UVAR,
66 R MASS, DX0, DY0, DZ0,
67 S RX0, RY0, RZ0, SLIPRING_STRAND,
68 T DFS, RING_SLIP, X02, LMIN,
69 U SLIPRING_ID, UPDATE_FLAG, RETRACTOR_ID, ADD_NODE1,
70 V ADD_NODE2, NC1, NC2, NC3,
71 W X1DP, X2DP, X3DP, VX3,
72 X VY3, VZ3, FLAG_SLIPRING_UPDATE, FLAG_RETRACTOR_UPDATE,
73 Y SENSOR_TAB, UINER, FR_ID, FRAM_FACTOR,
74 Z EPS_OLD, FX_B2, DPX_B2, YIELDX_B2,
75 1 XX_OLD_B2, FXEP_B2, POSX_B2, EPS_OLD_B2,
76 2 NFT , NSENSOR, STF, SANIN,
77 3 IRESP, SNPC)
78C-----------------------------------------------
79C M o d u l e s
80C-----------------------------------------------
81 USE python_funct_mod
82 USE redef3_mod
83 USE seatbelt_mod
84 USE redef_seatbelt_mod
85 USE sensor_mod
86C-----------------------------------------------
87C I m p l i c i t T y p e s
88C-----------------------------------------------
89#include "implicit_f.inc"
90#include "comlock.inc"
91C-----------------------------------------------
92C G l o b a l P a r a m e t e r s
93C-----------------------------------------------
94#include "mvsiz_p.inc"
95C-----------------------------------------------
96C C o m m o n B l o c k s
97C-----------------------------------------------
98#include "param_c.inc"
99#include "com04_c.inc"
100#include "com08_c.inc"
101#include "scr14_c.inc"
102#include "scr17_c.inc"
103#include "units_c.inc"
104#include "com01_c.inc"
105#include "impl1_c.inc"
106C-----------------------------------------------
107C D u m m y A r g u m e n t s
108C-----------------------------------------------
109 TYPE(python_) :: PYTHON
110 INTEGER, INTENT(IN) :: STF !< Size of TF
111 INTEGER, INTENT(IN) :: SANIN !< Size of ANIM
112 INTEGER, INTENT(IN) :: IRESP !< Single response flag
113 INTEGER, INTENT(IN) :: SNPC !< Size of NPF
114 INTEGER, INTENT(IN) :: NFT,NSENSOR
115 INTEGER NPF(*),IGEO(NPROPGI,*),NEL,NGL(*),PID(*),MID(*),NUVAR,
116 . IPM(NPROPMI,*),NC1(*),NC2(*),NC3(*),ADD_NODE1(*),ADD_NODE2(*),
117 . SLIPRING_ID(*),UPDATE_FLAG(*),RETRACTOR_ID(*),SLIPRING_STRAND(*),
118 . FLAG_SLIPRING_UPDATE,FLAG_RETRACTOR_UPDATE,FR_ID(*)
119C REAL
120 my_real
121 . skew(lskew,*), geo(npropg,*), fx(*), fy(*), fz(*), e(*), dx(*),
122 . dy(*), dz(*), tf(stf), off(*), dpx(*), dpy(*), dpz(*), fxep(*),
123 . fyep(*), fzep(*), x0(*), y0(*), z0(*), xmom(*), ymom(*),
124 . zmom(*), rx(*), ry(*), rz(*), rpx(*), rpy(*), rpz(*), xmep(*),
125 . ymep(*), zmep(*), dpx2(*), dpy2(*), dpz2(*),rpx2(*), rpy2(*),
126 . rpz2(*),anim(sanin),posx(*),posy(*),posz(*),posxx(*),
127 . posyy(*),poszz(*),fr_wave(*),e6(nel,6),
128 . exx2(mvsiz), eyx2(mvsiz), ezx2(mvsiz),
129 . exy2(mvsiz), eyy2(mvsiz), ezy2(mvsiz),
130 . exz2(mvsiz), eyz2(mvsiz), ezz2(mvsiz),
131 . crit_new(*), x0_err(mvsiz),yieldx(*),yieldy(*),
132 . yieldz(*),yieldx2(*),yieldy2(*),yieldz2(*),
133 . exx(mvsiz), eyx(mvsiz), ezx(mvsiz), exy(mvsiz),
134 . eyy(mvsiz), ezy(mvsiz), exz(mvsiz), eyz(mvsiz),
135 . ezz(mvsiz), xcr(mvsiz), rx1(mvsiz), rx2(mvsiz),
136 . ry1(mvsiz), ry2(mvsiz), rz1(mvsiz), rz2(mvsiz),
137 . xin(mvsiz),ak(mvsiz),xm(mvsiz),xkm(mvsiz),xcm(mvsiz),
138 . xkr(mvsiz),vx1(mvsiz),vx2(mvsiz),vy1(mvsiz),vy2(mvsiz),
139 . vz1(mvsiz),vz2(mvsiz),uvar(nuvar,*),uparam(*),mass(*),
140 . dx0(*),dy0(*),dz0(*),rx0(*),ry0(*),rz0(*),lmin(*),dfs(*),
141 . ring_slip(*),x02(*),vx3(mvsiz),vy3(mvsiz),vz3(mvsiz),
142 . uiner(*),fram_factor(*)
143 my_real ,INTENT(INOUT) :: eps_old(mvsiz),fx_b2(mvsiz),dpx_b2(mvsiz),yieldx_b2(mvsiz),
144 . xx_old_b2(mvsiz),fxep_b2(mvsiz),posx_b2(mvsiz),eps_old_b2(mvsiz)
145 TARGET :: uvar
146 DOUBLE PRECISION ALDP(MVSIZ),AL2DP(MVSIZ),X1DP(3,*),X2DP(3,*),X3DP(3,*)
147 TYPE (SENSOR_STR_) ,DIMENSION(NSENSOR) :: SENSOR_TAB
148C-----------------------------------------------
149C L o c a l V a r i a b l e s
150C-----------------------------------------------
151 INTEGER INDX(MVSIZ),
152 . iecrou(mvsiz), ifunc(mvsiz), ifv(mvsiz), ifunc2(mvsiz),
153 . i,k,ileng, j, kk, ifail(mvsiz),ifail2(mvsiz),adher(mvsiz),
154 . nindx,israte, ifunc3(mvsiz),i1,i2,i3,i4,i5,i6,i7,i8,
155 . i9,i10,i11,i12,i13,i14,if1,if2,if3,if4,iadbuf,nupar,nn1,nn2,ret,
156 . flag_retractor_update_old,flag,index1(mvsiz),index2(mvsiz),compt,
157 . dir,strd,ii,index_slip(mvsiz),compt2
158C REAL
159 my_real
160 . xk(mvsiz), yk(mvsiz),
161 . zk(mvsiz),xc(mvsiz), yc(mvsiz), zc(mvsiz),xh(mvsiz),
162 . xhr(mvsiz),dxold(mvsiz), dyold(mvsiz), dzold(mvsiz),
163 . b(mvsiz), d(mvsiz), epla(mvsiz),
164 . dv(mvsiz),vrt(mvsiz),vrr(mvsiz),
165 . dmn(mvsiz),dmx(mvsiz),xl0(mvsiz),crit(mvsiz),
166 . xn(mvsiz),ff(mvsiz),lscale(mvsiz),ee(mvsiz),gf3(mvsiz),
167 . hx(mvsiz), hy(mvsiz), hz(mvsiz),fx_max(mvsiz),e_offset(mvsiz),
168 . xk_comp(mvsiz),mx_max(mvsiz),dfs_old(mvsiz),not_used,xl02(mvsiz),
169 . xkp(mvsiz),dfx(mvsiz)
170 my_real
171 . at,dt05,xka,yka,zka,cc,cn,xa,xb,dlim,vfail,
172 . x21, y21, z21, epxy, epxz,
173 . vx21, vy21, vz21, ryavp, rzavp,eyzp,exzp,
174 . ryav, rzav,den, c, cp, exyp,
175 . x21phi, y21phi, z21phi, vx21phi, vy21phi, vz21phi,
176 . ryav1, rzav1, ryav1p, rzav1p,asrate,norm,gap,lmin_b2(mvsiz),
177 . eb(mvsiz),dxoldb(mvsiz),xkp_b2(mvsiz),dx_b2(mvsiz),a1
178 DOUBLE PRECISION X0DP(MVSIZ),X0DP_B2(MVSIZ),EX2DP(MVSIZ),EY2DP(MVSIZ),EZ2DP(MVSIZ),
179 . ALDP_B2(MVSIZ)
180 my_real ,DIMENSION(:), POINTER :: COORD_OLD
181C-----------------------------------------------
182C
183 not_used = zero
184C
185C-----------------------------------------------
186C-- FRAM_FACTOR - used for 2D seatbelts --
187C = 1.0 -> 1D Seatbelt spring
188C = 1/(Nb_fram -1) ----> 2D Seatbelt spring inside belt
189C = 0.5/(Nb_fram -1) --> 2D Seatbelt spring on belt edge
190C-----------------------------------------------
191C
192 nupar = 4
193 i1 = nupar
194 i2 = i1 + 6
195 i3 = i2 + 6
196 i4 = i3 + 6
197 i5 = i4 + 6
198 i6 = i5 + 6
199 i7 = i6 + 6
200 i8 = i7 + 6
201 i9 = i8 + 6
202 i10 = i9 + 6
203 i11 = i10 + 6
204 i12 = i11 + 6
205 i13 = i12 + 6
206 i14 = i13 + 6
207 nupar = nupar + 84
208 DO i=1,nel
209C
210 index1(i) = i
211 iadbuf= ipm(7,mid(i)) - 1
212 epla(i)=zero
213 dfx(i)=zero
214 xm(i)=mass(i)
215C
216 xk(i)=uparam(iadbuf + i11 + 1)
217 yk(i)=uparam(iadbuf + i11 + 2)
218 zk(i)=uparam(iadbuf + i11 + 3)
219C
220 xc(i)=uparam(iadbuf + i12 + 1) * fram_factor(i)
221 yc(i)=uparam(iadbuf + i12 + 2) * fram_factor(i)
222 zc(i)=uparam(iadbuf + i12 + 3) * fram_factor(i)
223C
224 xk_comp(i) = uparam(iadbuf + 117)*geo(1,pid(i))
225C
226 xka=xk(i)*uparam(iadbuf + i1 + 1)*fram_factor(i)
227 yka=yk(i)*uparam(iadbuf + i1 + 2)*fram_factor(i)
228 zka=zk(i)*uparam(iadbuf + i1 + 3)*fram_factor(i)
229 xkm(i)= max(xka,yka,zka,xk_comp(i))
230C
231 hx(i) = uparam(iadbuf + i14 + 1)
232 hy(i) = uparam(iadbuf + i14 + 2)
233 hz(i) = uparam(iadbuf + i14 + 3)
234C
235 xh(i)= max(hx(i),hy(i),hz(i))
236 xcm(i)= max(xc(i),yc(i),zc(i))
237 xcm(i)= xcm(i)+xh(i)
238
239 xkr(i)= max(yka,zka) * aldp(i) * aldp(i)
240 xcr(i)= (max(yc(i),zc(i)) + max(hy(i),hz(i))) * aldp(i) * aldp(i)
241 vrt(i) = uparam(iadbuf + nupar + 1)
242 vrr(i) = uparam(iadbuf + nupar + 2)
243 ifail(i) = nint(uparam(iadbuf + 1 ))
244 ifail2(i)= nint(uparam(iadbuf + 3 ))
245C
246 e_offset(i) = uparam(iadbuf + 118)
247 IF (tt == zero) lmin(i) = uparam(iadbuf + 119)
248 fx_max(i) = uparam(iadbuf + 120)
249 mx_max(i) = uparam(iadbuf + 121)
250 ENDDO
251C
252 IF (inispri /= 0 .and. tt == zero) THEN
253 DO i=1,nel
254 xl0(i)= x0(i)
255! if not initialized length
256 IF (xl0(i) == zero) xl0(i) = aldp(i)
257 ENDDO
258 ENDIF
259C
260 IF (tt == zero) THEN
261 DO i=1,nel
262 x0(i)= aldp(i) ! cast double vers My_real
263 ENDDO
264 ENDIF
265C
266 IF (scodver >= 101) THEN
267 IF (tt == zero) THEN
268 DO i=1,nel
269 x0_err(i)= aldp(i)-x0(i) ! difference between double and My_real
270 ENDDO
271 ENDIF
272 ENDIF
273C
274 IF ( inispri /= 0 .and. tt == zero) THEN
275 DO i=1,nel
276 x0(i)= xl0(i)
277 ENDDO
278 ENDIF
279C
280 DO i=1,nel
281 x0dp(i)= x0(i) ! cast double vers My_real
282 ENDDO
283C
284 IF (scodver >= 101) THEN
285 DO i=1,nel
286 x0dp(i)= x0(i) + x0_err(i) ! difference between double and My_real
287 ENDDO
288 ENDIF
289C---------------------
290C TRANSLATIONS
291C---------------------
292 DO i=1,nel
293 dxold(i)=dx(i)
294 dyold(i)=dy(i)
295 dzold(i)=dz(i)
296 ENDDO
297!
298 IF (inispri /= 0 .and. tt == zero) THEN
299 DO i=1,nel
300 dxold(i)=dx0(i)
301 dyold(i)=dy0(i)
302 dzold(i)=dz0(i)
303 ENDDO
304 ENDIF
305C
306 dt05=half*dt1
307 IF (ismdisp > 0) THEN
308 DO i=1,nel
309 vx21 = vx2(i)-vx1(i)
310 vy21 = vy2(i)-vy1(i)
311 vz21 = vz2(i)-vz1(i)
312 dx(i) = dxold(i)+(vx21*exx(i)+vy21*eyx(i)+vz21*ezx(i))*dt1
313 dy(i) = dyold(i)+(vx21*exy(i)+vy21*eyy(i)+vz21*ezy(i))*dt1
314 dz(i) = dzold(i)+(vx21*exz(i)+vy21*eyz(i)+vz21*ezz(i))*dt1
315C
316 x21 = (rx2(i)+rx1(i))
317 y21 = (ry2(i)+ry1(i))
318 z21 = (rz2(i)+rz1(i))
319C
320 ryav1 = (x21*exy2(i)+y21*eyy2(i)+z21*ezy2(i))
321 rzav1 = (x21*exz2(i)+y21*eyz2(i)+z21*ezz2(i))
322C
323 ryav = dt05 * ryav1
324 rzav = dt05 * rzav1
325C
326 dy(i) = dy(i) - rzav * al2dp(i)
327 dz(i) = dz(i) + ryav * al2dp(i)
328C
329 crit(i) = zero
330 ENDDO
331 ELSE
332 DO i=1,nel
333 vx21 = vx2(i)-vx1(i)
334 vy21 = vy2(i)-vy1(i)
335 vz21 = vz2(i)-vz1(i)
336C
337 epxy = (vx21*exy2(i)+vy21*eyy2(i)+vz21*ezy2(i))*dt05
338 epxz = (vx21*exz2(i)+vy21*eyz2(i)+vz21*ezz2(i))*dt05
339C
340 x21 = (rx2(i)+rx1(i))
341 y21 = (ry2(i)+ry1(i))
342 z21 = (rz2(i)+rz1(i))
343C
344 ryav1 = (x21*exy2(i)+y21*eyy2(i)+z21*ezy2(i))
345 rzav1 = (x21*exz2(i)+y21*eyz2(i)+z21*ezz2(i))
346C
347 at=epxz/max(al2dp(i),em30)
348 at=atan(at)
349 ryav = dt05 * (ryav1) + two * at
350 at=epxy/max(al2dp(i),em30)
351 at=atan(at)
352 rzav = dt05 * (rzav1) - two * at
353C
354 dx(i) = aldp(i) - x0dp(i)
355 dy(i) = dyold(i) - rzav * al2dp(i)
356 dz(i) = dzold(i) + ryav * al2dp(i)
357C
358 crit(i) = zero
359 ENDDO
360 ENDIF !(ISMDISP > 0) THEN
361C
362 DO i=1,nel
363 iadbuf = ipm(7,mid(i)) - 1
364 ileng = nint(uparam(iadbuf + 2))
365 IF (ileng /= 0) THEN
366 xl0(i)=max(x0dp(i),lmin(i))
367 ELSE
368 xl0(i)=one
369 ENDIF
370 ENDDO
371C-------------------------------
372 nindx = 0
373 if1 = 0
374 if2 = 6
375 if3 = 12
376 if4 = 18
377 compt = 0
378C
379 DO i=1,nel
380 iadbuf = ipm(7,mid(i)) - 1
381 ifunc(i) = ipm(10 + if1 + 1,mid(i))
382 ifv(i) = ipm(10 + if2 + 1,mid(i))
383C unloading curve activated only if cross point is passed
384 IF (yieldx(i) > uparam(iadbuf + 125)) THEN
385 ifunc2(i)= ipm(10 + if3 + 1,mid(i))
386 ELSE
387 ifunc2(i)= ifunc(i)
388 ENDIF
389 ifunc3(i)= ipm(10 + if4 + 1,mid(i))
390 iecrou(i)= nint(uparam(iadbuf + i13 + 1))
391C
392C-- parameters scaled by fram_factor for springs of 2D seatbelts
393 ak(i) = uparam(iadbuf + i1 + 1) * fram_factor(i)
394 b(i) = uparam(iadbuf + i2 + 1)
395 d(i) = uparam(iadbuf + i3 + 1)
396 ee(i) = uparam(iadbuf + i4 + 1)
397 gf3(i) = uparam(iadbuf + i5 + 1)
398 ff(i) = uparam(iadbuf + i6 + 1)
399 lscale(i)= uparam(iadbuf + i7 + 1)
400 dmn(i) = uparam(iadbuf + i8 + 1)
401 dmx(i) = uparam(iadbuf + i9 + 1)
402C
403 IF (update_flag(i) /= 0) THEN
404 compt = compt + 1
405 index2(compt) = i
406 ENDIF
407 ENDDO
408C
409 IF (nslipring + nretractor > 0) THEN
410C
411 IF (compt > 0) THEN
412C
413C---------------------
414C SLIPRING / RETRACTOR UPDATE - step3
415C---------------------
416C
417 DO ii=1,compt
418 i = index2(ii)
419C
420 IF (slipring_id(i) > 0) THEN
421C-- slipring update step3
422 j = slipring_id(i)
423 k = fr_id(i)
424 IF (slipring_strand(i) == 1) THEN
425 lmin(i) = uparam(iadbuf + 119)
426 xl0(i)=max(x0dp(i),lmin(i))
427C-- Loading of internal variables for computation of strand 2
428 x02(i) = slipring(j)%FRAM(k)%RESIDUAL_LENGTH(2)
429 fx_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR2(1)
430 dpx_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR2(2)
431 yieldx_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR2(3)
432 xx_old_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR2(4)
433 fxep_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR2(5)
434 posx_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR2(6)
435 eps_old_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR2(7)
436 IF (update_flag(i) <= 0) ring_slip(i) = x0(i)
437 ELSEIF (slipring_strand(i) == 2) THEN
438 lmin(i) = uparam(iadbuf + 119)
439 xl0(i)=max(x0dp(i),lmin(i))
440C-- Loading of internal variables for computation of strand 2
441 x02(i) = slipring(j)%FRAM(k)%RESIDUAL_LENGTH(1)
442 fx_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR1(1)
443 dpx_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR1(2)
444 yieldx_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR1(3)
445 xx_old_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR1(4)
446 fxep_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR1(5)
447 posx_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR1(6)
448 eps_old_b2(i) = slipring(j)%FRAM(k)%INTVAR_STR1(7)
449 ENDIF
450C
451 ELSEIF ((retractor_id(i) > 0).AND.(update_flag(i) == -1)) THEN
452C-- retractor update of mouth element step3
453 off(i) = one
454C
455 ELSEIF ((retractor_id(i) > 0).AND.(update_flag(i) == -2)) THEN
456C-- retractor update of new entered element step3
457 off(i) = zero
458 update_flag(i) = 0
459 x0(i) = zero
460C
461 ENDIF
462C
463 ENDDO
464C
465 ENDIF
466C
467 compt = 0
468 compt2 = 0
469 DO i=1,nel
470 adher(i) = 0
471 IF (slipring_strand(i) > 0) THEN
472 compt = compt + 1
473 compt2 = compt2 + 1
474 index2(compt) = i
475 index_slip(compt2) = i
476 IF (slipring(slipring_id(i))%NFRAM > 1) THEN
477 yc(i) = 0.1*xc(i)
478 zc(i) = 0.1*xc(i)
479 yk(i) = 0.01*xk(i)
480 zk(i) = 0.01*xk(i)
481 ELSE
482 yc(i) = zero
483 zc(i) = zero
484 ENDIF
485C- No damping in slipring
486 xc(i) = zero
487 xcm(i)= max(xc(i),yc(i),zc(i))+xh(i)
488 xk_comp(i) = xk(i)
489 fx_max(i) = ep20
490 ELSEIF (slipring_strand(i) == -1) THEN
491C- No damping in retractor
492 compt = compt + 1
493 index2(compt) = i
494 xc(i) = zero
495 yc(i) = zero
496 zc(i) = zero
497 xcm(i)= max(xc(i),yc(i),zc(i))+xh(i)
498 xk_comp(i) = xk(i)
499 fx_max(i) = ep20
500 ENDIF
501 ENDDO
502C
503 ENDIF
504C
505 IF (nuvar >=1) coord_old => uvar(1,1:nel)
506C
507 flag = 1
508 CALL redef_seatbelt(python,
509 1 fx, xk, dx, fxep,
510 2 dxold, dpx, tf, npf,
511 3 xc, off, e6(1,1), anim,
512 4 anim_fe(11),posx(1), xl0, dmn,
513 5 dmx, lscale, yieldx, ak,
514 6 iecrou, ifunc, ifunc2, coord_old,
515 7 fx_max, xk_comp, nel, index1,
516 8 flag, xkp, eps_old, fram_factor,
517 9 nft, snpc, stf, sanin,
518 a dt1, impl_s, idyna, nel)
519C
520
521C
522 IF (nslipring + nretractor > 0) THEN
523C---------------------
524C MATERIAL FLOW COMPUTATION FOR SLIPRING/RETRACTOR
525C---------------------
526C
527C- Computation of length of 2nd strand
528 DO j=1,compt2
529 i = index_slip(j)
530 strd = slipring_strand(i)
531 dir = slipring(slipring_id(i))%FRAM(fr_id(i))%STRAND_DIRECTION(strd)
532C--- Depending on the orientation of the spring the 2nd strand is defined by IXR(2)/IXR(3) or IXR(1)/IXR(3)
533 IF (((strd==1).AND.(dir ==1)).OR.((strd==2).AND.(dir==-1))) THEN
534 ex2dp(i)=x2dp(1,i)-x3dp(1,i)
535 ey2dp(i)=x2dp(2,i)-x3dp(2,i)
536 ez2dp(i)=x2dp(3,i)-x3dp(3,i)
537 ELSE
538 ex2dp(i)=x3dp(1,i)-x1dp(1,i)
539 ey2dp(i)=x3dp(2,i)-x1dp(2,i)
540 ez2dp(i)=x3dp(3,i)-x1dp(3,i)
541 ENDIF
542 norm= max(em15,sqrt(ex2dp(i)*ex2dp(i)+ey2dp(i)*ey2dp(i)+ez2dp(i)*ez2dp(i)))
543 ex2dp(i)= ex2dp(i)/norm
544 ey2dp(i)= ey2dp(i)/norm
545 ez2dp(i)= ez2dp(i)/norm
546 aldp_b2(i)= norm
547 ENDDO
548C
549 IF (tt == zero) THEN
550 DO j=1,compt2
551 i = index_slip(j)
552 x02(i)= aldp_b2(i)
553 ENDDO
554 ENDIF
555C
556 DO j=1,compt2
557 i = index_slip(j)
558 x0dp_b2(i)= x02(i) ! cast double vers My_real
559 xl02(i)=max(x0dp_b2(i),lmin(i))
560 dx_b2(i) = aldp_b2(i) - x0dp_b2(i)
561C for 2n strand - unloading curve activated only if cross point is passed
562 IF (yieldx_b2(i) > uparam(iadbuf + 125)) THEN
563 ifunc2(i)= ipm(10 + if3 + 1,mid(i))
564 ELSE
565 ifunc2(i)= ifunc(i)
566 ENDIF
567 ENDDO
568C
569C-- Computation of stress of 2nd strand for sliprings
570 flag = 2
571 eb(1:mvsiz) = zero
572 dxoldb(1:mvsiz) = zero
573 CALL redef_seatbelt(python,
574 1 fx_b2, xk, dx_b2, fxep_b2,
575 2 dxoldb, dpx_b2, tf, npf,
576 3 xc, off, eb, anim,
577 4 anim_fe(11),posx_b2, xl02, dmn,
578 5 dmx, lscale, yieldx_b2, ak,
579 6 iecrou, ifunc, ifunc2, xx_old_b2,
580 7 fx_max, xk_comp, compt2, index_slip,
581 8 flag, xkp_b2, eps_old_b2, fram_factor,
582 9 nft, snpc, stf, sanin,
583 a dt1, impl_s, idyna, mvsiz)
584C
585 DO j=1,compt2
586 i = index_slip(j)
587 IF (aldp(i) > aldp_b2(i)) xkp_b2(i) = xkp(i)
588 IF (aldp_b2(i) > aldp(i)) xkp(i) = xkp_b2(i)
589 ENDDO
590
591C-- Computation of material flow for slipring and retractor
592 CALL material_flow(dfs,dfs_old,aldp,slipring_strand,xkp,
593 1 off,x0,x02,lmin,update_flag,
594 2 ring_slip,slipring_id,xl0,dx,dxold,
595 3 exx,eyx,ezx,x1dp,x2dp,
596 4 x3dp,adher,nc1,nc2,nc3,
597 5 flag_slipring_update,add_node1,add_node2,vx1,vy1,
598 6 vz1,vx2,vy2,vz2,vx3,
599 7 vy3,vz3,xc,retractor_id,flag_retractor_update,
600 8 sensor_tab,x0dp,fr_id,dfx,fx,fx_b2,
601 9 aldp_b2,x0dp_b2,ex2dp,ey2dp,ez2dp,xl02,
602 a xkp_b2,compt,index2,nsensor)
603C
604 DO i=1,nel
605 IF (slipring_strand(i) > 0) THEN
606 fx(i) = fx(i) + dfx(i)
607 e6(i,1) = e6(i,1) + (dx(i)-dxold(i))*dfx(i)*half
608 dx(i) = aldp(i)-x0dp(i)
609 xkm(i) = max(xkm(i),fram_factor(i)*xk(i)*(one + xl0(i)/xl02(i)))
610 slipring(slipring_id(i))%FRAM(fr_id(i))%SLIP_FORCE(slipring_strand(i)) = fx(i)
611 ELSEIF (slipring_strand(i) < 0) THEN
612 fx(i) = fx(i) + dfx(i)
613 e6(i,1) = e6(i,1) + (dx(i)-dxold(i))*dfx(i)*half
614 dx(i) = aldp(i)-x0dp(i)
615 retractor(retractor_id(i))%RET_FORCE = fx(i)
616 ENDIF
617 END DO
618C
619 ENDIF
620C
621C---------------------
622C
623 DO i=1,nel
624 cc = uparam(iadbuf + nupar + 3)
625 cn = uparam(iadbuf + nupar + 9)
626 xa = uparam(iadbuf + nupar + 15)
627 xb = uparam(iadbuf + nupar + 21)
628 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i)/= zero) THEN
629 IF (ifail2(i) == 0) THEN
630 xa = one
631 xb = two
632 IF (dx(i) > zero) THEN
633 dlim = dx(i) / dmx(i)
634 ELSE
635 dlim = dx(i) / dmn(i)
636 ENDIF
637 ELSE
638 vfail = cc * (abs(dv(i)/vrt(i)))**cn
639 IF (ifail2(i) == 1) THEN
640 IF (dx(i) > zero) THEN
641 dlim = dx(i) / (dmx(i) + vfail)
642 ELSE
643 dlim = dx(i) / (dmn(i) - vfail)
644 ENDIF
645 ELSEIF (ifail2(i) == 2) THEN
646 IF (fx(i) > zero) THEN
647 dlim = fx(i) / (dmx(i) + vfail)
648 ELSE
649 dlim = fx(i) / (dmn(i) - vfail)
650 ENDIF
651 ELSEIF (ifail2(i) == 3) THEN
652 dlim = max(zero,e6(i,1)) / (dmx(i) + vfail)
653 ENDIF
654 ENDIF
655 IF (ifail(i) == 0) THEN
656! Uniaxial rupture
657 IF ((xa*dlim) > one) THEN
658 off(i) = zero
659 nindx = nindx + 1
660 indx(nindx) = i
661 idel7nok = 1
662 ENDIF
663 ELSE
664! Multiaxial rupture
665 crit(i)= crit(i) + xa *(dlim/xl0(i))**xb
666 ENDIF
667 ENDIF
668 ENDDO
669C
670 DO i=1,nel
671 iadbuf = ipm(7,mid(i)) - 1
672 ifunc(i) = 0
673 ifv(i) = ipm(10 + if2 + 2,mid(i))
674 ifunc2(i)= ipm(10 + if3 + 2,mid(i))
675 ifunc3(i)= ipm(10 + if4 + 2,mid(i))
676 iecrou(i)= nint(uparam(iadbuf + i13 + 2))
677 IF (iecrou(i) > 0) ifunc(i) = -1 ! But why????
678 ak(i) = uparam(iadbuf + i1 + 2)
679 b(i) = uparam(iadbuf + i2 + 2)
680 d(i) = uparam(iadbuf + i3 + 2)
681 ee(i) = uparam(iadbuf + i4 + 2)
682 gf3(i) = uparam(iadbuf + i5 + 2)
683 ff(i) = uparam(iadbuf + i6 + 2)
684 lscale(i)= uparam(iadbuf + i7 + 2)
685 dmn(i) = uparam(iadbuf + i8 + 2)
686 dmx(i) = uparam(iadbuf + i9 + 2)
687 ENDDO
688C
689 kk = 1 + numelr * anim_fe(11)
690 IF (nuvar >= 2) coord_old => uvar(2,1:nel)
691 CALL redef3(python,
692 1 fy, yk, dy, fyep,
693 2 dyold, dpy, tf, npf,
694 3 yc, off, e6(1,2), dpy2,
695 4 anim(kk), anim_fe(12),posy,
696 5 xl0, dmn, dmx, dv,
697 6 ff, lscale, ee, gf3,
698 7 ifunc3, yieldy, aldp, ak,
699 8 b, d, iecrou, ifunc,
700 9 ifv, ifunc2, epla, coord_old,
701 a nel, nft, stf, sanin, dt1,
702 b iresp, impl_s, idyna, snpc,
703 c fx_max=fx_max)
704C
705 DO i=1,nel
706 iadbuf = ipm(7,mid(i)) - 1
707 cc = uparam(iadbuf + nupar + 4)
708 cn = uparam(iadbuf + nupar + 10)
709 xa = uparam(iadbuf + nupar + 16)
710 xb = uparam(iadbuf + nupar + 22)
711 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
712 IF (ifail2(i) == 0) THEN
713 xa = one
714 xb = two
715 IF (dy(i) > zero) THEN
716 dlim = dy(i) / dmx(i)
717 ELSE
718 dlim = dy(i) / dmn(i)
719 ENDIF
720 ELSE
721 vfail = cc * (abs(dv(i)/vrt(i)))**cn
722 IF (ifail2(i) == 1) THEN
723 IF (dy(i) > zero) THEN
724 dlim = dy(i) / (dmx(i) + vfail)
725 ELSE
726 dlim = dy(i) / (dmn(i) - vfail)
727 ENDIF
728 ELSEIF (ifail2(i) == 2) THEN
729 IF (fy(i) > zero) THEN
730 dlim = fy(i) / (dmx(i) + vfail)
731 ELSE
732 dlim = fy(i) / (dmn(i) - vfail)
733 ENDIF
734 ELSEIF (ifail2(i) == 3) THEN
735 dlim = max(zero,e6(i,2)) / (dmx(i) + vfail)
736 ENDIF
737 ENDIF
738 IF (ifail(i) == 0) THEN
739! Uniaxial rupture
740 IF ((xa*dlim) > one) THEN
741 off(i) = zero
742 nindx = nindx + 1
743 indx(nindx) = i
744 idel7nok = 1
745 ENDIF
746 ELSE
747! Multiaxial rupture
748 crit(i)= crit(i) + xa *(dlim/xl0(i))**xb
749 ENDIF
750 ENDIF
751 ENDDO
752C
753 DO i=1,nel
754 iadbuf = ipm(7,mid(i)) - 1
755 ifunc(i) = 0
756 ifv(i) = ipm(10 + if2 + 3,mid(i))
757 ifunc2(i)= ipm(10 + if3 + 3,mid(i))
758 ifunc3(i)= ipm(10 + if4 + 3,mid(i))
759 iecrou(i)= nint(uparam(iadbuf + i13 + 3))
760 IF (iecrou(i) > 0) ifunc(i) = -1
761 ak(i) = uparam(iadbuf + i1 + 3)
762 b(i) = uparam(iadbuf + i2 + 3)
763 d(i) = uparam(iadbuf + i3 + 3)
764 ee(i) = uparam(iadbuf + i4 + 3)
765 gf3(i) = uparam(iadbuf + i5 + 3)
766 ff(i) = uparam(iadbuf + i6 + 3)
767 lscale(i)= uparam(iadbuf + i7 + 3)
768 dmn(i) = uparam(iadbuf + i8 + 3)
769 dmx(i) = uparam(iadbuf + i9 + 3)
770 ENDDO
771C
772 kk = 1 + numelr * (anim_fe(11)+anim_fe(12))
773 IF (nuvar >= 3) coord_old => uvar(3,1:nel)
774
775 CALL redef3(python,
776 1 fz, zk, dz, fzep,
777 2 dzold, dpz, tf, npf,
778 3 zc, off, e6(1,3), dpz2,
779 4 anim(kk), anim_fe(13),posz,
780 5 xl0, dmn, dmx, dv,
781 6 ff, lscale, ee, gf3,
782 7 ifunc3, yieldz, aldp, ak,
783 8 b, d, iecrou, ifunc,
784 9 ifv, ifunc2, epla, coord_old,
785 a nel, nft, stf, sanin, dt1,
786 b iresp, impl_s, idyna, snpc,
787 c fx_max=fx_max)
788C
789 DO i=1,nel
790 iadbuf = ipm(7,mid(i)) - 1
791 cc = uparam(iadbuf + nupar + 5)
792 cn = uparam(iadbuf + nupar + 11)
793 xa = uparam(iadbuf + nupar + 17)
794 xb = uparam(iadbuf + nupar + 23)
795 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
796 IF (ifail2(i) == 0) THEN
797 xa = one
798 xb = two
799 IF (dz(i) > zero)THEN
800 dlim = dz(i) / dmx(i)
801 ELSE
802 dlim = dz(i) / dmn(i)
803 ENDIF
804 ELSE
805 vfail = cc * (abs(dv(i)/vrt(i)))**cn
806 IF (ifail2(i) == 1) THEN
807 IF (dz(i) > zero) THEN
808 dlim = dz(i) / (dmx(i) + vfail)
809 ELSE
810 dlim = dz(i) / (dmn(i) - vfail)
811 ENDIF
812 ELSEIF (ifail2(i) == 2) THEN
813 IF (fz(i) > zero) THEN
814 dlim = fz(i) / (dmx(i) + vfail)
815 ELSE
816 dlim = fz(i) / (dmn(i) - vfail)
817 ENDIF
818 ELSEIF (ifail2(i) == 3) THEN
819 dlim = max(zero,e6(i,3)) / (dmx(i) + vfail)
820 ENDIF
821 ENDIF
822 IF (ifail(i) == 0) THEN
823! Uniaxial rupture
824 IF ((xa*dlim) > one) THEN
825 off(i) = zero
826 nindx = nindx + 1
827 indx(nindx) = i
828 idel7nok = 1
829 ENDIF
830 ELSE
831! Multiaxial rupture
832 crit(i)= crit(i) + xa *(dlim/xl0(i))**xb
833 ENDIF
834 ENDIF
835 ENDDO
836C---------------------
837C ROTATIONS
838C---------------------
839 DO i=1,nel
840 iadbuf= ipm(7,mid(i)) - 1
841 xin(i)= uiner(i)
842 xk(i) = uparam(iadbuf + i11 + 4)
843 xc(i) = uparam(iadbuf + i12 + 4)
844 yk(i) = uparam(iadbuf + i11 + 5)
845 yc(i) = uparam(iadbuf + i12 + 5)
846 zk(i) = uparam(iadbuf + i11 + 6)
847 zc(i) = uparam(iadbuf + i12 + 6)
848 hx(i) = uparam(iadbuf + i14 + 4)
849 hy(i) = uparam(iadbuf + i14 + 5)
850 hz(i) = uparam(iadbuf + i14 + 6)
851C
852 xhr(i)= max(hx(i),hy(i),hz(i))
853
854 xkr(i)= max(xk(i)*uparam(iadbuf + i1 + 4),
855 . yk(i)*uparam(iadbuf + i1 + 5),
856 . zk(i)*uparam(iadbuf + i1 + 6))+xkr(i)
857 xcr(i)= max(xc(i),yc(i),zc(i)) + xhr(i) +xcr(i)+xh(i)
858 ENDDO
859C
860 DO i=1,nel
861 dxold(i)=rx(i)
862 dyold(i)=ry(i)
863 dzold(i)=rz(i)
864 ENDDO
865!
866 IF ( inispri /= 0 .AND. tt == zero) THEN
867 DO i=1,nel
868 dxold(i)=rx0(i)
869 dyold(i)=ry0(i)
870 dzold(i)=rz0(i)
871 ENDDO
872 ENDIF
873C
874 DO i=1,nel
875 x21 = (rx2(i)-rx1(i))*dt1
876 y21 = (ry2(i)-ry1(i))*dt1
877 z21 = (rz2(i)-rz1(i))*dt1
878 rx(i) = dxold(i) + x21*exx2(i)+y21*eyx2(i)+z21*ezx2(i)
879 ry(i) = dyold(i) + x21*exy2(i)+y21*eyy2(i)+z21*ezy2(i)
880 rz(i) = dzold(i) + x21*exz2(i)+y21*eyz2(i)+z21*ezz2(i)
881 ENDDO
882C-------------------------------
883 DO i=1,nel
884 iadbuf = ipm(7,mid(i)) - 1
885 ifunc(i) = 0
886 ifv(i) = ipm(10 + if2 + 4,mid(i))
887 ifunc2(i)= ipm(10 + if3 + 4,mid(i))
888 ifunc3(i)= ipm(10 + if4 + 4,mid(i))
889 iecrou(i)= nint(uparam(iadbuf + i13 + 4))
890 IF (iecrou(i) > 0) ifunc(i) = -1
891 ak(i) = uparam(iadbuf + i1 + 4)
892 b(i) = uparam(iadbuf + i2 + 4)
893 d(i) = uparam(iadbuf + i3 + 4)
894 ee(i) = uparam(iadbuf + i4 + 4)
895 gf3(i) = uparam(iadbuf + i5 + 4)
896 ff(i) = uparam(iadbuf + i6 + 4)
897 lscale(i)= uparam(iadbuf + i7 + 4)
898 dmn(i) = uparam(iadbuf + i8 + 4)
899 dmx(i) = uparam(iadbuf + i9 + 4)
900 ENDDO
901 IF (nuvar >= 4) coord_old => uvar(4,1:nel)
902
903 CALL redef3(python,
904 1 xmom, xk, rx, xmep,
905 2 dxold, rpx, tf, npf,
906 3 xc, off, e6(1,4), rpx2,
907 4 anim, 0, posxx,
908 5 xl0, dmn, dmx, dv,
909 6 ff, lscale, ee, gf3,
910 7 ifunc3, yieldx2, aldp, ak,
911 8 b, d, iecrou, ifunc,
912 9 ifv, ifunc2, epla, coord_old,
913 a nel, nft, stf, sanin, dt1,
914 b iresp, impl_s, idyna, snpc,
915 c fx_max=mx_max)
916C
917 DO i=1,nel
918 iadbuf= ipm(7,mid(i)) - 1
919 cc = uparam(iadbuf + nupar + 6)
920 cn = uparam(iadbuf + nupar + 12)
921 xa = uparam(iadbuf + nupar + 18)
922 xb = uparam(iadbuf + nupar + 24)
923 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
924 IF (ifail2(i) == 0) THEN
925 xa = one
926 xb = two
927 IF (rx(i) > zero) THEN
928 dlim = rx(i) / dmx(i)
929 ELSE
930 dlim = rx(i) / dmn(i)
931 ENDIF
932 ELSE
933 vfail = cc * (abs(dv(i)/vrr(i)))**cn
934 IF (ifail2(i) == 1) THEN
935 IF (rx(i) > zero) THEN
936 dlim = rx(i) / (dmx(i) + vfail)
937 ELSE
938 dlim = rx(i) / (dmn(i) - vfail)
939 ENDIF
940 ELSEIF (ifail2(i) == 2) THEN
941 IF(xmom(i)>0.)THEN
942 dlim = xmom(i)/(dmx(i) + vfail)
943 ELSE
944 dlim = xmom(i)/(dmn(i) - vfail)
945 ENDIF
946 ELSEIF (ifail2(i) == 3) THEN
947 dlim = max(zero,e6(i,4)) / (dmx(i) + vfail)
948 ENDIF
949 ENDIF
950 IF (ifail(i) == 0) THEN
951! Uniaxial rupture
952 IF ((xa*dlim) > one) THEN
953 off(i)= zero
954 nindx = nindx + 1
955 indx(nindx) = i
956 idel7nok = 1
957 ENDIF
958 ELSE
959! Multiaxial rupture
960 crit(i)= crit(i) + xa *(dlim/xl0(i))**xb
961 ENDIF
962 ENDIF
963 ENDDO
964C
965 DO i=1,nel
966 iadbuf = ipm(7,mid(i)) - 1
967 ifunc(i) = 0
968 ifv(i) = ipm(10 + if2 + 5,mid(i))
969 ifunc2(i)= ipm(10 + if3 + 5,mid(i))
970 ifunc3(i)= ipm(10 + if4 + 5,mid(i))
971 iecrou(i)= nint(uparam(iadbuf + i13 + 5))
972 IF (iecrou(i) > 0) ifunc(i) = -1
973 ak(i) = uparam(iadbuf + i1 + 5)
974 b(i) = uparam(iadbuf + i2 + 5)
975 d(i) = uparam(iadbuf + i3 + 5)
976 ee(i) = uparam(iadbuf + i4 + 5)
977 gf3(i) = uparam(iadbuf + i5 + 5)
978 ff(i) = uparam(iadbuf + i6 + 5)
979 lscale(i)= uparam(iadbuf + i7 + 5)
980 dmn(i) = uparam(iadbuf + i8 + 5)
981 dmx(i) = uparam(iadbuf + i9 + 5)
982 ENDDO
983 IF (nuvar >= 5) coord_old => uvar(5,1:nel)
984
985 CALL redef3(python,
986 1 ymom, yk, ry, ymep,
987 2 dyold, rpy, tf, npf,
988 3 yc, off, e6(1,5), rpy2,
989 4 anim, 0, posyy,
990 5 xl0, dmn, dmx, dv,
991 6 ff, lscale, ee, gf3,
992 7 ifunc3, yieldy2, aldp, ak,
993 8 b, d, iecrou, ifunc,
994 9 ifv, ifunc2, epla, coord_old,
995 a nel, nft, stf, sanin, dt1,
996 b iresp, impl_s, idyna, snpc,
997 c fx_max=mx_max)
998C
999 DO i=1,nel
1000 iadbuf= ipm(7,mid(i)) - 1
1001 cc = uparam(iadbuf + nupar + 7)
1002 cn = uparam(iadbuf + nupar + 13)
1003 xa = uparam(iadbuf + nupar + 19)
1004 xb = uparam(iadbuf + nupar + 25)
1005 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
1006 IF (ifail2(i) == 0) THEN
1007 xa = one
1008 xb = two
1009 IF (ry(i) > zero) THEN
1010 dlim = ry(i) / dmx(i)
1011 ELSE
1012 dlim = ry(i) / dmn(i)
1013 ENDIF
1014 ELSE
1015 vfail = cc * (abs(dv(i)/vrr(i)))**cn
1016 IF (ifail2(i) == 1) THEN
1017 IF (ry(i) > zero) THEN
1018 dlim = ry(i) / (dmx(i) + vfail)
1019 ELSE
1020 dlim = ry(i) / (dmn(i) - vfail)
1021 ENDIF
1022 ELSEIF (ifail2(i) == 2) THEN
1023 IF (ymom(i) > zero)THEN
1024 dlim = ymom(i)/(dmx(i) + vfail)
1025 ELSE
1026 dlim = ymom(i)/(dmn(i) - vfail)
1027 ENDIF
1028 ELSEIF (ifail2(i) == 3) THEN
1029 dlim = max(zero,e6(i,5)) / (dmx(i) + vfail)
1030 ENDIF
1031 ENDIF
1032 IF (ifail(i) == 0) THEN
1033! Uniaxial rupture
1034 IF ((xa*dlim) > 1) THEN
1035 off(i) = zero
1036 nindx = nindx + 1
1037 indx(nindx) = i
1038 idel7nok = 1
1039 ENDIF
1040 ELSE
1041! Multiaxial rupture
1042 crit(i)= crit(i) + xa *(dlim/xl0(i))**xb
1043 ENDIF
1044 ENDIF
1045 ENDDO
1046C
1047 DO i=1,nel
1048 iadbuf = ipm(7,mid(i)) - 1
1049 ifunc(i) = 0
1050 ifv(i) = ipm(10 + if2 + 6,mid(i))
1051 ifunc2(i)= ipm(10 + if3 + 6,mid(i))
1052 ifunc3(i)= ipm(10 + if4 + 6,mid(i))
1053 iecrou(i)= nint(uparam(iadbuf + i13 + 6))
1054 IF (iecrou(i) > 0) ifunc(i) = -1
1055 ak(i) = uparam(iadbuf + i1 + 6)
1056 b(i) = uparam(iadbuf + i2 + 6)
1057 d(i) = uparam(iadbuf + i3 + 6)
1058 ee(i) = uparam(iadbuf + i4 + 6)
1059 gf3(i) = uparam(iadbuf + i5 + 6)
1060 ff(i) = uparam(iadbuf + i6 + 6)
1061 lscale(i)= uparam(iadbuf + i7 + 6)
1062 dmn(i) = uparam(iadbuf + i8 + 6)
1063 dmx(i) = uparam(iadbuf + i9 + 6)
1064 ENDDO
1065 IF (nuvar >= 6) coord_old => uvar(6,1:nel)
1066 CALL redef3(python,
1067 1 zmom, zk, rz, zmep,
1068 2 dzold, rpz, tf, npf,
1069 3 zc, off, e6(1,6), rpz2,
1070 4 anim, 0, poszz,
1071 5 xl0, dmn, dmx, dv,
1072 6 ff, lscale, ee, gf3,
1073 7 ifunc3, yieldz2, aldp, ak,
1074 8 b, d, iecrou, ifunc,
1075 9 ifv, ifunc2, epla, coord_old,
1076 a nel, nft, stf, sanin, dt1,
1077 b iresp, impl_s, idyna, snpc,
1078 c fx_max=mx_max)
1079C
1080 DO i=1,nel
1081 iadbuf= ipm(7,mid(i)) - 1
1082 cc = uparam(iadbuf + nupar + 8)
1083 cn = uparam(iadbuf + nupar + 14)
1084 xa = uparam(iadbuf + nupar + 20)
1085 xb = uparam(iadbuf + nupar + 26)
1086 IF (off(i) == one .AND. dmx(i) /= zero .AND. dmn(i) /= zero) THEN
1087 IF (ifail2(i) == 0) THEN
1088 xa = one
1089 xb = two
1090 IF (rz(i) > zero) THEN
1091 dlim = rz(i) / dmx(i)
1092 ELSE
1093 dlim = rz(i) / dmn(i)
1094 ENDIF
1095 ELSE
1096 vfail = cc * (abs(dv(i)/vrr(i)))**cn
1097 IF (ifail2(i) == 1) THEN
1098 IF (rz(i) > zero)THEN
1099 dlim = rz(i) / (dmx(i) + vfail)
1100 ELSE
1101 dlim = rz(i) / (dmn(i) - vfail)
1102 ENDIF
1103 ELSEIF (ifail2(i) == 2) THEN
1104 IF (zmom(i) > zero)THEN
1105 dlim = zmom(i)/(dmx(i) + vfail)
1106 ELSE
1107 dlim = zmom(i)/(dmn(i) - vfail)
1108 ENDIF
1109 ELSEIF (ifail2(i) == 3) THEN
1110 dlim = max(zero,e6(i,6)) / (dmx(i) + vfail)
1111 ENDIF
1112 ENDIF
1113 IF (ifail(i) == 0) THEN
1114! Uniaxial rupture
1115 IF ((xa*dlim) > 1) THEN
1116 off(i) = zero
1117 nindx = nindx + 1
1118 indx(nindx) = i
1119 idel7nok = 1
1120 ENDIF
1121 ELSE
1122! Multiaxial rupture
1123 crit(i)= crit(i) + xa *(dlim/xl0(i))**xb
1124 ENDIF
1125 ENDIF
1126 ENDDO
1127C
1128 DO i=1,nel
1129 e(i) = e6(i,1)+e6(i,2)+e6(i,3)+e6(i,4)+e6(i,5)+e6(i,6)
1130 ENDDO
1131C-------------------------------
1132C COUPLED FAILURE
1133C-------------------------------
1134 DO i=1,nel
1135 iadbuf = ipm(7,mid(i)) - 1
1136 israte = nint(uparam(iadbuf + nupar + 27))
1137C---- smoothing factor alpha = 2PI*fc*dt/(2PI*fc*dt+1) ---
1138 asrate = (2*pi*uparam(iadbuf + nupar + 28)*dt1)/(one+2*pi*uparam(iadbuf + nupar + 28)*dt1)
1139 IF (israte /= 0) THEN
1140 crit(i) = asrate*crit(i) + (one - asrate)*crit_new(i)
1141 crit_new(i) = crit(i)
1142 ENDIF
1143 IF (off(i) == one .AND. ifail(i) == 1) THEN
1144 IF (crit(i) > one) THEN
1145 off(i)=zero
1146 nindx = nindx + 1
1147 indx(nindx) = i
1148 idel7nok = 1
1149 ENDIF
1150 ENDIF
1151 ENDDO
1152C
1153 DO j=1,nindx
1154 i = indx(j)
1155#include "lockon.inc"
1156 WRITE(iout, 1000) ngl(i)
1157 WRITE(istdo,1100) ngl(i),tt
1158#include "lockoff.inc"
1159 ENDDO
1160C-------------------------------
1161C COUPLED PLASTICITY
1162C-------------------------------
1163 CALL repla3(
1164 1 xk, rpx, tf, npf,
1165 2 iecrou, ifunc, ifv, epla,
1166 3 nel)
1167 CALL repla3(
1168 1 yk, rpy, tf, npf,
1169 2 iecrou, ifunc, ifv, epla,
1170 3 nel)
1171 CALL repla3(
1172 1 zk, rpz, tf, npf,
1173 2 iecrou, ifunc, ifv, epla,
1174 3 nel)
1175C
1176 DO i=1,nel
1177 iadbuf= ipm(7,mid(i)) - 1
1178 xk(i)=uparam(iadbuf + i11 + 1)
1179 yk(i)=uparam(iadbuf + i11 + 2)
1180 zk(i)=uparam(iadbuf + i11 + 3)
1181 ENDDO
1182C
1183 CALL repla3(
1184 1 xk, dpx, tf, npf,
1185 2 iecrou, ifunc, ifv, epla,
1186 3 nel)
1187 CALL repla3(
1188 1 yk, dpy, tf, npf,
1189 2 iecrou, ifunc, ifv, epla,
1190 3 nel)
1191 CALL repla3(
1192 1 zk, dpz, tf, npf,
1193 2 iecrou, ifunc, ifv, epla,
1194 3 nel)
1195 DO i=1,nel
1196C XM(I)=XM(I)*XL0(I)
1197 xkm(i)=xkm(i)/xl0(i)
1198 xcm(i)=xcm(i)/xl0(i)
1199C XIN(I)=XIN(I)*XL0(I)
1200 xkr(i)=xkr(i)/xl0(i)
1201 xcr(i)=xcr(i)/xl0(i)
1202 ENDDO
1203C---
1204 1000 FORMAT(1x,'-- RUPTURE OF SPRING ELEMENT NUMBER ',i10)
1205 1100 FORMAT(1x,'-- RUPTURE OF SPRING ELEMENT :',i10,' AT TIME :',g11.4)
1206C---
1207 RETURN
1208 END
my_real
#define my_real
Definition cppsort.cpp:32
norm
norm(diag(diag(diag(inv(mat))) -id.SOL), 2) % destroy mumps instance id.JOB
max
#define max(a, b)
Definition macros.h:21
material_flow
subroutine material_flow(dfs, dfs_old, aldp, slipring_strand, xk, off, al0, al02, lmin, update_flag, ring_slip, slipring_id, xl0, dl, dlold, exdp, eydp, ezdp, x1dp, x2dp, x3dp, adher, nc1, nc2, nc3, flag_slipring_update, add_node1, add_node2, vx1, vy1, vz1, vx2, vy2, vz2, vx3, vy3, vz3, xc, retractor_id, flag_retractor_update, sensor_tab, al0dp, fr_id, ddf, fx, fx2, aldp2, al0dp2, ex2dp, ey2dp, ez2dp, xl02, xk2, compt, index2, nsensor)
Definition material_flow.F:52
seatbelt_mod
Definition seatbelt_mod.F:55
seatbelt_mod::retractor
type(retractor_struct), dimension(:), allocatable retractor
Definition seatbelt_mod.F:155
seatbelt_mod::slipring
type(slipring_struct), dimension(:), allocatable slipring
Definition seatbelt_mod.F:110
r23l114def3
subroutine r23l114def3(python, skew, ipm, igeo, mid, pid, geo, uparam, 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, ngl, crit_new, x0_err, aldp, yieldx, yieldy, yieldz, yieldx2, yieldy2, yieldz2, 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, mass, dx0, dy0, dz0, rx0, ry0, rz0, slipring_strand, dfs, ring_slip, x02, lmin, slipring_id, update_flag, retractor_id, add_node1, add_node2, nc1, nc2, nc3, x1dp, x2dp, x3dp, vx3, vy3, vz3, flag_slipring_update, flag_retractor_update, sensor_tab, uiner, fr_id, fram_factor, eps_old, fx_b2, dpx_b2, yieldx_b2, xx_old_b2, fxep_b2, posx_b2, eps_old_b2, nft, nsensor, stf, sanin, iresp, snpc)
Definition r23l114def3.F:78
repla3
subroutine repla3(xk, dpx, tf, npf, iecrou, ifunc, ifv, epla, nel)
Definition repla3.F:39