ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
SatAdj Class Reference

#include <ERF_SatAdj.H>

Inheritance diagram for SatAdj:
Collaboration diagram for SatAdj:

Public Member Functions

 SatAdj ()
 
virtual ~SatAdj ()=default
 
void AdvanceSatAdj (const SolverChoice &)
 
void Define (SolverChoice &sc) override
 
void Init (const amrex::MultiFab &cons_in, const amrex::BoxArray &, const amrex::Geometry &geom, const amrex::Real &dt_advance, std::unique_ptr< amrex::MultiFab > &, std::unique_ptr< amrex::MultiFab > &) override
 
void Copy_State_to_Micro (const amrex::MultiFab &cons_in) override
 
void Copy_Micro_to_State (amrex::MultiFab &cons_in) override
 
void Update_Micro_Vars (amrex::MultiFab &cons_in) override
 
void Update_State_Vars (amrex::MultiFab &cons_in) override
 
void Advance (const amrex::Real &dt_advance, const SolverChoice &solverChoice) override
 
amrex::MultiFab * Qmoist_Ptr (const int &varIdx) override
 
int Qmoist_Size () override
 
int Qstate_Size () override
 
void Qmoist_Restart_Vars (const SolverChoice &, std::vector< int > &a_idx, std::vector< std::string > &a_names) const override
 
- Public Member Functions inherited from NullMoist
 NullMoist ()
 
virtual ~NullMoist ()=default
 

Static Public Member Functions

AMREX_GPU_HOST_DEVICE static AMREX_FORCE_INLINE amrex::Real NewtonIterSat (int &i, int &j, int &k, const amrex::Real &fac_cond, const amrex::Array4< amrex::Real > &tabs_array, const amrex::Array4< amrex::Real > &pres_array, const amrex::Array4< amrex::Real > &qv_array, const amrex::Array4< amrex::Real > &qc_array)
 

Private Types

using FabPtr = std::shared_ptr< amrex::MultiFab >
 

Private Attributes

int m_qmoist_size = 2
 
int m_qstate_size = 2
 
amrex::Vector< int > MicVarMap
 
amrex::Geometry m_geom
 
amrex::Real dt
 
amrex::Real m_fac_cond
 
amrex::Real m_rdOcp
 
amrex::Array< FabPtr, MicVar_SatAdj::NumVarsmic_fab_vars
 

Member Typedef Documentation

◆ FabPtr

using SatAdj::FabPtr = std::shared_ptr<amrex::MultiFab>
private

Constructor & Destructor Documentation

◆ SatAdj()

SatAdj::SatAdj ( )
inline
47 {}

◆ ~SatAdj()

virtual SatAdj::~SatAdj ( )
virtualdefault

Member Function Documentation

◆ Advance()

void SatAdj::Advance ( const amrex::Real &  dt_advance,
const SolverChoice solverChoice 
)
inlineoverridevirtual

Reimplemented from NullMoist.

98  {
99  dt = dt_advance;
100 
101  this->AdvanceSatAdj(solverChoice);
102  }
void AdvanceSatAdj(const SolverChoice &)
Definition: ERF_SatAdj.cpp:8
amrex::Real dt
Definition: ERF_SatAdj.H:205
Here is the call graph for this function:

◆ AdvanceSatAdj()

void SatAdj::AdvanceSatAdj ( const SolverChoice )

Compute Precipitation-related Microphysics quantities.

