AMReX-Hydro
AMReX-based hydro routines for low Mach number flows
 
Loading...
Searching...
No Matches
hydro_godunov_plm.H
Go to the documentation of this file.
1/**
2 * \file hydro_godunov_plm.H
3 *
4 * \addtogroup Godunov
5 * @{
6 */
7
8#ifndef HYDRO_PLM_GODUNOV_H
9#define HYDRO_PLM_GODUNOV_H
10
11#include <AMReX_Slopes_K.H>
12#include <AMReX_Gpu.H>
13#include <AMReX_FArrayBox.H>
14#include <AMReX_BCRec.H>
15#include <AMReX_BC_TYPES.H>
16#include <AMReX_Array.H>
17#include <iomanip>
18#include <hydro_constants.H>
19
20namespace PLM {
21
22void PredictVelOnXFace ( amrex::Box const& xebox, int ncomp,
23 amrex::Array4<amrex::Real> const& Imx, amrex::Array4<amrex::Real> const& Ipx,
24 amrex::Array4<amrex::Real const> const& q,
25 amrex::Array4<amrex::Real const> const& vcc,
26 const amrex::Geometry& geom,
27 amrex::Real dt,
28 amrex::Vector<amrex::BCRec> const& h_bcrec,
29 amrex::BCRec const* pbc,
30 amrex::Array4<int const> const& bc_arr = {});
31
32void PredictVelOnYFace ( amrex::Box const& yebox, int ncomp,
33 amrex::Array4<amrex::Real> const& Imy, amrex::Array4<amrex::Real> const& Ipy,
34 amrex::Array4<amrex::Real const> const& q,
35 amrex::Array4<amrex::Real const> const& vcc,
36 const amrex::Geometry& geom,
37 amrex::Real dt,
38 amrex::Vector<amrex::BCRec> const& h_bcrec,
39 amrex::BCRec const* pbc,
40 amrex::Array4<int const> const& bc_arr= {});
41
42#if (AMREX_SPACEDIM==3)
43void PredictVelOnZFace ( amrex::Box const& zebox, int ncomp,
44 amrex::Array4<amrex::Real> const& Imz, amrex::Array4<amrex::Real> const& Ipz,
45 amrex::Array4<amrex::Real const> const& q,
46 amrex::Array4<amrex::Real const> const& vcc,
47 const amrex::Geometry& geom,
48 amrex::Real dt,
49 amrex::Vector<amrex::BCRec> const& h_bcrec,
50 amrex::BCRec const* pbc,
51 amrex::Array4<int const> const& bc_arr= {});
52#endif
53
54
55AMREX_GPU_DEVICE AMREX_FORCE_INLINE
56void PredictStateOnXFace ( const int i, const int j, const int k, const int n,
57 const amrex::Real dt, const amrex::Real dx,
58 amrex::Real& Im, amrex::Real& Ip,
59 const amrex::Array4<const amrex::Real> &S,
60 const amrex::Real& umac,
61 const amrex::BCRec bc,
62 const int domain_ilo, const int domain_ihi,
63 const bool is_velocity )
64{
65 using namespace amrex;
66 {
67 bool extdir_or_ho_ilo = (bc.lo(0) == BCType::ext_dir) ||
68 (bc.lo(0) == BCType::hoextrap);
69 bool extdir_or_ho_ihi = (bc.hi(0) == BCType::ext_dir) ||
70 (bc.hi(0) == BCType::hoextrap);
71
72 Real upls, umns;
73
74 int order = 4;
75
76 if (i == domain_ilo && ((bc.lo(0) == BCType::ext_dir) ||
77 (bc.lo(0) == BCType::direction_dependent && umac >= 0.0)) )
78 {
79 umns = S(i-1,j,k,n);
80
81 if ( n==XVEL && is_velocity )
82 {
83 upls = S(i-1,j,k,n);
84 }
85 else
86 {
87 upls = S(i ,j,k,n) + Real(0.5) * (-Real(1.0) - umac * dt/dx) *
88 amrex_calc_xslope_extdir(i ,j,k,n,order,S, extdir_or_ho_ilo, extdir_or_ho_ihi, domain_ilo, domain_ihi);
89 }
90
91 }
92 else if (i == domain_ihi+1 && ((bc.hi(0) == BCType::ext_dir) ||
93 (bc.hi(0) == BCType::direction_dependent && umac <= 0.0)) )
94 {
95 upls = S(i ,j,k,n);
96
97 if ( n==XVEL && is_velocity )
98 {
99 umns = S(i,j,k,n);
100 }
101 else
102 {
103 umns = S(i-1,j,k,n) + Real(0.5) * ( Real(1.0) - umac * dt/dx) *
104 amrex_calc_xslope_extdir(i-1,j,k,n,order,S, extdir_or_ho_ilo, extdir_or_ho_ihi, domain_ilo, domain_ihi);
105 }
106 }
107 else
108 {
109 // Note that we still call the "extdir version" here because interior cells one
110 // away from the boundary will still see the boundary condition in the 4th order
111 // slope
112 upls = S(i ,j,k,n) + Real(0.5) * (-Real(1.0) - umac * dt/dx) *
113 amrex_calc_xslope_extdir(i ,j,k,n,order,S, extdir_or_ho_ilo, extdir_or_ho_ihi, domain_ilo, domain_ihi);
114 umns = S(i-1,j,k,n) + Real(0.5) * ( Real(1.0) - umac * dt/dx) *
115 amrex_calc_xslope_extdir(i-1,j,k,n,order,S, extdir_or_ho_ilo, extdir_or_ho_ihi, domain_ilo, domain_ihi);
116 }
117
118 Ip = umns;
119 Im = upls;
120 }
121}
122
123
124AMREX_GPU_DEVICE AMREX_FORCE_INLINE
125void PredictStateOnYFace ( const int i, const int j, const int k, const int n,
126 const amrex::Real dt, const amrex::Real dy,
127 amrex::Real& Im, amrex::Real& Ip,
128 const amrex::Array4<const amrex::Real> &S,
129 const amrex::Real& vmac,
130 const amrex::BCRec bc,
131 const int domain_jlo, const int domain_jhi,
132 const bool is_velocity )
133{
134 using namespace amrex;
135 {
136 bool extdir_or_ho_jlo = (bc.lo(1) == BCType::ext_dir) ||
137 (bc.lo(1) == BCType::hoextrap);
138 bool extdir_or_ho_jhi = (bc.hi(1) == BCType::ext_dir) ||
139 (bc.hi(1) == BCType::hoextrap);
140
141 Real vpls, vmns;
142
143 int order = 4;
144
145 if (j == domain_jlo && ((bc.lo(1) == BCType::ext_dir) ||
146 (bc.lo(1) == BCType::direction_dependent && vmac >= 0.0)) )
147 {
148 vmns = S(i,j-1,k,n);
149 if ( n==YVEL && is_velocity )
150 {
151 vpls = S(i,j-1,k,n);
152 }
153 else
154 {
155 vpls = S(i,j ,k,n) + Real(0.5) * (-Real(1.0) - vmac * dt/dy) *
156 amrex_calc_yslope_extdir(i,j ,k,n,order,S, extdir_or_ho_jlo, extdir_or_ho_jhi, domain_jlo, domain_jhi);
157 }
158 }
159 else if (j == domain_jhi+1 && ((bc.hi(1) == BCType::ext_dir) ||
160 (bc.hi(1) == BCType::direction_dependent && vmac <= 0.0)) )
161 {
162 vpls = S(i,j ,k,n);
163 if ( n==YVEL && is_velocity )
164 {
165 vmns = S(i,j ,k,n);
166 }
167 else
168 {
169 vmns = S(i,j-1,k,n) + Real(0.5) * ( Real(1.0) - vmac * dt/dy) *
170 amrex_calc_yslope_extdir(i,j-1,k,n,order,S, extdir_or_ho_jlo, extdir_or_ho_jhi, domain_jlo, domain_jhi);
171 }
172 }
173 else
174 {
175 // Note that we still call the "extdir version" here because interior cells one
176 // away from the boundary will still see the boundary condition in the 4th order
177 // slope
178 vpls = S(i,j ,k,n) + Real(0.5) * (-Real(1.0) - vmac * dt/dy) *
179 amrex_calc_yslope_extdir(i,j ,k,n,order,S, extdir_or_ho_jlo, extdir_or_ho_jhi, domain_jlo, domain_jhi);
180 vmns = S(i,j-1,k,n) + Real(0.5) * ( Real(1.0) - vmac * dt/dy) *
181 amrex_calc_yslope_extdir(i,j-1,k,n,order,S, extdir_or_ho_jlo, extdir_or_ho_jhi, domain_jlo, domain_jhi);
182 }
183
184 Ip = vmns;
185 Im = vpls;
186 }
187}
188
189#if (AMREX_SPACEDIM==3)
190
191AMREX_GPU_DEVICE AMREX_FORCE_INLINE
192void PredictStateOnZFace ( const int i, const int j, const int k, const int n,
193 const amrex::Real dt, const amrex::Real dz,
194 amrex::Real& Im, amrex::Real& Ip,
195 const amrex::Array4<const amrex::Real> &S,
196 const amrex::Real& wmac,
197 const amrex::BCRec bc,
198 const int domain_klo, const int domain_khi,
199 const bool is_velocity )
200{
201 using namespace amrex;
202 {
203 bool extdir_or_ho_klo = (bc.lo(2) == BCType::ext_dir) ||
204 (bc.lo(2) == BCType::hoextrap);
205 bool extdir_or_ho_khi = (bc.hi(2) == BCType::ext_dir) ||
206 (bc.hi(2) == BCType::hoextrap);
207
208 Real wpls, wmns;
209
210 int order = 4;
211
212 if (k == domain_klo && ((bc.lo(2) == BCType::ext_dir) ||
213 (bc.lo(2) == BCType::direction_dependent && wmac >= 0.0)) )
214 {
215 wmns = S(i,j,k-1,n);
216 if ( n == ZVEL && is_velocity )
217 {
218 wpls = S(i,j,k-1,n);
219 }
220 else
221 {
222 wpls = S(i,j,k ,n) + Real(0.5) * (-Real(1.0) - wmac * dt/dz) *
223 amrex_calc_zslope_extdir(i,j,k ,n,order,S, extdir_or_ho_klo, extdir_or_ho_khi, domain_klo, domain_khi);
224 }
225 }
226 else if (k == domain_khi+1 && ((bc.hi(2) == BCType::ext_dir) ||
227 (bc.hi(2) == BCType::direction_dependent && wmac <= 0.0)) )
228 {
229 wpls = S(i,j,k ,n);
230 if ( n == ZVEL && is_velocity )
231 {
232 wmns = S(i,j,k ,n);
233 }
234 else
235 {
236 wmns = S(i,j,k-1,n) + Real(0.5) * ( Real(1.0) - wmac * dt/dz) *
237 amrex_calc_zslope_extdir(i,j,k-1,n,order,S, extdir_or_ho_klo, extdir_or_ho_khi, domain_klo, domain_khi);
238 }
239 }
240 else
241 {
242 // Note that we still call the "extdir version" here because interior cells one
243 // away from the boundary will still see the boundary condition in the 4th order
244 // slope
245 wpls = S(i,j,k ,n) + Real(0.5) * (-Real(1.0) - wmac * dt/dz) *
246 amrex_calc_zslope_extdir(i,j,k ,n,order,S, extdir_or_ho_klo, extdir_or_ho_khi, domain_klo, domain_khi);
247 wmns = S(i,j,k-1,n) + Real(0.5) * ( Real(1.0) - wmac * dt/dz) *
248 amrex_calc_zslope_extdir(i,j,k-1,n,order,S, extdir_or_ho_klo, extdir_or_ho_khi, domain_klo, domain_khi);
249 }
250
251 Ip = wmns;
252 Im = wpls;
253 }
254}
255#endif
256
257}
258#endif
259/** @} */
AMReX_Array.H
AMReX_BCRec.H
AMReX_BC_TYPES.H
AMREX_FORCE_INLINE
#define AMREX_FORCE_INLINE
AMReX_FArrayBox.H
AMREX_GPU_DEVICE
#define AMREX_GPU_DEVICE
AMReX_Gpu.H
AMReX_Slopes_K.H
amrex::BCRec::hi
__host__ __device__ const int * hi() const &noexcept
amrex::BCRec::lo
__host__ __device__ const int * lo() const &noexcept
YVEL
#define YVEL
Definition hydro_constants.H:29
XVEL
#define XVEL
Definition hydro_constants.H:28
ZVEL
#define ZVEL
Definition hydro_constants.H:30
hydro_constants.H
PLM
Definition hydro_godunov_plm.H:20
PLM::PredictStateOnYFace
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void PredictStateOnYFace(const int i, const int j, const int k, const int n, const amrex::Real dt, const amrex::Real dy, amrex::Real &Im, amrex::Real &Ip, const amrex::Array4< const amrex::Real > &S, const amrex::Real &vmac, const amrex::BCRec bc, const int domain_jlo, const int domain_jhi, const bool is_velocity)
Definition hydro_godunov_plm.H:125
PLM::PredictStateOnXFace
AMREX_GPU_DEVICE AMREX_FORCE_INLINE void PredictStateOnXFace(const int i, const int j, const int k, const int n, const amrex::Real dt, const amrex::Real dx, amrex::Real &Im, amrex::Real &Ip, const amrex::Array4< const amrex::Real > &S, const amrex::Real &umac, const amrex::BCRec bc, const int domain_ilo, const int domain_ihi, const bool is_velocity)
Definition hydro_godunov_plm.H:56
PLM::PredictVelOnYFace
void PredictVelOnYFace(amrex::Box const &yebox, int ncomp, amrex::Array4< amrex::Real > const &Imy, amrex::Array4< amrex::Real > const &Ipy, amrex::Array4< amrex::Real const > const &q, amrex::Array4< amrex::Real const > const &vcc, const amrex::Geometry &geom, amrex::Real dt, amrex::Vector< amrex::BCRec > const &h_bcrec, amrex::BCRec const *pbc, amrex::Array4< int const > const &bc_arr={})
PLM::PredictVelOnXFace
void PredictVelOnXFace(amrex::Box const &xebox, int ncomp, amrex::Array4< amrex::Real > const &Imx, amrex::Array4< amrex::Real > const &Ipx, amrex::Array4< amrex::Real const > const &q, amrex::Array4< amrex::Real const > const &vcc, const amrex::Geometry &geom, amrex::Real dt, amrex::Vector< amrex::BCRec > const &h_bcrec, amrex::BCRec const *pbc, amrex::Array4< int const > const &bc_arr={})
amrex::amrex_calc_zslope_extdir
__device__ Real amrex_calc_zslope_extdir(int i, int j, int k, int n, int order, amrex::Array4< Real const > const &q, bool edlo, bool edhi, int domlo, int domhi) noexcept