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