tpism_utilities.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
---
tpism_utilities.cc (12402B)
---
1 /* Copyright (C) 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
20 #include "pism_utilities.hh"
21
22 #include <cstdarg> // va_list, va_start(), va_end()
23 #include <sstream> // istringstream, ostringstream
24 #include <cstdio> // vsnprintf
25
26 #include <mpi.h> // MPI_Get_library_version
27 #include <fftw3.h> // fftw_version
28 #include <gsl/gsl_version.h>
29
30 #include "pism/pism_config.hh" // Pism_USE_XXX, version info
31
32 // The following is a stupid kludge necessary to make NetCDF 4.x work in
33 // serial mode in an MPI program:
34 #ifndef MPI_INCLUDED
35 #define MPI_INCLUDED 1
36 #endif
37 #include <netcdf.h> // nc_inq_libvers
38
39 #if (Pism_USE_PROJ==1)
40 #include "pism/util/Proj.hh" // pj_release
41 #endif
42
43 #if (Pism_USE_JANSSON==1)
44 #include <jansson.h> // JANSSON_VERSION
45 #endif
46
47 #include <petsctime.h> // PetscTime
48
49 #include "error_handling.hh"
50
51 namespace pism {
52
53 //! Returns true if `str` ends with `suffix` and false otherwise.
54 bool ends_with(const std::string &str, const std::string &suffix) {
55 if (suffix.size() > str.size()) {
56 return false;
57 }
58
59 if (str.rfind(suffix) + suffix.size() == str.size()) {
60 return true;
61 }
62
63 return false;
64 }
65
66 template <class T>
67 std::string join_impl(const T& input, const std::string& separator) {
68 auto j = input.begin();
69 std::string result = *j;
70 ++j;
71 while (j != input.end()) {
72 result += separator + *j;
73 ++j;
74 }
75 return result;
76 }
77
78 //! Concatenate `strings`, inserting `separator` between elements.
79 std::string join(const std::vector<std::string> &strings, const std::string &separator) {
80 return join_impl(strings, separator);
81 }
82
83 std::string set_join(const std::set<std::string> &input, const std::string& separator) {
84 return join_impl(input, separator);
85 }
86
87 //! Transform a `separator`-separated list (a string) into a vector of strings.
88 std::vector<std::string> split(const std::string &input, char separator) {
89 std::istringstream input_list(input);
90 std::string token;
91 std::vector<std::string> result;
92
93 while (getline(input_list, token, separator)) {
94 if (not token.empty()) {
95 result.push_back(token);
96 }
97 }
98 return result;
99 }
100
101 //! Transform a `separator`-separated list (a string) into a set of strings.
102 std::set<std::string> set_split(const std::string &input, char separator) {
103 std::istringstream input_list(input);
104 std::string token;
105 std::set<std::string> result;
106
107 while (getline(input_list, token, separator)) {
108 if (not token.empty()) {
109 result.insert(token);
110 }
111 }
112 return result;
113 }
114
115 //! Checks if a vector of doubles is strictly increasing.
116 bool is_increasing(const std::vector<double> &a) {
117 int len = (int)a.size();
118 for (int k = 0; k < len-1; k++) {
119 if (a[k] >= a[k+1]) {
120 return false;
121 }
122 }
123 return true;
124 }
125
126 bool member(const std::string &string, const std::set<std::string> &set) {
127 return (set.find(string) != set.end());
128 }
129
130 void GlobalReduce(MPI_Comm comm, double *local, double *result, int count, MPI_Op op) {
131 int err = MPI_Allreduce(local, result, count, MPI_DOUBLE, op, comm);
132 PISM_C_CHK(err, 0, "MPI_Allreduce");
133 }
134
135 void GlobalMin(MPI_Comm comm, double *local, double *result, int count) {
136 GlobalReduce(comm, local, result, count, MPI_MIN);
137 }
138
139 void GlobalMax(MPI_Comm comm, double *local, double *result, int count) {
140 GlobalReduce(comm, local, result, count, MPI_MAX);
141 }
142
143 void GlobalSum(MPI_Comm comm, double *local, double *result, int count) {
144 GlobalReduce(comm, local, result, count, MPI_SUM);
145 }
146
147 unsigned int GlobalSum(MPI_Comm comm, unsigned int input) {
148 unsigned int result;
149 int err = MPI_Allreduce(&input, &result, 1, MPI_UNSIGNED, MPI_SUM, comm);
150 PISM_C_CHK(err, 0, "MPI_Allreduce");
151 return result;
152 }
153
154 int GlobalSum(MPI_Comm comm, int input) {
155 int result;
156 int err = MPI_Allreduce(&input, &result, 1, MPI_INT, MPI_SUM, comm);
157 PISM_C_CHK(err, 0, "MPI_Allreduce");
158 return result;
159 }
160
161 double GlobalMin(MPI_Comm comm, double local) {
162 double result;
163 GlobalMin(comm, &local, &result, 1);
164 return result;
165 }
166
167 double GlobalMax(MPI_Comm comm, double local) {
168 double result;
169 GlobalMax(comm, &local, &result, 1);
170 return result;
171 }
172
173 double GlobalSum(MPI_Comm comm, double local) {
174 double result;
175 GlobalSum(comm, &local, &result, 1);
176 return result;
177 }
178
179 static const int TEMPORARY_STRING_LENGTH = 32768;
180
181 std::string version() {
182 char buffer[TEMPORARY_STRING_LENGTH];
183 std::string result;
184
185 snprintf(buffer, sizeof(buffer), "PISM (%s)\n", pism::revision);
186 result += buffer;
187
188 snprintf(buffer, sizeof(buffer), "CMake %s.\n", pism::cmake_version);
189 result += buffer;
190
191 PetscGetVersion(buffer, TEMPORARY_STRING_LENGTH);
192 result += buffer;
193 result += "\n";
194
195 snprintf(buffer, sizeof(buffer), "PETSc configure: %s\n",
196 pism::petsc_configure_flags);
197 result += buffer;
198
199 // OpenMPI added MPI_Get_library_version in version 1.7 (relatively recently).
200 #ifdef OPEN_MPI
201 snprintf(buffer, TEMPORARY_STRING_LENGTH, "OpenMPI %d.%d.%d\n",
202 OMPI_MAJOR_VERSION, OMPI_MINOR_VERSION, OMPI_RELEASE_VERSION);
203 #else
204 // Assume that other MPI libraries implement this part of the MPI-3 standard...
205 int string_length = TEMPORARY_STRING_LENGTH;
206 MPI_Get_library_version(buffer, &string_length);
207 #endif
208 result += buffer;
209
210 snprintf(buffer, sizeof(buffer), "NetCDF %s.\n", nc_inq_libvers());
211 result += buffer;
212
213 snprintf(buffer, sizeof(buffer), "FFTW %s.\n", fftw_version);
214 result += buffer;
215
216 snprintf(buffer, sizeof(buffer), "GSL %s.\n", GSL_VERSION);
217 result += buffer;
218
219 #if (Pism_USE_PROJ==1)
220 snprintf(buffer, sizeof(buffer), "PROJ %s.\n", pj_release);
221 result += buffer;
222 #endif
223
224 #if (Pism_USE_JANSSON==1)
225 snprintf(buffer, sizeof(buffer), "Jansson %s.\n", JANSSON_VERSION);
226 result += buffer;
227 #endif
228
229 #if (Pism_BUILD_PYTHON_BINDINGS==1)
230 snprintf(buffer, sizeof(buffer), "SWIG %s.\n", pism::swig_version);
231 result += buffer;
232
233 snprintf(buffer, sizeof(buffer), "petsc4py %s.\n", pism::petsc4py_version);
234 result += buffer;
235 #endif
236
237 return result;
238 }
239
240
241 //! Return time since the beginning of the run, in hours.
242 double wall_clock_hours(MPI_Comm com, double start_time) {
243 int rank = 0;
244 double result = 0.0;
245
246 MPI_Comm_rank(com, &rank);
247
248 ParallelSection rank0(com);
249 try {
250 if (rank == 0) {
251 result = (get_time() - start_time) / 3600.0;
252 }
253 } catch (...) {
254 rank0.failed();
255 }
256 rank0.check();
257
258 MPI_Bcast(&result, 1, MPI_DOUBLE, 0, com);
259
260 return result;
261 }
262
263 //! Creates a time-stamp used for the history NetCDF attribute.
264 std::string timestamp(MPI_Comm com) {
265 time_t now;
266 tm tm_now;
267 char date_str[50];
268 now = time(NULL);
269 localtime_r(&now, &tm_now);
270 // Format specifiers for strftime():
271 // %F = ISO date format, %T = Full 24 hour time, %Z = Time Zone name
272 strftime(date_str, sizeof(date_str), "%F %T %Z", &tm_now);
273
274 MPI_Bcast(date_str, 50, MPI_CHAR, 0, com);
275
276 return std::string(date_str);
277 }
278
279 //! Creates a string with the user name, hostname and the time-stamp (for history strings).
280 std::string username_prefix(MPI_Comm com) {
281 PetscErrorCode ierr;
282
283 char username[50];
284 ierr = PetscGetUserName(username, sizeof(username));
285 PISM_CHK(ierr, "PetscGetUserName");
286 if (ierr != 0) {
287 username[0] = '\0';
288 }
289 char hostname[100];
290 ierr = PetscGetHostName(hostname, sizeof(hostname));
291 PISM_CHK(ierr, "PetscGetHostName");
292 if (ierr != 0) {
293 hostname[0] = '\0';
294 }
295
296 std::ostringstream message;
297 message << username << "@" << hostname << " " << timestamp(com) << ": ";
298
299 std::string result = message.str();
300 unsigned int length = result.size();
301 MPI_Bcast(&length, 1, MPI_UNSIGNED, 0, com);
302
303 result.resize(length);
304 MPI_Bcast(&result[0], length, MPI_CHAR, 0, com);
305
306 return result;
307 }
308
309 //! \brief Uses argc and argv to create the string with current PISM
310 //! command-line arguments.
311 std::string args_string() {
312 int argc;
313 char **argv;
314 PetscErrorCode ierr = PetscGetArgs(&argc, &argv);
315 PISM_CHK(ierr, "PetscGetArgs");
316
317 std::string cmdstr, argument;
318 for (int j = 0; j < argc; j++) {
319 argument = argv[j];
320
321 // enclose arguments containing spaces with double quotes:
322 if (argument.find(" ") != std::string::npos) {
323 argument = "\"" + argument + "\"";
324 }
325
326 cmdstr += std::string(" ") + argument;
327 }
328 cmdstr += "\n";
329
330 return cmdstr;
331 }
332
333 //! \brief Adds a suffix to a filename.
334 /*!
335 * Returns filename + separator + suffix + .nc if the original filename had the
336 * .nc suffix, otherwise filename + separator. If the old filename had the form
337 * "name + separator + more stuff + .nc", then removes the string after the
338 * separator.
339 */
340 std::string filename_add_suffix(const std::string &filename,
341 const std::string &separator,
342 const std::string &suffix) {
343 std::string basename = filename, result;
344
345 // find where the separator begins:
346 std::string::size_type j = basename.rfind(separator);
347 if (j == std::string::npos) {
348 j = basename.rfind(".nc");
349 }
350
351 // if the separator was not found, find the .nc suffix:
352 if (j == std::string::npos) {
353 j = basename.size();
354 }
355
356 // cut off everything starting from the separator (or the .nc suffix):
357 basename.resize(static_cast<int>(j));
358
359 result = basename + separator + suffix;
360
361 if (ends_with(filename, ".nc")) {
362 result += ".nc";
363 }
364
365 return result;
366 }
367
368 double get_time() {
369 PetscLogDouble result;
370 PetscErrorCode ierr = PetscTime(&result); PISM_CHK(ierr, "PetscTime");
371 return result;
372 }
373
374 std::string printf(const char *format, ...) {
375 std::string result(1024, ' ');
376 va_list arglist;
377 size_t length;
378
379 va_start(arglist, format);
380 if((length = vsnprintf(&result[0], result.size(), format, arglist)) > result.size()) {
381 result.reserve(length);
382 vsnprintf(&result[0], result.size(), format, arglist);
383 }
384 va_end(arglist);
385 return result.substr(0, length);
386 }
387
388 /*!
389 * Validate a format string. In this application a format string should contain `%` exactly
390 * once, followed by `s` (i.e. `%s`).
391 *
392 * Throws RuntimeError if the provided string is invalid.
393 */
394 void validate_format_string(const std::string &format) {
395 if (format.find("%s") == std::string::npos) {
396 throw RuntimeError::formatted(PISM_ERROR_LOCATION, "format string %s does not contain %%s",
397 format.c_str());
398 }
399
400 if (format.find("%") != format.rfind("%")) {
401 throw RuntimeError::formatted(PISM_ERROR_LOCATION, "format string %s contains more than one %%",
402 format.c_str());
403 }
404 }
405
406 double vector_min(const std::vector<double> &input) {
407 double my_min = input[0];
408 for (auto x : input) {
409 my_min = std::min(x, my_min);
410 }
411 return my_min;
412 }
413
414 double vector_max(const std::vector<double> &input) {
415 double my_max = input[0];
416 for (auto x : input) {
417 my_max = std::max(x, my_max);
418 }
419 return my_max;
420 }
421
422
423 /*!
424 * Fletcher's checksum
425 *
426 * See https://en.wikipedia.org/wiki/Fletcher%27s_checksum#Optimizations
427 */
428 uint64_t fletcher64(const uint32_t *data, size_t length) {
429 // Accumulating a sum of block_size unsigned 32-bit integers in an unsigned 64-bit
430 // integer will not lead to an overflow.
431 //
432 // This constant is found by solving n * (n + 1) / 2 * (2^32 - 1) < (2^64 - 1).
433 static const size_t block_size = 92681;
434
435 uint64_t c0 = 0, c1 = 0;
436 while (length != 0) {
437 size_t block = std::min(block_size, length);
438
439 for (size_t i = 0; i < block; ++i) {
440 c0 = c0 + *data++;
441 c1 = c1 + c0;
442 }
443
444 c0 = c0 % UINT32_MAX;
445 c1 = c1 % UINT32_MAX;
446
447 length = length > block_size ? length - block_size : 0;
448 }
449 return (c1 << 32 | c0);
450 }
451
452 } // end of namespace pism