ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_MakeBuoyancy.cpp File Reference
#include <AMReX_MultiFab.H>
#include <AMReX_ArrayLim.H>
#include <AMReX_GpuContainers.H>
#include <ERF_Constants.H>
#include <ERF_EOS.H>
#include <ERF_IndexDefines.H>
#include <ERF_SrcHeaders.H>
#include <ERF_BuoyancyUtils.H>
#include <ERF_EB.H>
Include dependency graph for ERF_MakeBuoyancy.cpp:

Functions

void make_buoyancy (int lev, const Vector< MultiFab > &S_data, const MultiFab &S_prim, const MultiFab &qt, MultiFab &buoyancy, const Geometry geom, const SolverChoice &solverChoice, const MultiFab &base_state, const int n_qstate, const eb_ &ebfact, const int anelastic)
 

Function Documentation

◆ make_buoyancy()

void make_buoyancy ( int  lev,
const Vector< MultiFab > &  S_data,
const MultiFab &  S_prim,
const MultiFab &  qt,
MultiFab &  buoyancy,
const Geometry  geom,
const SolverChoice &  solverChoice,
const MultiFab &  base_state,
const int  n_qstate,
const eb_ &  ebfact,
const int  anelastic 
)

Compute the vertical perturbational gravitational force density [N/m^3] used by the dry-momentum equation. Compressible selector values 1–4 select density, temperature (2 and 3 are identical), or potential- temperature formulations. Anelastic flow dispatches to its active dry or moist kernel independently of the stored selector.

