tOrographicPrecipitation.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
---
tOrographicPrecipitation.cc (4146B)
---
1 // Copyright (C) 2011, 2012, 2013, 2014, 2015, 2016, 2017, 2018, 2019, 2020 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 "OrographicPrecipitation.hh"
20
21 #include "OrographicPrecipitationSerial.hh"
22 #include "pism/coupler/util/options.hh"
23 #include "pism/geometry/Geometry.hh"
24 #include "pism/util/ConfigInterface.hh"
25 #include "pism/util/IceGrid.hh"
26 #include "pism/util/Time.hh"
27
28 namespace pism {
29 namespace atmosphere {
30
31 OrographicPrecipitation::OrographicPrecipitation(IceGrid::ConstPtr grid,
32 std::shared_ptr<AtmosphereModel> in)
33 : AtmosphereModel(grid, in) {
34
35 m_precipitation = allocate_precipitation(grid);
36
37 m_work0 = m_precipitation->allocate_proc0_copy();
38
39 const int
40 Mx = m_grid->Mx(),
41 My = m_grid->My(),
42 Z = m_config->get_number("atmosphere.orographic_precipitation.grid_size_factor"),
43 Nx = Z * (Mx - 1) + 1,
44 Ny = Z * (My - 1) + 1;
45
46 ParallelSection rank0(m_grid->com);
47 try {
48 if (m_grid->rank() == 0) {
49 m_serial_model.reset(new OrographicPrecipitationSerial(*m_config,
50 Mx, My,
51 m_grid->dx(), m_grid->dy(),
52 Nx, Ny));
53 }
54 } catch (...) {
55 rank0.failed();
56 }
57 rank0.check();
58 }
59
60 OrographicPrecipitation::~OrographicPrecipitation() {
61 // empty
62 }
63
64 const IceModelVec2S &OrographicPrecipitation::mean_precipitation_impl() const {
65 return *m_precipitation;
66 }
67
68 void OrographicPrecipitation::init_impl(const Geometry &geometry) {
69 (void)geometry;
70
71 m_input_model->init(geometry);
72
73 m_log->message(2, "* Initializing the atmosphere model computing precipitation using the\n"
74 " Linear Theory of Orographic Precipitation model with scalar wind speeds...\n");
75
76 m_reference = "R. B. Smith and I. Barstad, 2004.\n"
77 "A Linear Theory of Orographic Precipitation. J. Atmos. Sci. 61, 1377-1391.";
78
79 m_precipitation->metadata().set_string("source", m_reference);
80 }
81
82
83 void OrographicPrecipitation::update_impl(const Geometry &geometry, double t, double dt) {
84 m_input_model->update(geometry, t, dt);
85
86 geometry.ice_surface_elevation.put_on_proc0(*m_work0);
87
88 ParallelSection rank0(m_grid->com);
89 try {
90 if (m_grid->rank() == 0) { // processor zero updates the precipitation
91 m_serial_model->update(*m_work0);
92
93 PetscErrorCode ierr = VecCopy(m_serial_model->precipitation(), *m_work0);
94 PISM_CHK(ierr, "VecCopy");
95 }
96 } catch (...) {
97 rank0.failed();
98 }
99 rank0.check();
100
101 m_precipitation->get_from_proc0(*m_work0);
102
103 // convert from mm/s to kg / (m^2 s):
104 double water_density = m_config->get_number("constants.fresh_water.density");
105 m_precipitation->scale(1e-3 * water_density);
106 }
107
108 void OrographicPrecipitation::precip_time_series_impl(int i, int j,
109 std::vector<double> &result) const {
110
111 for (unsigned int k = 0; k < m_ts_times.size(); k++) {
112 result[k] = (*m_precipitation)(i, j);
113 }
114 }
115
116 void OrographicPrecipitation::begin_pointwise_access_impl() const {
117 m_input_model->begin_pointwise_access();
118 m_precipitation->begin_access();
119 }
120
121 void OrographicPrecipitation::end_pointwise_access_impl() const {
122 m_precipitation->end_access();
123 m_input_model->end_pointwise_access();
124 }
125
126 } // end of namespace atmosphere
127 } // end of namespace pism