toutput.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
---
toutput.cc (9490B)
---
1 // Copyright (C) 2004-2019 Jed Brown, Ed Bueler and 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 <cstring> // strncpy
20 #include <cstdio> // snprintf
21
22 #include <algorithm>
23 #include <set>
24
25 #include "IceModel.hh"
26
27 #include "pism/util/IceGrid.hh"
28 #include "pism/util/ConfigInterface.hh"
29 #include "pism/util/Diagnostic.hh"
30 #include "pism/util/Time.hh"
31 #include "pism/util/error_handling.hh"
32 #include "pism/util/io/File.hh"
33 #include "pism/util/pism_options.hh"
34
35 #include "pism/util/Vars.hh"
36 #include "pism/util/io/io_helpers.hh"
37 #include "pism/util/Profiling.hh"
38 #include "pism/util/pism_utilities.hh"
39 #include "pism/util/projection.hh"
40 #include "pism/util/Component.hh"
41 #include "pism/energy/utilities.hh"
42
43 namespace pism {
44
45 MaxTimestep reporting_max_timestep(const std::vector<double> ×, double t,
46 const std::string &description) {
47
48 const size_t N = times.size();
49 if (t >= times.back()) {
50 return MaxTimestep();
51 }
52
53 size_t j = 0;
54 double dt = 0.0;
55 if (t < times[0]) {
56 j = -1;
57 } else {
58 j = gsl_interp_bsearch(×[0], t, 0, N - 1);
59 }
60
61 dt = times[j + 1] - t;
62
63 // now make sure that we don't end up taking a time-step of less than 1
64 // second long
65 if (dt < 1.0) {
66 if (j + 2 < N) {
67 return MaxTimestep(times[j + 2] - t, description);
68 } else {
69 return MaxTimestep(description);
70 }
71 } else {
72 return MaxTimestep(dt, description);
73 }
74 }
75
76 //! Write time-independent metadata to a file.
77 void IceModel::write_metadata(const File &file, MappingTreatment mapping_flag,
78 HistoryTreatment history_flag) {
79 if (mapping_flag == WRITE_MAPPING) {
80 write_mapping(file);
81 }
82
83 m_config->write(file);
84
85 if (history_flag == PREPEND_HISTORY) {
86 VariableMetadata tmp = m_output_global_attributes;
87
88 std::string old_history = file.read_text_attribute("PISM_GLOBAL", "history");
89
90 tmp.set_name("PISM_GLOBAL");
91 tmp.set_string("history", tmp.get_string("history") + old_history);
92
93 io::write_attributes(file, tmp, PISM_DOUBLE);
94 } else {
95 io::write_attributes(file, m_output_global_attributes, PISM_DOUBLE);
96 }
97 }
98
99 //! Save model state in NetCDF format.
100 /*!
101 Calls save_variables() to do the actual work.
102 */
103 void IceModel::save_results() {
104 {
105 update_run_stats();
106
107 auto str = pism::printf(
108 "PISM done. Performance stats: %.4f wall clock hours, %.4f proc.-hours, %.4f model years per proc.-hour.",
109 m_run_stats.get_number("wall_clock_hours"),
110 m_run_stats.get_number("processor_hours"),
111 m_run_stats.get_number("model_years_per_processor_hour"));
112
113 prepend_history(str);
114 }
115
116 std::string filename = m_config->get_string("output.file_name");
117
118 if (filename.empty()) {
119 m_log->message(2, "WARNING: output.file_name is empty. Using unnamed.nc instead.\n");
120 filename = "unnamed.nc";
121 }
122
123 if (not ends_with(filename, ".nc")) {
124 m_log->message(2,
125 "PISM WARNING: output file name does not have the '.nc' suffix!\n");
126 }
127
128 const Profiling &profiling = m_ctx->profiling();
129
130 profiling.begin("io.model_state");
131 if (m_config->get_string("output.size") != "none") {
132 m_log->message(2, "Writing model state to file `%s'...\n", filename.c_str());
133 File file(m_grid->com,
134 filename,
135 string_to_backend(m_config->get_string("output.format")),
136 PISM_READWRITE_MOVE,
137 m_ctx->pio_iosys_id());
138
139 write_metadata(file, WRITE_MAPPING, PREPEND_HISTORY);
140
141 write_run_stats(file);
142
143 save_variables(file, INCLUDE_MODEL_STATE, m_output_vars,
144 m_time->current());
145 }
146 profiling.end("io.model_state");
147 }
148
149 void IceModel::write_mapping(const File &file) {
150 // only write mapping if it is set.
151 const VariableMetadata &mapping = m_grid->get_mapping_info().mapping;
152 std::string name = mapping.get_name();
153 if (mapping.has_attributes()) {
154 if (not file.find_variable(name)) {
155 file.define_variable(name, PISM_DOUBLE, {});
156 }
157 io::write_attributes(file, mapping, PISM_DOUBLE);
158
159 // Write the PROJ string to mapping:proj_params (for CDO).
160 std::string proj = m_grid->get_mapping_info().proj;
161 if (not proj.empty()) {
162 file.write_attribute(name, "proj_params", proj);
163 }
164 }
165 }
166
167 void IceModel::write_run_stats(const File &file) {
168 update_run_stats();
169 if (not file.find_variable(m_run_stats.get_name())) {
170 file.define_variable(m_run_stats.get_name(), PISM_DOUBLE, {});
171 }
172 io::write_attributes(file, m_run_stats, PISM_DOUBLE);
173 }
174
175 void IceModel::save_variables(const File &file,
176 OutputKind kind,
177 const std::set<std::string> &variables,
178 double time,
179 IO_Type default_diagnostics_type) {
180
181 // define the time dimension if necessary (no-op if it is already defined)
182 io::define_time(file, *m_grid->ctx());
183 // define the "timestamp" (wall clock time since the beginning of the run)
184 // Note: it is time-dependent, so we need to define time first.
185 io::define_timeseries(m_timestamp, file, PISM_FLOAT);
186 // append to the time dimension
187 io::append_time(file, *m_config, time);
188
189 // Write metadata *before* everything else:
190 //
191 // FIXME: we should write this to variables instead of attributes because NetCDF-4 crashes after
192 // about 2^16 attribute modifications per variable. :-(
193 write_run_stats(file);
194
195 if (kind == INCLUDE_MODEL_STATE) {
196 define_model_state(file);
197 }
198 define_diagnostics(file, variables, default_diagnostics_type);
199
200 // Done defining variables
201
202 {
203 // Note: we don't use "variables" (an argument of this method) here because it
204 // contains PISM's names of diagnostic quantities which (in some cases) map to more
205 // than one NetCDF variable. Moreover, here we're concerned with file contents, not
206 // the list of requested variables.
207 std::set<std::string> var_names;
208 unsigned int n_vars = file.nvariables();
209 for (unsigned int k = 0; k < n_vars; ++k) {
210 var_names.insert(file.variable_name(k));
211 }
212
213 // If this output file contains variables lat and lon...
214 if (member("lat", var_names) and member("lon", var_names)) {
215
216 // add the coordinates attribute to all variables that use x and y dimensions
217 for (auto v : var_names) {
218 std::set<std::string> dims;
219 for (auto d : file.dimensions(v)) {
220 dims.insert(d);
221 }
222
223 if (not member(v, {"lat", "lon", "lat_bnds", "lon_bnds"}) and
224 member("x", dims) and member("y", dims)) {
225 file.write_attribute(v, "coordinates", "lat lon");
226 }
227 }
228
229 // and if it also contains lat_bnds and lon_bnds, add the bounds attribute to lat
230 // and lon.
231 if (member("lat_bnds", var_names) and member("lon_bnds", var_names)) {
232 file.write_attribute("lat", "bounds", "lat_bnds");
233 file.write_attribute("lon", "bounds", "lon_bnds");
234 }
235 }
236 }
237
238 if (kind == INCLUDE_MODEL_STATE) {
239 write_model_state(file);
240 }
241 write_diagnostics(file, variables);
242
243 // find out how much time passed since the beginning of the run and save it to the output file
244 {
245 unsigned int time_length = file.dimension_length(m_config->get_string("time.dimension_name"));
246 size_t start = time_length > 0 ? static_cast<size_t>(time_length - 1) : 0;
247 io::write_timeseries(file, m_timestamp, start,
248 wall_clock_hours(m_grid->com, m_start_time));
249 }
250 }
251
252 void IceModel::define_diagnostics(const File &file, const std::set<std::string> &variables,
253 IO_Type default_type) {
254 for (auto variable : variables) {
255 auto diag = m_diagnostics.find(variable);
256
257 if (diag != m_diagnostics.end()) {
258 diag->second->define(file, default_type);
259 }
260 }
261 }
262
263 //! \brief Writes variables listed in vars to filename, using nctype to write
264 //! fields stored in dedicated IceModelVecs.
265 void IceModel::write_diagnostics(const File &file, const std::set<std::string> &variables) {
266 for (auto variable : variables) {
267 auto diag = m_diagnostics.find(variable);
268
269 if (diag != m_diagnostics.end()) {
270 diag->second->compute()->write(file);
271 }
272 }
273 }
274
275 void IceModel::define_model_state(const File &file) {
276 for (auto v : m_model_state) {
277 v->define(file);
278 }
279
280 for (auto m : m_submodels) {
281 m.second->define_model_state(file);
282 }
283
284 for (auto d : m_diagnostics) {
285 d.second->define_state(file);
286 }
287 }
288
289 void IceModel::write_model_state(const File &file) {
290 for (auto v : m_model_state) {
291 v->write(file);
292 }
293
294 for (auto m : m_submodels) {
295 m.second->write_model_state(file);
296 }
297
298 for (auto d : m_diagnostics) {
299 d.second->write_state(file);
300 }
301 }
302
303 } // end of namespace pism