OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
alefvm_aflux3.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_aflux3 ../engine/source/ale/alefvm/alefvm_aflux3.F
25!||--- called by ------------------------------------------------------
26!|| aflux0 ../engine/source/ale/aflux0.F
27!||--- uses -----------------------------------------------------
28!|| ale_connectivity_mod ../common_source/modules/ale/ale_connectivity_mod.F
29!|| alefvm_mod ../common_source/modules/ale/alefvm_mod.F
30!|| element_mod ../common_source/modules/elements/element_mod.F90
31!|| i22tri_mod ../common_source/modules/interfaces/cut-cell-search_mod.F
32!||====================================================================
33 SUBROUTINE alefvm_aflux3(PM ,IXS ,W ,FLUX ,
34 2 FLU1 ,ALE_CONNECT ,
35 3 IPM ,NV46 ,X ,
36 4 NEL )
37C-----------------------------------------------
38C D e s c r i p t i o n
39C-----------------------------------------------
40C 'alefvm' is related to a collocated scheme (built from FVM and based on Godunov scheme)
41C which was temporarily introduced for experimental option /INTER/TYPE22 (FSI coupling with cut cell method)
42C This cut cell method is not completed, abandoned, and is not an official option.
43C There is no other use for this scheme which is automatically enabled when /INTER/TYPE22 is defined (INT22>0 => IALEFVM=1).
44C
45C This subroutine is treating an uncut cell.
46C-----------------------------------------------
47C M o d u l e s
48C-----------------------------------------------
49 USE i22tri_mod
50 USE alefvm_mod
51 USE ale_connectivity_mod
52 use element_mod , only : nixs
53C-----------------------------------------------
54C I m p l i c i t T y p e s
55C-----------------------------------------------
56#include "implicit_f.inc"
57C-----------------------------------------------
58C G l o b a l P a r a m e t e r s
59C-----------------------------------------------
60#include "mvsiz_p.inc"
61C-----------------------------------------------
62C C o m m o n B l o c k s
63C-----------------------------------------------
64#include "vect01_c.inc"
65#include "param_c.inc"
66#include "com01_c.inc"
67#include "com08_c.inc"
68C-----------------------------------------------
69C D u m m y A r g u m e n t s
70C-----------------------------------------------
71 INTEGER :: IXS(NIXS,*), NV46, IPM(NPROPMI,*),NEL
72 my_real :: PM(NPROPM,*), FLUX(6,*), FLU1(*), X(3,*),W(3,*)
73 TYPE(t_ale_connectivity), INTENT(IN) :: ALE_CONNECT
74C-----------------------------------------------
75C L o c a l V a r i a b l e s
76C-----------------------------------------------
77 INTEGER MAT(MVSIZ), NC1(MVSIZ), NC2(MVSIZ), NC3(MVSIZ), NC4(MVSIZ), NC5(MVSIZ), NC6(MVSIZ),
78 . NC7(MVSIZ), NC8(MVSIZ), I, II,J,IMAT,IALEFVM_FLG,IDV, ICF(4,6),
79 . IX1,IX2,IX3,IX4
80 my_real
81 . x1(mvsiz), x2(mvsiz), x3(mvsiz), x4(mvsiz), x5(mvsiz), x6(mvsiz), x7(mvsiz), x8(mvsiz),
82 . y1(mvsiz), y2(mvsiz), y3(mvsiz), y4(mvsiz), y5(mvsiz), y6(mvsiz), y7(mvsiz), y8(mvsiz), z1(mvsiz),
83 . z2(mvsiz), z3(mvsiz), z4(mvsiz), z5(mvsiz), z6(mvsiz), z7(mvsiz), z8(mvsiz), n1x(mvsiz),
84 . n2x(mvsiz), n3x(mvsiz), n4x(mvsiz), n5x(mvsiz), n6x(mvsiz), n1y(mvsiz), n2y(mvsiz), n3y(mvsiz),
85 . n4y(mvsiz), n5y(mvsiz), n6y(mvsiz), n1z(mvsiz), n2z(mvsiz), n3z(mvsiz), n4z(mvsiz), n5z(mvsiz),
86 . n6z(mvsiz), flux1(mvsiz), flux2(mvsiz), flux3(mvsiz), flux4(mvsiz), flux5(mvsiz),
87 . flux6(mvsiz), vx1(mvsiz), vx2(mvsiz), vx3(mvsiz), vx4(mvsiz), vx5(mvsiz), vx6(mvsiz),
88 . vy1(mvsiz), vy2(mvsiz), vy3(mvsiz), vy4(mvsiz), vy5(mvsiz), vy6(mvsiz), vz1(mvsiz), vz2(mvsiz),
89 . vz3(mvsiz), vz4(mvsiz), vz5(mvsiz), vz6(mvsiz),
90 . reduc,upwl(6,mvsiz),valvois(6,nv46,mvsiz), valel(6,mvsiz),vec(3,6),cf(3,mvsiz),
91 . z(3),zadj(3),zcf(3),zzadj(3),zcf_,zzadj_,cos, ps, denom, denom1,denom2,lambda,xsi,wface(3,6,mvsiz),sr1,sr2,
92 . surf(6), term2, n(3,6,mvsiz)
93 LOGICAL debug_outp
94 INTEGER :: idbf, idbl, NC(8,MVSIZ), IAD2, LGTH
95C-----------------------------------------------
96 DATA icf/1,4,3,2,3,4,8,7,5,6,7,8,1,2,6,5,2,3,7,6,1,5,8,4/
97C-----------------------------------------------
98C D e s c r i p t i o n
99C-----------------------------------------------
100C This subroutine enables to compute face velocities
101C from centroid values [rho*U] stored in Fcell vector
102C (later used for centroid forces)
103C Face velocities are stored in vectors F_Face (ALEFVM_MOD)
104C This same vector is used later to store forces on faces
105C (used to compute centroid force)
106 !IPM(251,MAT(1)) = I_ALE_SOLVER
107 ! 0: Default = 1 expect if /ALE/MAT or /EULER/MAT has IFROM flag defined.
108 ! 1 : FEM
109 ! 2 : FVM U average
110 ! 3 : FVM rho.U average
111 ! 4 : FVM rho.c.U average
112 ! 5 : Godunov Acoustic
113 ! 6 : experimental
114C-----------------------------------------------
115C P r e - C o n d i t i o n s
116C-----------------------------------------------
117 IF(alefvm_param%IEnabled == 0) RETURN
118 imat = ixs(1,1+nft)
119 ialefvm_flg = ipm(251,imat)
120 IF(ialefvm_flg <= 1)RETURN
121C-----------------------------------------------
122C S o u r c e L i n e s
123C-----------------------------------------------
124
125 DO i=1,nel
126 ii = i + nft
127 mat(i) = ixs(1,ii)
128 !---8 local node numbers NC1 TO NC8 for solid element I ---!
129 nc1(i)=ixs(2,ii)
130 nc2(i)=ixs(3,ii)
131 nc3(i)=ixs(4,ii)
132 nc4(i)=ixs(5,ii)
133 nc5(i)=ixs(6,ii)
134 nc6(i)=ixs(7,ii)
135 nc7(i)=ixs(8,ii)
136 nc8(i)=ixs(9,ii)
137 !
138 !---Coordinates of the 8 nodes
139 x1(i)=x(1,nc1(i))
140 y1(i)=x(2,nc1(i))
141 z1(i)=x(3,nc1(i))
142 !
143 x2(i)=x(1,nc2(i))
144 y2(i)=x(2,nc2(i))
145 z2(i)=x(3,nc2(i))
146 !
147 x3(i)=x(1,nc3(i))
148 y3(i)=x(2,nc3(i))
149 z3(i)=x(3,nc3(i))
150 !
151 x4(i)=x(1,nc4(i))
152 y4(i)=x(2,nc4(i))
153 z4(i)=x(3,nc4(i))
154 !
155 x5(i)=x(1,nc5(i))
156 y5(i)=x(2,nc5(i))
157 z5(i)=x(3,nc5(i))
158 !
159 x6(i)=x(1,nc6(i))
160 y6(i)=x(2,nc6(i))
161 z6(i)=x(3,nc6(i))
162 !
163 x7(i)=x(1,nc7(i))
164 y7(i)=x(2,nc7(i))
165 z7(i)=x(3,nc7(i))
166 !
167 x8(i)=x(1,nc8(i))
168 y8(i)=x(2,nc8(i))
169 z8(i)=x(3,nc8(i))
170 ENDDO
171
172 !======================================================!
173 ! FACE VELOCITIES !
174 !======================================================!
175 DO i=1,nel
176 wface(1,1,i) = fourth*(w(1,nc1(i))+w(1,nc2(i))+w(1,nc3(i))+w(1,nc4(i)))
177 wface(2,1,i) = fourth*(w(2,nc1(i))+w(2,nc2(i))+w(2,nc3(i))+w(2,nc4(i)))
178 wface(3,1,i) = fourth*(w(3,nc1(i))+w(3,nc2(i))+w(3,nc3(i))+w(3,nc4(i)))
179
180 wface(1,2,i) = fourth*(w(1,nc3(i))+w(1,nc4(i))+w(1,nc7(i))+w(1,nc8(i)))
181 wface(2,2,i) = fourth*(w(2,nc3(i))+w(2,nc4(i))+w(2,nc7(i))+w(2,nc8(i)))
182 wface(3,2,i) = fourth*(w(3,nc3(i))+w(3,nc4(i))+w(3,nc7(i))+w(3,nc8(i)))
183
184
185 wface(1,3,i) = fourth*(w(1,nc5(i))+w(1,nc6(i))+w(1,nc7(i))+w(1,nc8(i)))
186 wface(2,3,i) = fourth*(w(2,nc5(i))+w(2,nc6(i))+w(2,nc7(i))+w(2,nc8(i)))
187 wface(3,3,i) = fourth*(w(3,nc5(i))+w(3,nc6(i))+w(3,nc7(i))+w(3,nc8(i)))
188
189 wface(1,4,i) = fourth*(w(1,nc1(i))+w(1,nc2(i))+w(1,nc5(i))+w(1,nc6(i)))
190 wface(2,4,i) = fourth*(w(2,nc1(i))+w(2,nc2(i))+w(2,nc5(i))+w(2,nc6(i)))
191 wface(3,4,i) = fourth*(w(3,nc1(i))+w(3,nc2(i))+w(3,nc5(i))+w(3,nc6(i)))
192
193 wface(1,5,i) = fourth*(w(1,nc2(i))+w(1,nc3(i))+w(1,nc6(i))+w(1,nc7(i)))
194 wface(2,5,i) = fourth*(w(2,nc2(i))+w(2,nc3(i))+w(2,nc6(i))+w(2,nc7(i)))
195 wface(3,5,i) = fourth*(w(3,nc2(i))+w(3,nc3(i))+w(3,nc6(i))+w(3,nc7(i)))
196
197 wface(1,6,i) = fourth*(w(1,nc1(i))+w(1,nc4(i))+w(1,nc5(i))+w(1,nc8(i)))
198 wface(2,6,i) = fourth*(w(2,nc1(i))+w(2,nc4(i))+w(2,nc5(i))+w(2,nc8(i)))
199 wface(3,6,i) = fourth*(w(3,nc1(i))+w(3,nc4(i))+w(3,nc5(i))+w(3,nc8(i)))
200 ENDDO
201
202 !======================================================!
203 ! INITIALIZATION : ADJAENT VALUES !
204 !======================================================!
205 !FCELL = [rhoU]. this quantity must be know for all groups
206 DO i=1,nel
207 ii = i + nft
208 valel(1,i) = alefvm_buffer%FCELL(1,ii) !rho.ux
209 valel(2,i) = alefvm_buffer%FCELL(2,ii) !rho.uy
210 valel(3,i) = alefvm_buffer%FCELL(3,ii) !rho.uz
211 valel(4,i) = alefvm_buffer%FCELL(4,ii) !rho
212 valel(5,i) = alefvm_buffer%FCELL(5,ii) !ssp
213 valel(6,i) = alefvm_buffer%FCELL(6,ii) !P
214 iad2 = ale_connect%ee_connect%iad_connect(ii)
215 lgth = ale_connect%ee_connect%iad_connect(ii+1) - iad2
216 DO j=1,lgth
217 idv = ale_connect%ee_connect%connected(iad2 + j - 1)
218 IF(idv <= 0) THEN
219 valvois(1,j,i) = -alefvm_buffer%FCELL(1,ii)
220 valvois(2,j,i) = -alefvm_buffer%FCELL(2,ii)
221 valvois(3,j,i) = -alefvm_buffer%FCELL(3,ii)
222 valvois(4,j,i) = alefvm_buffer%FCELL(4,ii)
223 valvois(5,j,i) = alefvm_buffer%FCELL(5,ii)
224 valvois(6,j,i) = alefvm_buffer%FCELL(6,ii)
225 ELSE
226 valvois(1,j,i) = alefvm_buffer%FCELL(1,idv)
227 valvois(2,j,i) = alefvm_buffer%FCELL(2,idv)
228 valvois(3,j,i) = alefvm_buffer%FCELL(3,idv)
229 valvois(4,j,i) = alefvm_buffer%FCELL(4,idv)
230 valvois(5,j,i) = alefvm_buffer%FCELL(5,idv)
231 valvois(6,j,i) = alefvm_buffer%FCELL(6,idv)
232 ENDIF
233 enddo!next J
234 enddo!next I
235
236 !======================================================!
237 ! NORMAL VECTORS ON EACH DACE !
238 ! 2S[n] = [diag1] x [diag2] !
239 ! where !
240 ! [n] : unitary normal vector on face !
241 !======================================================!
242 DO i=1,nel
243 ! Face-1
244 n1x(i)=(y3(i)-y1(i))*(z2(i)-z4(i)) - (z3(i)-z1(i))*(y2(i)-y4(i))
245 n1y(i)=(z3(i)-z1(i))*(x2(i)-x4(i)) - (x3(i)-x1(i))*(z2(i)-z4(i))
246 n1z(i)=(x3(i)-x1(i))*(y2(i)-y4(i)) - (y3(i)-y1(i))*(x2(i)-x4(i))
247 ! Face-2
248 n2x(i)=(y7(i)-y4(i))*(z3(i)-z8(i)) - (z7(i)-z4(i))*(y3(i)-y8(i))
249 n2y(i)=(z7(i)-z4(i))*(x3(i)-x8(i)) - (x7(i)-x4(i))*(z3(i)-z8(i))
250 n2z(i)=(x7(i)-x4(i))*(y3(i)-y8(i)) - (y7(i)-y4(i))*(x3(i)-x8(i))
251 ! Face-3
252 n3x(i)=(y6(i)-y8(i))*(z7(i)-z5(i)) - (z6(i)-z8(i))*(y7(i)-y5(i))
253 n3y(i)=(z6(i)-z8(i))*(x7(i)-x5(i)) - (x6(i)-x8(i))*(z7(i)-z5(i))
254 n3z(i)=(x6(i)-x8(i))*(y7(i)-y5(i)) - (y6(i)-y8(i))*(x7(i)-x5(i))
255 ! Face-4
256 n4x(i)=(y2(i)-y5(i))*(z6(i)-z1(i)) - (z2(i)-z5(i))*(y6(i)-y1(i))
257 n4y(i)=(z2(i)-z5(i))*(x6(i)-x1(i)) - (x2(i)-x5(i))*(z6(i)-z1(i))
258 n4z(i)=(x2(i)-x5(i))*(y6(i)-y1(i)) - (y2(i)-y5(i))*(x6(i)-x1(i))
259 ! Face-5
260 n5x(i)=(y7(i)-y2(i))*(z6(i)-z3(i)) - (z7(i)-z2(i))*(y6(i)-y3(i))
261 n5y(i)=(z7(i)-z2(i))*(x6(i)-x3(i)) - (x7(i)-x2(i))*(z6(i)-z3(i))
262 n5z(i)=(x7(i)-x2(i))*(y6(i)-y3(i)) - (y7(i)-y2(i))*(x6(i)-x3(i))
263 ! face-6
264 n6x(i)=(y8(i)-y1(i))*(z4(i)-z5(i)) - (z8(i)-z1(i))*(y4(i)-y5(i))
265 n6y(i)=(z8(i)-z1(i))*(x4(i)-x5(i)) - (x8(i)-x1(i))*(z4(i)-z5(i))
266 n6z(i)=(x8(i)-x1(i))*(y4(i)-y5(i)) - (y8(i)-y1(i))*(x4(i)-x5(i))
267 ! Stack
268 n(1,1,i) = n1x(i)
269 n(2,1,i) = n1y(i)
270 n(3,1,i) = n1z(i)
271 !
272 n(1,2,i) = n2x(i)
273 n(2,2,i) = n2y(i)
274 n(3,2,i) = n2z(i)
275 !
276 n(1,3,i) = n3x(i)
277 n(2,3,i) = n3y(i)
278 n(3,3,i) = n3z(i)
279 !
280 n(1,4,i) = n4x(i)
281 n(2,4,i) = n4y(i)
282 n(3,4,i) = n4z(i)
283 !
284 n(1,5,i) = n5x(i)
285 n(2,5,i) = n5y(i)
286 n(3,5,i) = n5z(i)
287 !
288 n(1,6,i) = n6x(i)
289 n(2,6,i) = n6y(i)
290 n(3,6,i) = n6z(i)
291 END DO
292
293 !======================================================!
294 ! RELATIVE VELOCITIES ON EACH FACE !
295 !======================================================!
296 !UPDATE LATER WITH LINEAR INTERPOLATION INSTEAD OF MEAN VALUE
297 !MEAN VALUE OK FOR CARTESIAN REGULAR MESH
298
299 IF(alefvm_param%IFORM/=0)ialefvm_flg = alefvm_param%IFORM !for debug purpose if set in alefvm_init.F
300
301 SELECT CASE(alefvm_param%ISOLVER)
302
303
304 CASE(1)
305 !CENTERED FVM - U average
306 DO i=1,nel
307 ii = i + nft
308 DO j=1,nv46
309 vec(1,j) = half * (valel(1,i)/valel(4,i) + valvois(1,j,i)/valvois(4,j,i))
310 vec(2,j) = half * (valel(2,i)/valel(4,i) + valvois(2,j,i)/valvois(4,j,i))
311 vec(3,j) = half * (valel(3,i)/valel(4,i) + valvois(3,j,i)/valvois(4,j,i))
312 vec(1,j) = vec(1,j) - wface(1,j,i)
313 vec(2,j) = vec(2,j) - wface(2,j,i)
314 vec(3,j) = vec(3,j) - wface(3,j,i)
315 enddo!next J
316 DO j=1,nv46
317 alefvm_buffer%F_FACE(1,j,ii) = vec(1,j)
318 alefvm_buffer%F_FACE(2,j,ii) = vec(2,j)
319 alefvm_buffer%F_FACE(3,j,ii) = vec(3,j)
320 ENDDO
321 !---face-1
322 vx1(i) = vec(1,1)
323 vy1(i) = vec(2,1)
324 vz1(i) = vec(3,1)
325 !---face-2
326 vx2(i) = vec(1,2)
327 vy2(i) = vec(2,2)
328 vz2(i) = vec(3,2)
329 !---face-3
330 vx3(i) = vec(1,3)
331 vy3(i) = vec(2,3)
332 vz3(i) = vec(3,3)
333 !---face-4
334 vx4(i) = vec(1,4)
335 vy4(i) = vec(2,4)
336 vz4(i) = vec(3,4)
337 !---face-5
338 vx5(i) = vec(1,5)
339 vy5(i) = vec(2,5)
340 vz5(i) = vec(3,5)
341 !---face-6
342 vx6(i) = vec(1,6)
343 vy6(i) = vec(2,6)
344 vz6(i) = vec(3,6)
345 enddo!next I
346
347 CASE(2)
348 !CENTERED FVM - rho.U average
349 DO i=1,nel
350 ii = i + nft
351 DO j=1,nv46
352 vec(1,j) = (valel(1,i) + valvois(1,j,i)) / (valel(4,i)+valvois(4,j,i))
353 vec(2,j) = (valel(2,i) + valvois(2,j,i)) / (valel(4,i)+valvois(4,j,i))
354 vec(3,j) = (valel(3,i) + valvois(3,j,i)) / (valel(4,i)+valvois(4,j,i))
355 vec(1,j) = vec(1,j) - wface(1,j,i)
356 vec(2,j) = vec(2,j) - wface(2,j,i)
357 vec(3,j) = vec(3,j) - wface(3,j,i)
358 enddo!next J
359 DO j=1,nv46
360 alefvm_buffer%F_FACE(1,j,ii) = vec(1,j)
361 alefvm_buffer%F_FACE(2,j,ii) = vec(2,j)
362 alefvm_buffer%F_FACE(3,j,ii) = vec(3,j)
363 ENDDO
364 !---face-1
365 vx1(i) = vec(1,1)
366 vy1(i) = vec(2,1)
367 vz1(i) = vec(3,1)
368 !---face-2
369 vx2(i) = vec(1,2)
370 vy2(i) = vec(2,2)
371 vz2(i) = vec(3,2)
372 !---face-3
373 vx3(i) = vec(1,3)
374 vy3(i) = vec(2,3)
375 vz3(i) = vec(3,3)
376 !---face-4
377 vx4(i) = vec(1,4)
378 vy4(i) = vec(2,4)
379 vz4(i) = vec(3,4)
380 !---face-5
381 vx5(i) = vec(1,5)
382 vy5(i) = vec(2,5)
383 vz5(i) = vec(3,5)
384 !---face-6
385 vx6(i) = vec(1,6)
386 vy6(i) = vec(2,6)
387 vz6(i) = vec(3,6)
388 enddo!next I
389
390 CASE(3)
391 !CENTERED FVM - Roe-average
392 DO i=1,nel
393 ii = i + nft
394 sr1 = sqrt(valel(4,i))
395 DO j=1,nv46
396 sr2 = sqrt(valvois(4,j,i))
397 vec(1,j) = (valel(1,i)/sr1 + valvois(1,j,i)/sr2) / (sr1+sr2)
398 vec(2,j) = (valel(2,i)/sr1 + valvois(2,j,i)/sr2) / (sr1+sr2)
399 vec(3,j) = (valel(3,i)/sr1 + valvois(3,j,i)/sr2) / (sr1+sr2)
400 vec(1,j) = vec(1,j) - wface(1,j,i)
401 vec(2,j) = vec(2,j) - wface(2,j,i)
402 vec(3,j) = vec(3,j) - wface(3,j,i)
403 enddo!next J
404 DO j=1,nv46
405 alefvm_buffer%F_FACE(1,j,ii) = vec(1,j)
406 alefvm_buffer%F_FACE(2,j,ii) = vec(2,j)
407 alefvm_buffer%F_FACE(3,j,ii) = vec(3,j)
408 ENDDO
409 !---face-1
410 vx1(i) = vec(1,1)
411 vy1(i) = vec(2,1)
412 vz1(i) = vec(3,1)
413 !---face-2
414 vx2(i) = vec(1,2)
415 vy2(i) = vec(2,2)
416 vz2(i) = vec(3,2)
417 !---face-3
418 vx3(i) = vec(1,3)
419 vy3(i) = vec(2,3)
420 vz3(i) = vec(3,3)
421 !---face-4
422 vx4(i) = vec(1,4)
423 vy4(i) = vec(2,4)
424 vz4(i) = vec(3,4)
425 !---face-5
426 vx5(i) = vec(1,5)
427 vy5(i) = vec(2,5)
428 vz5(i) = vec(3,5)
429 !---face-6
430 vx6(i) = vec(1,6)
431 vy6(i) = vec(2,6)
432 vz6(i) = vec(3,6)
433 enddo!next I
434
435 CASE(4)
436 !CENTERED FVM - rho.c.U average
437 DO i=1,nel
438 ii = i + nft
439 DO j=1,nv46
440 IF(dt1==zero)THEN
441 vec(1,j) = (valel(1,i) + valvois(1,j,i) )
442 . / (valel(4,i) + valvois(4,j,i) )
443 vec(2,j) = (valel(2,i) + valvois(2,j,i) )
444 . / (valel(4,i) + valvois(4,j,i) )
445 vec(3,j) = (valel(3,i) + valvois(3,j,i) )
446 . / (valel(4,i) + valvois(4,j,i) )
447 ELSE
448 vec(1,j) = (valel(1,i)*valel(5,i) + valvois(1,j,i)*valvois(5,j,i))
449 . / (valel(4,i)*valel(5,i) + valvois(4,j,i)*valvois(5,j,i))
450 vec(2,j) = (valel(2,i)*valel(5,i) + valvois(2,j,i)*valvois(5,j,i))
451 . / (valel(4,i)*valel(5,i) + valvois(4,j,i)*valvois(5,j,i))
452 vec(3,j) = (valel(3,i)*valel(5,i) + valvois(3,j,i)*valvois(5,j,i))
453 . / (valel(4,i)*valel(5,i) + valvois(4,j,i)*valvois(5,j,i))
454 ENDIF
455 vec(1,j) = vec(1,j) - wface(1,j,i)
456 vec(2,j) = vec(2,j) - wface(2,j,i)
457 vec(3,j) = vec(3,j) - wface(3,j,i)
458 enddo!next J
459 DO j=1,nv46
460 alefvm_buffer%F_FACE(1,j,ii) = vec(1,j)
461 alefvm_buffer%F_FACE(2,j,ii) = vec(2,j)
462 alefvm_buffer%F_FACE(3,j,ii) = vec(3,j)
463 ENDDO
464 !---face-1
465 vx1(i) = vec(1,1)
466 vy1(i) = vec(2,1)
467 vz1(i) = vec(3,1)
468 !---face-2
469 vx2(i) = vec(1,2)
470 vy2(i) = vec(2,2)
471 vz2(i) = vec(3,2)
472 !---face-3
473 vx3(i) = vec(1,3)
474 vy3(i) = vec(2,3)
475 vz3(i) = vec(3,3)
476 !---face-4
477 vx4(i) = vec(1,4)
478 vy4(i) = vec(2,4)
479 vz4(i) = vec(3,4)
480 !---face-5
481 vx5(i) = vec(1,5)
482 vy5(i) = vec(2,5)
483 vz5(i) = vec(3,5)
484 !---face-6
485 vx6(i) = vec(1,6)
486 vy6(i) = vec(2,6)
487 vz6(i) = vec(3,6)
488 enddo!next I
489
490 CASE(5)
491 !Godunov Acoustic for Riemann problem
492 DO i=1,nel
493 ii = i + nft
494 IF(dt1==zero)THEN
495 DO j=1,nv46
496 vec(1,j) = (valel(1,i) + valvois(1,j,i))
497 . / (valel(4,i) + valvois(4,j,i))
498 vec(2,j) = (valel(2,i) + valvois(2,j,i))
499 . / (valel(4,i) + valvois(4,j,i))
500 vec(3,j) = (valel(3,i) + valvois(3,j,i))
501 . / (valel(4,i) + valvois(4,j,i))
502 vec(1,j) = vec(1,j) - wface(1,j,i)
503 vec(2,j) = vec(2,j) - wface(2,j,i)
504 vec(3,j) = vec(3,j) - wface(3,j,i)
505 enddo!next J
506 ELSE
507 surf(1) = one/sqrt(n1x(i)*n1x(i)+n1y(i)*n1y(i)+n1z(i)*n1z(i))
508 surf(2) = one/sqrt(n2x(i)*n2x(i)+n2y(i)*n2y(i)+n2z(i)*n2z(i))
509 surf(3) = one/sqrt(n3x(i)*n3x(i)+n3y(i)*n3y(i)+n3z(i)*n3z(i))
510 surf(4) = one/sqrt(n4x(i)*n4x(i)+n4y(i)*n4y(i)+n4z(i)*n4z(i))
511 surf(5) = one/sqrt(n5x(i)*n5x(i)+n5y(i)*n5y(i)+n5z(i)*n5z(i))
512 surf(6) = one/sqrt(n6x(i)*n6x(i)+n6y(i)*n6y(i)+n6z(i)*n6z(i))
513 DO j=1,nv46
514 term2 = surf(j) * ( valel(6,i)-valvois(6,j,i) ) / (valel(4,i)*valel(5,i)+valvois(4,j,i)*valvois(5,j,i))
515 !TERM2 = TERM2 * ZERO
516 vec(1,j) = (valel(1,i)*valel(5,i) + valvois(1,j,i)*valvois(5,j,i))
517 . / (valel(4,i)*valel(5,i)+valvois(4,j,i)*valvois(5,j,i)) + term2*n(1,j,i)
518 vec(2,j) = (valel(2,i)*valel(5,i) + valvois(2,j,i)*valvois(5,j,i))
519 . / (valel(4,i)*valel(5,i)+valvois(4,j,i)*valvois(5,j,i)) + term2*n(2,j,i)
520 vec(3,j) = (valel(3,i)*valel(5,i) + valvois(3,j,i)*valvois(5,j,i))
521 . / (valel(4,i)*valel(5,i)+valvois(4,j,i)*valvois(5,j,i)) + term2*n(3,j,i)
522 vec(1,j) = vec(1,j) - wface(1,j,i)
523 vec(2,j) = vec(2,j) - wface(2,j,i)
524 vec(3,j) = vec(3,j) - wface(3,j,i)
525 enddo!next J
526 ENDIF
527 DO j=1,nv46
528 alefvm_buffer%F_FACE(1,j,ii) = vec(1,j)
529 alefvm_buffer%F_FACE(2,j,ii) = vec(2,j)
530 alefvm_buffer%F_FACE(3,j,ii) = vec(3,j)
531 ENDDO
532 !---face-1
533 vx1(i) = vec(1,1)
534 vy1(i) = vec(2,1)
535 vz1(i) = vec(3,1)
536 !---face-2
537 vx2(i) = vec(1,2)
538 vy2(i) = vec(2,2)
539 vz2(i) = vec(3,2)
540 !---face-3
541 vx3(i) = vec(1,3)
542 vy3(i) = vec(2,3)
543 vz3(i) = vec(3,3)
544 !---face-4
545 vx4(i) = vec(1,4)
546 vy4(i) = vec(2,4)
547 vz4(i) = vec(3,4)
548 !---face-5
549 vx5(i) = vec(1,5)
550 vy5(i) = vec(2,5)
551 vz5(i) = vec(3,5)
552 !---face-6
553 vx6(i) = vec(1,6)
554 vy6(i) = vec(2,6)
555 vz6(i) = vec(3,6)
556 enddo!next I
557
558 CASE(6)
559 !INTERPOLATED
560 !WARNING : ADD SYMMETRY LATER IF FORMULATION IS KEPT (mean value for identical value on both sides)
561 DO i=1,nel
562 ii = i + nft
563 nc(1,i) = ixs(2,ii)
564 nc(2,i) = ixs(3,ii)
565 nc(3,i) = ixs(4,ii)
566 nc(4,i) = ixs(5,ii)
567 nc(5,i) = ixs(6,ii)
568 nc(6,i) = ixs(7,ii)
569 nc(7,i) = ixs(8,ii)
570 nc(8,i) = ixs(9,ii)
571 z(1) = one_over_8*sum(x(1,nc(1:8,i)))
572 z(2) = one_over_8*sum(x(2,nc(1:8,i)))
573 z(3) = one_over_8*sum(x(3,nc(1:8,i)))
574 iad2 = ale_connect%ee_connect%iad_connect(ii)
575 lgth = ale_connect%ee_connect%iad_connect(ii + 1) - iad2
576 DO j=1,lgth
577 idv = ale_connect%ee_connect%connected(iad2 + j - 1)
578 IF(idv<=0)THEN
579 vec(1:3,j) = zero
580 cycle
581 ENDIF
582 ix1=ixs(icf(1,j)+1,ii)
583 ix2=ixs(icf(2,j)+1,ii)
584 ix3=ixs(icf(3,j)+1,ii)
585 ix4=ixs(icf(4,j)+1,ii)
586 cf(1,i) = fourth*(x(1,ix1)+x(1,ix2)+x(1,ix3)+x(1,ix4))
587 cf(2,i) = fourth*(x(2,ix1)+x(2,ix2)+x(2,ix3)+x(2,ix4))
588 cf(3,i) = fourth*(x(3,ix1)+x(3,ix2)+x(3,ix3)+x(3,ix4))
589 nc(1,i)=ixs(2,idv)
590 nc(2,i)=ixs(3,idv)
591 nc(3,i)=ixs(4,idv)
592 nc(4,i)=ixs(5,idv)
593 nc(5,i)=ixs(6,idv)
594 nc(6,i)=ixs(7,idv)
595 nc(7,i)=ixs(8,idv)
596 nc(8,i)=ixs(9,idv)
597 zadj(1) = one_over_8*sum(x(1,nc(1:8,i)))
598 zadj(2) = one_over_8*sum(x(2,nc(1:8,i)))
599 zadj(3) = one_over_8*sum(x(3,nc(1:8,i)))
600 zzadj(1) = zadj(1)-z(1)
601 zzadj(2) = zadj(2)-z(2)
602 zzadj(3) = zadj(3)-z(3)
603 zcf(1) = cf(1,i) - z(1)
604 zcf(2) = cf(2,i) - z(2)
605 zcf(3) = cf(3,i) - z(3)
606 ps = zcf(1)*zzadj(1) + zcf(2)*zzadj(2) + zcf(3)*zzadj(3)
607 zcf_ = zcf(1)**2 + zcf(2)**2 + zcf(3)**2
608 zzadj_ = zzadj(1)**2 + zzadj(2)**2 + zzadj(3)**2
609 denom1 = sqrt(zcf_)
610 denom2 = sqrt(zzadj_)
611 denom = denom1*denom2
612 cos = ps/denom
613 xsi = sqrt(zcf_) * cos
614 lambda = xsi / denom2
615 !interpoler VALEL -> VALVOIS utilisant lambda
616 vec(1,j) = valel(1,i) + lambda*(valvois(1,j,i)/valvois(4,j,i)-valel(1,i)/valel(4,i)) - wface(1,j,i)
617 vec(2,j) = valel(2,i) + lambda*(valvois(2,j,i)/valvois(4,j,i)-valel(2,i)/valel(4,i)) - wface(2,j,i)
618 vec(3,j) = valel(3,i) + lambda*(valvois(3,j,i)/valvois(4,j,i)-valel(3,i)/valel(4,i)) - wface(3,j,i)
619 enddo!next J
620 DO j=1,nv46
621 alefvm_buffer%F_FACE(1,j,ii) = vec(1,j)
622 alefvm_buffer%F_FACE(2,j,ii) = vec(2,j)
623 alefvm_buffer%F_FACE(3,j,ii) = vec(3,j)
624 ENDDO
625 !---face-1
626 vx1(i) = vec(1,1)
627 vy1(i) = vec(2,1)
628 vz1(i) = vec(3,1)
629 !---face-2
630 vx2(i) = vec(1,2)
631 vy2(i) = vec(2,2)
632 vz2(i) = vec(3,2)
633 !---face-3
634 vx3(i) = vec(1,3)
635 vy3(i) = vec(2,3)
636 vz3(i) = vec(3,3)
637 !---face-4
638 vx4(i) = vec(1,4)
639 vy4(i) = vec(2,4)
640 vz4(i) = vec(3,4)
641 !---face-5
642 vx5(i) = vec(1,5)
643 vy5(i) = vec(2,5)
644 vz5(i) = vec(3,5)
645 !---face-6
646 vx6(i) = vec(1,6)
647 vy6(i) = vec(2,6)
648 vz6(i) = vec(3,6)
649 enddo!next I
650
651 END SELECT !I_ALE_SOLVER
652
653 !======================================================!
654 ! FLUXES CALCULATION ON EACH FACE !
655 ! FLUX_face = [0.5*V_face] . [2S*n] !
656 !======================================================!
657 DO i=1,nel
658 flux1(i)=half*(vx1(i)*n1x(i)+vy1(i)*n1y(i)+vz1(i)*n1z(i))
659 flux2(i)=half*(vx2(i)*n2x(i)+vy2(i)*n2y(i)+vz2(i)*n2z(i))
660 flux3(i)=half*(vx3(i)*n3x(i)+vy3(i)*n3y(i)+vz3(i)*n3z(i))
661 flux4(i)=half*(vx4(i)*n4x(i)+vy4(i)*n4y(i)+vz4(i)*n4z(i))
662 flux5(i)=half*(vx5(i)*n5x(i)+vy5(i)*n5y(i)+vz5(i)*n5z(i))
663 flux6(i)=half*(vx6(i)*n6x(i)+vy6(i)*n6y(i)+vz6(i)*n6z(i))
664 ENDDO
665
666 !DEBUG-OUTPUT---------------!
667 if(alefvm_param%IOUTP_FLUX /= 0)then
668 debug_outp = .false.
669 if(alefvm_param%IOUTP_FLUX>0)then
670 do i=lft,llt
671 ii = nft + i
672 if(ixs(11,ii)==alefvm_param%IOUTP_FLUX)THEN
673 debug_outp = .true.
674 idbf = i
675 idbl = i
676 EXIT
677 endif
678 enddo
679 elseif(alefvm_param%IOUTP_FLUX==-1)then
680 debug_outp=.true.
681 idbf = lft
682 idbl = llt
683 endif
684!#!include "lockon.inc"
685 if(debug_outp)then
686!#!include "lockon.inc"
687 print *, " |----alefvm_aflux3.F-----|"
688 print *, " | THREAD INFORMATION |"
689 print *, " |------------------------|"
690 print *, " NCYCLE =", ncycle
691 do i=idbf,idbl
692 ii = nft + i
693 print *, " brique=", ixs(11,nft+i)
694 write (*,fmt='(A,6E26.14)') " Flux(1:6) =", flux1(i),flux2(i),flux3(i),flux4(i),flux5(i),flux6(i)
695 write(*,fmt='(A24,1A20)') " ",
696 . "#--------- cell-----------"
697 write (*,fmt='(A,1E26.14)') " V_cell-X =", alefvm_buffer%FCELL(1,ii)
698 write (*,fmt='(A,1E26.14)') " V_cell-Y =", alefvm_buffer%FCELL(2,ii)
699 write (*,fmt='(A,1E26.14)') " V_cell-Z =", alefvm_buffer%FCELL(3,ii)
700 write(*,fmt='(A24,6A26)') " ",
701 . "#--------- adj_1 ---------","#--------- adj_2 ---------",
702 . "#--------- adj_3 ---------","#--------- adj_4 ---------",
703 . "#--------- adj_5 ---------","#--------- adj_6 --------#"
704 write (*,fmt='(A,6E26.14)') " V_faces-X =", alefvm_buffer%F_FACE(1,1:6,ii)
705 write (*,fmt='(A,6E26.14)') " V_faces-Y =", alefvm_buffer%F_FACE(2,1:6,ii)
706 write (*,fmt='(A,6E26.14)') " V_faces-Z =", alefvm_buffer%F_FACE(3,1:6,ii)
707 print *, " "
708 enddo
709!#!include "lockoff.inc"
710 endif
711 endif
712 !-----------------------------------------!
713
714
715 !======================================================!
716 ! TRIMATERIAL CASE INITIALIZATION (LAW51) !
717 ! -->RETURN !
718 !======================================================!
719 IF(nint(pm(19,mat(1)))==51)THEN
720 DO i=1,nel
721 flux(1,i)=flux1(i)
722 flux(2,i)=flux2(i)
723 flux(3,i)=flux3(i)
724 flux(4,i)=flux4(i)
725 flux(5,i)=flux5(i)
726 flux(6,i)=flux6(i)
727 ENDDO !next I
728 RETURN
729 ENDIF
730
731 !======================================================!
732 ! BOUNDARY FACE : no volume flux by default !
733 ! slip wall bc !
734 !======================================================!
735 DO j=1,6
736 DO i=1,nel
737 upwl(j,i)=pm(16,mat(i))
738 ENDDO !next I
739 ENDDO !next J
740
741 DO i=1,nel
742 reduc=pm(92,mat(i))
743 iad2 = ale_connect%ee_connect%iad_connect(i + nft)
744 !---face1---!
745 ii=ale_connect%ee_connect%connected(iad2 + 1 - 1)
746 IF(ii==0)THEN
747 flux1(i)=flux1(i)*reduc
748 ENDIF
749 !---face2---!
750 ii=ale_connect%ee_connect%connected(iad2 + 2 - 1)
751 IF(ii==0)THEN
752 flux2(i)=flux2(i)*reduc
753 ENDIF
754 !---face3---!
755 ii=ale_connect%ee_connect%connected(iad2 + 3 - 1)
756 IF(ii==0)THEN
757 flux3(i)=flux3(i)*reduc
758 ENDIF
759 !---face4---!
760 ii=ale_connect%ee_connect%connected(iad2 + 4 - 1)
761 IF(ii==0)THEN
762 flux4(i)=flux4(i)*reduc
763 ENDIF
764 !---face5---!
765 ii=ale_connect%ee_connect%connected(iad2 + 5 - 1)
766 IF(ii==0)THEN
767 flux5(i)=flux5(i)*reduc
768 ENDIF
769 !---face6---!
770 ii=ale_connect%ee_connect%connected(iad2 + 6 - 1)
771 IF(ii==0)THEN
772 flux6(i)=flux6(i)*reduc
773 ENDIF
774 ENDDO !next I
775
776 DO i=1,nel
777 flux(1,i)=flux1(i)-upwl(1,i)*abs(flux1(i))
778 flux(2,i)=flux2(i)-upwl(2,i)*abs(flux2(i))
779 flux(3,i)=flux3(i)-upwl(3,i)*abs(flux3(i))
780 flux(4,i)=flux4(i)-upwl(4,i)*abs(flux4(i))
781 flux(5,i)=flux5(i)-upwl(5,i)*abs(flux5(i))
782 flux(6,i)=flux6(i)-upwl(6,i)*abs(flux6(i))
783
784 flu1(i) =flux1(i)+upwl(1,i)*abs(flux1(i))
785 . +flux2(i)+upwl(2,i)*abs(flux2(i))
786 . +flux3(i)+upwl(3,i)*abs(flux3(i))
787 . +flux4(i)+upwl(4,i)*abs(flux4(i))
788 . +flux5(i)+upwl(5,i)*abs(flux5(i))
789 . +flux6(i)+upwl(6,i)*abs(flux6(i))
790 ENDDO !next I
791
792
793
794
795
796 RETURN
797 END
alefvm_aflux3
subroutine alefvm_aflux3(pm, ixs, w, flux, flu1, ale_connect, ipm, nv46, x, nel)
Definition alefvm_aflux3.F:37
ale_connectivity_mod
Definition ale_connectivity_mod.F:225
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
i22tri_mod
Definition cut-cell-search_mod.F:212