9 {
11 
12  // Expose for GPU
13  Real d_fac_cond = m_fac_cond;
14  Real rdOcp = m_rdOcp;
15 
16  // get the temperature, dentisy, theta, qt and qc from input
17  for ( MFIter mfi(*tabs,TilingIfNotGPU()); mfi.isValid(); ++mfi) {
18 
19  const auto& tbx = mfi.tilebox();
20 
21  auto qv_array = mic_fab_vars[MicVar_SatAdj::qv]->array(mfi);
22  auto qc_array = mic_fab_vars[MicVar_SatAdj::qc]->array(mfi);
23  auto tabs_array = mic_fab_vars[MicVar_SatAdj::tabs]->array(mfi);
24  auto theta_array = mic_fab_vars[MicVar_SatAdj::theta]->array(mfi);
25  auto pres_array = mic_fab_vars[MicVar_SatAdj::pres]->array(mfi);
26 
27  ParallelFor(tbx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
28  {
29  qc_array(i,j,k) = std::max(0.0, qc_array(i,j,k));
30 
31  //------- Evaporation/condensation
32  Real qsat;
33  erf_qsatw(tabs_array(i,j,k), pres_array(i,j,k), qsat);
34 
35  // There is enough moisutre to drive to equilibrium
36  if ((qv_array(i,j,k)+qc_array(i,j,k)) > qsat) {
37 
38  // Update temperature
39  tabs_array(i,j,k) = NewtonIterSat(i, j, k ,
40  d_fac_cond, tabs_array, pres_array,
41  qv_array , qc_array );
42 
43  // Update theta (constant pressure)
44  theta_array(i,j,k) = getThgivenPandT(tabs_array(i,j,k), 100.0*pres_array(i,j,k), rdOcp);
45 
46  //
47  // We cannot blindly relax to qsat, but we can convert qc/qi -> qv.
48  // The concept here is that if we put all the moisture into qv and modify
49  // the temperature, we can then check if qv > qsat occurs (for final T/P/qv).
50  // If the reduction in T/qsat and increase in qv does trigger the
51  // aforementioned condition, we can do Newton iteration to drive qv = qsat.
52  //
53  } else {
54  // Changes in each component
55  Real delta_qc = qc_array(i,j,k);
56 
57  // Partition the change in non-precipitating q
58  qv_array(i,j,k) += qc_array(i,j,k);
59  qc_array(i,j,k) = 0.0;
60 
61  // Update temperature (endothermic since we evap/sublime)
62  tabs_array(i,j,k) -= d_fac_cond * delta_qc;
63 
64  // Update theta
65  theta_array(i,j,k) = getThgivenPandT(tabs_array(i,j,k), 100.0*pres_array(i,j,k), rdOcp);
66 
67  // Verify assumption that qv > qsat does not occur
68  erf_qsatw(tabs_array(i,j,k), pres_array(i,j,k), qsat);
69  if (qv_array(i,j,k) > qsat) {
70 
71  // Update temperature
72  tabs_array(i,j,k) = NewtonIterSat(i, j, k ,
73  d_fac_cond , tabs_array, pres_array,
74  qv_array , qc_array );
75 
76  // Update theta
77  theta_array(i,j,k) = getThgivenPandT(tabs_array(i,j,k), 100.0*pres_array(i,j,k), rdOcp);
78 
79  }
80  }
81  });
82  }
83 }
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getThgivenPandT(const amrex::Real T, const amrex::Real P, const amrex::Real rdOcp)
Definition: ERF_EOS.H:18
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void erf_qsatw(amrex::Real t, amrex::Real p, amrex::Real &qsatw)
Definition: ERF_MicrophysicsUtils.H:158
amrex::Real m_fac_cond
Definition: ERF_SatAdj.H:208
amrex::Array< FabPtr, MicVar_SatAdj::NumVars > mic_fab_vars
Definition: ERF_SatAdj.H:212
amrex::Real m_rdOcp
Definition: ERF_SatAdj.H:209
AMREX_GPU_HOST_DEVICE static AMREX_FORCE_INLINE amrex::Real NewtonIterSat(int &i, int &j, int &k, const amrex::Real &fac_cond, const amrex::Array4< amrex::Real > &tabs_array, const amrex::Array4< amrex::Real > &pres_array, const amrex::Array4< amrex::Real > &qv_array, const amrex::Array4< amrex::Real > &qc_array)
Definition: ERF_SatAdj.H:129
@ tabs
Definition: ERF_Kessler.H:32
@ theta
Definition: ERF_SatAdj.H:31
@ pres
Definition: ERF_SatAdj.H:33
@ qv
Definition: ERF_SatAdj.H:35
@ tabs
Definition: ERF_SatAdj.H:32
@ qc
Definition: ERF_SatAdj.H:36

Referenced by Advance().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ Copy_Micro_to_State()

