PISM, A Parallel Ice Sheet Model  stable v2.0.4 committed by Constantine Khrulev on 2022-05-25 12:02:27 -0800
EnthalpyConverter.hh
Go to the documentation of this file.
1 // Copyright (C) 2009-2011, 2013, 2014, 2015, 2016, 2021 Andreas Aschwanden, Ed Bueler and Constantine Khroulev
2 //
3 // This file is part of PISM.
4 //
5 // PISM is free software; you can redistribute it and/or modify it under the
6 // terms of the GNU General Public License as published by the Free Software
7 // Foundation; either version 3 of the License, or (at your option) any later
8 // version.
9 //
10 // PISM is distributed in the hope that it will be useful, but WITHOUT ANY
11 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
12 // FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
13 // details.
14 //
15 // You should have received a copy of the GNU General Public License
16 // along with PISM; if not, write to the Free Software
17 // Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
18 
19 #ifndef __enthalpyConverter_hh
20 #define __enthalpyConverter_hh
21 
22 #include <vector>
23 #include <memory>
24 
25 namespace pism {
26 
27 class Config;
28 
29 //! Converts between specific enthalpy and temperature or liquid content.
30 /*!
31  Use this way, for example within IceModel with Config config member:
32  \code
33  #include "enthalpyConverter.hh"
34 
35  EnthalpyConverter EC(&config); // runs constructor; do after initialization of Config config
36  ...
37  for (...) {
38  ...
39  E_s = EC.enthalpy_cts(p);
40  ... etc ...
41  }
42  \endcode
43 
44  The three methods that get the enthalpy from temperatures and liquid
45  fractions, namely enthalpy(), enthalpy_permissive() are more strict
46  about error checking. They throw RuntimeError if their arguments are
47  invalid.
48 
49  This class is documented by [\ref AschwandenBuelerKhroulevBlatter].
50 */
52 public:
53  EnthalpyConverter(const Config &config);
54  virtual ~EnthalpyConverter() = default;
55 
56  typedef std::shared_ptr<EnthalpyConverter> Ptr;
57 
58  bool is_temperate(double E, double P) const;
59  bool is_temperate_relaxed(double E, double P) const;
60 
61  double temperature(double E, double P) const;
62  double melting_temperature(double P) const;
63  double pressure_adjusted_temperature(double E, double P) const;
64 
65  double water_fraction(double E, double P) const;
66 
67  double enthalpy(double T, double omega, double P) const;
68  double enthalpy_cts(double P) const;
69  double enthalpy_liquid(double P) const;
70  double enthalpy_permissive(double T, double omega, double P) const;
71 
72  double c() const;
73  double L(double T_pm) const;
74 
75  double pressure(double depth) const;
76  void pressure(const std::vector<double> &depth,
77  unsigned int ks, std::vector<double> &result) const;
78 protected:
79  void validate_E_P(double E, double P) const;
80  void validate_T_omega_P(double T, double omega, double P) const;
81 
82  double temperature_cold(double E) const;
83  double enthalpy_cold(double T) const;
84 
85  //! melting temperature of pure water at atmospheric pressure
86  double m_T_melting;
87  //! latent heat of fusion of water at atmospheric pressure
88  double m_L;
89  //! specific heat capacity of ice
90  double m_c_i;
91  //! specific heat capacity of pure water
92  double m_c_w;
93  //! density of ice
94  double m_rho_i;
95  //! acceleration due to gravity
96  double m_g;
97  //! atmospheric pressure
98  double m_p_air;
99  //! beta in the Clausius-Clapeyron relation (@f$ \diff{T_m}{p} = - \beta @f$).
100  double m_beta;
101  //! temperature tolerance used in `is_temperate()` in cold ice mode
103  //! reference temperature in the definition of ice enthalpy
104  double m_T_0;
105  //! @brief if cold ice methods are selected, use `is_temperate()`
106  //! check based on temperature, not enthalpy
108 };
109 
110 
111 //! An EnthalpyConverter for use in verification tests.
112 
113 /*! Treats ice at any temperature below 10^6 Kelvin as cold (= zero liquid fraction).
114 
115  The pressure dependence of the pressure-melting temperature is neglected.c;
116 
117  Note: Any instance of FlowLaw uses an EnthalpyConverter; this is
118  the one used in cold mode verification code.
119 
120 
121  This is the special enthalpy converter that is used in
122  temperature-based verification tests only.
123 
124  In these tests ice temperatures in an exact solution may exceed the
125  pressure-melting temperature, but we still want to pretend that this
126  ice is "cold" to ensure that the map from enthalpy to temperature is
127  one-to-one. (Normally enthalpy is mapped to the (temperature, water
128  fraction) pair; here water fraction is zero, so enthalpy <-->
129  (temperature, 0.)
130 
131  So, I had to pick a threshold (melting) temperature that is above
132  all ice temperatures. 10^6K was chosen.
133 */
135 public:
136  ColdEnthalpyConverter(const Config &config);
137  virtual ~ColdEnthalpyConverter() = default;
138 };
139 
140 } // end of namespace pism
141 
142 #endif // __enthalpyConverter_hh
143 
virtual ~ColdEnthalpyConverter()=default
ColdEnthalpyConverter(const Config &config)
An EnthalpyConverter for use in verification tests.
A class for storing and accessing PISM configuration flags and parameters.
double enthalpy(double T, double omega, double P) const
Compute enthalpy from absolute temperature, liquid water fraction, and pressure.
EnthalpyConverter(const Config &config)
void validate_E_P(double E, double P) const
double m_rho_i
density of ice
bool is_temperate(double E, double P) const
double melting_temperature(double P) const
Get melting temperature from pressure p.
std::shared_ptr< EnthalpyConverter > Ptr
double m_T_tolerance
temperature tolerance used in is_temperate() in cold ice mode
double temperature(double E, double P) const
Get absolute (not pressure-adjusted) ice temperature (K) from enthalpy and pressure.
double m_beta
beta in the Clausius-Clapeyron relation ( ).
virtual ~EnthalpyConverter()=default
bool m_do_cold_ice_methods
if cold ice methods are selected, use is_temperate() check based on temperature, not enthalpy
double enthalpy_cold(double T) const
Convert temperature into enthalpy (cold case).
double enthalpy_liquid(double P) const
Compute the maximum allowed value of ice enthalpy (corresponds to ).
double m_c_i
specific heat capacity of ice
double c() const
Specific heat capacity of ice.
double m_g
acceleration due to gravity
double enthalpy_permissive(double T, double omega, double P) const
Compute enthalpy more permissively than enthalpy().
double m_L
latent heat of fusion of water at atmospheric pressure
double pressure_adjusted_temperature(double E, double P) const
Get pressure-adjusted ice temperature, in Kelvin, from enthalpy and pressure.
double enthalpy_cts(double P) const
Get enthalpy E_s(p) at cold-temperate transition point from pressure p.
bool is_temperate_relaxed(double E, double P) const
A relaxed version of is_temperate().
double temperature_cold(double E) const
Convert enthalpy into temperature (cold case).
double pressure(double depth) const
Get pressure in ice from depth below surface using the hydrostatic assumption.
double m_T_melting
melting temperature of pure water at atmospheric pressure
double water_fraction(double E, double P) const
Get liquid water fraction from enthalpy and pressure.
double m_c_w
specific heat capacity of pure water
double m_p_air
atmospheric pressure
double L(double T_pm) const
void validate_T_omega_P(double T, double omega, double P) const
double m_T_0
reference temperature in the definition of ice enthalpy
Converts between specific enthalpy and temperature or liquid content.