tPointwiseIsostasy.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
---
tPointwiseIsostasy.cc (4048B)
---
1 // Copyright (C) 2010, 2011, 2013, 2014, 2015, 2016, 2017, 2018, 2019 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 #include "BedDef.hh"
20 #include "pism/util/IceGrid.hh"
21 #include "pism/util/Time.hh"
22 #include "pism/util/ConfigInterface.hh"
23 #include "pism/util/Vars.hh"
24 #include "pism/util/MaxTimestep.hh"
25
26 namespace pism {
27 namespace bed {
28
29 PointwiseIsostasy::PointwiseIsostasy(IceGrid::ConstPtr g)
30 : BedDef(g) {
31 m_load_last.create(m_grid, "load_last", WITH_GHOSTS, m_config->get_number("grid.max_stencil_width"));
32 }
33
34 PointwiseIsostasy::~PointwiseIsostasy() {
35 // empty
36 }
37
38 void PointwiseIsostasy::init_impl(const InputOptions &opts, const IceModelVec2S &ice_thickness,
39 const IceModelVec2S &sea_level_elevation) {
40
41 m_log->message(2,
42 "* Initializing the pointwise isostasy bed deformation model...\n");
43
44 BedDef::init_impl(opts, ice_thickness, sea_level_elevation);
45
46 // store the initial load
47 compute_load(m_topg, ice_thickness, sea_level_elevation, m_load_last);
48 }
49
50
51 void PointwiseIsostasy::bootstrap_impl(const IceModelVec2S &bed_elevation,
52 const IceModelVec2S &bed_uplift,
53 const IceModelVec2S &ice_thickness,
54 const IceModelVec2S &sea_level_elevation) {
55 BedDef::bootstrap_impl(bed_elevation, bed_uplift, ice_thickness, sea_level_elevation);
56
57 // store initial load and bed elevation
58 compute_load(bed_elevation, ice_thickness, sea_level_elevation, m_load_last);
59 m_topg_last.copy_from(bed_elevation);
60 }
61
62 MaxTimestep PointwiseIsostasy::max_timestep_impl(double t) const {
63 (void) t;
64 return MaxTimestep("bed_def iso");
65 }
66
67 //! Updates the pointwise isostasy model.
68 void PointwiseIsostasy::update_impl(const IceModelVec2S &ice_thickness,
69 const IceModelVec2S &sea_level_elevation,
70 double t, double dt) {
71 (void) t;
72
73 const double
74 mantle_density = m_config->get_number("bed_deformation.mantle_density"),
75 load_density = m_config->get_number("constants.ice.density"),
76 ocean_density = m_config->get_number("constants.sea_water.density"),
77 ice_density = m_config->get_number("constants.ice.density"),
78 f = load_density / mantle_density;
79
80 //! Our goal: topg = topg_last - f*(load - load_last)
81
82 IceModelVec::AccessList list{&m_topg, &m_topg_last,
83 &ice_thickness, &sea_level_elevation, &m_load_last};
84
85 ParallelSection loop(m_grid->com);
86 try {
87 for (Points p(*m_grid); p; p.next()) {
88 const int i = p.i(), j = p.j();
89
90 double load = compute_load(m_topg(i, j),
91 ice_thickness(i, j),
92 sea_level_elevation(i, j),
93 ice_density, ocean_density);
94
95 m_topg(i, j) = m_topg_last(i, j) - f * (load - m_load_last(i, j));
96 m_load_last(i, j) = load;
97 }
98 } catch (...) {
99 loop.failed();
100 }
101 loop.check();
102
103 //! Finally, we need to update bed uplift, topg_last and load_last.
104 compute_uplift(m_topg, m_topg_last, dt, m_uplift);
105
106 m_topg_last.copy_from(m_topg);
107
108 //! Increment the topg state counter. SIAFD relies on this!
109 m_topg.inc_state_counter();
110 }
111
112 } // end of namespace bed
113 } // end of namespace pism