void SatAdj::Copy_Micro_to_State ( amrex::MultiFab &  cons_in)
overridevirtual

Updates conserved and microphysics variables in the provided MultiFabs from the internal MultiFabs that store Microphysics module data.

Parameters
[out]consConserved variables
[out]qmoistqv, qc

Reimplemented from NullMoist.

13 {
14  // Get the temperature, density, theta, qt and qp from input
15  for (MFIter mfi(cons,TilingIfNotGPU()); mfi.isValid(); ++mfi) {
16  const auto& tbx = mfi.tilebox();
17 
18  auto states_arr = cons.array(mfi);
19 
20  auto rho_arr = mic_fab_vars[MicVar_SatAdj::rho]->array(mfi);
21  auto theta_arr = mic_fab_vars[MicVar_SatAdj::theta]->array(mfi);
22  auto qv_arr = mic_fab_vars[MicVar_SatAdj::qv]->array(mfi);
23  auto qc_arr = mic_fab_vars[MicVar_SatAdj::qc]->array(mfi);
24 
25  // get potential total density, temperature, qt, qp
26  ParallelFor(tbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
27  {
28  states_arr(i,j,k,RhoTheta_comp) = rho_arr(i,j,k)*theta_arr(i,j,k);
29  states_arr(i,j,k,RhoQ1_comp) = rho_arr(i,j,k)*qv_arr(i,j,k);
30  states_arr(i,j,k,RhoQ2_comp) = rho_arr(i,j,k)*qc_arr(i,j,k);
31  });
32  }
33 
34  // Fill interior ghost cells and periodic boundaries
35  cons.FillBoundary(m_geom.periodicity());
36 }
#define RhoTheta_comp
Definition: ERF_IndexDefines.H:37
#define RhoQ2_comp
Definition: ERF_IndexDefines.H:43
#define RhoQ1_comp
Definition: ERF_IndexDefines.H:42
amrex::Geometry m_geom
Definition: ERF_SatAdj.H:202
@ rho
Definition: ERF_SatAdj.H:30
@ cons
Definition: ERF_IndexDefines.H:129

Referenced by Update_State_Vars().

Here is the caller graph for this function:

◆ Copy_State_to_Micro()

void SatAdj::Copy_State_to_Micro ( const amrex::MultiFab &  cons_in)
overridevirtual

Initializes the Microphysics module.

Parameters
[in]cons_inConserved variables input

Reimplemented from NullMoist.

44 {
45  // Get the temperature, density, theta, qt and qp from input
46  for (MFIter mfi(cons_in); mfi.isValid(); ++mfi) {
47  const auto& tbx = mfi.tilebox();
48 
49  auto states_array = cons_in.array(mfi);
50 
51  auto qv_array = mic_fab_vars[MicVar_SatAdj::qv]->array(mfi);
52  auto qc_array = mic_fab_vars[MicVar_SatAdj::qc]->array(mfi);
53 
54  auto rho_array = mic_fab_vars[MicVar_SatAdj::rho]->array(mfi);
55  auto theta_array = mic_fab_vars[MicVar_SatAdj::theta]->array(mfi);
56  auto tabs_array = mic_fab_vars[MicVar_SatAdj::tabs]->array(mfi);
57  auto pres_array = mic_fab_vars[MicVar_SatAdj::pres]->array(mfi);
58 
59  // Get pressure, theta, temperature, density, and qt, qp
60  ParallelFor(tbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
61  {
62  rho_array(i,j,k) = states_array(i,j,k,Rho_comp);
63  theta_array(i,j,k) = states_array(i,j,k,RhoTheta_comp)/states_array(i,j,k,Rho_comp);
64  qv_array(i,j,k) = states_array(i,j,k,RhoQ1_comp)/states_array(i,j,k,Rho_comp);
65  qc_array(i,j,k) = states_array(i,j,k,RhoQ2_comp)/states_array(i,j,k,Rho_comp);
66 
67  tabs_array(i,j,k) = getTgivenRandRTh(states_array(i,j,k,Rho_comp),
68  states_array(i,j,k,RhoTheta_comp),
69  qv_array(i,j,k));
70 
71  // Pressure in [mbar] for qsat evaluation
72  pres_array(i,j,k) = getPgivenRTh(states_array(i,j,k,RhoTheta_comp), qv_array(i,j,k)) * 0.01;
73  });
74  }
75 }
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getPgivenRTh(const amrex::Real rhotheta, const amrex::Real qv=0.)
Definition: ERF_EOS.H:84
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getTgivenRandRTh(const amrex::Real rho, const amrex::Real rhotheta, const amrex::Real qv=0.0)
Definition: ERF_EOS.H:46
#define Rho_comp
Definition: ERF_IndexDefines.H:36

