tPatersonBuddCold.cc - 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
---
tPatersonBuddCold.cc (2536B)
---
1 /* Copyright (C) 2015, 2016, 2018 PISM Authors
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
20 #include <cmath>
21
22 #include "PatersonBuddCold.hh"
23
24 namespace pism {
25 namespace rheology {
26
27 PatersonBuddCold::PatersonBuddCold(const std::string &prefix,
28 const Config &config,
29 EnthalpyConverter::Ptr ec)
30 : PatersonBudd(prefix, config, ec) {
31 m_name = "Paterson-Budd (cold case)";
32 }
33
34 double PatersonBuddCold::tempFromSoftness(double A) const {
35 return - m_Q_cold / (m_ideal_gas_constant * (log(A) - log(m_A_cold)));
36 }
37
38 PatersonBuddCold::~PatersonBuddCold() {
39 // empty
40 }
41
42 // takes care of hardness...
43 double PatersonBuddCold::softness_from_temp(double T_pa) const {
44 return m_A_cold * exp(-m_Q_cold / (m_ideal_gas_constant * T_pa));
45 }
46
47 // ignores pressure and uses non-pressure-adjusted temperature
48 double PatersonBuddCold::flow_from_temp(double stress, double temp,
49 double , double) const {
50 return softness_from_temp(temp) * pow(stress,m_n-1);
51 }
52
53
54 // Rather than make this part of the base class, we just check at some reference values.
55 bool FlowLawIsPatersonBuddCold(const FlowLaw &flow_law, const Config &config,
56 EnthalpyConverter::Ptr EC) {
57 static const struct {double s, E, p, gs;} v[] = {
58 {1e3, 223, 1e6, 1e-3}, {450000, 475000, 500000, 525000}, {5e4, 268, 5e6, 3e-3}, {1e5, 273, 8e6, 5e-3}};
59 PatersonBuddCold cpb("stress_balance.sia.", config, EC); // This is unmodified cold Paterson-Budd
60 for (int i=0; i<4; i++) {
61 const double left = flow_law.flow(v[i].s, v[i].E, v[i].p, v[i].gs),
62 right = cpb.flow(v[i].s, v[i].E, v[i].p, v[i].gs);
63 if (fabs((left - right)/left)>1.0e-15) {
64 return false;
65 }
66 }
67 return true;
68 }
69
70 } // end of namespace rheology
71 } // end of namespace pism