Parameters
[in]levlevel
[in]S_datacurrent solution
[in]S_primprimitive variables (i.e. conserved variables divided by density)
[in]qttotal water in the state
[out]buoyancybuoyancy term computed here
[in]geomContainer for geometric information
[in]solverChoiceContainer for solver parameters
[in]base_statebase state
[in]n_qstateNumber of moist variables used by the current model
[in]ebfactContainer of EB information
[in]anelasticAre we solving the anelastic equations (1 if yes, 0 if no)
45 {
46  BL_PROFILE("make_buoyancy()");
47 
48  const Array<Real,AMREX_SPACEDIM> grav{zero, zero, -solverChoice.gravity};
49  const GpuArray<Real,AMREX_SPACEDIM> grav_gpu{grav[0], grav[1], grav[2]};
50 
51  const int klo = geom.Domain().smallEnd()[2];
52  const int khi = geom.Domain().bigEnd()[2] + 1;
53 
54  MultiFab r0 (base_state, make_alias, BaseState::r0_comp , 1);
55  MultiFab p0 (base_state, make_alias, BaseState::p0_comp , 1);
56  MultiFab th0(base_state, make_alias, BaseState::th0_comp, 1);
57  MultiFab qv0(base_state, make_alias, BaseState::qv0_comp, 1);
58 
59 #ifdef _OPENMP
60 #pragma omp parallel if (amrex::Gpu::notInLaunchRegion())
61 #endif
62  for ( MFIter mfi(buoyancy,TilingIfNotGPU()); mfi.isValid(); ++mfi)
63  {
64  Box tbz = mfi.tilebox();
65 
66  // We don't compute a source term for z-momentum on the bottom or top boundary
67  if (tbz.smallEnd(2) == klo) tbz.growLo(2,-1);
68  if (tbz.bigEnd(2) == khi) tbz.growHi(2,-1);
69 
70  const Array4<const Real> & cell_data = S_data[IntVars::cons].array(mfi);
71  const Array4<const Real> & cell_prim = S_prim.array(mfi);
72  const Array4<const Real> & qt_arr = qt.array(mfi);
73  const Array4< Real> & buoyancy_fab = buoyancy.array(mfi);
74 
75  // Base state density and pressure
76  const Array4<const Real>& r0_arr = r0.const_array(mfi);
77  const Array4<const Real>& p0_arr = p0.const_array(mfi);
78  const Array4<const Real>& th0_arr = th0.const_array(mfi);
79  const Array4<const Real>& qv0_arr = qv0.const_array(mfi);
80 
81  if (solverChoice.terrain_type != TerrainType::EB) {
82 
83  if ( anelastic && (solverChoice.moisture_type == MoistureType::None) )
84  {
85  // ******************************************************************************************
86  // Dry anelastic
87  // ******************************************************************************************
88  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
89  {
90  //
91  // Return the dry anelastic potential-temperature force density.
92  //
93  buoyancy_fab(i, j, k) = buoyancy_dry_anelastic(i,j,k,grav_gpu[2],
94  r0_arr,th0_arr,cell_data);
95  });
96  }
97  else if ( anelastic && (solverChoice.moisture_type != MoistureType::None) )
98  {
99  // ******************************************************************************************
100  // Moist anelastic
101  // ******************************************************************************************
102  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
103  {
104  //
105  // Use the active moist potential-temperature force density.
106  buoyancy_fab(i, j, k) = buoyancy_moist_Thpert(i,j,k,n_qstate,grav_gpu[2],
107  r0_arr,th0_arr,qv0_arr,cell_prim,qt_arr);
108  });
109  }
110  else if ( !anelastic && (solverChoice.moisture_type == MoistureType::None) )
111  {
112  // ******************************************************************************************
113  // Dry compressible
114  // ******************************************************************************************
115  if (solverChoice.buoyancy_type[lev] == 1) {
116 
117  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
118  {
119  //
120  // Exact density-perturbation force density.
121  //
122  buoyancy_fab(i, j, k) = buoyancy_rhopert(i,j,k,grav_gpu[2],
123  r0_arr,qv0_arr,cell_data,qt_arr);
124  });
125  }
126  else if (solverChoice.buoyancy_type[lev] == 2 || solverChoice.buoyancy_type[lev] == 3)
127  {
128  // Legacy selectors 2 and 3 intentionally share this approximation.
129  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
130  {
131  //
132  // Approximate temperature-perturbation force density.
133  //
134  buoyancy_fab(i, j, k) = buoyancy_dry_Tpert(i,j,k,grav_gpu[2],
135  r0_arr,p0_arr,th0_arr,cell_data);
136  });
137  }
138  else if (solverChoice.buoyancy_type[lev] == 4)
139  {
140  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
141  {
142  //
143  // Approximate potential-temperature force density.
144  //
145  buoyancy_fab(i, j, k) = buoyancy_dry_Thpert(i,j,k,grav_gpu[2],
146  r0_arr,th0_arr,cell_prim);
147  });
148  } // buoyancy_type for dry compressible
149  }
150  else // if ( !anelastic && (solverChoice.moisture_type != MoistureType::None) )
151  {
152  // ******************************************************************************************
153  // Moist compressible
154  // ******************************************************************************************
155 
156  if ( (solverChoice.moisture_type == MoistureType::Kessler_NoRain) ||
157  (solverChoice.moisture_type == MoistureType::SAM) ||
158  (solverChoice.moisture_type == MoistureType::SAM_NoPrecip_NoIce) )
159  {
160  AMREX_ALWAYS_ASSERT(solverChoice.buoyancy_type[lev] == 1);
161  }
162 
163  if (solverChoice.buoyancy_type[lev] == 1)
164  {
165  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
166  {
167  buoyancy_fab(i, j, k) = buoyancy_rhopert(i,j,k,grav_gpu[2],
168  r0_arr,qv0_arr,cell_data,qt_arr);
169  });
170  }
171  else if (solverChoice.buoyancy_type[lev] == 2 || solverChoice.buoyancy_type[lev] == 3)
172  {
173  // Legacy selectors 2 and 3 intentionally share this approximation.
174  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
175  {
176  buoyancy_fab(i, j, k) = buoyancy_moist_Tpert(i,j,k,n_qstate,grav_gpu[2],
177  r0_arr,th0_arr,qv0_arr,p0_arr,
178  cell_prim,cell_data,qt_arr);
179  });
180  }
181  else if (solverChoice.buoyancy_type[lev] == 4)
182  {
183  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
184  {
185  buoyancy_fab(i, j, k) = buoyancy_moist_Thpert(i,j,k,n_qstate,grav_gpu[2],
186  r0_arr,th0_arr,qv0_arr,cell_prim,qt_arr);
187  });
188  }
189  } // moist compressible
190 
191  } else {
192 
193  if ( anelastic && (solverChoice.moisture_type == MoistureType::None) ) {
194 
195  if (grav_gpu[2]==0) {
196  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
197  {
198  buoyancy_fab(i, j, k) = zero;
199  });
200  } else {
201 
202  Array4<const EBCellFlag> cellflg = (ebfact.get_const_factory())->getMultiEBCellFlagFab()[mfi].const_array();
203 
204  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
205  {
206  buoyancy_fab(i, j, k) = buoyancy_dry_anelastic_eb(i,j,k,grav_gpu[2],
207  r0_arr,th0_arr,cell_data,cellflg);
208  });
209  }
210  }
211  else
212  {
213  if (grav_gpu[2]==0) {
214  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
215  {
216  buoyancy_fab(i, j, k) = zero;
217  });
218  } else {
219  // Currently, only dry compressible is supported
220  AMREX_ASSERT( !anelastic && (solverChoice.moisture_type == MoistureType::None) && solverChoice.buoyancy_type[lev] == 1 );
221 
222  Array4<const EBCellFlag> cellflg = (ebfact.get_const_factory())->getMultiEBCellFlagFab()[mfi].const_array();
223 
224  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
225  {
226  buoyancy_fab(i, j, k) = buoyancy_rhopert_eb(i,j,k,grav_gpu[2],
227  r0_arr,qv0_arr,cell_data,qt_arr,cellflg);
228  });
229  }
230  }
231  } // TerrainType::EB
232  } // mfi
233 }
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_dry_Thpert(int &i, int &j, int &k, const amrex::Real &grav_gpu, const amrex::Array4< const amrex::Real > &r0_arr, const amrex::Array4< const amrex::Real > &th0_arr, const amrex::Array4< const amrex::Real > &cell_prim)
Definition: ERF_BuoyancyUtils.H:239
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_rhopert(int &i, int &j, int &k, const amrex::Real &grav_gpu, const amrex::Array4< const amrex::Real > &r0_arr, const amrex::Array4< const amrex::Real > &qv0_arr, const amrex::Array4< const amrex::Real > &cell_data, const amrex::Array4< const amrex::Real > &qt_arr)
Definition: ERF_BuoyancyUtils.H:113
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_moist_Tpert(int &i, int &j, int &k, const int &n_qstate, const amrex::Real &grav_gpu, const amrex::Array4< const amrex::Real > &r0_arr, const amrex::Array4< const amrex::Real > &th0_arr, const amrex::Array4< const amrex::Real > &qv0_arr, const amrex::Array4< const amrex::Real > &p0_arr, const amrex::Array4< const amrex::Real > &cell_prim, const amrex::Array4< const amrex::Real > &cell_data, const amrex::Array4< const amrex::Real > &qt_arr)
Definition: ERF_BuoyancyUtils.H:275
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_moist_Thpert(int &i, int &j, int &k, const int &n_qstate, const amrex::Real &grav_gpu, const amrex::Array4< const amrex::Real > &r0_arr, const amrex::Array4< const amrex::Real > &th0_arr, const amrex::Array4< const amrex::Real > &qv0_arr, const amrex::Array4< const amrex::Real > &cell_prim, const amrex::Array4< const amrex::Real > &qt_arr)
Definition: ERF_BuoyancyUtils.H:326
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_dry_anelastic_eb(int &i, int &j, int &k, const amrex::Real &grav_gpu, const amrex::Array4< const amrex::Real > &r0_arr, const amrex::Array4< const amrex::Real > &th0_arr, const amrex::Array4< const amrex::Real > &cell_data, amrex::Array4< amrex::EBCellFlag const > const &flag)
Definition: ERF_BuoyancyUtils.H:55
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_rhopert_eb(int &i, int &j, int &k, const amrex::Real &grav_gpu, const amrex::Array4< const amrex::Real > &r0_arr, const amrex::Array4< const amrex::Real > &qv0_arr, const amrex::Array4< const amrex::Real > &cell_data, const amrex::Array4< const amrex::Real > &qt_arr, amrex::Array4< amrex::EBCellFlag const > const &flag)
Definition: ERF_BuoyancyUtils.H:144
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_dry_anelastic(int &i, int &j, int &k, const amrex::Real &grav_gpu, const amrex::Array4< const amrex::Real > &r0_arr, const amrex::Array4< const amrex::Real > &th0_arr, const amrex::Array4< const amrex::Real > &cell_data)
Definition: ERF_BuoyancyUtils.H:22
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_dry_Tpert(int &i, int &j, int &k, const amrex::Real &grav_gpu, const amrex::Array4< const amrex::Real > &r0_arr, const amrex::Array4< const amrex::Real > &p0_arr, const amrex::Array4< const amrex::Real > &th0_arr, const amrex::Array4< const amrex::Real > &cell_data)
Definition: ERF_BuoyancyUtils.H:202
const int klo
Definition: ERF_InitCustomPert_ABL.H:75
const int khi
Definition: ERF_InitCustomPert_Bubble.H:21
AMREX_ALWAYS_ASSERT(bx.length()[2]==khi+1)
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);})
constexpr amrex::Real zero
Definition: ERF_NumericalConstants.H:36
const std::unique_ptr< amrex::EBFArrayBoxFactory > & get_const_factory() const noexcept
Return the cell-centered EB factory.
Definition: ERF_EB.H:102
@ qv0_comp
Definition: ERF_IndexDefines.H:83
@ p0_comp
Definition: ERF_IndexDefines.H:80
@ th0_comp
Definition: ERF_IndexDefines.H:82
@ r0_comp
Definition: ERF_IndexDefines.H:79
@ cons
Definition: ERF_IndexDefines.H:235
@ None
Definition: ERF_LSM_Types.H:6
@ qt
Definition: ERF_Kessler.H:30
real(c_double), parameter p0
Definition: ERF_module_model_constants.F90:40
real(kind=kind_phys), parameter, private r0
Definition: ERF_module_mp_wdm6.F90:75
amrex::Real gravity
Effective gravitational acceleration.
Definition: ERF_DataStruct.H:2343
MoistureType moisture_type
Moisture or microphysics model.
Definition: ERF_DataStruct.H:2524
amrex::Vector< int > buoyancy_type
Per-level buoyancy formulation selector.
Definition: ERF_DataStruct.H:2263
static TerrainType terrain_type
Terrain or immersed-boundary representation.
Definition: ERF_DataStruct.H:2224

Referenced by ERF::compute_max_buoyancy_gradp_diagnostic(), and ERF::FillPlot3DVars().

Here is the call graph for this function:
Here is the caller graph for this function: