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