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