ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_SLM.H
Go to the documentation of this file.
1 #ifndef ERF_SLM_H
2 #define ERF_SLM_H
3 
4 #include <string>
5 #include <vector>
6 #include <memory>
7 
8 #include <AMReX_FArrayBox.H>
9 #include <AMReX_Geometry.H>
10 #include <AMReX_MultiFabUtil.H>
11 
12 #include <ERF_NullSurf.H>
13 #include <ERF_Constants.H>
14 #include <ERF_IndexDefines.H>
15 #include <ERF_DataStruct.H>
16 
17 namespace LsmVar_SLM {
18  enum {
19  // independent variables
20  theta = 0,
21  NumVars
22  };
23 }
24 
25 /* See Chen & Dudhia (2001) https://doi.org/10.1175/1520-0493(2001)129<0569:CAALSH>2.0.CO;2 */
26 class SLM : public NullSurf {
27 
28  using FabPtr = std::shared_ptr<amrex::MultiFab>;
29 
30 public:
31  // Constructor
32  SLM () {}
33 
34  // Destructor
35  virtual ~SLM () = default;
36 
37  // Set thermo and grid properties
38  void
39  Define (SolverChoice& /*sc*/) override
40  {
41  // NOTE: We should parse things from sc here,
42  // but they are hard coded because this
43  // is a demonstration for now.
44  }
45 
46  // Initialize data structures
47  void
48  Init (const int& /*lev*/,
49  const amrex::MultiFab& cons_in,
50  const amrex::Geometry& geom,
51  const amrex::Real& dt) override;
52 
53  // Wrapper to do all the updating
54  void
55  Advance (const amrex::Real& dt) override
56  {
57  m_dt = dt;
58  this->ComputeFluxes();
59  this->AdvanceSLM();
60  this->ComputeTsurf();
61  }
62 
63  // Compute surface temperature
64  void
65  ComputeTsurf ();
66 
67  // Compute diffusive fluxes
68  void
69  ComputeFluxes ();
70 
71  // Advance the lsm state vars
72  void
73  AdvanceSLM ();
74 
75  // Get state vars from lsm class
76  amrex::MultiFab*
77  Lsm_Data_Ptr (const int& varIdx) override
78  {
79  int lsmIdx = LsmVarMap[varIdx];
80  AMREX_ALWAYS_ASSERT(lsmIdx < SLM::m_lsm_size && lsmIdx>=0);
81  return lsm_fab_vars[lsmIdx].get();
82  }
83 
84  // Get flux vars from lsm class
85  amrex::MultiFab*
86  Lsm_Flux_Ptr (const int& varIdx) override
87  {
88  int lsmIdx = LsmVarMap[varIdx];
89  AMREX_ALWAYS_ASSERT(lsmIdx < SLM::m_lsm_size && lsmIdx>=0);
90  return lsm_fab_flux[lsmIdx].get();
91  }
92 
93  // Get lsm geometry
94  amrex::Geometry
95  Lsm_Geom ( ) override { return m_lsm_geom; }
96 
97  // Get number of vars lsm class contains
98  int
99  Lsm_Data_Size () override { return SLM::m_lsm_size; }
100 
101  // Get variable names
102  std::string
103  Lsm_VarName (const int& varIdx) override
104  {
105  int lsmIdx = LsmVarMap[varIdx];
106  AMREX_ALWAYS_ASSERT(lsmIdx < SLM::m_lsm_size && lsmIdx>=0);
107  return LsmVarName[lsmIdx];
108  }
109 
110 private:
111  // number of lsm variables (theta)
112  int m_lsm_size = 1;
113 
114  // LsmVar map (state indices -> LsmVar enum)
115  amrex::Vector<int> LsmVarMap;
116 
117  // Lsm varnames
118  amrex::Vector<std::string> LsmVarName;
119 
120  // geometry for atmosphere
121  amrex::Geometry m_geom;
122 
123  // geometry for lsm
124  amrex::Geometry m_lsm_geom;
125 
126  // timestep
128 
129  // domain klo-1 or lsm khi
130  int khi_lsm;
131 
132  // constants
133  // amrex::Real m_fac_cond;
134  // amrex::Real m_fac_fus;
135  // amrex::Real m_fac_sub;
136  // amrex::Real m_gOcp;
137 
138  // independent variables
139  amrex::Array<FabPtr, LsmVar_SLM::NumVars> lsm_fab_vars;
140 
141  // flux array for conjugate transfer
142  amrex::Array<FabPtr, LsmVar_SLM::NumVars> lsm_fab_flux;
143 
144  // Vars that should be parsed
145  // ==========================
146  // Number of grid points in z
147  int m_nz_lsm = 30;
148 
149  // Size of grid spacing in z
150  amrex::Real m_dz_lsm = 0.1; // 3 [m] below surface
151 
152  // Specific heat
153  amrex::Real m_cp_soil = 1.26e6; // [J/m^3 K]
154 
155  // Conductivity
156  amrex::Real m_k_soil = 0.2; // dry soil [W/m K]
157 
158  // Theta Dirichlet value at lowest point below surface
159  amrex::Real m_theta_dir = 400.0; // [K]
160 
161  // Thermal diffusivity
163 };
164 #endif
amrex::Real Real
Definition: ERF_ShocInterface.H:16
Definition: ERF_NullSurf.H:8
Definition: ERF_SLM.H:26
void ComputeTsurf()
Definition: ERF_SLM.cpp:64
amrex::MultiFab * Lsm_Flux_Ptr(const int &varIdx) override
Definition: ERF_SLM.H:86
SLM()
Definition: ERF_SLM.H:32
void Init(const int &, const amrex::MultiFab &cons_in, const amrex::Geometry &geom, const amrex::Real &dt) override
Definition: ERF_SLM.cpp:7
int khi_lsm
Definition: ERF_SLM.H:130
virtual ~SLM()=default
void Advance(const amrex::Real &dt) override
Definition: ERF_SLM.H:55
amrex::Array< FabPtr, LsmVar_SLM::NumVars > lsm_fab_flux
Definition: ERF_SLM.H:142
amrex::Real m_cp_soil
Definition: ERF_SLM.H:153
std::string Lsm_VarName(const int &varIdx) override
Definition: ERF_SLM.H:103
amrex::Vector< std::string > LsmVarName
Definition: ERF_SLM.H:118
amrex::Real m_dt
Definition: ERF_SLM.H:127
int m_lsm_size
Definition: ERF_SLM.H:112
amrex::Real m_d_soil
Definition: ERF_SLM.H:162
amrex::Real m_dz_lsm
Definition: ERF_SLM.H:150
amrex::Array< FabPtr, LsmVar_SLM::NumVars > lsm_fab_vars
Definition: ERF_SLM.H:139
amrex::Geometry m_geom
Definition: ERF_SLM.H:121
amrex::Geometry Lsm_Geom() override
Definition: ERF_SLM.H:95
amrex::MultiFab * Lsm_Data_Ptr(const int &varIdx) override
Definition: ERF_SLM.H:77
amrex::Real m_k_soil
Definition: ERF_SLM.H:156
void Define(SolverChoice &) override
Definition: ERF_SLM.H:39
void ComputeFluxes()
Definition: ERF_SLM.cpp:83
std::shared_ptr< amrex::MultiFab > FabPtr
Definition: ERF_SLM.H:28
void AdvanceSLM()
Definition: ERF_SLM.cpp:108
amrex::Geometry m_lsm_geom
Definition: ERF_SLM.H:124
int m_nz_lsm
Definition: ERF_SLM.H:147
int Lsm_Data_Size() override
Definition: ERF_SLM.H:99
amrex::Vector< int > LsmVarMap
Definition: ERF_SLM.H:115
amrex::Real m_theta_dir
Definition: ERF_SLM.H:159
Definition: ERF_SLM.H:17
@ NumVars
Definition: ERF_SLM.H:21
@ theta
Definition: ERF_SLM.H:20
Definition: ERF_DataStruct.H:123