ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
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 (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 ( 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 for computing the buoyancy term to be used in the evolution equation for the z-component of momentum in the slow integrator. There are three options for how buoyancy is computed (two are the same in the absence of moisture).

Parameters
[in]S_datacurrent solution
[in]S_primprimitive variables (i.e. conserved variables divided by density)
[out]buoyancybuoyancy term computed here
[in]qmoistmoisture variables (in order: qv, qc, qi, ...)
[in]qv_dlateral average of cloud vapor
[in]qc_dlateral average of cloud vapor
[in]qd_dlateral average of cloud vapor
[in]geomContainer for geometric information
[in]solverChoiceContainer for solver parameters
[in]r0Reference (hydrostatically stratified) density
[in]n_qstateNumber of moist variables used by the current model
42 {
43  BL_PROFILE("make_buoyancy()");
44 
45  const Array<Real,AMREX_SPACEDIM> grav{0.0, 0.0, -solverChoice.gravity};
46  const GpuArray<Real,AMREX_SPACEDIM> grav_gpu{grav[0], grav[1], grav[2]};
47 
48  const int klo = geom.Domain().smallEnd()[2];
49  const int khi = geom.Domain().bigEnd()[2] + 1;
50 
51  Real rd_over_cp = solverChoice.rdOcp;
52  //Real rv_over_rd = R_v/R_d;
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 -rho0 g (thetaprime / theta0)
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  // Return -rho0 g (thetaprime / theta0)
106  //
107  //buoyancy_fab(i, j, k) = buoyancy_moist_anelastic(i,j,k,grav_gpu[2],rv_over_rd,
108  // r0_arr,th0_arr,qv0_arr,cell_data,qt_arr);
109 
110  // NOTE: Using the type 4, which we formally derived.
111  // The above has errors and needs rederiving.
112  buoyancy_fab(i, j, k) = buoyancy_moist_Thpert(i,j,k,n_qstate,grav_gpu[2],
113  r0_arr,th0_arr,qv0_arr,cell_prim,qt_arr);
114  });
115  }
116  else if ( !anelastic && (solverChoice.moisture_type == MoistureType::None) )
117  {
118  // ******************************************************************************************
119  // Dry compressible
120  // ******************************************************************************************
121  if (solverChoice.buoyancy_type == 1) {
122 
123  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
124  {
125  //
126  // Return -rho0 g (thetaprime / theta0)
127  //
128  buoyancy_fab(i, j, k) = buoyancy_rhopert(i,j,k,grav_gpu[2],
129  r0_arr,cell_data,qt_arr);
130  });
131  }
132  else if (solverChoice.buoyancy_type == 2 || solverChoice.buoyancy_type == 3)
133  {
134  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
135  {
136  //
137  // Return -rho0 g (Tprime / T0)
138  //
139  buoyancy_fab(i, j, k) = buoyancy_dry_Tpert(i,j,k,grav_gpu[2],rd_over_cp,
140  r0_arr,p0_arr,th0_arr,cell_data);
141  });
142  }
143  else if (solverChoice.buoyancy_type == 4)
144  {
145  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
146  {
147  //
148  // Return -rho0 g (Theta_prime / Theta_0)
149  //
150  buoyancy_fab(i, j, k) = buoyancy_dry_Thpert(i,j,k,grav_gpu[2],
151  r0_arr,th0_arr,cell_data);
152  });
153  } // buoyancy_type for dry compressible
154  }
155  else // if ( !anelastic && (solverChoice.moisture_type != MoistureType::None) )
156  {
157  // ******************************************************************************************
158  // Moist compressible
159  // ******************************************************************************************
160 
161  if ( (solverChoice.moisture_type == MoistureType::Kessler_NoRain) ||
162  (solverChoice.moisture_type == MoistureType::SAM) ||
163  (solverChoice.moisture_type == MoistureType::SAM_NoPrecip_NoIce) )
164  {
165  AMREX_ALWAYS_ASSERT(solverChoice.buoyancy_type == 1);
166  }
167 
168  if (solverChoice.buoyancy_type == 1)
169  {
170  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
171  {
172  buoyancy_fab(i, j, k) = buoyancy_rhopert(i,j,k,grav_gpu[2],
173  r0_arr,cell_data,qt_arr);
174  });
175  }
176  else if (solverChoice.buoyancy_type == 2 || solverChoice.buoyancy_type == 3)
177  {
178 
179  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
180  {
181  buoyancy_fab(i, j, k) = buoyancy_moist_Tpert(i,j,k,n_qstate,grav_gpu[2],rd_over_cp,
182  r0_arr,th0_arr,qv0_arr,p0_arr,
183  cell_prim,cell_data,qt_arr);
184  });
185  }
186  else if (solverChoice.buoyancy_type == 4)
187  {
188  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
189  {
190  buoyancy_fab(i, j, k) = buoyancy_moist_Thpert(i,j,k,n_qstate,grav_gpu[2],
191  r0_arr,th0_arr,qv0_arr,cell_prim,qt_arr);
192  });
193  }
194  } // moist compressible
195 
196  } else {
197 
198  // Currently, only dry compressible is supported
199  AMREX_ASSERT( !anelastic && (solverChoice.moisture_type == MoistureType::None) && solverChoice.buoyancy_type == 1 );
200 
201  Array4<const EBCellFlag> cellflg = (ebfact.get_const_factory())->getMultiEBCellFlagFab()[mfi].const_array();
202 
203  ParallelFor(tbz, [=] AMREX_GPU_DEVICE (int i, int j, int k)
204  {
205  buoyancy_fab(i, j, k) = buoyancy_rhopert_eb(i,j,k,grav_gpu[2],r0_arr,cell_data,qt_arr,cellflg);
206  });
207  } // TerrainType::EB
208  } // mfi
209 }
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_moist_Tpert(int &i, int &j, int &k, const int &n_qstate, amrex::Real const &grav_gpu, amrex::Real const &rd_over_cp, 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:165
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_moist_Thpert(int &i, int &j, int &k, const int &n_qstate, amrex::Real const &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:202
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_dry_anelastic(int &i, int &j, int &k, amrex::Real const &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:10
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_dry_Thpert(int &i, int &j, int &k, amrex::Real const &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:147
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_rhopert_eb(int &i, int &j, int &k, amrex::Real const &grav_gpu, const amrex::Array4< const amrex::Real > &r0_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:86
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_dry_Tpert(int &i, int &j, int &k, amrex::Real const &grav_gpu, amrex::Real const &rd_over_cp, 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:121
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real buoyancy_rhopert(int &i, int &j, int &k, amrex::Real const &grav_gpu, const amrex::Array4< const amrex::Real > &r0_arr, const amrex::Array4< const amrex::Real > &cell_data, const amrex::Array4< const amrex::Real > &qt_arr)
Definition: ERF_BuoyancyUtils.H:72
const std::unique_ptr< amrex::EBFArrayBoxFactory > & get_const_factory() const noexcept
Definition: ERF_EB.H:46
@ qv0_comp
Definition: ERF_IndexDefines.H:67
@ p0_comp
Definition: ERF_IndexDefines.H:64
@ th0_comp
Definition: ERF_IndexDefines.H:66
@ r0_comp
Definition: ERF_IndexDefines.H:63
@ cons
Definition: ERF_IndexDefines.H:158
@ qt
Definition: ERF_Kessler.H:27
real(c_double), parameter p0
Definition: ERF_module_model_constants.F90:40
amrex::Real rdOcp
Definition: ERF_DataStruct.H:808
amrex::Real gravity
Definition: ERF_DataStruct.H:806
MoistureType moisture_type
Definition: ERF_DataStruct.H:850
static TerrainType terrain_type
Definition: ERF_DataStruct.H:751
int buoyancy_type
Definition: ERF_DataStruct.H:782
Here is the call graph for this function: