OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
s4deri3.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/.
23CCgw|============================================================
24!||====================================================================
25!|| s4deri3 ../starter/source/elements/solid/solide4/s4deri3.F
26!||--- called by ------------------------------------------------------
27!|| inirig_mat ../starter/source/elements/initia/inirig_mat.F
28!|| inivoid ../starter/source/elements/initia/inivoid.F
29!|| multifluid_init3t ../starter/source/multifluid/multifluid_init3t.F
30!|| s4init3 ../starter/source/elements/solid/solide4/s4init3.F
31!|| s4refsta3 ../starter/source/elements/solid/solide4/s4refsta3.F
32!||--- calls -----------------------------------------------------
33!|| ancmsg ../starter/source/output/message/message.F
34!||--- uses -----------------------------------------------------
35!|| message_mod ../starter/share/message_module/message_mod.F
36!||====================================================================
37 SUBROUTINE s4deri3(VOL ,VEUL ,GEO ,IGEO ,RX ,
38 . RY ,RZ ,SX ,SY ,SZ ,
39 . TX ,TY ,TZ ,
40 . X1 ,X2 ,X3 ,X4 ,Y1 ,Y2 ,
41 . Y3 ,Y4 ,Z1 ,Z2 ,Z3 ,Z4 ,
42 . PX1 ,PX2 ,PX3 ,PX4 ,
43 . PY1 ,PY2 ,PY3 ,PY4 ,
44 . PZ1 ,PZ2 ,PZ3 ,PZ4 ,JAC_I,
45 . DELTAX,DET ,NGL ,NGEO ,MXT ,
46 . PM ,VOLDP )
47 USE message_mod
48C-----------------------------------------------
49C I m p l i c i t T y p e s
50C-----------------------------------------------
51#include "implicit_f.inc"
52C-----------------------------------------------
53C G l o b a l P a r a m e t e r s
54C-----------------------------------------------
55#include "mvsiz_p.inc"
56C-----------------------------------------------
57C C o m m o n B l o c k s
58C-----------------------------------------------
59#include "scr03_c.inc"
60#include "vect01_c.inc"
61#include "param_c.inc"
62#include "com01_c.inc"
63#include "scr17_c.inc"
64C-----------------------------------------------
65C D u m m y A r g u m e n t s
66C-----------------------------------------------
67 INTEGER IGEO(NPROPGI,*), MXT(MVSIZ)
68 DOUBLE PRECISION
69 . X1(MVSIZ), X2(MVSIZ), X3(MVSIZ), X4(MVSIZ),
70 . Y1(MVSIZ), Y2(MVSIZ), Y3(MVSIZ), Y4(MVSIZ),
71 . Z1(MVSIZ), Z2(MVSIZ), Z3(MVSIZ), Z4(MVSIZ),VOLDP(*)
72C REAL
73 my_real
74 . VOL(*), VEUL(LVEUL,*),GEO(NPROPG,*), PM(NPROPM,*),
75 . RX(*), RY(*), RZ(*), SX(*), SY(*), SZ(*), TX(*), TY(*), TZ(*), DET(*),
76 . px1(mvsiz), px2(mvsiz), px3(mvsiz), px4(mvsiz),
77 . py1(mvsiz), py2(mvsiz), py3(mvsiz), py4(mvsiz),
78 . pz1(mvsiz), pz2(mvsiz), pz3(mvsiz), pz4(mvsiz),
79 . deltax(mvsiz),jac_i(10,mvsiz)
80 INTEGER NGL(*), NGEO(*)
81C-----------------------------------------------
82C L o c a l V a r i a b l e s
83C-----------------------------------------------
84 INTEGER I
85 my_real
86 . B1(MVSIZ), B2(MVSIZ), B3(MVSIZ), B4(MVSIZ),
87 . C1(MVSIZ), C2(MVSIZ), C3(MVSIZ), C4(MVSIZ),
88 . D1(MVSIZ), D2(MVSIZ), D3(MVSIZ), D4(MVSIZ)
89 double precision
90 . x41, y41, z41, x42, y42, z42, x43, y43, z43,b1dp,c1dp,d1dp
91 my_real
92 . d ,pxx, pyy, pzz, pxy, pyz, pxz, gfac, aa, bb, p, ld
93C=======================================================================
94 DO i=lft,llt
95 x43 = x4(i) - x3(i)
96 y43 = y4(i) - y3(i)
97 z43 = z4(i) - z3(i)
98 x41 = x4(i) - x1(i)
99 y41 = y4(i) - y1(i)
100 z41 = z4(i) - z1(i)
101 x42 = x4(i) - x2(i)
102 y42 = y4(i) - y2(i)
103 z42 = z4(i) - z2(i)
104 rx(i) = -x41
105 ry(i) = -y41
106 rz(i) = -z41
107 sx(i) = -x42
108 sy(i) = -y42
109 sz(i) = -z42
110 tx(i) = -x43
111 ty(i) = -y43
112 tz(i) = -z43
113C
114 b1dp = y43*z42 - y42*z43
115 b1(i) = b1dp
116 b2(i) = y41*z43 - y43*z41
117 b3(i) = y42*z41 - y41*z42
118 b4(i) = -(b1(i) + b2(i) + b3(i))
119C
120 c1dp = z43*x42 - z42*x43
121 c1(i) = c1dp
122 c2(i) = z41*x43 - z43*x41
123 c3(i) = z42*x41 - z41*x42
124 c4(i) = -(c1(i) + c2(i) + c3(i))
125C
126 d1dp = x43*y42 - x42*y43
127 d1(i) = d1dp
128 d2(i) = x41*y43 - x43*y41
129 d3(i) = x42*y41 - x41*y42
130 d4(i) = -(d1(i) + d2(i) + d3(i))
131C
132 voldp(i) = (x41*b1dp + y41*c1dp + z41*d1dp)*one_over_6
133 det(i) = voldp(i)
134 vol(i) = det(i)
135 ENDDO
136C
137 DO i=lft,llt
138 IF (det(i) <= zero) THEN
139 IF (igeo(11,ngeo(i)) /= 0) THEN
140 CALL ancmsg(msgid=245,
141 . msgtype=msgerror,
142 . anmode=aninfo,
143 . i1=ngl(i))
144 vol(i) = em20 ! to prevent crash
145 det(i) = em20
146 ELSE
147 CALL ancmsg(msgid=635,
148 . msgtype=msgwarning,
149 . anmode=aninfo,
150 . i1=ngl(i))
151 ENDIF
152 ENDIF
153 ENDDO
154C
155 DO i=lft,llt
156 d = one/det(i)/six
157 px1(i)=-b1(i)*d
158 py1(i)=-c1(i)*d
159 pz1(i)=-d1(i)*d
160 px2(i)=-b2(i)*d
161 py2(i)=-c2(i)*d
162 pz2(i)=-d2(i)*d
163 px3(i)=-b3(i)*d
164 py3(i)=-c3(i)*d
165 pz3(i)=-d3(i)*d
166 px4(i)=-b4(i)*d
167 py4(i)=-c4(i)*d
168 pz4(i)=-d4(i)*d
169 ENDDO
170
171
172 IF(idt1sol==0)THEN
173
174 DO i=lft,llt
175 d = max(px1(i)*px1(i)+py1(i)*py1(i)+pz1(i)*pz1(i),
176 . px2(i)*px2(i)+py2(i)*py2(i)+pz2(i)*pz2(i),
177 . px3(i)*px3(i)+py3(i)*py3(i)+pz3(i)*pz3(i),
178 . px4(i)*px4(i)+py4(i)*py4(i)+pz4(i)*pz4(i))
179 deltax(i) = one / sqrt(d)
180 END DO
181
182 ELSEIF(iformdt==0)THEN
183 DO i=lft,llt
184 d = px1(i)*px1(i)+py1(i)*py1(i)+pz1(i)*pz1(i)
185 . + px2(i)*px2(i)+py2(i)*py2(i)+pz2(i)*pz2(i)
186 . + px3(i)*px3(i)+py3(i)*py3(i)+pz3(i)*pz3(i)
187 . + px4(i)*px4(i)+py4(i)*py4(i)+pz4(i)*pz4(i)
188 deltax(i) = one / sqrt(d)
189 END DO
190
191 ELSEIF(iformdt==1)THEN
192
193 gfac=pm(105,mxt(1))
194 ld =two*sqrt(max(one-gfac,zero))+one
195 DO i=lft,llt
196 pxx=px1(i)*px1(i)+px2(i)*px2(i)+px3(i)*px3(i)+px4(i)*px4(i)
197 pyy=py1(i)*py1(i)+py2(i)*py2(i)+py3(i)*py3(i)+py4(i)*py4(i)
198 pzz=pz1(i)*pz1(i)+pz2(i)*pz2(i)+pz3(i)*pz3(i)+pz4(i)*pz4(i)
199 pxy=px1(i)*py1(i)+px2(i)*py2(i)+px3(i)*py3(i)+px4(i)*py4(i)
200 pxz=px1(i)*pz1(i)+px2(i)*pz2(i)+px3(i)*pz3(i)+px4(i)*pz4(i)
201 pyz=py1(i)*pz1(i)+py2(i)*pz2(i)+py3(i)*pz3(i)+py4(i)*pz4(i)
202C
203 aa = -(pxx+pyy+pzz)
204 bb = (pxx*pyy+pxx*pzz+pyy*pzz-pxy**2-pxz**2-pyz**2)
205 p = bb-third*aa*aa
206 d = two*sqrt(third*max(-p,zero))-third*aa
207C
208 d = ld*d
209C
210 deltax(i) = one / sqrt(d)
211 END DO
212
213 ELSEIF(iformdt==2)THEN
214
215 gfac=pm(105,mxt(1))
216 DO i=lft,llt
217 pxx=px1(i)*px1(i)+px2(i)*px2(i)+px3(i)*px3(i)+px4(i)*px4(i)
218 pyy=py1(i)*py1(i)+py2(i)*py2(i)+py3(i)*py3(i)+py4(i)*py4(i)
219 pzz=pz1(i)*pz1(i)+pz2(i)*pz2(i)+pz3(i)*pz3(i)+pz4(i)*pz4(i)
220 pxy=px1(i)*py1(i)+px2(i)*py2(i)+px3(i)*py3(i)+px4(i)*py4(i)
221 pxz=px1(i)*pz1(i)+px2(i)*pz2(i)+px3(i)*pz3(i)+px4(i)*pz4(i)
222 pyz=py1(i)*pz1(i)+py2(i)*pz2(i)+py3(i)*pz3(i)+py4(i)*pz4(i)
223C
224 aa = -(pxx+pyy+pzz)
225 bb = gfac*(pxx*pyy+pxx*pzz+pyy*pzz-pxy**2-pxz**2-pyz**2)
226 p = bb-third*aa*aa
227 d = two*sqrt(third*max(-p,zero))-third*aa
228C
229 deltax(i) = one / sqrt(d)
230 END DO
231
232 END IF
233
234 IF (ismstr==10.OR.ismstr==12) THEN
235 DO i=lft,llt
236 jac_i(1,i) = px1(i)
237 jac_i(2,i) = px2(i)
238 jac_i(3,i) = px3(i)
239 jac_i(4,i) = py1(i)
240 jac_i(5,i) = py2(i)
241 jac_i(6,i) = py3(i)
242 jac_i(7,i) = pz1(i)
243 jac_i(8,i) = pz2(i)
244 jac_i(9,i) = pz3(i)
245 jac_i(10,i) = vol(i)
246 ENDDO
247 END IF
248C
249 IF (jeul * (1 - imulti_fvm) /= 0 .AND. nxref==0) THEN
250 DO i=lft,llt
251 veul(1,i) =-px1(i)
252 veul(2,i) = px2(i)
253 veul(3,i) =-px3(i)
254 veul(4,i) = px4(i)
255 veul(5,i) =-py1(i)
256 veul(6,i) = py2(i)
257 veul(7,i) =-py3(i)
258 veul(8,i) = py4(i)
259 veul(9,i) =-pz1(i)
260 veul(10,i)= pz2(i)
261 veul(11,i)=-pz3(i)
262 veul(12,i)= pz4(i)
263 veul(13,i)= deltax(i)
264C------------------------------------------
265C CALCUL DE LA NORMALE A CHAQUE FACE
266C------------------------------------------
267 veul(14,i) = zero
268 veul(15,i) = half*b1(i)
269 veul(16,i) = zero
270 veul(17,i) = half*b2(i)
271 veul(18,i) = half*b3(i)
272 veul(19,i) = half*b4(i)
273C
274 veul(20,i) = zero
275 veul(21,i) = half*c1(i)
276 veul(22,i) = zero
277 veul(23,i) = half*c2(i)
278 veul(24,i) = half*c3(i)
279 veul(25,i) = half*c4(i)
280C
281 veul(26,i) = zero
282 veul(27,i) = half*d1(i)
283 veul(28,i) = zero
284 veul(29,i) = half*d2(i)
285 veul(30,i) = half*d3(i)
286 veul(31,i) = half*d4(i)
287 ENDDO
288C
289 IF (igeo(11,ngeo(lft)) == 15) THEN
290 DO i=lft,llt
291 vol(i)=vol(i)*geo(1,ngeo(i))
292 voldp(i)=voldp(i)*geo(1,ngeo(i))
293 ENDDO
294 ENDIF
295 ENDIF
296C-----------
297 RETURN
298 END
max
#define max(a, b)
Definition macros.h:21
message_mod
Definition message_mod.F:1249
s4deri3
subroutine s4deri3(vol, veul, geo, igeo, rx, ry, rz, sx, sy, sz, tx, ty, tz, x1, x2, x3, x4, y1, y2, y3, y4, z1, z2, z3, z4, px1, px2, px3, px4, py1, py2, py3, py4, pz1, pz2, pz3, pz4, jac_i, deltax, det, ngl, ngeo, mxt, pm, voldp)
Definition s4deri3.F:47
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:889