ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
ERF_EulerianMicrophysics.H
Go to the documentation of this file.
1 /*! @file ERF_EulerianMicrophysics.H
2  * \brief Contains the Eulerian microphysics class
3  */
4 
5 #ifndef ERF_EULERIANMICROPHYSICS_H
6 #define ERF_EULERIANMICROPHYSICS_H
7 
8 #include "ERF_NullMoist.H"
9 #include "ERF_SAM.H"
10 #include "ERF_Kessler.H"
11 #include "ERF_Morrison.H"
12 #include "ERF_SatAdj.H"
13 #include "ERF_Microphysics.H"
14 
15 /*! \brief Eulerian microphysics interface */
17 
18 public:
19 
20  /*! \brief Null constructor */
22 
23  /*! \brief default destructor */
24  ~EulerianMicrophysics () = default;
25 
26  /*! \brief Constructor: create the moisture model */
27  EulerianMicrophysics (const int& nlev, /*!< Number of AMR levels */
28  const MoistureType& a_model_type /*!< moisture model */)
29  {
30  AMREX_ASSERT( Microphysics::modelType(a_model_type) == MoistureModelType::Eulerian );
31  m_moist_model.resize(nlev);
32  if (a_model_type == MoistureType::SAM ||
33  a_model_type == MoistureType::SAM_NoIce ||
34  a_model_type == MoistureType::SAM_NoPrecip_NoIce) {
35  SetModel<SAM>();
36  } else if (a_model_type == MoistureType::Kessler ||
37  a_model_type == MoistureType::Kessler_NoRain) {
38  SetModel<Kessler>();
39  } else if (a_model_type == MoistureType::Morrison ||
40  a_model_type == MoistureType::Morrison_NoIce) {
41  SetModel<Morrison>();
42  } else if (a_model_type == MoistureType::SatAdj) {
43  SetModel<SatAdj>();
44  } else if (a_model_type == MoistureType::None) {
45  SetModel<NullMoist>();
46  } else {
47  amrex::Abort("EulerianMicrophysics: Dont know this moisture_type!") ;
48  }
49  }
50 
51  /*! \brief Define the moisture model */
52  void Define (const int& lev, /*!< AMR level */
53  SolverChoice& sc /*!< Solver choice object */) override
54  {
55  m_moist_model[lev]->Define(sc);
56  }
57 
58  /*! \brief Initialize the moisture model */
59  void Init (const int& lev, /*!< AMR level */
60  const amrex::MultiFab& cons_in, /*!< Conserved state variables */
61  const amrex::BoxArray& grids, /*!< Grids */
62  const amrex::Geometry& geom, /*!< Geometry */
63  const amrex::Real& dt_advance, /*!< Time step */
64  std::unique_ptr<amrex::MultiFab>& z_phys_nd, /*< Nodal z heights */
65  std::unique_ptr<amrex::MultiFab>& detJ_cc /*< CC Jacobian determinants */) override
66  {
67  m_moist_model[lev]->Init(cons_in, grids, geom, dt_advance,
68  z_phys_nd, detJ_cc);
69  }
70 
71  /*! \brief Advance the moisture model for one time step */
72  void Advance (const int& lev, /*!< AMR level */
73  const amrex::Real& dt_advance, /*!< Time step */
74  const int&, /*!< iteration number */
75  const amrex::Real&, /*!< current time */
76  const SolverChoice &solverChoice, /*!< Solver choice object */
77  amrex::Vector<amrex::Vector<amrex::MultiFab>>&, /*!< Dycore state variables */
78  const amrex::Vector<std::unique_ptr<amrex::MultiFab>>&, /*!< terrain */
79  const amrex::GpuArray<ERF_BC, AMREX_SPACEDIM*2>& ) override
80  {
81  m_moist_model[lev]->Advance(dt_advance, solverChoice);
82  }
83 
84  /*! \brief update microphysics variables from ERF state variables */
85  void Update_Micro_Vars_Lev (const int& lev, /*! AMR level */
86  amrex::MultiFab& cons_in /*!< Conserved state variables */) override
87  {
88  m_moist_model[lev]->Update_Micro_Vars(cons_in);
89  }
90 
91  /*! \brief update ERF state variables from microphysics variables */
92  void Update_State_Vars_Lev (const int& lev, /*!< AMR level */
93  amrex::MultiFab& cons_in /*!< Conserved state variables */) override
94  {
95  m_moist_model[lev]->Update_State_Vars(cons_in);
96  }
97 
98  /*! \brief get pointer to a moisture variable */
99  amrex::MultiFab* Get_Qmoist_Ptr (const int& lev, /*!< AMR level */
100  const int& varIdx /*!< moisture variable index */) override
101  {
102  return m_moist_model[lev]->Qmoist_Ptr(varIdx);
103  }
104 
105  /*! \brief get the number of moisture model variables */
106  int Get_Qmoist_Size (const int& /* lev */) override
107  {
108  return m_moist_model[0]->Qmoist_Size();
109  }
110 
111  /*! \brief get the number of microphysics conserved moist (water-related) state variables */
112  int Get_Qstate_Moist_Size () override
113  {
114  return m_moist_model[0]->Qstate_Moist_Size();
115  }
116 
117  /*! \brief get the number of microphysics conserved non-moist (non-water, i.e., other vapor/condensed species) state variables */
118  int Get_Qstate_NonMoist_Size () override
119  {
120  return m_moist_model[0]->Qstate_NonMoist_Size();
121  }
122 
123  /*! \brief get the indices and names of moisture model variables for restart
124  at a given level */
125  void Get_Qmoist_Restart_Vars (const int a_lev, /*!< level */
126  const SolverChoice& a_sc, /*!< Solver choice object */
127  std::vector<int>& a_idx, /*!< indices */
128  std::vector<std::string>& a_names /*!< names */ ) const override
129  {
130  m_moist_model[a_lev]->Qmoist_Restart_Vars(a_sc, a_idx, a_names);
131  }
132 
133  /*! \brief Returns a list of additional plot variable names */
134  virtual void GetPlotVarNames (amrex::Vector<std::string>& a_vec) const override
135  {
136  m_moist_model[0]->GetPlotVarNames(a_vec);
137  }
138 
139  /*! \brief Fills in a MultiFab for plotting */
140  virtual void GetPlotVar (const std::string& a_name,
141  amrex::MultiFab& a_mf,
142  const int a_lev) const override
143  {
144  m_moist_model[a_lev]->GetPlotVar(a_name, a_mf);
145  }
146 
147 protected:
148 
149  /*! \brief Create and set the specified moisture model */
150  template<class NewMoistModel>
151  void SetModel ()
152  {
153  for (int lev(0); lev<m_moist_model.size(); ++lev) {
154  m_moist_model[lev] = std::make_unique<NewMoistModel>();
155  }
156  }
157 
158 private:
159  amrex::Vector<std::unique_ptr<NullMoist>> m_moist_model; /*!< moisture model */
160 };
161 #endif
Contains the base class for microphysics.
Eulerian microphysics interface.
Definition: ERF_EulerianMicrophysics.H:16
void Get_Qmoist_Restart_Vars(const int a_lev, const SolverChoice &a_sc, std::vector< int > &a_idx, std::vector< std::string > &a_names) const override
get the indices and names of moisture model variables for restart at a given level
Definition: ERF_EulerianMicrophysics.H:125
void Update_Micro_Vars_Lev(const int &lev, amrex::MultiFab &cons_in) override
update microphysics variables from ERF state variables
Definition: ERF_EulerianMicrophysics.H:85
EulerianMicrophysics(const int &nlev, const MoistureType &a_model_type)
Constructor: create the moisture model.
Definition: ERF_EulerianMicrophysics.H:27
~EulerianMicrophysics()=default
default destructor
amrex::Vector< std::unique_ptr< NullMoist > > m_moist_model
Definition: ERF_EulerianMicrophysics.H:159
amrex::MultiFab * Get_Qmoist_Ptr(const int &lev, const int &varIdx) override
get pointer to a moisture variable
Definition: ERF_EulerianMicrophysics.H:99
void Define(const int &lev, SolverChoice &sc) override
Define the moisture model.
Definition: ERF_EulerianMicrophysics.H:52
int Get_Qstate_NonMoist_Size() override
get the number of microphysics conserved non-moist (non-water, i.e., other vapor/condensed species) s...
Definition: ERF_EulerianMicrophysics.H:118
void Update_State_Vars_Lev(const int &lev, amrex::MultiFab &cons_in) override
update ERF state variables from microphysics variables
Definition: ERF_EulerianMicrophysics.H:92
virtual void GetPlotVar(const std::string &a_name, amrex::MultiFab &a_mf, const int a_lev) const override
Fills in a MultiFab for plotting.
Definition: ERF_EulerianMicrophysics.H:140
void SetModel()
Create and set the specified moisture model.
Definition: ERF_EulerianMicrophysics.H:151
virtual void GetPlotVarNames(amrex::Vector< std::string > &a_vec) const override
Returns a list of additional plot variable names.
Definition: ERF_EulerianMicrophysics.H:134
void Advance(const int &lev, const amrex::Real &dt_advance, const int &, const amrex::Real &, const SolverChoice &solverChoice, amrex::Vector< amrex::Vector< amrex::MultiFab >> &, const amrex::Vector< std::unique_ptr< amrex::MultiFab >> &, const amrex::GpuArray< ERF_BC, AMREX_SPACEDIM *2 > &) override
Advance the moisture model for one time step.
Definition: ERF_EulerianMicrophysics.H:72
void Init(const int &lev, const amrex::MultiFab &cons_in, const amrex::BoxArray &grids, const amrex::Geometry &geom, const amrex::Real &dt_advance, std::unique_ptr< amrex::MultiFab > &z_phys_nd, std::unique_ptr< amrex::MultiFab > &detJ_cc) override
Initialize the moisture model.
Definition: ERF_EulerianMicrophysics.H:59
int Get_Qmoist_Size(const int &) override
get the number of moisture model variables
Definition: ERF_EulerianMicrophysics.H:106
int Get_Qstate_Moist_Size() override
get the number of microphysics conserved moist (water-related) state variables
Definition: ERF_EulerianMicrophysics.H:112
EulerianMicrophysics()
Null constructor.
Definition: ERF_EulerianMicrophysics.H:21
Base class for microphysics interface.
Definition: ERF_Microphysics.H:14
static MoistureModelType modelType(const MoistureType a_moisture_type)
query if a specified moisture model is Eulerian or Lagrangian
Definition: ERF_Microphysics.H:86
Definition: ERF_DataStruct.H:91