OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
lag_rby_cond.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!|| lag_rby_cond ../engine/source/tools/lagmul/lag_rby_cond.F
25!||--- called by ------------------------------------------------------
26!|| lag_rby ../engine/source/tools/lagmul/lag_rby.F
27!||====================================================================
28 SUBROUTINE lag_rby_cond(NPBYL ,LPBYL ,RBYL ,MASS ,INER ,
29 2 X ,V ,VR ,A ,AR ,
30 3 IADLL ,LLL ,COMNTAG,NN ,NC )
31C-----------------------------------------------
32C I m p l i c i t T y p e s
33C-----------------------------------------------
34#include "implicit_f.inc"
35C-----------------------------------------------
36C C o m m o n B l o c k s
37C-----------------------------------------------
38#include "param_c.inc"
39#include "com08_c.inc"
40C-----------------------------------------------
41C D u m m y A r g u m e n t s
42C-----------------------------------------------
43 INTEGER NN, NC
44 INTEGER LLL(*),IADLL(*),NPBYL(NNPBY,*),LPBYL(*),COMNTAG(*)
45C REAL
47 . rbyl(nrby,*),x(3,*),v(3,*),vr(3,*),a(3,*),ar(3,*),
48 . mass(*),iner(*)
49C-----------------------------------------------
50C L o c a l V a r i a b l e s
51C-----------------------------------------------
52 INTEGER I,J,IC,IK,N,NS,MSL,MSL2,M,IFX,IFR,NFIX,NFRE
53 my_real XX,YY,ZZ,XY,YZ,XZ,IXX,IYY,IZZ,IXY,IXZ,IYZ,
54 . jxx,jyy,jzz,jxy,jxz,jyz,jxy2,jxz2,jyz2,det,
55 . b1,b2,b3,c1,c2,c3,vx1,vx2,vx3,wa1,wa2,wa3,mmas,usdt,ddt,
56 . vi(3),vg(3),ag(3),xm(3),vtm(3),vrm(3),atm(3),arm(3),
57 . xrot(9),ipr(3)
58C======================================================================|
59 msl = npbyl(2,nn)
60 msl2= msl*2
61 nfix= 0
62 nfre= 0
63 ifx = 0
64 ifr = 0
65C--- main
66 ns = lpbyl(msl)
67 IF (comntag(ns)==1) THEN
68 nfix = nfix + 1
69 lpbyl(msl+nfix) = ns
70 ELSE
71 nfre = nfre + 1
72 lpbyl(msl2+nfre) = ns
73 ENDIF
74C--- and secnds
75 DO n=1,msl-1
76 ns = lpbyl(n)
77 IF (comntag(ns)==1) THEN
78 nfix = nfix + 1
79 lpbyl(msl+nfix) = ns
80 ELSE
81 nfre = nfre + 1
82 lpbyl(msl2+nfre) = ns
83 ENDIF
84 ENDDO
85 IF (nfix<=1) nfix = 0
86C--------------------------------------
87 IF (nfix/=0) THEN
88 ifx = lpbyl(msl+1)
89 jxx = zero
90 jyy = zero
91 jzz = zero
92 jxy = zero
93 jyz = zero
94 jxz = zero
95 mmas = zero
96 DO j=1,3
97 xm(j) = zero
98 vg(j) = zero
99 ag(j) = zero
100 vtm(j) = zero
101 atm(j) = zero
102 vrm(j) = zero
103 arm(j) = zero
104 ENDDO
105C--- CONDENSATION: CDG, masse, transl velocity, accel
106 DO i=1,nfix
107 n = lpbyl(msl+i)
108 mmas = mmas + mass(n)
109 DO j=1,3
110 xm(j)= xm(j) + x(j,n)*mass(n)
111 vtm(j)=vtm(j) + v(j,n)*mass(n)
112 atm(j)=atm(j) + a(j,n)*mass(n)
113 ENDDO
114 ENDDO
115 DO j=1,3
116 xm(j) = xm(j) / mmas
117 vtm(j) = vtm(j) / mmas
118 atm(j) = atm(j) / mmas
119 ENDDO
120 DO i=1,nfix
121 n = lpbyl(msl+i)
122 xx = x(1,n) - xm(1)
123 yy = x(2,n) - xm(2)
124 zz = x(3,n) - xm(3)
125c
126c VG(1)=VG(1) + VR(1,N)*INER(N)+MASS(N)*(YY*V(3,N)-ZZ*V(2,N))
127c VG(2)=VG(2) + VR(2,N)*INER(N)+MASS(N)*(ZZ*V(1,N)-XX*V(3,N))
128c VG(3)=VG(3) + VR(3,N)*INER(N)+MASS(N)*(XX*V(2,N)-YY*V(1,N))
129c
130c sum of moments
131 ag(1)=ag(1) + ar(1,n)*iner(n)+mass(n)*(yy*a(3,n)-zz*a(2,n))
132 ag(2)=ag(2) + ar(2,n)*iner(n)+mass(n)*(zz*a(1,n)-xx*a(3,n))
133 ag(3)=ag(3) + ar(3,n)*iner(n)+mass(n)*(xx*a(2,n)-yy*a(1,n))
134c---
135 xy=xx*yy
136 yz=yy*zz
137 xz=xx*zz
138c
139 xx=xx*xx
140 yy=yy*yy
141 zz=zz*zz
142 ixx = iner(n)+(yy+zz)*mass(n)
143 iyy = iner(n)+(xx+zz)*mass(n)
144 izz = iner(n)+(xx+yy)*mass(n)
145 ixy = xy*mass(n)
146 iyz = yz*mass(n)
147 ixz = xz*mass(n)
148 jxx = jxx + ixx
149 jyy = jyy + iyy
150 jzz = jzz + izz
151 jxy = jxy - ixy
152 jyz = jyz - iyz
153 jxz = jxz - ixz
154 ENDDO
155C-----------------------------
156 jxy2 = jxy*jxy
157 jyz2 = jyz*jyz
158 jxz2 = jxz*jxz
159 det = jxx*jyy*jzz-jxx*jyz2-jyy*jxz2-jzz*jxy2-two*jxy*jyz*jxz
160 det = one / det
161 b1 = det*(jzz*jyy-jyz2)
162 b2 = det*(jxx*jzz-jxz2)
163 b3 = det*(jyy*jxx-jxy2)
164 c1 = det*(jxx*jyz+jxz*jxy)
165 c2 = det*(jyy*jxz+jxy*jyz)
166 c3 = det*(jzz*jxy+jyz*jxz)
167C-----------------------------
168 vrm(1) = vr(1,ifx)
169 vrm(2) = vr(2,ifx)
170 vrm(3) = vr(3,ifx)
171C
172 vg(1) = vrm(1)*jxx + vrm(2)*jxy + vrm(3)*jxz
173 vg(2) = vrm(1)*jxy + vrm(2)*jyy + vrm(3)*jyz
174 vg(3) = vrm(1)*jxz + vrm(2)*jyz + vrm(3)*jzz
175C
176c VRM(1)= VG(1)*B1 + VG(2)*C3 + VG(3)*C2
177c VRM(2)= VG(1)*C3 + VG(2)*B2 + VG(3)*C1
178c VRM(3)= VG(1)*C2 + VG(2)*C1 + VG(3)*B3
179
180C
181 ag(1) = ag(1) - vrm(2)*vg(3) + vrm(3)*vg(2)
182 ag(2) = ag(2) - vrm(3)*vg(1) + vrm(1)*vg(3)
183 ag(3) = ag(3) - vrm(1)*vg(2) + vrm(2)*vg(1)
184C
185 arm(1)= ag(1)*b1 + ag(2)*c3 + ag(3)*c2
186 arm(2)= ag(1)*c3 + ag(2)*b2 + ag(3)*c1
187 arm(3)= ag(1)*c2 + ag(2)*c1 + ag(3)*b3
188C-----------------------------
189 IF (nfre==0) THEN
190C--- Total condensation => direct solution + decondensation
191 usdt = one / dt12
192 ddt = half * dt12
193 DO j=1,3
194 vrm(j) = vrm(j) + arm(j)*dt12
195 ENDDO
196 DO n=1,msl
197 ns = lpbyl(n)
198 DO j=1,3
199 ar(j,ns) = (vrm(j)-vr(j,ns)) * usdt
200 ENDDO
201 xx = x(1,ns)-xm(1)
202 yy = x(2,ns)-xm(2)
203 zz = x(3,ns)-xm(3)
204 vx1 = vrm(2)*zz - vrm(3)*yy
205 vx2 = vrm(3)*xx - vrm(1)*zz
206 vx3 = vrm(1)*yy - vrm(2)*xx
207 a(1,ns) = atm(1) + usdt*
208 . (vtm(1)-v(1,ns)+vx1+ddt*(vrm(2)*vx3-vrm(3)*vx2))
209 a(2,ns) = atm(2) + usdt*
210 . (vtm(2)-v(2,ns)+vx2+ddt*(vrm(3)*vx1-vrm(1)*vx3))
211 a(3,ns) = atm(3) + usdt*
212 . (vtm(3)-v(3,ns)+vx3+ddt*(vrm(1)*vx2-vrm(2)*vx1))
213
214c A(1,NS) = ATM(1) + USDT*(VTM(1) - V(1,NS) + VX1)
215c A(2,NS) = ATM(2) + USDT*(VTM(2) - V(2,NS) + VX2)
216c A(3,NS) = ATM(3) + USDT*(VTM(3) - V(3,NS) + VX3)
217 ENDDO
218 ELSE
219C--- partial condensation : save condensed values for further treatment
220 ifr = lpbyl(2*msl+1)
221 rbyl(10,nn) = mass(ifx)
222 rbyl(14,nn) = v(1,ifx)
223 rbyl(15,nn) = v(2,ifx)
224 rbyl(16,nn) = v(3,ifx)
225 rbyl(17,nn) = vr(1,ifx)
226 rbyl(18,nn) = vr(2,ifx)
227 rbyl(19,nn) = vr(3,ifx)
228 rbyl(20,nn) = a(1,ifx)
229 rbyl(21,nn) = a(2,ifx)
230 rbyl(22,nn) = a(3,ifx)
231 rbyl(23,nn) = ar(1,ifx)
232 rbyl(24,nn) = ar(2,ifx)
233 rbyl(25,nn) = ar(3,ifx)
234C
235 rbyl(1,nn) = b1
236 rbyl(2,nn) = b2
237 rbyl(3,nn) = b3
238 rbyl(4,nn) = c1
239 rbyl(5,nn) = c2
240 rbyl(6,nn) = c3
241 rbyl(11,nn) = xm(1)
242 rbyl(12,nn) = xm(2)
243 rbyl(13,nn) = xm(3)
244C
245 mass(ifx) = mmas
246 v(1,ifx) = vtm(1)
247 v(2,ifx) = vtm(2)
248 v(3,ifx) = vtm(3)
249 vr(1,ifx) = vrm(1)
250 vr(2,ifx) = vrm(2)
251 vr(3,ifx) = vrm(3)
252 a(1,ifx) = atm(1)
253 a(2,ifx) = atm(2)
254 a(3,ifx) = atm(3)
255 ar(1,ifx) = arm(1)
256 ar(2,ifx) = arm(2)
257 ar(3,ifx) = arm(3)
258 ENDIF
259 ENDIF
260 npbyl(4,nn) = nfix
261 npbyl(5,nn) = nfre
262 npbyl(7,nn) = ifx
263 npbyl(8,nn) = ifr
264C---
265 RETURN
266 END
267!||====================================================================
268!|| rby_decond ../engine/source/tools/lagmul/lag_rby_cond.F
269!||--- called by ------------------------------------------------------
270!|| lag_mult ../engine/source/tools/lagmul/lag_mult.F
271!|| lag_multp ../engine/source/tools/lagmul/lag_mult.F
272!||--- calls -----------------------------------------------------
273!|| ancmsg ../engine/source/output/message/message.F
274!|| arret ../engine/source/system/arret.F
275!||--- uses -----------------------------------------------------
276!|| message_mod ../engine/share/message_module/message_mod.F
277!||====================================================================
278 SUBROUTINE rby_decond(X ,V ,VR ,A ,AR ,
279 2 IADLL ,LLL ,JLL ,XLL ,LAMBDA ,
280 3 MASS ,INER ,RBYL ,NPBYL ,LPBYL ,
281 4 NC ,NCR )
282C-----------------------------------------------
283C M o d u l e s
284C-----------------------------------------------
285 USE message_mod
286C-----------------------------------------------
287C I m p l i c i t T y p e s
288C-----------------------------------------------
289#include "implicit_f.inc"
290C-----------------------------------------------
291C C o m m o n B l o c k s
292C-----------------------------------------------
293#include "param_c.inc"
294#include "lagmult.inc"
295#include "com08_c.inc"
296
297
298C-----------------------------------------------
299C D u m m y A r g u m e n t s
300C-----------------------------------------------
301 INTEGER NC, NCR
302 INTEGER IADLL(*),LLL(*),JLL(*),NPBYL(NNPBY,*),LPBYL(*)
303C REAL
304 my_real
305 . RBYL(NRBY,*),XLL(*),X(3,*),V(3,*),VR(3,*),A(3,*),AR(3,*),
306 . mass(*),iner(*),lambda(*)
307C-----------------------------------------------
308C L o c a l V a r i a b l e s
309C-----------------------------------------------
310 INTEGER I,J,K,JF,IC,IK,IR,IFX,IFR,N,NS,NFIX,NFRE,MSL,TNSL
311 my_real
312 . XX,YY,ZZ,VX1,VX2,VX3,USDT,DDT,XM(3),VTM(3),VRM(3),ATM(3),ARM(3)
313C======================================================================|
314 USDT = one / dt12
315 ddt = half *dt12
316 tnsl = 0
317 ic = ncr
318 DO ir = 1,nrbylag
319 msl = npbyl(2,ir)
320 ifx = npbyl(7,ir)
321 nfix = npbyl(4,ir)
322 nfre = npbyl(5,ir)
323 IF (nfix>0.AND.nfre>0) THEN
324C--- calculate acceleration of condensed node: a = ao + [M]-1[L]t la
325 DO k = 1,3
326 ic = ic + 1
327 DO ik=iadll(ic),iadll(ic+1)-1
328 i = lll(ik)
329 j = jll(ik)
330 xll(ik) = xll(ik)*lambda(ic)
331 IF (j<=3) THEN
332c A(J,I) = A(J,I) - XLL(IK)/MASS(I)
333 CALL ancmsg(msgid=117,anmode=aninfo,
334 . i1=i)
335 CALL arret(2)
336 ELSEIF(i/=ifx) THEN
337 j = j-3
338 ar(j,i) = ar(j,i) - xll(ik)/iner(i)
339 ELSEIF (xll(ik)/=0.) THEN
340 IF(j==4) THEN
341 ar(1,ifx) = ar(1,ifx) - xll(ik)*rbyl(1,ir)
342 ar(2,ifx) = ar(2,ifx) - xll(ik)*rbyl(6,ir)
343 ar(3,ifx) = ar(3,ifx) - xll(ik)*rbyl(5,ir)
344 ELSEIF(j==5) THEN
345 ar(1,ifx) = ar(1,ifx) - xll(ik)*rbyl(6,ir)
346 ar(2,ifx) = ar(2,ifx) - xll(ik)*rbyl(2,ir)
347 ar(3,ifx) = ar(3,ifx) - xll(ik)*rbyl(4,ir)
348 ELSE
349 ar(1,ifx) = ar(1,ifx) - xll(ik)*rbyl(5,ir)
350 ar(2,ifx) = ar(2,ifx) - xll(ik)*rbyl(4,ir)
351 ar(3,ifx) = ar(3,ifx) - xll(ik)*rbyl(3,ir)
352 ENDIF
353 ENDIF
354 ENDDO
355 ENDDO
356C--- decondensation
357 mass(ifx) = rbyl(10,ir)
358 DO j=1,3
359 vtm(j) = v(j,ifx)
360 vrm(j) = vr(j,ifx)
361 atm(j) = a(j,ifx)
362 arm(j) = ar(j,ifx)
363 xm(j) = rbyl(10+j,ir)
364 v(j,ifx) = rbyl(13+j,ir)
365 vr(j,ifx) = rbyl(16+j,ir)
366 a(j,ifx) = rbyl(19+j,ir)
367 ar(j,ifx) = rbyl(22+j,ir)
368 ENDDO
369 k = tnsl+msl
370 DO j=1,3
371 vrm(j) = vrm(j) + arm(j)*dt12
372 ENDDO
373 DO n=1,nfix
374 ns = lpbyl(k+n)
375 DO j=1,3
376 ar(j,ns) = (vrm(j)-vr(j,ns)) * usdt
377 ENDDO
378 xx = x(1,ns) - xm(1)
379 yy = x(2,ns) - xm(2)
380 zz = x(3,ns) - xm(3)
381c
382 vx1 = vrm(2)*zz - vrm(3)*yy
383 vx2 = vrm(3)*xx - vrm(1)*zz
384 vx3 = vrm(1)*yy - vrm(2)*xx
385c pas de second ordre
386c A(1,NS) = ATM(1) + USDT*(VTM(1) + VX1 - V(1,NS))
387c A(2,NS) = ATM(2) + USDT*(VTM(2) + VX2 - V(2,NS))
388c A(3,NS) = ATM(3) + USDT*(VTM(3) + VX3 - V(3,NS))
389c calcul du second ordre
390 a(1,ns) = atm(1) + usdt*
391 . (vtm(1)-v(1,ns)+vx1+ddt*(vrm(2)*vx3-vrm(3)*vx2))
392 a(2,ns) = atm(2) + usdt*
393 . (vtm(2)-v(2,ns)+vx2+ddt*(vrm(3)*vx1-vrm(1)*vx3))
394 a(3,ns) = atm(3) + usdt*
395 . (vtm(3)-v(3,ns)+vx3+ddt*(vrm(1)*vx2-vrm(2)*vx1))
396 ENDDO
397 ENDIF
398 tnsl = tnsl + 3*msl
399 ENDDO
400C---
401 RETURN
402 END
#define my_real
Definition cppsort.cpp:32
subroutine lag_rby_cond(npbyl, lpbyl, rbyl, mass, iner, x, v, vr, a, ar, iadll, lll, comntag, nn, nc)
subroutine rby_decond(x, v, vr, a, ar, iadll, lll, jll, xll, lambda, mass, iner, rbyl, npbyl, lpbyl, nc, ncr)
subroutine ancmsg(msgid, msgtype, anmode, i1, i2, i3, i4, i5, i6, i7, i8, i9, i10, i11, i12, i13, i14, i15, i16, i17, i18, i19, i20, r1, r2, r3, r4, r5, r6, r7, r8, r9, c1, c2, c3, c4, c5, c6, c7, c8, c9, prmode)
Definition message.F:889
subroutine arret(nn)
Definition arret.F:87