/home/runner/work/amr-wind/amr-wind/amr-wind/equation_systems/vof/split_advection.H Source File

AMR-Wind API: /home/runner/work/amr-wind/amr-wind/amr-wind/equation_systems/vof/split_advection.H Source File
AMR-Wind API v0.1.0
CFD solver for wind plant simulations
Loading...
Searching...
No Matches
split_advection.H
Go to the documentation of this file.
1#ifndef SPLIT_ADVECTION_H_
2#define SPLIT_ADVECTION_H_
3
4#include <cmath>
6#include "AMReX_FArrayBox.H"
7#include "AMReX_BCRec.H"
8#include "AMReX_REAL.H"
9
10using namespace amrex::literals;
11
12namespace amr_wind::multiphase {
13
14AMREX_GPU_DEVICE AMREX_FORCE_INLINE void eulerian_implicit(
15 const int i,
16 const int j,
17 const int k,
18 const int dir,
19 const amrex::Real dtdx,
20 const amrex::Real velL,
21 const amrex::Real velR,
22 amrex::Array4<amrex::Real const> const& volfrac,
23 amrex::Array4<amrex::Real> const& vofL,
24 amrex::Array4<amrex::Real> const& vofR)
25{
26 constexpr amrex::Real tiny =
27 std::numeric_limits<amrex::Real>::epsilon() * 1.0e4_rt;
28 amrex::Real mx = 0.0_rt, my = 0.0_rt, mz = 0.0_rt, alpha = 0.0_rt;
29 amrex::Real x0, deltax;
30 amrex::Real aL = velL * dtdx;
31 amrex::Real aR = velR * dtdx;
32
33 vofL(i, j, k) = 0.0_rt;
34 vofR(i, j, k) = 0.0_rt;
35
36 if (std::abs(volfrac(i, j, k) - 1.0_rt) <= tiny) {
37 if (aL < 0.0_rt) {
38 vofL(i, j, k) = 1.0_rt;
39 }
40 if (aR > 0.0_rt) {
41 vofR(i, j, k) = 1.0_rt;
42 }
43 } else if (volfrac(i, j, k) > tiny) {
44 fit_plane(i, j, k, volfrac, mx, my, mz, alpha);
45 // Eulerian advection
46 x0 = 0.0_rt;
47 if (dir == 0) {
48 if (aL < 0.0_rt) {
49 deltax = -aL;
50 vofL(i, j, k) = cut_volume(mx, my, mz, alpha, x0, deltax);
51 }
52 if (aR > 0.0_rt) {
53 x0 = 1.0_rt - aR;
54 deltax = aR;
55 vofR(i, j, k) = cut_volume(mx, my, mz, alpha, x0, deltax);
56 }
57 } else if (dir == 1) {
58 if (aL < 0.0_rt) {
59 deltax = -aL;
60 vofL(i, j, k) = cut_volume(my, mz, mx, alpha, x0, deltax);
61 }
62 if (aR > 0.0_rt) {
63 x0 = 1.0_rt - aR;
64 deltax = aR;
65 vofR(i, j, k) = cut_volume(my, mz, mx, alpha, x0, deltax);
66 }
67 } else if (dir == 2) {
68 if (aL < 0.0_rt) {
69 deltax = -aL;
70 vofL(i, j, k) = cut_volume(mz, mx, my, alpha, x0, deltax);
71 }
72 if (aR > 0.0_rt) {
73 x0 = 1.0_rt - aR;
74 deltax = aR;
75 vofR(i, j, k) = cut_volume(mz, mx, my, alpha, x0, deltax);
76 }
77 }
78 }
79}
80
81AMREX_GPU_DEVICE AMREX_FORCE_INLINE void balance_eulerian_fluxes(
82 const int i,
83 const int j,
84 const int k,
85 const int dir,
86 const amrex::Real dxi,
87 const amrex::Real velL,
88 const amrex::Real velR,
89 amrex::Array4<amrex::Real> const& volfrac,
90 amrex::Array4<amrex::Real const> const& fluxF,
91 amrex::Array4<amrex::Real const> const& fluxC)
92{
93
94 if (dir == 0) {
95 volfrac(i, j, k) += (fluxF(i, j, k) - fluxF(i + 1, j, k)) * dxi +
96 fluxC(i, j, k) * (velR - velL);
97 } else if (dir == 1) {
98 volfrac(i, j, k) += (fluxF(i, j, k) - fluxF(i, j + 1, k)) * dxi +
99 fluxC(i, j, k) * (velR - velL);
100 } else if (dir == 2) {
101 volfrac(i, j, k) += (fluxF(i, j, k) - fluxF(i, j, k + 1)) * dxi +
102 fluxC(i, j, k) * (velR - velL);
103 }
104 // Do clipping
105 volfrac(i, j, k) = amrex::max<amrex::Real>(
106 0.0_rt, amrex::min<amrex::Real>(1.0_rt, volfrac(i, j, k)));
107}
108
109AMREX_GPU_DEVICE AMREX_FORCE_INLINE void fluxes_bc_save(
110 const int i,
111 const int j,
112 const int k,
113 const int dir,
114 const amrex::Real disp,
115 amrex::Array4<amrex::Real> const& f_f,
116 amrex::Array4<amrex::Real> const& vofL,
117 amrex::Array4<amrex::Real> const& vofR,
118 amrex::Array4<amrex::Real> const& advalpha_f,
119 amrex::GpuArray<BC, AMREX_SPACEDIM * 2> BCs,
120 const int domlo,
121 const int domhi)
122{
123 auto bclo = BCs[amrex::Orientation(dir, amrex::Orientation::low)];
124 auto bchi = BCs[amrex::Orientation(dir, amrex::Orientation::high)];
125
126 if (dir == 0) {
127 // For wall BCs, do not allow flow into or out of domain
128 if (bclo == BC::no_slip_wall || bclo == BC::slip_wall ||
129 bclo == BC::wall_model || bclo == BC::symmetric_wall) {
130 if (i == domlo) {
131 vofR(domlo - 1, j, k) = 0.0_rt;
132 vofL(domlo, j, k) = 0.0_rt;
133 }
134 }
135
136 if (bchi == BC::no_slip_wall || bchi == BC::slip_wall ||
137 bchi == BC::wall_model || bchi == BC::symmetric_wall) {
138 if (i == domhi + 1) {
139 vofL(domhi + 1, j, k) = 0.0_rt;
140 vofR(domhi, j, k) = 0.0_rt;
141 }
142 }
143 // Store volume fraction of each flux and multiply by integral
144 // displacement to get the flux term
145 advalpha_f(i, j, k) = vofR(i - 1, j, k) + vofL(i, j, k);
146 f_f(i, j, k) = advalpha_f(i, j, k) * disp;
147 // Note: only one of vofR or vofL can be nonzero, due to the nature of
148 // eulerian_implicit(), so this sum is not for the sake of accumulation
149 } else if (dir == 1) {
150 if (bclo == BC::no_slip_wall || bclo == BC::slip_wall ||
151 bclo == BC::wall_model || bclo == BC::symmetric_wall) {
152 if (j == domlo) {
153 vofR(i, domlo - 1, k) = 0.0_rt;
154 vofL(i, domlo, k) = 0.0_rt;
155 }
156 }
157 if (bchi == BC::no_slip_wall || bchi == BC::slip_wall ||
158 bchi == BC::wall_model || bchi == BC::symmetric_wall) {
159 if (j == domhi + 1) {
160 vofL(i, domhi + 1, k) = 0.0_rt;
161 vofR(i, domhi, k) = 0.0_rt;
162 }
163 }
164 advalpha_f(i, j, k) = vofR(i, j - 1, k) + vofL(i, j, k);
165 f_f(i, j, k) = advalpha_f(i, j, k) * disp;
166 } else if (dir == 2) {
167 if (bclo == BC::no_slip_wall || bclo == BC::slip_wall ||
168 bclo == BC::wall_model || bclo == BC::symmetric_wall) {
169 if (k == domlo) {
170 vofR(i, j, domlo - 1) = 0.0_rt;
171 vofL(i, j, domlo) = 0.0_rt;
172 }
173 }
174 if (bchi == BC::no_slip_wall || bchi == BC::slip_wall ||
175 bchi == BC::wall_model || bchi == BC::symmetric_wall) {
176 if (k == domhi + 1) {
177 vofL(i, j, domhi + 1) = 0.0_rt;
178 vofR(i, j, domhi) = 0.0_rt;
179 }
180 }
181 advalpha_f(i, j, k) = vofR(i, j, k - 1) + vofL(i, j, k);
182 f_f(i, j, k) = advalpha_f(i, j, k) * disp;
183 }
184}
185
186AMREX_GPU_DEVICE AMREX_FORCE_INLINE void c_mask(
187 const int i,
188 const int j,
189 const int k,
190 amrex::Array4<amrex::Real> const& volfrac,
191 amrex::Array4<amrex::Real> const& volfrac_masked)
192{
193 if (volfrac(i, j, k) > 0.5_rt) {
194 volfrac_masked(i, j, k) = 1.0_rt;
195 } else {
196 volfrac_masked(i, j, k) = 0.0_rt;
197 }
198}
199
200AMREX_GPU_DEVICE AMREX_FORCE_INLINE void remove_vof_debris(
201 const int i,
202 const int j,
203 const int k,
204 amrex::Array4<amrex::Real> const& volfrac)
205{
206 amrex::Real small_vof = 1.0e-6_rt;
207 amrex::Real volFxL = volfrac(i - 1, j, k);
208 amrex::Real volFxR = volfrac(i + 1, j, k);
209 amrex::Real volFyL = volfrac(i, j - 1, k);
210 amrex::Real volFyR = volfrac(i, j + 1, k);
211 amrex::Real volFzL = volfrac(i, j, k - 1);
212 amrex::Real volFzR = volfrac(i, j, k + 1);
213
214 if (volfrac(i, j, k) > 0.0_rt && volFxL < small_vof && volFxR < small_vof &&
215 volFyL < small_vof && volFyR < small_vof && volFzL < small_vof &&
216 volFzR < small_vof) {
217 volfrac(i, j, k) = 0.0_rt;
218 }
219}
220
221} // namespace amr_wind::multiphase
222#endif // SPLIT_ADVECTION.H
@ slip_wall
Definition incflo_enums.H:12
@ wall_model
Definition incflo_enums.H:13
@ symmetric_wall
Definition incflo_enums.H:17
@ no_slip_wall
Definition incflo_enums.H:11
Definition height_functions.H:8
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void fit_plane(const int i, const int j, const int k, amrex::Array4< amrex::Real const > const &volfrac, amrex::Real &mx, amrex::Real &my, amrex::Real &mz, amrex::Real &alpha)
Definition volume_fractions.H:457
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void balance_eulerian_fluxes(const int i, const int j, const int k, const int dir, const amrex::Real dxi, const amrex::Real velL, const amrex::Real velR, amrex::Array4< amrex::Real > const &volfrac, amrex::Array4< amrex::Real const > const &fluxF, amrex::Array4< amrex::Real const > const &fluxC)
Definition split_advection.H:81
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void eulerian_implicit(const int i, const int j, const int k, const int dir, const amrex::Real dtdx, const amrex::Real velL, const amrex::Real velR, amrex::Array4< amrex::Real const > const &volfrac, amrex::Array4< amrex::Real > const &vofL, amrex::Array4< amrex::Real > const &vofR)
Definition split_advection.H:14
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void c_mask(const int i, const int j, const int k, amrex::Array4< amrex::Real > const &volfrac, amrex::Array4< amrex::Real > const &volfrac_masked)
Definition split_advection.H:186
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real cut_volume(const amrex::Real m1, const amrex::Real m2, const amrex::Real m3, const amrex::Real alpha, const amrex::Real r0, const amrex::Real dr0) noexcept
Definition volume_fractions.H:368
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void fluxes_bc_save(const int i, const int j, const int k, const int dir, const amrex::Real disp, amrex::Array4< amrex::Real > const &f_f, amrex::Array4< amrex::Real > const &vofL, amrex::Array4< amrex::Real > const &vofR, amrex::Array4< amrex::Real > const &advalpha_f, amrex::GpuArray< BC, AMREX_SPACEDIM *2 > BCs, const int domlo, const int domhi)
Definition split_advection.H:109
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void remove_vof_debris(const int i, const int j, const int k, amrex::Array4< amrex::Real > const &volfrac)
Definition split_advection.H:200