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