OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
alemuscl_upwind.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!|| alemuscl_upwind ../engine/source/ale/alemuscl/alemuscl_upwind.F
25!||--- called by ------------------------------------------------------
26!|| afluxt ../engine/source/ale/ale51/afluxt.F
27!||--- uses -----------------------------------------------------
28!|| ale_connectivity_mod ../common_source/modules/ale/ale_connectivity_mod.F
29!|| alemuscl_mod ../common_source/modules/ale/alemuscl_mod.F
30!|| element_mod ../common_source/modules/elements/element_mod.F90
31!|| i22bufbric_mod ../common_source/modules/interfaces/cut-cell-search_mod.F
32!|| i22tri_mod ../common_source/modules/interfaces/cut-cell-search_mod.F
33!|| segvar_mod ../engine/share/modules/segvar_mod.F
34!||====================================================================
35 SUBROUTINE alemuscl_upwind(FLUX, ALE_CONNECT, X, IXS, FLUX_VOIS,
36 . N4_VOIS, ITAB, NV46, ITRIMAT, SEGVAR)
37C-----------------------------------------------
38C D e s c r i p t i o n
39C This subroutines performs the following steps:
40C 1 - compute a gradient for volume fraction ALPH
41C (calls GRADIENT_RECONSTRUCTION)
42C 2 - reconstruct a value for volume fraction on each edge of the mesh
43C based on an affine approximation
44C 3 - upwind this value on the edge and store it in the flux
45C-----------------------------------------------
46 USE i22bufbric_mod
47 USE i22tri_mod
48 USE alemuscl_mod
49 USE segvar_mod
50 USE ale_connectivity_mod
51 use element_mod , only :nixs
52C-----------------------------------------------
53C I m p l i c i t T y p e s
54C-----------------------------------------------
55#include "implicit_f.inc"
56C-----------------------------------------------
57C C o m m o n B l o c k s
58C-----------------------------------------------
59#include "spmd_c.inc"
60#include "vect01_c.inc"
61#include "com04_c.inc"
62C-----------------------------------------------
63C D u m m y A r g u m e n t s
64C-----------------------------------------------
65 INTEGER, INTENT(IN) :: NV46
66 my_real, INTENT(OUT) :: flux(nv46, *)
67 my_real, INTENT(IN) :: x(3,numnod)
68 INTEGER, INTENT(IN) :: IXS(NIXS,NUMELS)
69 my_real, INTENT(OUT) :: flux_vois(numels+nsvois, nv46)
70 INTEGER, INTENT(OUT) :: N4_VOIS(NUMELS+NSVOIS,8)
71 INTEGER, INTENT(IN) :: ITAB(NUMNOD)
72 INTEGER, INTENT(IN) :: ITRIMAT
73 TYPE(t_segvar),INTENT(IN) :: SEGVAR
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 :: I, II, KK, JJ, IAD2, IAD3
79 INTEGER :: NEIGHBOOR_LIST(NV46), FACE_NEIGHBOOR(NV46)
80 my_real :: alphak
81 my_real :: xk, yk, zk
82 my_real :: xf, yf, zf
83 INTEGER :: FACE_TO_NODE_LOCAL_ID(6, 4)
84 my_real :: norm(3), a(3), b(3), c(3), surf, surf1, surf2
85C-----------------------------------------------
86C S o u r c e L i n e s
87C-----------------------------------------------
88!!! Once for all, associate node local id to a face number
89!!! Face 1
90 face_to_node_local_id(1, 1) = 1 ; face_to_node_local_id(1, 2) = 4
91 face_to_node_local_id(1, 3) = 3 ; face_to_node_local_id(1, 4) = 2
92!!! Face 2
93 face_to_node_local_id(2, 1) = 3 ; face_to_node_local_id(2, 2) = 4
94 face_to_node_local_id(2, 3) = 8 ; face_to_node_local_id(2, 4) = 7
95!!! Face 3
96 face_to_node_local_id(3, 1) = 5 ; face_to_node_local_id(3, 2) = 6
97 face_to_node_local_id(3, 3) = 7 ; face_to_node_local_id(3, 4) = 8
98!!! Face 4
99 face_to_node_local_id(4, 1) = 1 ; face_to_node_local_id(4, 2) = 2
100 face_to_node_local_id(4, 3) = 6 ; face_to_node_local_id(4, 4) = 5
101!!! Face 5
102 face_to_node_local_id(5, 1) = 2 ; face_to_node_local_id(5, 2) = 3
103 face_to_node_local_id(5, 3) = 7 ; face_to_node_local_id(5, 4) = 6
104!!! Face 6
105 face_to_node_local_id(6, 1) = 1 ; face_to_node_local_id(6, 2) = 5
106 face_to_node_local_id(6, 3) = 8 ; face_to_node_local_id(6, 4) = 4
107
108!!! First of all, compute gradient for alpha
109 DO i = lft, llt
110 ii = i + nft
111 iad2 = ale_connect%ee_connect%iad_connect(ii)
112 !!!centroid element
113 xk = alemuscl_buffer%ELCENTER(ii,1) ;
114 yk = alemuscl_buffer%ELCENTER(ii,2) ;
115 zk = alemuscl_buffer%ELCENTER(ii,3)
116 !!! Neighbors
117 DO kk = 1, nv46
118 !!! Only for outgoing fluxes
119 IF (flux(kk, ii) > zero) THEN
120 !!! Storing neighbor indexes
121 neighboor_list(kk) = ale_connect%ee_connect%connected(iad2 + kk - 1)
122 face_neighboor(kk) = kk
123 IF (neighboor_list(kk) <= 0) THEN
124 IF(neighboor_list(kk)==0)neighboor_list(kk) = ii
125 !case <0 is for eBCS. -NEIGHBOR_LIST is then the segment number
126 ELSEIF (neighboor_list(kk) <= numels) THEN
127 iad3 = ale_connect%ee_connect%iad_connect(neighboor_list(kk))
128 !!! Store the face number to which II and NEIGHBOR_LIST(KK) are adjacent
129 DO jj = 1, nv46
130 IF (ale_connect%ee_connect%connected(iad3 + jj - 1) == ii) THEN
131 face_neighboor(kk) = jj
132 ENDIF
133 ENDDO ! JJ = 1, NV46
134 ENDIF
135
136 !!! Face centroid
137 xf = zero
138 yf = zero
139 zf = zero
140
141 a(1:3) = x(1:3, ixs(face_to_node_local_id(kk, 1) + 1, ii))
142 b(1:3) = x(1:3, ixs(face_to_node_local_id(kk, 2) + 1, ii))
143 c(1:3) = x(1:3, ixs(face_to_node_local_id(kk, 3) + 1, ii))
144
145 norm(1) = (b(2) - a(2)) * (c(3) - a(3)) - (b(3) - a(3)) * (c(2) - a(2))
146 norm(2) = (b(3) - a(3)) * (c(1) - a(1)) - (b(1) - a(1)) * (c(3) - a(3))
147 norm(3) = (b(1) - a(1)) * (c(2) - a(2)) - (b(2) - a(2)) * (c(1) - a(1))
148
149 surf1 = half * abs(sqrt(norm(1) * norm(1) + norm(2) * norm(2) + norm(3) * norm(3)))
150 xf = surf1 * third * (a(1) + b(1) + c(1))
151 yf = surf1 * third * (a(2) + b(2) + c(2))
152 zf = surf1 * third * (a(3) + b(3) + c(3))
153
154 a(1:3) = x(1:3, ixs(face_to_node_local_id(kk, 1) + 1, ii))
155 b(1:3) = x(1:3, ixs(face_to_node_local_id(kk, 3) + 1, ii))
156 c(1:3) = x(1:3, ixs(face_to_node_local_id(kk, 4) + 1, ii))
157
158 norm(1) = (b(2) - a(2)) * (c(3) - a(3)) - (b(3) - a(3)) * (c(2) - a(2))
159 norm(2) = (b(3) - a(3)) * (c(1) - a(1)) - (b(1) - a(1)) * (c(3) - a(3))
160 norm(3) = (b(1) - a(1)) * (c(2) - a(2)) - (b(2) - a(2)) * (c(1) - a(1))
161
162 surf2 = half * abs(sqrt(norm(1) * norm(1) + norm(2) * norm(2) + norm(3) * norm(3)))
163 xf = xf + surf2 * third * (a(1) + b(1) + c(1))
164 yf = yf + surf2 * third * (a(2) + b(2) + c(2))
165 zf = zf + surf2 * third * (a(3) + b(3) + c(3))
166
167 surf = surf1 + surf2
168 xf = xf / surf
169 yf = yf / surf
170 zf = zf / surf
171
172 !!! Reconstruct second order value for ALPHA(II) on the face
173 alphak = alemuscl_buffer%VOLUME_FRACTION(ii,itrimat)
174 . + alemuscl_buffer%GRAD(ii,1,itrimat) * (xf - xk)
175 . + alemuscl_buffer%GRAD(ii,2,itrimat) * (yf - yk)
176 . + alemuscl_buffer%GRAD(ii,3,itrimat) * (zf - zk)
177 !!! Partial volume flux is then computed as:
178 flux(kk, ii) = alphak * flux(kk, ii)
179 IF (neighboor_list(kk) > 0)THEN
180 IF (neighboor_list(kk) <= numels .AND. neighboor_list(kk) > 0) THEN
181 !!! The opposite of the flux goes to the neighboord
182 flux(face_neighboor(kk), neighboor_list(kk)) = -flux(kk, ii)
183 ELSE
184 !!! ALE51_ANTIDIFF3
185 flux_vois(ii, kk) = flux(kk, ii)
186 n4_vois(ii, 1) = itab(ixs(2, ii))
187 n4_vois(ii, 2) = itab(ixs(3, ii))
188 n4_vois(ii, 3) = itab(ixs(4, ii))
189 n4_vois(ii, 4) = itab(ixs(5, ii))
190 n4_vois(ii, 5) = itab(ixs(6, ii))
191 n4_vois(ii, 6) = itab(ixs(7, ii))
192 n4_vois(ii, 7) = itab(ixs(8, ii))
193 n4_vois(ii, 8) = itab(ixs(9, ii))
194 ENDIF
195 ENDIF
196 ENDIF ! (FLUX(KK, II) > ZERO)
197 ENDDO ! KK = 1, NV46
198 ENDDO ! I = LFT, LLT
199
200C-----------------------------------------------
201C incoming volume fluxes from EBCS
202C-----------------------------------------------
203 IF(nsegflu > 0)THEN
204 DO i = lft, llt
205 ii = i + nft
206 iad2 = ale_connect%ee_connect%iad_connect(ii)
207 DO kk=1,nv46
208 IF(flux(kk,ii) < zero .AND. ale_connect%ee_connect%connected(iad2 + kk - 1) < 0)THEN
209 flux(kk,ii) = segvar%PHASE_ALPHA(itrimat,-ale_connect%ee_connect%connected(iad2 + kk - 1))*flux(kk,ii)
210 ENDIF
211 ENDDO
212 ENDDO
213 ENDIF
214
215C-----------------------------------------------
216 END SUBROUTINE alemuscl_upwind
alemuscl_upwind
subroutine alemuscl_upwind(flux, ale_connect, x, ixs, flux_vois, n4_vois, itab, nv46, itrimat, segvar)
Definition alemuscl_upwind.F:37
my_real
#define my_real
Definition cppsort.cpp:32
norm
norm(diag(diag(diag(inv(mat))) -id.SOL), 2) % destroy mumps instance id.JOB
ale_connectivity_mod
Definition ale_connectivity_mod.F:225
alemuscl_mod
Definition alemuscl_mod.F:55
alemuscl_mod::alemuscl_buffer
type(alemuscl_buffer_) alemuscl_buffer
Definition alemuscl_mod.F:75
i22bufbric_mod
Definition cut-cell-search_mod.F:95
i22tri_mod
Definition cut-cell-search_mod.F:212
segvar_mod
Definition segvar_mod.F:63
ale_connectivity_mod::t_ale_connectivity
Definition ale_connectivity_mod.F:256
segvar_mod::t_segvar
Definition segvar_mod.F:71