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