tEnthalpyConverter.hh - pism - [fork] customized build of PISM, the parallel ice sheet model (tillflux branch)
HTML git clone git://src.adamsgaard.dk/pism
DIR Log
DIR Files
DIR Refs
DIR LICENSE
---
tEnthalpyConverter.hh (4706B)
---
1 // Copyright (C) 2009-2011, 2013, 2014, 2015, 2016 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 */
51 class EnthalpyConverter {
52 public:
53 EnthalpyConverter(const Config &config);
54 virtual ~EnthalpyConverter();
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_m) 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
102 double m_T_tolerance;
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
107 bool m_do_cold_ice_methods;
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 */
134 class ColdEnthalpyConverter : public EnthalpyConverter {
135 public:
136 ColdEnthalpyConverter(const Config &config);
137 virtual ~ColdEnthalpyConverter();
138 };
139
140 } // end of namespace pism
141
142 #endif // __enthalpyConverter_hh
143