Referenced by Update_Micro_Vars().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ Define()

void SatAdj::Define ( SolverChoice sc)
inlineoverridevirtual

Reimplemented from NullMoist.

58  {
59  m_fac_cond = lcond / sc.c_p;
60  m_rdOcp = sc.rdOcp;
61  }
constexpr amrex::Real lcond
Definition: ERF_Constants.H:66
amrex::Real rdOcp
Definition: ERF_DataStruct.H:628
amrex::Real c_p
Definition: ERF_DataStruct.H:627

◆ Init()

void SatAdj::Init ( const amrex::MultiFab &  cons_in,
const amrex::BoxArray &  ,
const amrex::Geometry &  geom,
const amrex::Real &  dt_advance,
std::unique_ptr< amrex::MultiFab > &  ,
std::unique_ptr< amrex::MultiFab > &   
)
overridevirtual

Initializes the Microphysics module.

Parameters
[in]cons_inConserved variables input
[in]qc_inCloud variables input
[in,out]qv_inVapor variables input
[in]qi_inIce variables input
[in]gridsThe boxes on which we will evolve the solution
[in]geomGeometry associated with these MultiFabs and grids
[in]dt_advanceTimestep for the advance

Reimplemented from NullMoist.

23 {
24  dt = dt_advance;
25  m_geom = geom;
26 
27  MicVarMap.resize(m_qmoist_size);
29 
30  // initialize microphysics variables
31  for (auto ivar = 0; ivar < MicVar_SatAdj::NumVars; ++ivar) {
32  mic_fab_vars[ivar] = std::make_shared<MultiFab>(cons_in.boxArray(), cons_in.DistributionMap(),
33  1, cons_in.nGrowVect());
34  mic_fab_vars[ivar]->setVal(0.);
35  }
36 }
amrex::Vector< int > MicVarMap
Definition: ERF_SatAdj.H:199
int m_qmoist_size
Definition: ERF_SatAdj.H:193
@ NumVars
Definition: ERF_SatAdj.H:37

◆ NewtonIterSat()

AMREX_GPU_HOST_DEVICE static AMREX_FORCE_INLINE amrex::Real SatAdj::NewtonIterSat ( int &  i,
int &  j,
int &  k,
const amrex::Real &  fac_cond,
const amrex::Array4< amrex::Real > &  tabs_array,
const amrex::Array4< amrex::Real > &  pres_array,
const amrex::Array4< amrex::Real > &  qv_array,
const amrex::Array4< amrex::Real > &  qc_array 
)
inlinestatic
135  {
136  // Solution tolerance
137  amrex::Real tol = 1.0e-8;
138 
139  // Saturation moisture fractions
140  amrex::Real qsat, dqsat;
141 
142  // Newton iteration vars
143  int niter;
144  amrex::Real fff, dfff;
145  amrex::Real dtabs, delta_qv;
146 
147  // Initial guess for temperature & pressure
148  amrex::Real tabs = tabs_array(i,j,k);
149  amrex::Real pres = pres_array(i,j,k);
150 
151  niter = 0;
152  dtabs = 1;
153 
154  //==================================================
155  // Newton iteration to qv=qsat (cloud phase only)
156  //==================================================
157  do {
158  // Saturation moisture fractions
159  erf_qsatw(tabs, pres, qsat);
160  erf_dtqsatw(tabs, pres, dqsat);
161 
162  // Function for root finding:
163  // 0 = -T_new + T_old + L_eff/C_p * (qv - qsat)
164  fff = -tabs + tabs_array(i,j,k) + fac_cond*(qv_array(i,j,k) - qsat);
165 
166  // Derivative of function (T_new iterated on)
167  dfff = -1.0 - fac_cond*dqsat;
168 
169  // Update the temperature
170  dtabs = -fff/dfff;
171  tabs += dtabs;
172 
173  // Update iteration
174  niter = niter+1;
175  } while(std::abs(dtabs) > tol && niter < 20);
176 
177  // Update qsat from last iteration (dq = dq/dt * dt)
178  qsat += dqsat*dtabs;
179 
180  // Changes in each component
181  delta_qv = qv_array(i,j,k) - qsat;
182 
183  // Partition the change in non-precipitating q
184  qv_array(i,j,k) = qsat;
185  qc_array(i,j,k) += delta_qv;
186 
187  // Return to temperature
188  return tabs;
189  }
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void erf_dtqsatw(amrex::Real t, amrex::Real p, amrex::Real &dtqsatw)
Definition: ERF_MicrophysicsUtils.H:170
Here is the call graph for this function:

