ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
ERF_LandSurface.H
Go to the documentation of this file.
1 #ifndef ERF_LANDSURFACE_H
2 #define ERF_LANDSURFACE_H
3 
4 #include <AMReX_PlotFileUtil.H>
5 #include <AMReX_MultiFabUtil.H>
6 
7 #include <ERF_NullSurf.H>
8 #include <ERF_SLM.H>
9 #include <ERF_MM5.H>
10 
11 #if ERF_USE_NOAH
12 #include <ERF_NOAH.H>
13 #endif
14 
15 class LandSurface {
16 
17 public:
18 
19  LandSurface () { }
20 
21  ~LandSurface () = default;
22 
23  void
24  ReSize (const int& nlev) { m_lsm_model.resize(nlev); }
25 
26  template<class NewSurfModel>
27  void
29  {
30  for (int lev(0); lev<m_lsm_model.size(); ++lev) {
31  m_lsm_model[lev] = std::make_unique<NewSurfModel>();
32  }
33  }
34 
35  void
36  Define (const int& lev,
37  SolverChoice& sc)
38  {
39  m_lsm_model[lev]->Define(sc);
40  }
41 
42  void
43  Init (const int& lev,
44  const amrex::MultiFab& cons_in,
45  const amrex::Geometry& geom,
46  const amrex::Real& dt_advance)
47  {
48  m_lsm_model[lev]->Init(cons_in, geom, dt_advance);
49  }
50 
51  void
52  Advance (const int& lev,
53  amrex::MultiFab& cons_in,
54  amrex::MultiFab& xvel_in,
55  amrex::MultiFab& yvel_in,
56  amrex::MultiFab* hfx3_out,
57  amrex::MultiFab* qfx3_out,
58  const amrex::Real& dt_advance,
59  const int& nstep)
60  {
61  m_lsm_model[lev]->Advance(lev, cons_in, xvel_in, yvel_in, hfx3_out, qfx3_out, dt_advance, nstep);
62  }
63 
64  void
65  Advance (const int& lev, const amrex::Real& dt_advance)
66  {
67  m_lsm_model[lev]->Advance(dt_advance);
68  }
69 
70  void
71  Update_Micro_Vars_Lev (const int& lev, amrex::MultiFab& cons_in)
72  {
73  m_lsm_model[lev]->Update_Micro_Vars(cons_in);
74  }
75 
76  void
77  Update_State_Vars_Lev (const int& lev, amrex::MultiFab& cons_in)
78  {
79  m_lsm_model[lev]->Update_State_Vars(cons_in);
80  }
81 
82  amrex::MultiFab*
83  Get_Data_Ptr (const int& lev, const int& varIdx) { return m_lsm_model[lev]->Lsm_Data_Ptr(varIdx); }
84 
85  amrex::MultiFab*
86  Get_Flux_Ptr (const int& lev, const int& varIdx) { return m_lsm_model[lev]->Lsm_Flux_Ptr(varIdx); }
87 
88  amrex::Geometry
89  Get_Lsm_Geom (const int& lev ) { return m_lsm_model[lev]->Lsm_Geom(); }
90 
91  int
92  Get_Data_Size () { return m_lsm_model[0]->Lsm_Data_Size(); }
93 
94  std::string
95  Get_VarName (const int& varIdx) { return m_lsm_model[0]->Lsm_VarName(varIdx); }
96 
97  void
98  Plot_Lsm_Data (amrex::Real time,
99  const amrex::Vector<int> &level_steps,
100  const amrex::Vector<amrex::IntVect> &ref_ratio)
101  {
102  int nlev = m_lsm_model.size();
103  int nvar = this->Get_Data_Size();
104  amrex::Vector<std::string> varnames;
105 
106  // Only write if we have valid pointers
107  if (this->Get_Data_Ptr(0,0)) {
108  varnames.resize(nvar);
109  m_lsm_geom_lev.resize(nlev);
110  m_lsm_data_lev.resize(nlev);
111  std::string plotfilename = amrex::Concatenate(plot_file_lsm, level_steps[0], 5);
112  for (int lev(0); lev<nlev; ++lev) {
113  m_lsm_geom_lev[lev] = this->Get_Lsm_Geom(lev);
114 
115  amrex::MultiFab* mf_lsm = this->Get_Data_Ptr(lev,0);
116  amrex::IntVect ng(0,0,0);
117  amrex::BoxArray ba = mf_lsm->boxArray();
118  amrex::DistributionMapping dm = mf_lsm->DistributionMap();
119  m_lsm_data_lev[lev].define(ba, dm, nvar, ng);
120  for (int n(0); n<nvar; ++n) {
121  mf_lsm = this->Get_Data_Ptr(lev,n);
122  amrex::MultiFab::Copy(m_lsm_data_lev[lev],*(mf_lsm),0,n,1,0);
123  if (lev==0) varnames[n] = this->Get_VarName(n);
124  }
125  }
126  WriteMultiLevelPlotfile (plotfilename, nlev, GetVecOfConstPtrs(m_lsm_data_lev),
127  varnames, m_lsm_geom_lev, time, level_steps, ref_ratio);
128  m_lsm_geom_lev.clear();
129  m_lsm_data_lev.clear();
130  }
131  }
132 
133 private:
134  // lsm model at each level
135  amrex::Vector<std::unique_ptr<NullSurf>> m_lsm_model;
136 
137  // plotfile prefix
138  std::string plot_file_lsm {"plt_lsm_"};
139 
140  // Vector of geometry
141  amrex::Vector<amrex::Geometry> m_lsm_geom_lev;
142 
143  // Vector of data pointers
144  amrex::Vector<amrex::MultiFab> m_lsm_data_lev;
145 };
146 #endif
Definition: ERF_LandSurface.H:15
amrex::Vector< amrex::Geometry > m_lsm_geom_lev
Definition: ERF_LandSurface.H:141
std::string Get_VarName(const int &varIdx)
Definition: ERF_LandSurface.H:95
void Advance(const int &lev, amrex::MultiFab &cons_in, amrex::MultiFab &xvel_in, amrex::MultiFab &yvel_in, amrex::MultiFab *hfx3_out, amrex::MultiFab *qfx3_out, const amrex::Real &dt_advance, const int &nstep)
Definition: ERF_LandSurface.H:52
void Advance(const int &lev, const amrex::Real &dt_advance)
Definition: ERF_LandSurface.H:65
int Get_Data_Size()
Definition: ERF_LandSurface.H:92
amrex::Geometry Get_Lsm_Geom(const int &lev)
Definition: ERF_LandSurface.H:89
~LandSurface()=default
void ReSize(const int &nlev)
Definition: ERF_LandSurface.H:24
amrex::Vector< amrex::MultiFab > m_lsm_data_lev
Definition: ERF_LandSurface.H:144
void Plot_Lsm_Data(amrex::Real time, const amrex::Vector< int > &level_steps, const amrex::Vector< amrex::IntVect > &ref_ratio)
Definition: ERF_LandSurface.H:98
amrex::Vector< std::unique_ptr< NullSurf > > m_lsm_model
Definition: ERF_LandSurface.H:135
void SetModel()
Definition: ERF_LandSurface.H:28
amrex::MultiFab * Get_Flux_Ptr(const int &lev, const int &varIdx)
Definition: ERF_LandSurface.H:86
void Init(const int &lev, const amrex::MultiFab &cons_in, const amrex::Geometry &geom, const amrex::Real &dt_advance)
Definition: ERF_LandSurface.H:43
LandSurface()
Definition: ERF_LandSurface.H:19
void Update_State_Vars_Lev(const int &lev, amrex::MultiFab &cons_in)
Definition: ERF_LandSurface.H:77
void Define(const int &lev, SolverChoice &sc)
Definition: ERF_LandSurface.H:36
amrex::MultiFab * Get_Data_Ptr(const int &lev, const int &varIdx)
Definition: ERF_LandSurface.H:83
std::string plot_file_lsm
Definition: ERF_LandSurface.H:138
void Update_Micro_Vars_Lev(const int &lev, amrex::MultiFab &cons_in)
Definition: ERF_LandSurface.H:71
Definition: ERF_DataStruct.H:86