tConstant.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
---
tConstant.cc (3148B)
---
1 // Copyright (C) 2011, 2012, 2013, 2014, 2015, 2016, 2017, 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 #include "Constant.hh"
20
21 #include "pism/util/Vars.hh"
22 #include "pism/util/ConfigInterface.hh"
23 #include "pism/util/IceGrid.hh"
24 #include "pism/util/iceModelVec.hh"
25 #include "pism/util/MaxTimestep.hh"
26 #include "pism/util/pism_utilities.hh"
27 #include "pism/geometry/Geometry.hh"
28
29 namespace pism {
30 namespace ocean {
31
32 Constant::Constant(IceGrid::ConstPtr g)
33 : CompleteOceanModel(g) {
34 // empty
35 }
36
37 Constant::~Constant() {
38 // empty
39 }
40
41 void Constant::update_impl(const Geometry &geometry, double t, double dt) {
42 (void) t;
43 (void) dt;
44
45 const IceModelVec2S &ice_thickness = geometry.ice_thickness;
46
47 const double
48 melt_rate = m_config->get_number("ocean.constant.melt_rate", "m second-1"),
49 ice_density = m_config->get_number("constants.ice.density"),
50 mass_flux = melt_rate * ice_density;
51
52 melting_point_temperature(ice_thickness, *m_shelf_base_temperature);
53
54 m_shelf_base_mass_flux->set(mass_flux);
55
56 m_melange_back_pressure_fraction->set(m_config->get_number("ocean.melange_back_pressure_fraction"));
57 }
58
59 void Constant::init_impl(const Geometry &geometry) {
60 (void) geometry;
61
62 if (not m_config->get_flag("ocean.always_grounded")) {
63 m_log->message(2, "* Initializing the constant ocean model...\n");
64 m_log->message(2, " Sub-shelf melt rate set to %f m/year.\n",
65 m_config->get_number("ocean.constant.melt_rate", "m year-1"));
66
67 }
68 }
69
70 MaxTimestep Constant::max_timestep_impl(double t) const {
71 (void) t;
72 return MaxTimestep("ocean constant");
73 }
74
75 /*!
76 * Compute melting point temperature of ice at the depth `depth`.
77 */
78 void Constant::melting_point_temperature(const IceModelVec2S& depth,
79 IceModelVec2S &result) const {
80 const double
81 T0 = m_config->get_number("constants.fresh_water.melting_point_temperature"),
82 beta_CC = m_config->get_number("constants.ice.beta_Clausius_Clapeyron"),
83 g = m_config->get_number("constants.standard_gravity"),
84 ice_density = m_config->get_number("constants.ice.density");
85
86 IceModelVec::AccessList list{&depth, &result};
87
88 for (Points p(*m_grid); p; p.next()) {
89 const int i = p.i(), j = p.j();
90 const double pressure = ice_density * g * depth(i, j); // FIXME issue #15
91
92 result(i, j) = T0 - beta_CC * pressure;
93 }
94 }
95
96 } // end of namespape ocean
97 } // end of namespace pism