ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_LagrangianMicrophysics.H
Go to the documentation of this file.
1 /*! @file ERF_LagrangianMicrophysics.H
2  * \brief Contains the Lagrangian microphysics class
3  */
4 
5 #ifndef ERF_LAGRANGIANMICROPHYSICS_H
6 #define ERF_LAGRANGIANMICROPHYSICS_H
7 
8 #ifdef ERF_USE_PARTICLES
9 
10 #include <utility>
11 #include <string>
12 
14 #include "ERF_Microphysics.H"
15 
16 /* forward declaration */
17 class ERFPC;
18 
19 /*! \brief Eulerian microphysics interface
20  *
21  * One key difference from #EulerianMicrophysics is that only the base AMR
22  * level has the moisture model. Thus, for higher AMR levels, the microphysics
23  * interface need not do anything. The number of conserved state variables for
24  * moisture (RhoQ1, RhoQ2, ...) are the same for all AMR levels; however, for
25  * levels other than the base, they just contain and evolve zeros. */
26 class LagrangianMicrophysics : public Microphysics {
27 
28 public:
29 
30  /*! \brief Null constructor */
31  LagrangianMicrophysics () { }
32 
33  /*! \brief default destructor */
34  ~LagrangianMicrophysics () = default;
35 
36  /*! \brief Constructor: create the moisture model */
37  LagrangianMicrophysics (const int& /* nlev */, /*!< Number of AMR levels */
38  const MoistureType& a_model_type /*!< moisture model */ )
39  {
40  AMREX_ASSERT( Microphysics::modelType(a_model_type) == MoistureModelType::Lagrangian );
41  amrex::Abort("No Lagrangian moisture model implemented yet!") ;
42  }
43 
44  /*! \brief Define the moisture model */
45  void Define (const int& lev, /*!< AMR level */
46  SolverChoice& sc /*!< Solver choice object */) override
47  {
48  if (lev > 0) return;
49  m_moist_model->Define(sc);
50  }
51 
52  /*! \brief Initialize the moisture model */
53  void Init (const int& lev, /*!< AMR level */
54  const amrex::MultiFab& cons_in, /*!< Conserved state variables */
55  const amrex::BoxArray& grids, /*!< Grids */
56  const amrex::Geometry& geom, /*!< Geometry */
57  const amrex::Real& dt_advance, /*!< Time step */
58  std::unique_ptr<amrex::MultiFab>& z_phys_nd,/*< Nodal z heights */
59  std::unique_ptr<amrex::MultiFab>& detJ_cc /*< CC Jacobian determinants */) override
60  {
61  if (lev > 0) return;
62  m_moist_model->Init(cons_in, grids, geom, dt_advance,
63  z_phys_nd, detJ_cc);
64  }
65 
66  /*! \brief Advance the moisture model for one time step */
67  void Advance (const int& lev, /*!< AMR level */
68  const amrex::Real& dt_advance, /*!< Time step */
69  const int& iter, /*!< iteration number */
70  const amrex::Real& time, /*!< current time */
71  const SolverChoice& /*solverChoice*/, /*!< Solver choice object */
72  amrex::Vector<amrex::Vector<amrex::MultiFab>>& a_vars, /*!< Dycore state variables */
73  const amrex::Vector<std::unique_ptr<amrex::MultiFab>>& a_z/*!< terrain */) override
74  {
75  if (lev > 0) return;
76  m_moist_model->Advance(dt_advance, iter, time, a_vars, a_z);
77  }
78 
79  /*! \brief update microphysics variables from ERF state variables */
80  void Update_Micro_Vars_Lev (const int& lev, /*! AMR level */
81  amrex::MultiFab& cons_in /*!< Conserved state variables */) override
82  {
83  if (lev > 0) return;
84  m_moist_model->Update_Micro_Vars(cons_in);
85  }
86 
87  /*! \brief update ERF state variables from microphysics variables */
88  void Update_State_Vars_Lev (const int& lev, /*!< AMR level */
89  amrex::MultiFab& cons_in /*!< Conserved state variables */) override
90  {
91  if (lev > 0) return;
92  m_moist_model->Update_State_Vars(cons_in);
93  }
94 
95  /*! \brief get pointer to a moisture variable */
96  amrex::MultiFab* Get_Qmoist_Ptr (const int& lev, /*!< AMR level */
97  const int& varIdx /*!< moisture variable index */) override
98  {
99  return (lev > 0 ? nullptr : m_moist_model->Qmoist_Ptr(varIdx));
100  }
101 
102  /*! \brief get the number of moisture model variables */
103  int Get_Qmoist_Size (const int& a_lev /*!< AMR level */) override
104  {
105  return (a_lev > 0 ? 0 : m_moist_model->Qmoist_Size());
106  }
107 
108  /*! \brief get the number of moisture-model-related conserved state variables */
109  int Get_Qstate_Size () override
110  {
111  return m_moist_model->Qstate_Size();
112  }
113 
114  /*! \brief initialize particles in particle container */
115  inline void initParticles ( std::unique_ptr<amrex::MultiFab>& z_phys_nd /*!< Nodal z heights */ )
116  {
117  m_moist_model->InitParticles(z_phys_nd);
118  }
119 
120  /*! \brief restart particles in particle container */
121  inline void restartParticles ( amrex::ParGDBBase* a_gdb, const std::string& a_fname)
122  {
123  m_moist_model->RestartParticles(a_gdb, a_fname);
124  }
125 
126  /*! \brief get the particle container from the moisture model */
127  inline ERFPC* getParticleContainer () const
128  {
129  return m_moist_model->getParticleContainer();
130  }
131 
132  /*! \brief get the name of the moisture model's particle container */
133  inline const std::string& getName () const
134  {
135  return m_moist_model->getName();
136  }
137 
138  /*! \brief get the indices and names of moisture model variables for restart
139  at a given level */
140  void Get_Qmoist_Restart_Vars ( const int a_lev, /*!< level */
141  const SolverChoice& a_sc, /*!< Solver choice object */
142  std::vector<int>& a_idx, /*!< indices */
143  std::vector<std::string>& a_names /*!< names */ ) const override
144  {
145  if (a_lev == 0) {
146  m_moist_model->Qmoist_Restart_Vars( a_sc, a_idx, a_names );
147  } else {
148  a_idx.clear();
149  a_names.clear();
150  }
151  }
152 
153 protected:
154 
155  /*! \brief Create and set the specified moisture model */
156  template<class NewMoistModel>
157  void SetModel ()
158  {
159  m_moist_model = std::make_unique<NewMoistModel>();
160  }
161 
162  std::unique_ptr<NullMoistLagrangian> m_moist_model; /*!< moisture model */
163 };
164 
165 #endif
166 #endif
Contains the base class for microphysics.
Contains the Lagrangian moisture model base class.
Base class for microphysics interface.
Definition: ERF_Microphysics.H:13
virtual int Get_Qmoist_Size(const int &)=0
get the number of moisture model variables
virtual void Define(const int &, SolverChoice &)=0
define the microphysics object
virtual void Init(const int &, const amrex::MultiFab &, const amrex::BoxArray &, const amrex::Geometry &, const amrex::Real &, std::unique_ptr< amrex::MultiFab > &, std::unique_ptr< amrex::MultiFab > &)=0
initialize the microphysics object
virtual void Update_State_Vars_Lev(const int &, amrex::MultiFab &)=0
update ERF state variables from microphysics variables
virtual void Update_Micro_Vars_Lev(const int &, amrex::MultiFab &)=0
update microphysics variables from ERF state variables
virtual int Get_Qstate_Size()=0
get the number of moisture-model-related conserved state variables
virtual void Advance(const int &, const amrex::Real &, const int &, const amrex::Real &, const SolverChoice &, amrex::Vector< amrex::Vector< amrex::MultiFab >> &, const amrex::Vector< std::unique_ptr< amrex::MultiFab >> &)=0
advance microphysics for one time step
virtual amrex::MultiFab * Get_Qmoist_Ptr(const int &, const int &)=0
get pointer to a moisture variable
virtual void Get_Qmoist_Restart_Vars(int, const SolverChoice &, std::vector< int > &, std::vector< std::string > &) const =0
get the indices and names of moisture model variables for restart at a given level
static MoistureModelType modelType(const MoistureType a_moisture_type)
query if a specified moisture model is Eulerian or Lagrangian
Definition: ERF_Microphysics.H:64
Definition: ERF_DataStruct.H:78