OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
alefvm_stress.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!|| alefvm_stress ../engine/source/ale/alefvm/alefvm_stress.F
25!||--- called by ------------------------------------------------------
26!|| sforc3 ../engine/source/elements/solid/solide/sforc3.F
27!||--- uses -----------------------------------------------------
28!|| alefvm_mod ../common_source/modules/ale/alefvm_mod.F
29!|| element_mod ../common_source/modules/elements/element_mod.F90
30!||====================================================================
31 SUBROUTINE alefvm_stress (
32 1 IXS, SIG, QVIS,
33 2 N1X, N2X, N3X, N4X, N5X, N6X,
34 3 N1Y, N2Y, N3Y, N4Y, N5Y, N6Y,
35 4 N1Z, N2Z, N3Z, N4Z, N5Z, N6Z,
36 5 X , IPM, RHO, VOL,
37 . NEL, MOM, SSP, SVIS)
38C-----------------------------------------------
39C D e s c r i p t i o n
40C-----------------------------------------------
41C 'alefvm' is related to a collocated scheme (built from FVM and based on Godunov scheme)
42C which was temporarily introduced for experimental option /INTER/TYPE22 (FSI coupling with cut cell method)
43C This cut cell method is not completed, abandoned, and is not an official option.
44C There is no other use for this scheme which is automatically enabled when /INTER/TYPE22 is defined (INT22>0 => IALEFVM=1).
45C
46C This subroutine is treating an uncut cell.
47C-----------------------------------------------
48C M o d u l e s
49C-----------------------------------------------
50 USE alefvm_mod
51 use element_mod , only : nixs
52C-----------------------------------------------
53C I m p l i c i t T y p e s
54C-----------------------------------------------
55#include "implicit_f.inc"
56C-----------------------------------------------
57C G l o b a l P a r a m e t e r s
58C-----------------------------------------------
59#include "mvsiz_p.inc"
60C-----------------------------------------------
61C C o m m o n B l o c k s
62C-----------------------------------------------
63#include "vect01_c.inc"
64#include "com01_c.inc"
65#include "param_c.inc"
66C-----------------------------------------------
67C D e s c r i p t i o n
68C-----------------------------------------------
69C This subroutines computes internal forces for
70C finite volume scheme (IALEFVM==1)
71C
72C If option is not detected in input file then
73C subroutine is unplugged
74C-----------------------------------------------
75C D u m m y A r g u m e n t s
76C-----------------------------------------------
77 INTEGER :: NEL
78 INTEGER :: IXS(NIXS,*),IPM(NPROPMI,*)
79 my_real :: SIG(NEL,6),QVIS(*),X(3,*),RHO(*), VOL(*)
80 my_real,INTENT(IN) :: MOM(NEL,3), SSP(*)
81 my_real :: N1X(*), N2X(*), N3X(*), N4X(*), N5X(*), N6X(*),
82 . N1Y(*), N2Y(*), N3Y(*), N4Y(*), N5Y(*), N6Y(*),
83 . n1z(*), n2z(*), n3z(*), n4z(*), n5z(*), n6z(*)
84 my_real, DIMENSION(MVSIZ,6), INTENT(INOUT) :: svis
85C-----------------------------------------------
86C L o c a l V a r i a b l e s
87C-----------------------------------------------
88 INTEGER :: I, II, IMAT, IALEFVM_FLG
89 INTEGER :: NC1(MVSIZ),NC2(MVSIZ),NC3(MVSIZ),NC4(MVSIZ),NC5(MVSIZ),NC6(MVSIZ),NC7(MVSIZ),NC8(MVSIZ)
90
91 my_real :: p(mvsiz) ,
92 . s1(mvsiz), s2(mvsiz), s3(mvsiz),
93 . s4(mvsiz), s5(mvsiz), s6(mvsiz),
94 . mass
95 my_real x1(mvsiz), x2(mvsiz), x3(mvsiz), x4(mvsiz), x5(mvsiz), x6(mvsiz), x7(mvsiz), x8(mvsiz),
96 . y1(mvsiz), y2(mvsiz), y3(mvsiz), y4(mvsiz), y5(mvsiz), y6(mvsiz), y7(mvsiz), y8(mvsiz),
97 . z1(mvsiz), z2(mvsiz), z3(mvsiz), z4(mvsiz), z5(mvsiz), z6(mvsiz), z7(mvsiz), z8(mvsiz),
98 . u_n(6), surf(6), norm(6)
99
100 INTEGER :: idbf,idbl
101 LOGICAL :: debug_outp
102
103C-----------------------------------------------
104C P r e - C o n d i t i o n s
105C-----------------------------------------------
106 IF(alefvm_param%IEnabled==0)RETURN
107 imat = ixs(1,1+nft)
108 ialefvm_flg = ipm(251,imat)
109 IF(ialefvm_flg <= 1)RETURN
110C-----------------------------------------------
111C S o u r c e L i n e s
112C-----------------------------------------------
113
114 ! attention: mmain must be processed for each group
115 ! in order to know the neighboring tensors (SIG,SVIS,QVIS) !
116 ! the assembly is then carried out after all calls to
117 ! ALEMAIN>SFORC3() pour chaque groupe.
118
119 !-------------------------------------------------------------!
120 ! TOTAL TENSOR !
121 !-------------------------------------------------------------!
122 DO i=1,nel
123 s1(i) = sig(i,1) + svis(i,1) - qvis(i)
124 s2(i) = sig(i,2) + svis(i,2) - qvis(i)
125 s3(i) = sig(i,3) + svis(i,3) - qvis(i)
126 s4(i) = sig(i,4) + svis(i,4)
127 s5(i) = sig(i,5) + svis(i,5)
128 s6(i) = sig(i,6) + svis(i,6)
129 ENDDO
130
131 DO i=1,nel
132 p(i) = -third*(s1(i)+s2(i)+s3(i))
133 ENDDO
134
135
136 !-------------------------------------------------------------!
137 ! NORMAL VECTOR FOR ALE !
138 !-------------------------------------------------------------!
139 DO i=1,nel
140 ii = i + nft
141 !---8 local node numbers NC1 TO NC8 for solid element I ---!
142 nc1(i)=ixs(2,ii)
143 nc2(i)=ixs(3,ii)
144 nc3(i)=ixs(4,ii)
145 nc4(i)=ixs(5,ii)
146 nc5(i)=ixs(6,ii)
147 nc6(i)=ixs(7,ii)
148 nc7(i)=ixs(8,ii)
149 nc8(i)=ixs(9,ii)
150 !
151 !---Coordinates of the 8 nodes
152 x1(i)=x(1,nc1(i))
153 y1(i)=x(2,nc1(i))
154 z1(i)=x(3,nc1(i))
155 !
156 x2(i)=x(1,nc2(i))
157 y2(i)=x(2,nc2(i))
158 z2(i)=x(3,nc2(i))
159 !
160 x3(i)=x(1,nc3(i))
161 y3(i)=x(2,nc3(i))
162 z3(i)=x(3,nc3(i))
163 !
164 x4(i)=x(1,nc4(i))
165 y4(i)=x(2,nc4(i))
166 z4(i)=x(3,nc4(i))
167 !
168 x5(i)=x(1,nc5(i))
169 y5(i)=x(2,nc5(i))
170 z5(i)=x(3,nc5(i))
171 !
172 x6(i)=x(1,nc6(i))
173 y6(i)=x(2,nc6(i))
174 z6(i)=x(3,nc6(i))
175 !
176 x7(i)=x(1,nc7(i))
177 y7(i)=x(2,nc7(i))
178 z7(i)=x(3,nc7(i))
179 !
180 x8(i)=x(1,nc8(i))
181 y8(i)=x(2,nc8(i))
182 z8(i)=x(3,nc8(i))
183 ENDDO
184 DO i=1,nel
185 ! Face-1
186 n1x(i)=(y3(i)-y1(i))*(z2(i)-z4(i)) - (z3(i)-z1(i))*(y2(i)-y4(i))
187 n1y(i)=(z3(i)-z1(i))*(x2(i)-x4(i)) - (x3(i)-x1(i))*(z2(i)-z4(i))
188 n1z(i)=(x3(i)-x1(i))*(y2(i)-y4(i)) - (y3(i)-y1(i))*(x2(i)-x4(i))
189 ! Face-2
190 n2x(i)=(y7(i)-y4(i))*(z3(i)-z8(i)) - (z7(i)-z4(i))*(y3(i)-y8(i))
191 n2y(i)=(z7(i)-z4(i))*(x3(i)-x8(i)) - (x7(i)-x4(i))*(z3(i)-z8(i))
192 n2z(i)=(x7(i)-x4(i))*(y3(i)-y8(i)) - (y7(i)-y4(i))*(x3(i)-x8(i))
193 ! Face-3
194 n3x(i)=(y6(i)-y8(i))*(z7(i)-z5(i)) - (z6(i)-z8(i))*(y7(i)-y5(i))
195 n3y(i)=(z6(i)-z8(i))*(x7(i)-x5(i)) - (x6(i)-x8(i))*(z7(i)-z5(i))
196 n3z(i)=(x6(i)-x8(i))*(y7(i)-y5(i)) - (y6(i)-y8(i))*(x7(i)-x5(i))
197 ! Face-4
198 n4x(i)=(y2(i)-y5(i))*(z6(i)-z1(i)) - (z2(i)-z5(i))*(y6(i)-y1(i))
199 n4y(i)=(z2(i)-z5(i))*(x6(i)-x1(i)) - (x2(i)-x5(i))*(z6(i)-z1(i))
200 n4z(i)=(x2(i)-x5(i))*(y6(i)-y1(i)) - (y2(i)-y5(i))*(x6(i)-x1(i))
201 ! Face-5
202 n5x(i)=(y7(i)-y2(i))*(z6(i)-z3(i)) - (z7(i)-z2(i))*(y6(i)-y3(i))
203 n5y(i)=(z7(i)-z2(i))*(x6(i)-x3(i)) - (x7(i)-x2(i))*(z6(i)-z3(i))
204 n5z(i)=(x7(i)-x2(i))*(y6(i)-y3(i)) - (y7(i)-y2(i))*(x6(i)-x3(i))
205 ! Face-6
206 n6x(i)=(y8(i)-y1(i))*(z4(i)-z5(i)) - (z8(i)-z1(i))*(y4(i)-y5(i))
207 n6y(i)=(z8(i)-z1(i))*(x4(i)-x5(i)) - (x8(i)-x1(i))*(z4(i)-z5(i))
208 n6z(i)=(x8(i)-x1(i))*(y4(i)-y5(i)) - (y8(i)-y1(i))*(x4(i)-x5(i))
209 ENDDO
210
211 !--------------------------------------------------------------!
212 ! PREPARING BUFFER FOR ALEFVM_SFINT3 !
213 ! 0 : cell centroid data !
214 ! 1:3 : INTEGRAL ON EACH FACE from Integral(DIV(SIGMA),Volume)!
215 ! 4 : max(u+x,u-c) !
216 !--------------------------------------------------------------!
217 DO i=1, nel
218 ii = i + nft
219 mass = rho(i)*vol(i)
220 norm(1) = sqrt(n1x(i)*n1x(i) + n1y(i)*n1y(i) + n1z(i)*n1z(i))
221 norm(2) = sqrt(n2x(i)*n2x(i) + n2y(i)*n2y(i) + n2z(i)*n2z(i))
222 norm(3) = sqrt(n3x(i)*n3x(i) + n3y(i)*n3y(i) + n3z(i)*n3z(i))
223 norm(4) = sqrt(n4x(i)*n4x(i) + n4y(i)*n4y(i) + n4z(i)*n4z(i))
224 norm(5) = sqrt(n5x(i)*n5x(i) + n5y(i)*n5y(i) + n5z(i)*n5z(i))
225 norm(6) = sqrt(n6x(i)*n6x(i) + n6y(i)*n6y(i) + n6z(i)*n6z(i))
226 u_n(1) = (mom(i,1)*n1x(i) + mom(i,2)*n1y(i) + mom(i,3)*n1z(i)) / (mass*norm(1))
227 u_n(2) = (mom(i,1)*n2x(i) + mom(i,2)*n2y(i) + mom(i,3)*n2z(i)) / (mass*norm(2))
228 u_n(3) = (mom(i,1)*n3x(i) + mom(i,2)*n3y(i) + mom(i,3)*n3z(i)) / (mass*norm(3))
229 u_n(4) = (mom(i,1)*n4x(i) + mom(i,2)*n4y(i) + mom(i,3)*n4z(i)) / (mass*norm(4))
230 u_n(5) = (mom(i,1)*n5x(i) + mom(i,2)*n5y(i) + mom(i,3)*n5z(i)) / (mass*norm(5))
231 u_n(6) = (mom(i,1)*n6x(i) + mom(i,2)*n6y(i) + mom(i,3)*n6z(i)) / (mass*norm(6))
232 surf(1) = half*norm(1)
233 surf(2) = half*norm(2)
234 surf(3) = half*norm(3)
235 surf(4) = half*norm(4)
236 surf(5) = half*norm(5)
237 surf(6) = half*norm(6)
238 !---storage
239 ! index1=1 centroid data
240 alefvm_buffer%F_FACE(1 ,1 ,ii) = rho(i)
241 alefvm_buffer%F_FACE(1 ,2 ,ii) = ssp(i)
242 alefvm_buffer%F_FACE(1 ,3 ,ii) = rho(i)*ssp(i)
243 alefvm_buffer%F_FACE(1 ,4 ,ii) = p(i)
244 alefvm_buffer%F_FACE(1 ,5 ,ii) = sqrt(mom(i,1)*mom(i,1)+mom(i,2)*mom(i,2)+mom(i,3)*mom(i,3))/rho(i)/vol(i) / ssp(i) ! material velocity (rho.U/rho) divided by sound speed is mach number : M=U/ssps
245 alefvm_buffer%F_FACE(1 ,6 ,ii) = zero !available
246 ! index1=2 Surface
247 alefvm_buffer%F_FACE(2 ,1:6,ii) = surf(1:6)
248 ! index1=3 Normal Velocity
249 alefvm_buffer%F_FACE(3 ,1:6,ii) = u_n(1:6)
250 ENDDO
251
252 !DEBUG-OUTPUT---------------!
253 if(alefvm_param%IOUTP_STRESS /= 0)then
254 debug_outp = .false.
255 if(alefvm_param%IOUTP_STRESS>0)then
256 do i=lft,llt
257 ii = nft + i
258 if(ixs(11,ii)==alefvm_param%IOUTP_STRESS)THEN
259 debug_outp = .true.
260 idbf = i
261 idbl = i
262 EXIT
263 endif
264 enddo
265 elseif(alefvm_param%IOUTP_STRESS==-1)then
266 debug_outp=.true.
267 idbf = lft
268 idbl = llt
269 endif
270 if(debug_outp)then
271!#!include "lockon.inc"
272 print *, " |----alefvm_stress.F-----|"
273 print *, " | THREAD INFORMATION |"
274 print *, " |------------------------|"
275 print *, " NCYCLE =", ncycle
276 do i=idbf,idbl
277 ii = nft + i
278 print *, " brique=", ixs(11,nft+i)
279 write(*,fmt='(A24,1A26)') " ",
280 . "#-stress Tensor (P+VIS+Q)#"
281
282 write (*,fmt='(A,3E26.14,A)') " | ", sig(i,1),sig(i,4),sig(i,6), " |"
283 write (*,fmt='(A,3E26.14,A)') " P =| ", sig(i,4),sig(i,2),sig(i,5), " |"
284 write (*,fmt='(A,3E26.14,A)') " |_", sig(i,6),sig(i,5),sig(i,3), "_|"
285 write (*,fmt='(A,3E26.14,A)') " | ", svis(i,1),svis(i,4),svis(i,6)," |"
286 write (*,fmt='(A,3E26.14,A)') " VIS =| ", svis(i,4),svis(i,2),svis(i,5)," |"
287 write (*,fmt='(A,3E26.14,A)') " |_", svis(i,6),svis(i,5),svis(i,3),"_|"
288 write (*,fmt='(A,3E26.14,A)') " | ", s1(i),s4(i),s6(i)," |"
289 write (*,fmt='(A,3E26.14,A)') " SIGMA = P+VIS+Q =| ", s4(i),s2(i),s5(i)," |"
290 write (*,fmt='(A,3E26.14,A)') " |_", s6(i),s5(i),s3(i),"_|"
291 write (*,fmt='(A,2E26.14)') " p = ",-third*sum(sig(i,1:3)),p(i)
292 write (*,fmt='(A,1E26.14)') " q = ",qvis(i)
293 write (*,fmt='(A,1E26.14)') " M = ",alefvm_buffer%F_FACE(1 ,5 ,ii)
294 write (*,fmt='(A,1E26.14)') " rho = ",alefvm_buffer%F_FACE(1 ,1 ,ii)
295 write (*,fmt='(A,1E26.14)') " u = ",alefvm_buffer%F_FACE(1 ,5 ,ii)*ssp(i)
296 write (*,fmt='(A,1E26.14)') " ssp = ",alefvm_buffer%F_FACE(1 ,2 ,ii)
297 write (*,fmt='(A,1E26.14)') " z = ",alefvm_buffer%F_FACE(1 ,3 ,ii)
298 write (*,fmt='(A,3E26.14)') " rho.U = ",mom(i,1:3)
299 write (*,fmt='(A,3E26.14)') " V = ",vol(i)
300 write(*,fmt='(A34,6A26)') " ",
301 . "#-------- face_1 ---------","#-------- face_2 ---------",
302 . "#-------- face_3 ---------","#-------- face_4 ---------",
303 . "#-------- face_5 ---------","#-------- face_6 --------#"
304 write (*,fmt='(A,8E26.14)') " SURF =", alefvm_buffer%F_FACE(1,1:6,ii)
305 write (*,fmt='(A,8E26.14)') " <U,n> =", alefvm_buffer%F_FACE(2,1:6,ii)
306 print *, " "
307 enddo
308!#!include "lockoff.inc"
309 endif
310 endif
311 !-----------------------------------------!
312
313 RETURN
314 END
alefvm_stress
subroutine alefvm_stress(ixs, sig, qvis, n1x, n2x, n3x, n4x, n5x, n6x, n1y, n2y, n3y, n4y, n5y, n6y, n1z, n2z, n3z, n4z, n5z, n6z, x, ipm, rho, vol, nel, mom, ssp, svis)
Definition alefvm_stress.F:38
norm
norm(diag(diag(diag(inv(mat))) -id.SOL), 2) % destroy mumps instance id.JOB
alefvm_mod
Definition alefvm_mod.F:88
alefvm_mod::alefvm_buffer
type(alefvm_buffer_), target alefvm_buffer
Definition alefvm_mod.F:120
alefvm_mod::alefvm_param
type(alefvm_param_), target alefvm_param
Definition alefvm_mod.F:121