OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
ruser36.F File Reference
#include "implicit_f.inc"
#include "mvsiz_p.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine ruser36 (nel, iprop, uvar, nuvar, fr_wave, fx, fy, fz, xmom, ymom, zmom, e, off, stifm, stifr, viscm, viscr, mass, xiner, dt, xl, vx, ry1, rz1, rx, ry2, rz2)

Function/Subroutine Documentation

◆ ruser36()

subroutine ruser36 ( integer nel,
integer iprop,
uvar,
integer nuvar,
fr_wave,
fx,
fy,
fz,
xmom,
ymom,
zmom,
e,
off,
stifm,
stifr,
viscm,
viscr,
mass,
xiner,
dt,
xl,
vx,
ry1,
rz1,
rx,
ry2,
rz2 )

Definition at line 34 of file ruser36.F.

41C-------------------------------------------------------------------------
42C I m p l i c i t T y p e s
43C-----------------------------------------------
44#include "implicit_f.inc"
45#include "mvsiz_p.inc"
46C----------------------------------------------------------
47C d u m m y a r g u m e n t s a n d f u n c t i o n
48C----------------------------------------------------------
49 INTEGER NEL,NUVAR,IPROP,
50 . GET_U_PNU,GET_U_PID,GET_U_MID,GET_U_MNU,
51 . KFUNC,KMAT,KPROP
52 my_real
53 . uvar(nuvar,*),dt ,
54 . fx(*), fy(*), fz(*), e(*), vx(*),mass(*) ,xiner(*),
55 . ry1(*), rz1(*), off(*), xmom(*), ymom(*),
56 . zmom(*), rx(*), ry2(*), rz2(*),xl(*),
57 . stifm(*) ,stifr(*) , viscm(*) ,viscr(*) ,fr_wave(*) ,
58 . get_u_mat, get_u_geo, get_u_func
59 EXTERNAL get_u_mnu,get_u_pnu,get_u_mid,get_u_pid,
60 . get_u_mat,get_u_geo, get_u_func
61 parameter(kfunc=29)
62 parameter(kmat=31)
63 parameter(kprop=33)
64C-----------------------------------------------
65C L O C A L V A R I A B L E S
66C-----------------------------------------------
67 INTEGER I,KCST,N0,
68 . IMAT1,IPROP1,
69 . IMAT2,IPROP2
70 INTEGER J,N,NINDX,NMAX,INDEX(MVSIZ),IFUNC
71 my_real
72 . fac,rho,area,ixx,iyy,izz,imyz,young,g,
73 . area1,ixx1,iyy1,izz1,rho1,young1,g1,
74 . area2,ixx2,iyy2,izz2,rho2,young2,g2
75 my_real r_y1,r_y2,r_z1,r_z2,r_y,r_z,r_x,facx,facy,facz,
76 . ay1,ay2,az1,az2,ay,az,by1,by2,bz1,bz2,by,bz,g3,
77 . cx1,cx2,cx,area_y(mvsiz),area_z(mvsiz),
78 . yld(mvsiz),dre,drg,dpla_i,dpla_j(mvsiz),dfe,dfg,
79 . err,f,df,yld_i,c1,fn,fm,y1,y2,
80 . pnx,pnx2,pny,pny2,pnz,pnz2,pmx,pmx2,pmy,pmy2,pmz,pmz2,
81 . nx2,ny2,nz2,mx2,my2,mz2,tempy,tempz,ktran,krot,
82 . nx(mvsiz),ny(mvsiz),nz(mvsiz),youngm(mvsiz),gm(mvsiz),
83 . mx(mvsiz),my(mvsiz),mz(mvsiz),svm(mvsiz),h(mvsiz),
84 . epsvxx(mvsiz),epsvxy(mvsiz),epsvxz(mvsiz),
85 . epscbx(mvsiz),epscby(mvsiz),epscbz(mvsiz),
86 . signxx(mvsiz),signxy(mvsiz),signxz(mvsiz),
87 . momnxx(mvsiz),momnyy(mvsiz),momnzz(mvsiz),
88 . gama(mvsiz),gamac(mvsiz),xl2(mvsiz),devol,
89 . pc1,dc1,pr1,ps1,dc2,pr2,ps2,dc,pr,ps,ueq,dtemp,
90 . sig01,sig02, sig0,hpla,h1,h2,m1,m2,m,dsigy,cc1,
91 . fac1,fac2
92 DATA nmax/10/
93C-----------------------------------------------
94 kcst = 0
95C
96 iprop1 = get_u_pnu(1,iprop,kprop)
97 area1 = get_u_geo(2,iprop1)
98 ixx1 = get_u_geo(3,iprop1)
99 iyy1 = get_u_geo(4,iprop1)
100 izz1 = get_u_geo(5,iprop1)
101 r_y1= get_u_geo(6,iprop1)
102 r_z1= get_u_geo(7,iprop1)
103 iprop2 = get_u_pnu(2,iprop,kprop)
104 area2 = get_u_geo(2,iprop2)
105 ixx2 = get_u_geo(3,iprop2)
106 iyy2 = get_u_geo(4,iprop2)
107 izz2 = get_u_geo(5,iprop2)
108 imat1 = get_u_pnu(1,iprop1,kmat)
109 ifunc = get_u_mnu(1,imat1,kfunc)
110 young1 = get_u_mat(7,imat1)
111 g1 = get_u_mat(6,imat1)
112 rho1 = get_u_mat(0,imat1)
113 ay1 = get_u_mat(8,imat1)
114 az1 = get_u_mat(9,imat1)
115 by1 = get_u_mat(10,imat1)
116 bz1 = get_u_mat(11,imat1)
117 cx1 = get_u_mat(12,imat1)
118 dc1 = get_u_mat(13,imat1)
119 pr1 = get_u_mat(14,imat1)
120 ps1 = get_u_mat(15,imat1)
121 sig01 = get_u_mat(16,imat1)
122 h1 = get_u_mat(17,imat1)
123 m1 = get_u_mat(18,imat1)
124 fac1 = get_u_mat(19,imat1)
125 imat2 = get_u_pnu(1,iprop2,kmat)
126 young2 = get_u_mat(7,imat2)
127 g2 = get_u_mat(6,imat2)
128 rho2 = get_u_mat(0,imat2)
129 ay2 = get_u_mat(8,imat2)
130 az2 = get_u_mat(9,imat2)
131 by2 = get_u_mat(10,imat2)
132 bz2 = get_u_mat(11,imat2)
133 cx2 = get_u_mat(12,imat2)
134 dc2 = get_u_mat(13,imat2)
135 pr2 = get_u_mat(14,imat2)
136 ps2 = get_u_mat(15,imat2)
137 sig02 = get_u_mat(16,imat2)
138 h2 = get_u_mat(17,imat2)
139 m2 = get_u_mat(18,imat2)
140 r_y2 = get_u_geo(6,iprop2)
141 r_z2 = get_u_geo(7,iprop2)
142 rho = half*(rho1+rho2)
143 area = half*(area1+area2)
144 fac = rho*area
145 ixx = half*(ixx1+ixx2)
146 iyy = half*(iyy1+iyy2)
147 izz = half*(izz1+izz2)
148 r_y = half*(r_y1+r_y2)
149 r_z = half*(r_z1+r_z2)
150 r_x = half*(r_y+r_z)
151 cx = half*(cx1+cx2)
152 ay = half*(ay1+ay2)
153 az = half*(az1+az2)
154 by = half*(by1+by2)
155 bz = half*(bz1+bz2)
156 imyz = max(iyy,izz)
157 young = half*(young1+young2)
158 g = half*(g1+g2)
159 facx = r_x/max(ixx,em20)
160 facy = r_y/max(iyy,em20)
161 facz = r_z/max(izz,em20)
162 tempy = g*area/max(twelve*young*iyy,em20)
163 tempz = g*area/max(twelve*young*izz,em20)
164 g3 = three*g
165 n0 = 11
166 dc = half*(dc1+dc2)
167 pr = half*(pr1+pr2)
168 ps = half*(ps1+ps2)
169 sig0 = half*(sig01+sig02)
170 hpla = half*(h1+h2)
171 m = half*(m1+m2)
172C
173C---- UVAR(N0,I) -> PLA , N0+1 -> SY0 , N0+2 -> EXX_P , N0+3 -> EYY_P
174C---- UVAR(N0+4,I) -> ENDOMMAGEMENT
175C
176 DO i=1,nel
177 youngm(i) = young*(1-uvar(n0+4,i))
178 gm(i) = g*(1-uvar(n0+4,i))
179C== WE FIRST CALCULATE THE VELOCITIES OF GENERALIZED DEFORMATIONS ====
180 dtemp = one/max(xl(i),em20)
181 epsvxx(i) = vx(i)*dtemp
182 epsvxy(i) = -half*(rz1(i) + rz2(i))
183 epsvxz(i) = half*(ry1(i) + ry2(i))
184 epscbx(i) = r_x*rx(i)*dtemp
185 epscby(i) = r_y*(ry2(i) - ry1(i))*dtemp
186 epscbz(i) = r_z*(rz2(i) - rz1(i))*dtemp
187C== ON CALCULE CONTRAINTES GENERALISEES ====
188 xl2(i) = xl(i)*xl(i)
189 area_y(i) = area/(one +tempy*xl2(i))
190 area_z(i) = area/(one +tempz*xl2(i))
191 signxx(i) = fx(i)/max(area,em20)+youngm(i)*epsvxx(i)*dt
192 signxy(i) = fy(i)/max(area_y(i),em20)+gm(i)*epsvxy(i)*dt
193 signxz(i) = fz(i)/max(area_z(i),em20)+gm(i)*epsvxz(i)*dt
194 momnxx(i) = xmom(i)*facx+gm(i)*epscbx(i)*dt
195 momnyy(i) = ymom(i)*facy+youngm(i)*epscby(i)*dt
196 momnzz(i) = zmom(i)*facz+youngm(i)*epscbz(i)*dt
197 ENDDO
198C
199C--- UVAR(11,I)-> DEF PLASTIQUE EQUIVALENTE -
200C--- UVAR(N0,I)-> PLA, N0+1->SY0; N0+2->EXX_P; N0+3->EYY_P;
201C--- UVAR(N0+3,I)-> COURBURE PLASTIQUE Y,Z UVAR(N0+2,I)-> COURBURE PLASTIQUE X -
202C
203 DO i=1,nel
204 IF (ifunc /= 0) THEN
205 y1= fac1*get_u_func(ifunc,uvar(n0,i),y2)
206 y2= fac1*y2
207 yld(i) = max(y1,em20)
208 m=1
209 ELSE
210 y1 = sig0+hpla*(uvar(n0,i))**m
211 y2 = hpla
212 ENDIF
213 yld(i) = max(y1,em20)
214 h(i) = max(y2,zero)
215 uvar(n0+1,i)=min(yld(i),uvar(n0+1,i))
216C
217 pc1 = abs(six*youngm(i)*(uvar(n0+2,i)/uvar(n0+1,i)))
218 c1=one-fourth*exp(-pc1)
219 gamac(i)=cx*three_over_4/max(c1,em20)
220C
221 pc1 = six*youngm(i)*(uvar(n0+3,i)/uvar(n0+1,i))
222 cc1 = (three*pi/sixteen)**2
223 c1=one - (one - cc1)*exp(-pc1)
224 gama(i)=cc1/max(c1,em20)
225C
226 nx(i)=signxx(i)
227 ny(i)=ay*signxy(i)
228 nz(i)=az*signxz(i)
229 mx(i)=gamac(i)*momnxx(i)
230 my(i)=by*gama(i)*momnyy(i)
231 mz(i)=bz*gama(i)*momnzz(i)
232 svm(i)=nx(i)*nx(i)+three*(ny(i)*ny(i)+nz(i)*nz(i)+mx(i)*mx(i))+
233 1 my(i)*my(i)+mz(i)*mz(i)
234 svm(i)=sqrt(svm(i))
235 ENDDO
236 fac2 = g/max(em20,young)
237 DO i=1,nel
238C
239C TIME STEP
240C
241 dtemp = one/max(xl(i),em20)
242 ktran = max(area,fac2*area_y(i),fac2*area_z(i))*dtemp
243 krot = four *max(iyy*dtemp,izz*dtemp)
244 stifm(i) = youngm(i)*ktran
245 stifr(i) = max( gm(i)*ixx*dtemp,youngm(i)*krot)
246 viscm(i) = zero
247 viscr(i) = zero
248 mass(i) = xl(i)*fac
249 xiner(i) = xl(i)*rho*max(ixx,imyz+area*xl2(i)/twelve)
250 ENDDO
251C
252 nindx=0
253 DO i=1,nel
254 IF (svm(i) > yld(i) .AND. off(i) == one) THEN
255 nindx=nindx+1
256 index(nindx)=i
257 ENDIF
258 ENDDO
259C
260 IF (nindx /= 0) THEN
261 DO j=1,nindx
262 i=index(j)
263 dpla_j(i)=(svm(i)-yld(i))/(g3+h(i))
264 ENDDO
265C
266C WARNING, THE VERSION OF THE SPRING MODEL (Q.ZENG/A.COMBESCURE/C.BENOIT)
267C NE CONSIDERE PAS D'ERREUR CUMULEE SUR PLUSIEURS ELEMENTS...
268C
269C --> EXISTE VERSION AVEC TEST MULTI-ELEMENTS (ERR_MULTI/RUSER29.F)
270C
271C
272C
273 n=0
274 err =2*em4
275 DO WHILE (err > em4 .AND. n < nmax)
276 DO j=1,nindx
277 i=index(j)
278 dpla_i=dpla_j(i)
279 IF (ifunc /= 0) THEN
280 yld_i =yld(i)+h(i)*dpla_i
281 ELSEIF (uvar(n0,i) > em20) THEN
282 yld_i =yld(i)+h(i)*m*dpla_i*(uvar(n0,i)+dpla_i)**(m-1)
283 ELSE
284 yld_i =yld(i)
285 ENDIF
286 dre =youngm(i)*dpla_i/yld_i
287 drg =g*dpla_i/yld_i
288 pnx =one/(one+dre)
289 pny =one/(one+three*drg*ay**2)
290 pnz =one/(one+three*drg*az**2)
291 pmx =one/(one+three*drg*gamac(i)**2)
292 pmy =one/(one+dre*(by*gama(i))**2)
293 pmz =one/(one+dre*(bz*gama(i))**2)
294 pnx2=pnx*pnx
295 pny2=pny*pny
296 pnz2=pnz*pnz
297 pmx2=pmx*pmx
298 pmy2=pmy*pmy
299 pmz2=pmz*pmz
300 nx2=nx(i)*nx(i)
301 ny2=ny(i)*ny(i)
302 nz2=nz(i)*nz(i)
303 mx2=mx(i)*mx(i)
304 my2=my(i)*my(i)
305 mz2=mz(i)*mz(i)
306 fn=nx2*pnx2+3.*(ny2*pny2+nz2*pnz2)
307 fm=3.*mx2*pmx2+my2*pmy2+mz2*pmz2
308 f = fn+fm-yld_i*yld_i
309 dfe = nx2*pnx2*pnx+my2*pmy2*pmy*(by*gama(i))**2
310 . +mz2*pmz2*pmz*(bz*gama(i))**2
311 dfg = ny2*pny2*pny*ay**2+nz2*pnz2*pnz*az**2
312 . +mx2*pmx2*pmx*gamac(i)**2
313 IF (ifunc /= 0) THEN
314 dsigy = -h(i)*yld_i*two
315 ELSEIF (uvar(n0,i)>em20) THEN
316 dsigy = -(h(i)*yld_i*m*(dpla_i+uvar(n0,i))**(m-1))*two
317 ELSE
318 dsigy = zero
319 ENDIF
320 df =(-dfe*(youngm(i)-dre*h(i))/yld_i
321 . -nine*dfg*(g-drg*h(i))/yld_i)*two
322 df=df+dsigy
323 err=abs(f/df)
324 dpla_j(i)=max(zero,dpla_i-f/df)
325 ENDDO ! DO J=1,NINDX
326 n=n+1
327 ENDDO ! DO WHILE (ERR > EM4 .AND. N < NMAX)
328C
329C CONTRAINTES PLASTIQUEMENT ADMISSIBLES
330C
331 DO j=1,nindx
332 i=index(j)
333 uvar(n0,i) = uvar(n0,i) + dpla_j(i)
334 dpla_i=dpla_j(i)
335C INCREMENT N+1
336 IF (ifunc /= 0) THEN
337 yld_i =yld(i)+h(i)*dpla_i
338 ELSEIF (uvar(n0,i) > em20) THEN
339 yld_i =yld(i)+h(i)*m*dpla_i*(uvar(n0,i))**(m-1)
340 ELSE
341 yld_i =yld(i)
342 ENDIF
343 c1 = dpla_i/yld_i
344 dre =youngm(i)*c1
345 drg =g*c1
346 pnx =one/(one+dre)
347 pny =one/(one+three*drg*ay**2)
348 pnz =one/(one+three*drg*az**2)
349 pmx =one/(one+three*drg*gamac(i)**2)
350 pmy =one/(one+dre*(by*gama(i))**2)
351 pmz =one/(one+dre*(bz*gama(i))**2)
352C
353 signxx(i) = signxx(i)*pnx
354 signxy(i) = signxy(i)*pny
355 signxz(i) = signxz(i)*pnz
356 momnxx(i) = momnxx(i)*pmx
357 momnyy(i) = momnyy(i)*pmy
358 momnzz(i) = momnzz(i)*pmz
359C
360 uvar(n0+2,i) = uvar(n0+2,i)+three*c1*gamac(i)**2*abs(momnxx(i))
361 ueq = c1*sqrt( signxx(i)**2
362 1 + gama(i)**4
363 1 *(by**4*momnyy(i)**2 + bz**4*momnzz(i)**2) )
364 uvar(n0+3,i) = uvar(n0+3,i) + ueq
365C
366C ENDOMAGEMENT
367C
368 IF (uvar(n0,i) > ps) THEN
369 devol = dc*dpla_i/(pr - ps)
370 ELSE
371 devol = zero
372 ENDIF
373C
374 uvar(n0+4,i) = uvar(n0+4,i)+devol
375C
376 IF (uvar(n0,i) >= pr .OR. uvar(n0+4,i) >= dc) THEN
377 off(i) = zero
378 uvar(n0+4,i) = dc
379 ENDIF
380 ENDDO ! DO J=1,NINDX
381 ENDIF ! IF (NINDX /= 0)
382C
383C TRANSLATE FORCES AND MOMENTS :
384C
385 DO i=1,nel
386 fx(i) = signxx(i)*area*off(i)
387 fy(i) = signxy(i)*area_y(i)*off(i)
388 fz(i) = signxz(i)*area_z(i)*off(i)
389 xmom(i) = momnxx(i)*off(i)/facx
390 ymom(i) = momnyy(i)*off(i)/facy
391 zmom(i) = momnzz(i)*off(i)/facz
392 ENDDO
393C---
394 RETURN
my_real
#define my_real
Definition cppsort.cpp:32
area
subroutine area(d1, x, x2, y, y2, eint, stif0)
Definition hm_read_prop32.F:567
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
get_u_pid
integer function get_u_pid(ip)
Definition uaccess.F:625
get_u_pnu
integer function get_u_pnu(ivar, ip, k)
Definition uaccess.F:481
get_u_mid
integer function get_u_mid(im)
Definition uaccess.F:667
get_u_mnu
integer function get_u_mnu(ivar, im, k)
Definition uaccess.F:564