tConstantPIK.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
---
tConstantPIK.cc (4029B)
---
1 // Copyright (C) 2008-2019 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 "ConstantPIK.hh"
20 #include "pism/util/io/File.hh"
21 #include "pism/util/Vars.hh"
22 #include "pism/util/IceGrid.hh"
23 #include "pism/util/pism_utilities.hh"
24 #include "pism/util/MaxTimestep.hh"
25 #include "pism/geometry/Geometry.hh"
26
27 namespace pism {
28 namespace surface {
29
30 ///// Constant-in-time surface model for accumulation,
31 ///// ice surface temperature parameterized as in PISM-PIK dependent on latitude and surface elevation
32
33
34 PIK::PIK(IceGrid::ConstPtr grid, std::shared_ptr<atmosphere::AtmosphereModel> atmosphere)
35 : SurfaceModel(grid) {
36 (void) atmosphere;
37
38 m_mass_flux = allocate_mass_flux(grid);
39 m_temperature = allocate_temperature(grid);
40 }
41
42 void PIK::init_impl(const Geometry &geometry) {
43 (void) geometry;
44
45 m_log->message(2,
46 "* Initializing the constant-in-time surface processes model PIK.\n"
47 " It reads surface mass balance directly from the file and holds it constant.\n"
48 " Ice upper-surface temperature is parameterized as in Martin et al. 2011, equation (1).\n"
49 " Any choice of atmosphere coupler (option '-atmosphere') is ignored.\n");
50
51 InputOptions opts = process_input_options(m_grid->com, m_config);
52
53 // read snow precipitation rate from file
54 m_log->message(2,
55 " reading surface mass balance rate 'climatic_mass_balance' from %s ... \n",
56 opts.filename.c_str());
57 if (opts.type == INIT_BOOTSTRAP) {
58 m_mass_flux->regrid(opts.filename, CRITICAL); // fails if not found!
59 } else {
60 m_mass_flux->read(opts.filename, opts.record); // fails if not found!
61 }
62
63 // parameterizing the ice surface temperature 'ice_surface_temp'
64 m_log->message(2,
65 " parameterizing the ice surface temperature 'ice_surface_temp' ... \n");
66 }
67
68 MaxTimestep PIK::max_timestep_impl(double t) const {
69 (void) t;
70 return MaxTimestep("surface PIK");
71 }
72
73 void PIK::update_impl(const Geometry &geometry, double t, double dt) {
74 (void) t;
75 (void) dt;
76
77 const IceModelVec2S
78 &surface_elevation = geometry.ice_surface_elevation,
79 &latitude = geometry.latitude;
80
81 IceModelVec::AccessList list{ m_temperature.get(), &surface_elevation, &latitude };
82
83 for (Points p(*m_grid); p; p.next()) {
84 const int i = p.i(), j = p.j();
85 (*m_temperature)(i, j) = 273.15 + 30 - 0.0075 * surface_elevation(i, j) - 0.68775 * latitude(i, j) * (-1.0);
86 }
87
88 dummy_accumulation(*m_mass_flux, *m_accumulation);
89 dummy_melt(*m_mass_flux, *m_melt);
90 dummy_runoff(*m_mass_flux, *m_runoff);
91
92 }
93
94 const IceModelVec2S &PIK::mass_flux_impl() const {
95 return *m_mass_flux;
96 }
97
98 const IceModelVec2S &PIK::temperature_impl() const {
99 return *m_temperature;
100 }
101
102 const IceModelVec2S &PIK::accumulation_impl() const {
103 return *m_accumulation;
104 }
105
106 const IceModelVec2S &PIK::melt_impl() const {
107 return *m_melt;
108 }
109
110 const IceModelVec2S &PIK::runoff_impl() const {
111 return *m_runoff;
112 }
113
114 void PIK::define_model_state_impl(const File &output) const {
115 m_mass_flux->define(output);
116 SurfaceModel::define_model_state_impl(output);
117 }
118
119 void PIK::write_model_state_impl(const File &output) const {
120 m_mass_flux->write(output);
121 SurfaceModel::write_model_state_impl(output);
122 }
123
124 } // end of namespace surface
125 } // end of namespace pism