◆ Qmoist_Ptr()

amrex::MultiFab* SatAdj::Qmoist_Ptr ( const int &  varIdx)
inlineoverridevirtual

Reimplemented from NullMoist.

106  {
107  AMREX_ALWAYS_ASSERT(varIdx < m_qmoist_size);
108  return mic_fab_vars[MicVarMap[varIdx]].get();
109  }

◆ Qmoist_Restart_Vars()

void SatAdj::Qmoist_Restart_Vars ( const SolverChoice ,
std::vector< int > &  a_idx,
std::vector< std::string > &  a_names 
) const
inlineoverridevirtual

Reimplemented from NullMoist.

121  {
122  a_idx.clear();
123  a_names.clear();
124  }

◆ Qmoist_Size()

int SatAdj::Qmoist_Size ( )
inlineoverridevirtual

Reimplemented from NullMoist.

112 { return SatAdj::m_qmoist_size; }

◆ Qstate_Size()

int SatAdj::Qstate_Size ( )
inlineoverridevirtual

Reimplemented from NullMoist.

115 { return SatAdj::m_qstate_size; }
int m_qstate_size
Definition: ERF_SatAdj.H:196

◆ Update_Micro_Vars()

void SatAdj::Update_Micro_Vars ( amrex::MultiFab &  cons_in)
inlineoverridevirtual

Reimplemented from NullMoist.

83  {
84  this->Copy_State_to_Micro(cons_in);
85  }
void Copy_State_to_Micro(const amrex::MultiFab &cons_in) override
Definition: ERF_InitSatAdj.cpp:43
Here is the call graph for this function:

◆ Update_State_Vars()

void SatAdj::Update_State_Vars ( amrex::MultiFab &  cons_in)
inlineoverridevirtual

Reimplemented from NullMoist.

90  {
91  this->Copy_Micro_to_State(cons_in);
92  }
void Copy_Micro_to_State(amrex::MultiFab &cons_in) override
Definition: ERF_UpdateSatAdj.cpp:12
Here is the call graph for this function:

Member Data Documentation

◆ dt

amrex::Real SatAdj::dt
private

Referenced by Advance().

◆ m_fac_cond

amrex::Real SatAdj::m_fac_cond
private

Referenced by Define().

◆ m_geom

amrex::Geometry SatAdj::m_geom
private

◆ m_qmoist_size

int SatAdj::m_qmoist_size = 2
private

Referenced by Qmoist_Ptr(), and Qmoist_Size().

◆ m_qstate_size

int SatAdj::m_qstate_size = 2
private

Referenced by Qstate_Size().

◆ m_rdOcp

amrex::Real SatAdj::m_rdOcp
private

Referenced by Define().

◆ mic_fab_vars

amrex::Array<FabPtr, MicVar_SatAdj::NumVars> SatAdj::mic_fab_vars
private

Referenced by Qmoist_Ptr().

◆ MicVarMap

amrex::Vector<int> SatAdj::MicVarMap
private

Referenced by Qmoist_Ptr().


The documentation for this class was generated from the following files: