OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
eflux3_int22_fvm.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!|| eflux3_int22_fvm ../engine/source/ale/alefvm/cut_cells/eflux3_int22_fvm.F
25!||--- called by ------------------------------------------------------
26!|| aflux0 ../engine/source/ale/aflux0.F
27!||--- calls -----------------------------------------------------
28!|| my_barrier ../engine/source/system/machine.F
29!||--- uses -----------------------------------------------------
30!|| alefvm_mod ../common_source/modules/ale/alefvm_mod.F
31!|| elbufdef_mod ../common_source/modules/mat_elem/elbufdef_mod.F90
32!|| i22bufbric_mod ../common_source/modules/interfaces/cut-cell-search_mod.F
33!|| i22tri_mod ../common_source/modules/interfaces/cut-cell-search_mod.f
34!||====================================================================
35 SUBROUTINE eflux3_int22_fvm(
36 . PM , IXS , FLUX , FLU1 ,
37 . IPARG , ELBUF_TAB, ITASK ,
38 . NV46 , IPM , X)
39C-----------------------------------------------
40C D e s c r i p t i o n
41C-----------------------------------------------
42C 'alefvm' is related to a collocated scheme (built from FVM and based on Godunov scheme)
43C which was temporarily introduced for experimental option /INTER/TYPE22 (FSI coupling with cut cell method)
44C This cut cell method is not completed, abandoned, and is not an official option.
45C There is no other use for this scheme which is automatically enabled when /INTER/TYPE22 is defined (INT22>0 => IALEFVM=1).
46C
47C This subroutine is computing fluxes for cut cells.
48C Result are stored in cut cell buffer (currently called BRICK_LIST)
49C Main Cell flux is already computed with the usual
50C subroutine (EFLUX3)
51C All Cell not in cut cell buffer only have uncut adjacent cells. Since
52C a cut cell has all its adjacent cells in the cut cell buffer.
53C Cells in cut cell buffer might have adjacent cut cells.
54 !IPM(251,MAT(1)) = I_ALE_SOLVER
55 ! 0: Default = 1 expect if /ALE/MAT or /EULER/MAT has IFROM flag defined.
56 ! 1 : FEM
57 ! 2 : FVM U average
58 ! 3 : FVM rho.U average
59 ! 4 : FVM rho.c.U average
60 ! 5 : Godunov Acoustic
61 ! 6 : experimental
62C-----------------------------------------------
63C M o d u l e s
64C-----------------------------------------------
65 USE i22bufbric_mod
66 USE i22tri_mod
67 USE elbufdef_mod
68 USE alefvm_mod
69C-----------------------------------------------
70C I m p l i c i t T y p e s
71C-----------------------------------------------
72#include "implicit_f.inc"
73C-----------------------------------------------
74C G l o b a l P a r a m e t e r s
75C-----------------------------------------------
76#include "mvsiz_p.inc"
77C-----------------------------------------------
78C C o m m o n B l o c k s
79C-----------------------------------------------
80#include "param_c.inc"
81#include "inter22.inc"
82#include "task_c.inc"
83#include "com01_c.inc"
84#include "com08_c.inc"
85C-----------------------------------------------
86C D u m m y A r g u m e n t s
87C-----------------------------------------------
88 INTEGER :: IXS(NIXS,*), IPARG(NPARG,*),ISILENT, NV46,IPM(NPROPMI,*)
89 my_real :: PM(NPROPM,*),FLUX(6,*), FLU1(*),X(3,*)
90 TYPE (ELBUF_STRUCT_), DIMENSION(NGROUP),TARGET :: ELBUF_TAB
91C-----------------------------------------------
92C L o c a l V a r i a b l e s
93C-----------------------------------------------
94 INTEGER MAT, MLW, NC(8), I, IE,J, K, Kv,L,ITASK, IDm, IDV
95 INTEGER IDvm, IBvm, NGvm, IEvm, IEV, Jm, IMAT, IALEFVM_FLG
96 INTEGER IB,IBv,NIN, NBCUT,ICELL,NCELL,NGm
97 my_real cellflux(6,9,nb,5),reduc,upwl(6)
98 my_real :: vf(3), norm(3,6), lnorm(6),term2
99 INTEGER :: NBF,NBL, MCELL,iNOD,ICELLv,numnod, numnodV
100 TYPE(g_bufel_) ,POINTER :: GBUF
101 my_real :: face , z(3), zadj(3), zzadj_, cf(3), zcf(3),zzadj(3)
102 my_real :: ps, lambda
103 INTEGER :: NUM, NADJ, IADJ, JV, NG, IE_M, IEm, IBm
104 my_real :: facev, ddvol, valel(6), valvois(6,6,5), sr1, sr2, srf
105 LOGICAL debug_outp
106C-----------------------------------------------
107C P r e - C o n d i t i o n
108C-----------------------------------------------
109 IF(int22 == 0)RETURN
110C-----------------------------------------------
111C S o u r c e L i n e s
112C-----------------------------------------------
113
114 !======================================================!
115 ! INITIALIZATIONS !
116 !======================================================!
117 ibvm = 0
118 vf = zero
119 nin = 1
120 nbf = 1+itask*nb/nthread
121 nbl = (itask+1)*nb/nthread
122 nbl = min(nbl,nb)
123
124 !======================================================!
125 ! DEBUG OUTPUT !
126 !======================================================!
127 !INTERFACE 22 ONLY - OUTPUT---------------!
128 debug_outp = .false.
129 if(ibug22_flux22/=0)then
130 debug_outp = .false.
131 if(ibug22_flux22>0)then
132 do ib=nbf,nbl
133 ie = brick_list(nin,ib)%id
134 if(ixs(11,ie)==ibug22_flux22)debug_outp=.true.
135 enddo
136 elseif(ibug22_flux22==-1)then
137 debug_outp = .true.
138 endif
139 endif
140 if(debug_outp .AND. itask==0)then
141 print *, " |---------eflux3_int22_fvm.F---------|"
142 print *, " | THREAD INFORMATION |"
143 print *, " |------------------------------------|"
144 print *, " NCYCLE =", ncycle
145 endif
146 !INTERFACE 22 ONLY - OUTPUT---------------!
147
148
149
150 !======================================================!
151 ! STEP A : CUT CELL FLUXES !
152 ! A.1 : main & SECND CELLS FLUXES !
153 !======================================================!
154 DO ib=nbf,nbl
155 ie = brick_list(nin,ib)%ID
156 mlw = brick_list(nin,ib)%MLW
157 ncell = brick_list(nin,ib)%NBCUT
158 mcell = brick_list(nin,ib)%MainID
159 icell = 0
160
161 !======================================================!
162 ! RESET FLUXES !
163 !======================================================!
164 flux(1:6,ie) = zero
165 flu1(ie) = zero
166
167 !======================================================!
168 ! NORMAL VECTORS ON EACH DACE !
169 ! 2S[n] = [diag1] x [diag2] !
170 ! where !
171 ! [n] : unitary normal vector on face !
172 !======================================================!
173 norm(1,1:6) = brick_list(nin,ib)%N(1:6,1) !VEUL(14:19,IE)
174 norm(2,1:6) = brick_list(nin,ib)%N(1:6,2) !VEUL(20:25,IE)
175 norm(3,1:6) = brick_list(nin,ib)%N(1:6,3) !VEUL(26:31,IE)
176 DO j=1,nv46
177 lnorm(j) = sqrt( norm(1,j)**2 + norm(2,j)**2 + norm(3,j)**2 )
178 norm(1:3,j) = norm(1:3,j) / lnorm(j)
179 ENDDO
180
181 !======================================================!
182 !COMPUTES LOCAL(VALEL) AND ADJACENT (VALVOIS) MOMENTUMS!
183 !======================================================!
184 DO WHILE (icell<=ncell) ! loop on polyhedron {1:NCELL} U {9}
185 icell = icell +1
186 IF (icell>ncell .AND. ncell/=0)icell=9
187 nbcut = brick_list(nin,ib)%NBCUT
188 jm = brick_list(nin,ib)%POLY(icell)%WhereIsMain(1)
189 iem = brick_list(nin,ib)%POLY(icell)%WhereIsMain(3)
190 ibm = brick_list(nin,ib)%POLY(icell)%WhereIsMain(4)
191 idm = brick_list(nin,ibm)%ID
192 ngm = brick_list(nin,ibm)%NG
193 gbuf =>elbuf_tab(ngm)%GBUF
194 valel(1) = alefvm_buffer%FCELL(1,iem) !rho.ux
195 valel(2) = alefvm_buffer%FCELL(2,iem) !rho.uy
196 valel(3) = alefvm_buffer%FCELL(3,iem) !rho.uz
197 valel(4) = alefvm_buffer%FCELL(4,iem) !rho
198 valel(5) = alefvm_buffer%FCELL(5,iem) !ssp
199 valel(6) = alefvm_buffer%FCELL(6,iem) !P
200 sr1 = sqrt(valel(4))
201 DO j=1,nv46
202 nadj = brick_list(nin,ib)%POLY(icell)%FACE(j)%NAdjCell
203 iev = brick_list(nin,ib)%Adjacent_Brick(j,1)
204 ibv = brick_list(nin,ib)%Adjacent_Brick(j,4)
205 jv = brick_list(nin,ib)%Adjacent_Brick(j,5)
206 cellflux(j,icell,ib,1:5) = zero
207 DO iadj=1,nadj
208 icellv = brick_list(nin,ib)%POLY(icell)%FACE(j)%Adjacent_Cell(iadj)
209 IF(icellv==0)THEN
210 valvois(j,1:3,iadj) = -valel(1:3)
211 valvois(j,4,iadj) = valel(4)
212 ELSE
213 !IF INSIDE CUT CELL
214 IF(ibv/=0)THEN
215 ievm = brick_list(nin,ibv)%POLY(icellv)%WhereIsMain(3)
216 ibvm = brick_list(nin,ibv)%POLY(icellv)%WhereIsMain(4)
217 ngvm = brick_list(nin,ibvm)%NG
218 idvm = brick_list(nin,ibvm)%IDLOC
219 valvois(j,1,iadj) = alefvm_buffer%FCELL(1,ievm)
220 valvois(j,2,iadj) = alefvm_buffer%FCELL(2,ievm)
221 valvois(j,3,iadj) = alefvm_buffer%FCELL(3,ievm)
222 valvois(j,4,iadj) = alefvm_buffer%FCELL(4,ievm)
223 valvois(j,5,iadj) = alefvm_buffer%FCELL(5,ievm)
224 valvois(j,6,iadj) = alefvm_buffer%FCELL(6,ievm)
225 ELSE
226 !NOT IN CUT CELL (FRONTIER TO USUAL DOMAIN)
227 valvois(j,1,iadj) = alefvm_buffer%FCELL(1,iev)
228 valvois(j,2,iadj) = alefvm_buffer%FCELL(2,iev)
229 valvois(j,3,iadj) = alefvm_buffer%FCELL(3,iev)
230 valvois(j,4,iadj) = alefvm_buffer%FCELL(4,iev)
231 valvois(j,5,iadj) = alefvm_buffer%FCELL(5,iev)
232 valvois(j,6,iadj) = alefvm_buffer%FCELL(6,iev)
233 ENDIF
234 ENDIF
235 sr2 = sqrt(valvois(j,4,iadj))
236 !======================================================!
237 ! FLUXES CALCULATION ON EACH FACE(J) AND EACH ADJ CELL !
238 !======================================================!
239 imat = ixs(1,ie)
240 ialefvm_flg = ipm(251,imat)
241 IF(alefvm_param%IFORM/=0)ialefvm_flg = alefvm_param%IFORM !for debug purpose if set in alefvm_init.F
242
243
244 SELECT CASE(alefvm_param%ISOLVER)
245 CASE(1)
246 !CENTERED FVM - U average
247 vf(1) = half * (valel(1)/valel(4)+valvois(j,1,iadj)/valvois(j,4,iadj))
248 vf(2) = half * (valel(2)/valel(4)+valvois(j,2,iadj)/valvois(j,4,iadj))
249 vf(3) = half * (valel(3)/valel(4)+valvois(j,3,iadj)/valvois(j,4,iadj))
250 CASE(2)
251 !CENTERED FVM - rho.U average
252 vf(1) = (valel(1)+valvois(j,1,iadj))/(valel(4)+valvois(j,4,iadj))
253 vf(2) = (valel(2)+valvois(j,2,iadj))/(valel(4)+valvois(j,4,iadj))
254 vf(3) = (valel(3)+valvois(j,3,iadj))/(valel(4)+valvois(j,4,iadj))
255 CASE(3)
256 !CENTERED FVM - Roe-average
257 vf(1) = (valel(1)/sr1+valvois(j,1,iadj)/sr2)/(sr1+sr2)
258 vf(2) = (valel(2)/sr1+valvois(j,2,iadj)/sr2)/(sr1+sr2)
259 vf(3) = (valel(3)/sr1+valvois(j,3,iadj)/sr2)/(sr1+sr2)
260 CASE(4)
261 IF(dt1==zero)THEN
262 vf(1) = (valel(1) +valvois(j,1,iadj) )
263 . /(valel(4) +valvois(j,4,iadj) )
264 vf(2) = (valel(2) +valvois(j,2,iadj) )
265 . /(valel(4) +valvois(j,4,iadj) )
266 vf(3) = (valel(3) +valvois(j,3,iadj) )
267 . /(valel(4) +valvois(j,4,iadj) )
268 ELSE
269
270 !CENTERED FVM - rho.c.U average
271 vf(1) = (valel(1)*valel(5)+valvois(j,1,iadj)*valvois(j,5,iadj))
272 . /(valel(4)*valel(5)+valvois(j,4,iadj)*valvois(j,5,iadj))
273 vf(2) = (valel(2)*valel(5)+valvois(j,2,iadj)*valvois(j,5,iadj))
274 . /(valel(4)*valel(5)+valvois(j,4,iadj)*valvois(j,5,iadj))
275 vf(3) = (valel(3)*valel(5)+valvois(j,3,iadj)*valvois(j,5,iadj))
276 . /(valel(4)*valel(5)+valvois(j,4,iadj)*valvois(j,5,iadj))
277
278 ENDIF
279
280 CASE(5)
281 !Godunoc/Riemann Acoustic
282 IF(dt1==zero)THEN
283 vf(1) = (valel(1) +valvois(j,1,iadj) )
284 . /(valel(4) +valvois(j,4,iadj) )
285 vf(2) = (valel(2) +valvois(j,2,iadj) )
286 . /(valel(4) +valvois(j,4,iadj) )
287 vf(3) = (valel(3) +valvois(j,3,iadj) )
288 . /(valel(4) +valvois(j,4,iadj) )
289 ELSE
290 term2 = ( valel(6)-valvois(j,6,iadj) )/ (valel(4)*valel(5)+valvois(j,4,iadj)*valvois(j,5,iadj))
291 vf(1) = (valel(1)*valel(5)+valvois(j,1,iadj)*valvois(j,5,iadj))
292 . /(valel(4)*valel(5)+valvois(j,4,iadj)*valvois(j,5,iadj)) + term2 *norm(1,j)
293 vf(2) = (valel(2)*valel(5)+valvois(j,2,iadj)*valvois(j,5,iadj))
294 . /(valel(4)*valel(5)+valvois(j,4,iadj)*valvois(j,5,iadj)) + term2 *norm(2,j)
295 vf(3) = (valel(3)*valel(5)+valvois(j,3,iadj)*valvois(j,5,iadj))
296 . /(valel(4)*valel(5)+valvois(j,4,iadj)*valvois(j,5,iadj)) + term2 *norm(3,j)
297 ENDIF
298 CASE(6)
299 z(1:3) = brick_list(nin,ibm)%SCellCenter(1:3)
300 IF(ibv /=0)THEN
301 zadj(1:3) = brick_list(nin,ibvm)%ScellCenter(1:3)
302 ELSE
303 nc(1) = ixs(2,iev)
304 nc(2) = ixs(3,iev)
305 nc(3) = ixs(4,iev)
306 nc(4) = ixs(5,iev)
307 nc(5) = ixs(6,iev)
308 nc(6) = ixs(7,iev)
309 nc(7) = ixs(8,iev)
310 nc(8) = ixs(9,iev)
311 zadj(1) = one_over_8*sum(x(1,nc(1:8)))
312 zadj(2) = one_over_8*sum(x(2,nc(1:8)))
313 zadj(3) = one_over_8*sum(x(3,nc(1:8)))
314 ENDIF
315
316 cf(1:3) = brick_list(nin,ib)%POLY(icell)%FACE(j)%Center(1:3)
317 IF(ibv /= 0)THEN
318 face = brick_list(nin,ib)%POLY(icell)%FACE(j)%Surf
319 facev = brick_list(nin,ibv)%POLY(icellv)%FACE(jv)%Surf
320 IF(facev<face) cf(1:3) = brick_list(nin,ibv)%POLY(icellv)%FACE(jv)%Center(1:3)
321 ENDIF
322 zzadj(1) = zadj(1)-z(1)
323 zzadj(2) = zadj(2)-z(2)
324 zzadj(3) = zadj(3)-z(3)
325 zcf(1) = cf(1) - z(1)
326 zcf(2) = cf(2) - z(2)
327 zcf(3) = cf(3) - z(3)
328 ps = zcf(1)*zzadj(1) + zcf(2)*zzadj(2) + zcf(3)*zzadj(3)
329 zzadj_ = zzadj(1)**2 + zzadj(2)**2 + zzadj(3)**2
330 lambda = ps / max(em20,zzadj_)
331 IF(lambda<zero .OR. lambda >one)then
332 print *, "labmda=", lambda
333c pause
334 endif
335 lambda = min(max(zero,lambda) , one)
336 lambda = sin(half*3.14159265358979d00*lambda)
337 lambda = lambda * lambda
338 !face density
339 sr1 = valel(4)
340 sr2 = valvois(j,4,iadj)
341 srf = sr1 + lambda*(sr2-sr1)
342 !face momentum density
343 vf(1) = valel(1) + lambda*(valvois(j,1,iadj)-valel(1))
344 vf(2) = valel(2) + lambda*(valvois(j,2,iadj)-valel(2))
345 vf(3) = valel(3) + lambda*(valvois(j,3,iadj)-valel(3))
346 !face velocity
347 vf(1) = vf(1) / srf
348 vf(2) = vf(2) / srf
349 vf(3) = vf(3) / srf
350 !CASE(2) gives
351 !VF(1) = (VALEL(1)+VALVOIS(J,1,IADj))/(VALEL(4)+VALVOIS(J,4,IADj))
352 !VF(2) = (VALEL(2)+VALVOIS(J,2,IADj))/(VALEL(4)+VALVOIS(J,4,IADj))
353 !VF(3) = (VALEL(3)+VALVOIS(J,3,IADj))/(VALEL(4)+VALVOIS(J,4,IADj))
354 END SELECT
355
356
357 cellflux(j,icell,ib,iadj) = (vf(1)*norm(1,j) + vf(2)*norm(2,j) + vf(3)*norm(3,j))
358 !if(ncycle>=1 )then
359 !if(ixs(11,ie)==1802081 )then
360 ! print *, "brick =", ixs(11,ie)
361 ! print *, "brickV =", ixs(11,brick_list(nin,ibv)%id), brick_list(nin,ibv)%id
362 ! print *, "eflux3 - brickID=", ixs(11,ie), icell
363 ! print *, "VF =", VF(1:3)
364 ! print *, "VALEL =", VALEL(1:3)
365 ! print *, "LAMBDA =", LAMBDA
366 ! print *, "VALVOIS =", VALVOIS(J,1:3,IADj)
367 ! print *, "Norm =", Norm(1:3,J)
368 ! print *, "ICELLv =", ICELLv
369 ! print *, "IBv =", IBV
370 ! print *, "FCELL =", FCELL(1:4,IEvm)
371 !endif
372 !endif
373
374 brick_list(nin,ib)%POLY(icell)%FACE(j)%Vel(1) = vf(1)
375 brick_list(nin,ib)%POLY(icell)%FACE(j)%Vel(2) = vf(2)
376 brick_list(nin,ib)%POLY(icell)%FACE(j)%Vel(3) = vf(3)
377
378
379 !======================================================!
380 ! REDUCTION FACTOR FOR POLYHEDRA FACES !
381 !======================================================!
382 !facette non communicante
383 numnod = brick_list(nin,ib)%POLY(icell)%NumNOD
384 k = 0
385 DO l=1,numnod
386 inod = brick_list(nin,ib)%POLY(icell)%ListNodID(l)
387 IF(int22_buf%IsNodeOnFace(inod,j))k = k +1
388 ENDDO
389 face = brick_list(nin,ib)%POLY(icell)%FACE(j)%Surf
390 IF(ibv==0)THEN
391 !calcul de facette
392 facev = face
393 numnodv = 8
394 kv = 1
395 ELSE
396 facev = brick_list(nin,ibv)%POLY(icellv)%FACE(jv)%Surf
397 numnod = brick_list(nin,ibv)%POLY(icellv)%NumNOD
398 kv = 0
399 DO l=1,numnod
400 inod = brick_list(nin,ibv)%POLY(icellv)%ListNodID(l)
401 IF(int22_buf%IsNodeOnFace(inod,jv))kv= kv +1
402 ENDDO
403 ENDIF
404 face = min(face,facev)
405 If(k==0 .OR. kv==0) face = zero !facette non communicante si pas de noeud de brique.
406 IF(ibv/=0 .AND. ibm==ibvm) face=zero !blocage flux entre le main et son escalve
407 cellflux(j,icell,ib,iadj) = face * cellflux(j,icell,ib,iadj)
408! write (*,FMT='(5I10,E30.16)') ixs(11,IE), J,ICELL,IB,IADj, Face
409 brick_list(nin,ib)%POLY(icell)%FACE(j)%Adjacent_FLUX(iadj) = cellflux(j,icell,ib,iadj)
410 enddo!next IADj
411 enddo!next J
412 enddo!next ICELL
413 enddo!next IB
414
415 !======================================================!
416 ! debug output !
417 !======================================================!
418 if(debug_outp)then
419 call my_barrier
420 if(itask==0)then
421 if(ibug22_flux22>0 .OR. ibug22_flux22==-1)then
422 print *, " |------e22flux3_int22_fvm.F------|"
423 do ib=1,nb
424 ie = brick_list(nin,ib)%Id
425 if (ibug22_flux22>0 .and. ibug22_flux22/= ixs(11,ie))cycle
426 icell = 0
427 ncell = brick_list(nin,ib)%nbcut
428 mcell = brick_list(nin,ib)%mainid
429 nbcut = ncell
430 DO WHILE (icell<=ncell) ! loop on polyhedron {1:NCELL} U {9}
431 icell = icell +1
432 IF (icell>ncell .AND. ncell/=0)icell=9
433 !print *, " brique =", ixs(11,ie)
434 !print *, " id =", brick_list(nin,ib)%id
435 !print *, " NCYCLE =", NCYCLE
436 !print *, " icell =", ICELL
437 !print *, " _mcell =", mcell
438 !print *, " nbcut =", nbcut
439 do i=1,6
440 nadj = brick_list(nin,ib)%POLY(icell)%FACE(i)%NAdjCell
441 !do k=1,NADJ
442 !IF(ABS(cellFLUX(i,ICELL,IB,K))<=EM20) cellFLUX(i,ICELL,IB,K) = ZERO
443 !enddo
444
445! write (*,FMT='(A,I10,A,5E26.14)') " phi(1:5,",i,")=", brick_list(nin,ib)%POLY(icell)%FACE(1)%Adjacent_FLUX(1:Nadj)
446! write (*,FMT='(A,I10,A,5E26.14)') " phi(1:5,",i,")=", brick_list(nin,ib)%POLY(icell)%FACE(2)%Adjacent_FLUX(1:Nadj)
447! write (*,FMT='(A,I10,A,5E26.14)') " phi(1:5,",i,")=", brick_list(nin,ib)%POLY(icell)%FACE(3)%Adjacent_FLUX(1:Nadj)
448! write (*,FMT='(A,I10,A,5E26.14)') " phi(1:5,",i,")=", brick_list(nin,ib)%POLY(icell)%FACE(4)%Adjacent_FLUX(1:Nadj)
449! write (*,FMT='(A,I10,A,5E26.14)') " phi(1:5,",i,")=", brick_list(nin,ib)%POLY(icell)%FACE(5)%Adjacent_FLUX(1:Nadj)
450! write (*,FMT='(A,I10,A,5E26.14)') " phi(1:5,",i,")=", brick_list(nin,ib)%POLY(icell)%FACE(6)%Adjacent_FLUX(1:Nadj)
451
452! write (*,FMT='(A,6E26.14)') " vel-x =", BRICK_LIST(NIN,IB)%POLY(ICELL)%FACE(1:6)%Vel(1)
453! write (*,FMT='(A,6E26.14)') " vel-y =", BRICK_LIST(NIN,IB)%POLY(ICELL)%FACE(1:6)%Vel(2)
454! write (*,FMT='(A,6E26.14)') " vel-z =", BRICK_LIST(NIN,IB)%POLY(ICELL)%FACE(1:6)%Vel(3)
455! write (*,FMT='(A,6E26.14)') " N-x =", BRICK_LIST(NIN,IB)%N(1:6,1)
456! write (*,FMT='(A,6E26.14)') " N-y =", BRICK_LIST(NIN,IB)%N(1:6,2)
457! write (*,FMT='(A,6E26.14)') " N-z =", BRICK_LIST(NIN,IB)%N(1:6,3)
458 enddo
459 enddo! next ICELL
460 enddo!next IB
461 print *, " ------------------------"
462 endif
463 endif
464 endif
465
466 !======================================================!
467
468
469 !==============!
470 CALL my_barrier
471 !==============!
472
473 !======================================================!
474 ! STEP B : NON CONFORM MESH !
475 ! USE CONSISTENT FLUX !
476 !======================================================!
477 DO ib=nbf,nbl
478 ie = brick_list(nin,ib)%ID
479 mlw = brick_list(nin,ib)%MLW
480 ncell = brick_list(nin,ib)%NBCUT
481 mcell = brick_list(nin,ib)%MainID
482 icell = 0
483 DO WHILE (icell<=ncell) ! loop on polyhedron {1:NCELL} U {9}
484 icell = icell +1
485 IF (icell>ncell .AND. ncell/=0)icell=9
486 DO j=1,6
487 brick_list(nin,ib)%POLY(icell)%FACE(j)%Adjacent_upwFLUX (1) = zero
488 brick_list(nin,ib)%POLY(icell)%FACE(j)%Adjacent_upwFLUX (2) = zero
489 brick_list(nin,ib)%POLY(icell)%FACE(j)%Adjacent_upwFLUX (3) = zero
490 brick_list(nin,ib)%POLY(icell)%FACE(j)%Adjacent_upwFLUX (4) = zero
491 brick_list(nin,ib)%POLY(icell)%FACE(j)%Adjacent_upwFLUX (5) = zero
492 enddo!next J
493
494 brick_list(nin,ib)%POLY(icell)%Adjacent_FLU1 = zero
495 !======================================================!
496 ! monomaterial upwind treatment !
497 !======================================================!
498 ie_m = brick_list(nin,ib)%POLY(icell)%WhereIsMain(3)
499 mat = ixs(1,ie_m)
500 upwl(1:6) = pm(16,mat)
501 reduc = pm(92,mat)
502 ddvol = zero
503 DO j=1,6
504 nadj = brick_list(nin,ib)%POLY(icell)%FACE(j)%NAdjCell
505 iclose = brick_list(nin,ib)%ClosedSurf(j) !pour l'hypothese de surface fermee : double coupe sur la facette et pas de partitionnement de l'autre cote => ON FERME.
506 IF(iclose==1) nadj=0
507 DO iadj = 1,nadj
508 idv = brick_list(nin,ib)%Adjacent_Brick(j,1)
509 ibv = brick_list(nin,ib)%Adjacent_Brick(j,4)
510 jv = brick_list(nin,ib)%Adjacent_Brick(j,5)
511 icellv = brick_list(nin,ib)%POLY(icell)%FACE(j)%Adjacent_Cell(iadj)
512 IF(idv==0)THEN
513 cellflux(j,icell,ib,iadj)=cellflux(j,icell,ib,iadj)*reduc
514 ELSEIF(idv>0)THEN
515 ng = brick_list(nin,ib)%NG
516 isilent = iparg(64,ng)
517 IF(isilent==1)THEN
518 upwl(j)=one
519 cellflux(j,icell,ib,iadj)=cellflux(j,icell,ib,iadj)*pm(92,ixs(1,idv))
520 ENDIF
521 ENDIF
522 brick_list(nin,ib)%POLY(icell)%FACE(j)%Adjacent_upwFLUX(iadj) =
523 . cellflux(j,icell,ib,iadj)-upwl(j)*abs(cellflux(j,icell,ib,iadj))
524 brick_list(nin,ib)%POLY(icell)%Adjacent_FLU1 =
525 . brick_list(nin,ib)%POLY(icell)%Adjacent_FLU1 + cellflux(j,icell,ib,iadj)+upwl(j)*abs(cellflux(j,icell,ib,iadj))
526 !-----DDVOL
527 ddvol = ddvol + cellflux(j,icell,ib,iadj)
528 !-----DDVOL*DT IS SUM FOR INCOMING AND OUTCOMING VOLUMES
529 enddo!next IADJ
530 enddo!next J
531
532 brick_list(nin,ib)%POLY(icell)%DDVOL = ddvol
533
534 !======================================================!
535 enddo!next ICELL
536 enddo!next IB
537
538 !==============!
539 CALL my_barrier
540 !==============!
541
542 !---------------------------------------------------------!
543 ! SECND CELLS STACK !
544 !---------------------------------------------------------!
545 !STACK Secnd cells values from ones connected to current main cell
546 IF(int22>0)THEN
547 nin = 1
548 DO ib=nbf,nbl
549 mlw = brick_list(nin,ib)%MLW
550 num = brick_list(nin,ib)%SecndList%Num
551 mcell = brick_list(nin,ib)%mainID
552 ddvol = zero
553 DO k=1,num
554 ibv = brick_list(nin,ib)%SecndList%IBV(k)
555 icellv = brick_list(nin,ib)%SecndList%ICELLv(k)
556 ddvol = ddvol + brick_list(nin,ibv)%POLY(icellv)%DDVOL
557 ENDDO
558 ddvol = ddvol + brick_list(nin,ib)%POLY(mcell)%DDVOL
559 !updating ddvol cumulative stack
560 brick_list(nin,ib)%POLY(mcell)%DDVOL = ddvol
561 enddo!next IB
562 ENDIF
563
564 !==============!
565 CALL my_barrier
566 !==============!
567
568 RETURN
569 END
norm
norm(diag(diag(diag(inv(mat))) -id.SOL), 2) % destroy mumps instance id.JOB
eflux3_int22_fvm
subroutine eflux3_int22_fvm(pm, ixs, flux, flu1, iparg, elbuf_tab, itask, nv46, ipm, x)
Definition eflux3_int22_fvm.F:39
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
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
i22bufbric_mod
Definition cut-cell-search_mod.F:95
i22bufbric_mod::brick_list
type(brick_entity), dimension(:,:), allocatable, target brick_list
Definition cut-cell-search_mod.F:105
i22bufbric_mod::int22_buf
type(int22_buf_t) int22_buf
Definition cut-cell-search_mod.F:130
i22tri_mod
Definition cut-cell-search_mod.F:212
i22tri_mod::nb
integer nb
Definition cut-cell-search_mod.F:487
i22tri_mod::ibug22_flux22
integer ibug22_flux22
Definition cut-cell-search_mod.F:296
my_barrier
subroutine my_barrier
Definition machine.F:31
upwind
subroutine upwind(rho, vis, vdx, vdy, vdz, r, s, t, deltax, gam, nel)
Definition upwind.F:35