ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
Functions
ERF_EBAdvectionSrcForScalars.H File Reference

Defines EB scalar-advection flux builders for higher-order stencils. More...

#include <ERF_IndexDefines.H>
#include <ERF_Interpolation.H>
#include <iostream>
#include <fstream>
#include <AMReX_Vector.H>
Include dependency graph for ERF_EBAdvectionSrcForScalars.H:
This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Functions

template<typename InterpType_H , typename InterpType_V >
void EBAdvectionSrcForScalarsWrapper (const amrex::Box &bx, const int &ncomp, const int &icomp, const amrex::GpuArray< const amrex::Array4< amrex::Real >, AMREX_SPACEDIM > flx_arr, const amrex::Array4< const amrex::Real > &cell_prim, const amrex::Array4< const amrex::Real > &avg_xmom, const amrex::Array4< const amrex::Real > &avg_ymom, const amrex::Array4< const amrex::Real > &avg_zmom, const amrex::Array4< const amrex::EBCellFlag > &cellflag, const amrex::Array4< const amrex::Real > &ax_arr, const amrex::Array4< const amrex::Real > &ay_arr, const amrex::Array4< const amrex::Real > &az_arr, const amrex::Real horiz_upw_frac, const amrex::Real vert_upw_frac)
 Build EB scalar fluxes with selected horizontal and vertical interpolation. More...
 
template<typename InterpType_H >
void EBAdvectionSrcForScalarsVert (const amrex::Box &bx, const int &ncomp, const int &icomp, const amrex::GpuArray< const amrex::Array4< amrex::Real >, AMREX_SPACEDIM > flx_arr, const amrex::Array4< const amrex::Real > &cell_prim, const amrex::Array4< const amrex::Real > &avg_xmom, const amrex::Array4< const amrex::Real > &avg_ymom, const amrex::Array4< const amrex::Real > &avg_zmom, const amrex::Array4< const amrex::EBCellFlag > &cellflag, const amrex::Array4< const amrex::Real > &ax_arr, const amrex::Array4< const amrex::Real > &ay_arr, const amrex::Array4< const amrex::Real > &az_arr, const amrex::Real horiz_upw_frac, const amrex::Real vert_upw_frac, const AdvType vert_adv_type)
 Dispatch EB scalar flux construction on the vertical advection scheme. More...
 

Detailed Description

Defines EB scalar-advection flux builders for higher-order stencils.

Function Documentation

◆ EBAdvectionSrcForScalarsVert()

template<typename InterpType_H >
void EBAdvectionSrcForScalarsVert ( const amrex::Box &  bx,
const int &  ncomp,
const int &  icomp,
const amrex::GpuArray< const amrex::Array4< amrex::Real >, AMREX_SPACEDIM >  flx_arr,
const amrex::Array4< const amrex::Real > &  cell_prim,
const amrex::Array4< const amrex::Real > &  avg_xmom,
const amrex::Array4< const amrex::Real > &  avg_ymom,
const amrex::Array4< const amrex::Real > &  avg_zmom,
const amrex::Array4< const amrex::EBCellFlag > &  cellflag,
const amrex::Array4< const amrex::Real > &  ax_arr,
const amrex::Array4< const amrex::Real > &  ay_arr,
const amrex::Array4< const amrex::Real > &  az_arr,
const amrex::Real  horiz_upw_frac,
const amrex::Real  vert_upw_frac,
const AdvType  vert_adv_type 
)

Dispatch EB scalar flux construction on the vertical advection scheme.

Template Parameters
InterpType_HHorizontal interpolation functor selected by the caller.
229 {
230  switch(vert_adv_type) {
231  case AdvType::Centered_2nd:
232  EBAdvectionSrcForScalarsWrapper<InterpType_H,CENTERED2>(bx, ncomp, icomp,
233  flx_arr, cell_prim,
234  avg_xmom, avg_ymom, avg_zmom,
235  cellflag, ax_arr, ay_arr, az_arr,
236  horiz_upw_frac, vert_upw_frac);
237  break;
238  case AdvType::Upwind_3rd:
239  EBAdvectionSrcForScalarsWrapper<InterpType_H,UPWIND3>(bx, ncomp, icomp,
240  flx_arr, cell_prim,
241  avg_xmom, avg_ymom, avg_zmom,
242  cellflag, ax_arr, ay_arr, az_arr,
243  horiz_upw_frac, vert_upw_frac);
244  break;
245  case AdvType::Centered_4th:
246  EBAdvectionSrcForScalarsWrapper<InterpType_H,CENTERED4>(bx, ncomp, icomp,
247  flx_arr, cell_prim,
248  avg_xmom, avg_ymom, avg_zmom,
249  cellflag, ax_arr, ay_arr, az_arr,
250  horiz_upw_frac, vert_upw_frac);
251  break;
252  case AdvType::Upwind_5th:
253  EBAdvectionSrcForScalarsWrapper<InterpType_H,UPWIND5>(bx, ncomp, icomp,
254  flx_arr, cell_prim,
255  avg_xmom, avg_ymom, avg_zmom,
256  cellflag, ax_arr, ay_arr, az_arr,
257  horiz_upw_frac, vert_upw_frac);
258  break;
259  case AdvType::Centered_6th:
260  EBAdvectionSrcForScalarsWrapper<InterpType_H,CENTERED6>(bx, ncomp, icomp,
261  flx_arr, cell_prim,
262  avg_xmom, avg_ymom, avg_zmom,
263  cellflag, ax_arr, ay_arr, az_arr,
264  horiz_upw_frac, vert_upw_frac);
265  break;
266  default:
267  AMREX_ASSERT_WITH_MESSAGE(false, "Unknown vertical advection scheme!");
268  }
269 }
AdvType::Upwind_3rd
@ Upwind_3rd
AdvType::Centered_4th
@ Centered_4th
AdvType::Upwind_5th
@ Upwind_5th
AdvType::Centered_6th
@ Centered_6th
AdvType::Centered_2nd
@ Centered_2nd
AMREX_ASSERT_WITH_MESSAGE
AMREX_ASSERT_WITH_MESSAGE(wbar_cutoff_min > wbar_cutoff_max, "ERROR: wbar_cutoff_min < wbar_cutoff_max")
Here is the call graph for this function:

◆ EBAdvectionSrcForScalarsWrapper()

template<typename InterpType_H , typename InterpType_V >
void EBAdvectionSrcForScalarsWrapper ( const amrex::Box &  bx,
const int &  ncomp,
const int &  icomp,
const amrex::GpuArray< const amrex::Array4< amrex::Real >, AMREX_SPACEDIM >  flx_arr,
const amrex::Array4< const amrex::Real > &  cell_prim,
const amrex::Array4< const amrex::Real > &  avg_xmom,
const amrex::Array4< const amrex::Real > &  avg_ymom,
const amrex::Array4< const amrex::Real > &  avg_zmom,
const amrex::Array4< const amrex::EBCellFlag > &  cellflag,
const amrex::Array4< const amrex::Real > &  ax_arr,
const amrex::Array4< const amrex::Real > &  ay_arr,
const amrex::Array4< const amrex::Real > &  az_arr,
const amrex::Real  horiz_upw_frac,
const amrex::Real  vert_upw_frac 
)

Build EB scalar fluxes with selected horizontal and vertical interpolation.

The wrapper reduces interpolation order locally when covered EB cells truncate the upwind stencil.

Template Parameters
InterpType_HHorizontal interpolation functor.
InterpType_VVertical interpolation functor.
36 {
37  // Instantiate structs for vert/horiz interp
38  InterpType_H interp_prim_h(cell_prim, horiz_upw_frac);
39  InterpType_V interp_prim_v(cell_prim, vert_upw_frac);
40 
41  const int ncells_h = interp_prim_h.GetUpwindCellNumber();
42  const int ncells_v = interp_prim_v.GetUpwindCellNumber();
43 
44  // Instantiate structs for lower-order vert/horiz interp
45  CENTERED2 interp_prim_CEN2(cell_prim, 0);
46  UPWIND3 interp_prim_UPW3(cell_prim, horiz_upw_frac);
47 
48  const amrex::Box xbx = amrex::surroundingNodes(bx,0).grow(amrex::IntVect(0, 1, 1));
49  const amrex::Box ybx = amrex::surroundingNodes(bx,1).grow(amrex::IntVect(1, 0, 1));
50  const amrex::Box zbx = amrex::surroundingNodes(bx,2).grow(amrex::IntVect(1, 1, 0));
51 
52  amrex::ParallelFor(xbx, ncomp,[=] AMREX_GPU_DEVICE (int i, int j, int k, int n) noexcept
53  {
54  const int cons_index = icomp + n;
55 
56  if ( ax_arr(i,j,k) > zero ) // Do we need a threshold value for the area fraction to avoid infinitesimal fractions?
57  {
58  // Find the highest upwind order based on the cell flag
59 
60  int icell = 0;
61 
62  if (avg_xmom(i,j,k)>zero)
63  {
64  for (int ii=0; ii<ncells_h; ++ii) {
65  if (cellflag(i-ii-1,j,k).isCovered()) {
66  break;
67  }
68  icell++;
69  }
70  }
71  else if (avg_xmom(i,j,k)<zero)
72  {
73  for (int ii=0; ii<ncells_h; ++ii) {
74  if (cellflag(i+ii,j,k).isCovered()) {
75  break;
76  }
77  icell++;
78  }
79  }
80 
81  // Interpolate the scalar variable (cell_prim) using the highest-order scheme
82 
83  const int prim_index = cons_index - 1;
84  amrex::Real interpx(zero);
85 
86  if (icell==ncells_h) {
87  interp_prim_h.InterpolateInX(i,j,k,prim_index,interpx,avg_xmom(i,j,k));
88  } else {
89  if (icell==1) {
90  interp_prim_CEN2.InterpolateInX(i,j,k,prim_index,interpx,avg_xmom(i,j,k));
91  } else if (icell==2) {
92  interp_prim_UPW3.InterpolateInX(i,j,k,prim_index,interpx,avg_xmom(i,j,k));
93  }
94  }
95  (flx_arr[0])(i,j,k,cons_index) = avg_xmom(i,j,k) * interpx;
96  }
97  else
98  {
99  (flx_arr[0])(i,j,k,cons_index) = zero;
100  }
101  });
102 
103  amrex::ParallelFor(ybx, ncomp,[=] AMREX_GPU_DEVICE (int i, int j, int k, int n) noexcept
104  {
105  const int cons_index = icomp + n;
106 
107  if ( ay_arr(i,j,k) > zero )
108  {
109  // Find the highest upwind order based on the cell flag
110 
111  int jcell = 0;
112 
113  if (avg_ymom(i,j,k)>zero)
114  {
115  for (int jj=0; jj<ncells_h; ++jj) {
116  if (cellflag(i,j-jj-1,k).isCovered()) {
117  break;
118  }
119  jcell++;
120  }
121  }
122  else if (avg_ymom(i,j,k)<zero)
123  {
124  for (int jj=0; jj<ncells_h; ++jj) {
125  if (cellflag(i,j+jj,k).isCovered()) {
126  break;
127  }
128  jcell++;
129  }
130  }
131 
132  // Interpolate the scalar variable (cell_prim) using the highest-order scheme
133 
134  const int prim_index = cons_index - 1;
135  amrex::Real interpy(zero);
136 
137  if (jcell==ncells_h) {
138  interp_prim_h.InterpolateInY(i,j,k,prim_index,interpy,avg_ymom(i,j,k));
139  } else {
140  if (jcell==1) {
141  interp_prim_CEN2.InterpolateInY(i,j,k,prim_index,interpy,avg_ymom(i,j,k));
142  } else if (jcell==2) {
143  interp_prim_UPW3.InterpolateInY(i,j,k,prim_index,interpy,avg_ymom(i,j,k));
144  }
145  }
146  (flx_arr[1])(i,j,k,cons_index) = avg_ymom(i,j,k) * interpy;
147  }
148  else
149  {
150  (flx_arr[1])(i,j,k,cons_index) = zero;
151  }
152  });
153 
154  interp_prim_UPW3.SetUpwinding(vert_upw_frac);
155 
156  amrex::ParallelFor(zbx, ncomp,[=] AMREX_GPU_DEVICE (int i, int j, int k, int n) noexcept
157  {
158  const int cons_index = icomp + n;
159 
160  if ( az_arr(i,j,k) > zero )
161  {
162  // Find the highest upwind order based on the cell flag
163 
164  int kcell = 0;
165 
166  if (avg_zmom(i,j,k)>zero)
167  {
168  for (int kk=0; kk<ncells_v; ++kk) {
169  if (cellflag(i,j,k-kk-1).isCovered()) {
170  break;
171  }
172  kcell++;
173  }
174  }
175  else if (avg_zmom(i,j,k)<zero)
176  {
177  for (int kk=0; kk<ncells_v; ++kk) {
178  if (cellflag(i,j,k+kk).isCovered()) {
179  break;
180  }
181  kcell++;
182  }
183  }
184 
185  // Interpolate the scalar variable (cell_prim) using the highest-order scheme
186 
187  const int prim_index = cons_index - 1;
188  amrex::Real interpz(zero);
189 
190  if (kcell==ncells_v) {
191  interp_prim_v.InterpolateInZ(i,j,k,prim_index,interpz,avg_zmom(i,j,k));
192  } else {
193  if (kcell==1) {
194  interp_prim_CEN2.InterpolateInZ(i,j,k,prim_index,interpz,avg_zmom(i,j,k));
195  } else if (kcell==2) {
196  interp_prim_UPW3.InterpolateInZ(i,j,k,prim_index,interpz,avg_zmom(i,j,k));
197  }
198  }
199  (flx_arr[2])(i,j,k,cons_index) = avg_zmom(i,j,k) * interpz;
200  }
201  else
202  {
203  (flx_arr[2])(i,j,k,cons_index) = zero;
204  }
205  });
206 }
zero
constexpr amrex::Real zero
Definition: ERF_Constants.H:8
ParallelFor
ParallelFor(fab_box, [=] AMREX_GPU_DEVICE(int i, int j, int k) { qrcuten_arr(i, j, k)=Real(0);qscuten_arr(i, j, k)=Real(0);qicuten_arr(i, j, k)=Real(0);})
Real
amrex::Real Real
Definition: ERF_ShocInterface.H:19
CENTERED2
Definition: ERF_Interpolation_UPW.H:10
UPWIND3
Definition: ERF_Interpolation_UPW.H:95
Here is the call graph for this function: