Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
utilities.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2005 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16#include <deal.II/base/config.h>
17
18// It's necessary to include winsock2.h before thread_local_storage.h,
19// because Intel implementation of TBB includes winsock.h,
20// and we'll get a conflict between winsock.h and winsock2.h otherwise.
21#ifdef DEAL_II_MSVC
22# include <winsock2.h>
23#endif
24
26#include <deal.II/base/mpi.h>
27#include <deal.II/base/point.h>
30
31#define BOOST_BIND_GLOBAL_PLACEHOLDERS
32#include <boost/archive/iterators/base64_from_binary.hpp>
33#include <boost/archive/iterators/binary_from_base64.hpp>
34#include <boost/archive/iterators/transform_width.hpp>
35#include <boost/iostreams/copy.hpp>
36#include <boost/lexical_cast.hpp>
37#include <boost/random.hpp>
38#undef BOOST_BIND_GLOBAL_PLACEHOLDERS
39
40#include <algorithm>
41#include <bitset>
42#include <cctype>
43#include <cerrno>
44#include <cmath>
45#include <cstddef>
46#include <cstdio>
47#include <ctime>
48#include <fstream>
49#include <iomanip>
50#include <iostream>
51#include <limits>
52#include <sstream>
53#include <string>
54
55#if defined(DEAL_II_HAVE_UNISTD_H) && defined(DEAL_II_HAVE_GETHOSTNAME)
56# include <unistd.h>
57#endif
58
59#ifndef DEAL_II_MSVC
60# include <cstdlib>
61#endif
62
63
64#ifdef DEAL_II_WITH_TRILINOS
65# ifdef DEAL_II_WITH_MPI
69
70# include <Epetra_MpiComm.h>
71# include <Teuchos_DefaultComm.hpp>
72# endif
73# include <Epetra_SerialComm.h>
74# include <Teuchos_RCP.hpp>
75#endif
76
78
79
80namespace Utilities
81{
83 unsigned int,
84 unsigned int,
85 << "When trying to convert " << arg1 << " to a string with "
86 << arg2 << " digits");
87 DeclException1(ExcInvalidNumber, unsigned int, << "Invalid number " << arg1);
89 std::string,
90 << "Can't convert the string " << arg1
91 << " to the desired type");
92
93
94 std::string
96 {
98 }
99
100
101 namespace
102 {
103 template <int dim,
104 typename Number,
105 int effective_dim,
106 typename LongDouble,
107 typename Integer>
108 std::vector<std::array<std::uint64_t, effective_dim>>
109 inverse_Hilbert_space_filling_curve_effective(
110 const std::vector<Point<dim, Number>> &points,
111 const Point<dim, Number> & bl,
112 const std::array<LongDouble, dim> & extents,
113 const std::bitset<dim> & valid_extents,
114 const int min_bits,
115 const Integer max_int)
116 {
117 std::vector<std::array<Integer, effective_dim>> int_points(points.size());
118
119 for (unsigned int i = 0; i < points.size(); ++i)
120 {
121 // convert into integers:
122 unsigned int eff_d = 0;
123 for (unsigned int d = 0; d < dim; ++d)
124 if (valid_extents[d])
125 {
126 Assert(extents[d] > 0, ExcInternalError());
127 const LongDouble v = (static_cast<LongDouble>(points[i][d]) -
128 static_cast<LongDouble>(bl[d])) /
129 extents[d];
130 Assert(v >= 0. && v <= 1., ExcInternalError());
131 AssertIndexRange(eff_d, effective_dim);
132 int_points[i][eff_d] =
133 static_cast<Integer>(v * static_cast<LongDouble>(max_int));
134 ++eff_d;
135 }
136 }
137
138 // note that we call this with "min_bits"
139 return inverse_Hilbert_space_filling_curve<effective_dim>(int_points,
140 min_bits);
141 }
142 } // namespace
143
144 template <int dim, typename Number>
145 std::vector<std::array<std::uint64_t, dim>>
147 const std::vector<Point<dim, Number>> &points,
148 const int bits_per_dim)
149 {
150 using Integer = std::uint64_t;
151 // take floating point number hopefully with mantissa >= 64bit
152 using LongDouble = long double;
153
154 // return if there is nothing to do
155 if (points.size() == 0)
156 return std::vector<std::array<std::uint64_t, dim>>();
157
158 // get bounding box:
159 Point<dim, Number> bl = points[0], tr = points[0];
160 for (const auto &p : points)
161 for (unsigned int d = 0; d < dim; ++d)
162 {
163 const double cid = p[d];
164 bl[d] = std::min(cid, bl[d]);
165 tr[d] = std::max(cid, tr[d]);
166 }
167
168 std::array<LongDouble, dim> extents;
169 std::bitset<dim> valid_extents;
170 for (unsigned int i = 0; i < dim; ++i)
171 {
172 extents[i] =
173 static_cast<LongDouble>(tr[i]) - static_cast<LongDouble>(bl[i]);
174 valid_extents[i] = (extents[i] > 0.);
175 }
176
177 // make sure our conversion from fractional coordinates to
178 // Integers work as expected, namely our cast (LongDouble)max_int
179 const int min_bits =
180 std::min(bits_per_dim,
181 std::min(std::numeric_limits<Integer>::digits,
182 std::numeric_limits<LongDouble>::digits));
183
184 // based on that get the maximum integer:
185 const Integer max_int = (min_bits == std::numeric_limits<Integer>::digits ?
186 std::numeric_limits<Integer>::max() :
187 (Integer(1) << min_bits) - 1);
188
189 const unsigned int effective_dim = valid_extents.count();
190 if (effective_dim == dim)
191 {
192 return inverse_Hilbert_space_filling_curve_effective<dim,
193 Number,
194 dim,
195 LongDouble,
196 Integer>(
197 points, bl, extents, valid_extents, min_bits, max_int);
198 }
199
200 // various degenerate cases
201 std::array<std::uint64_t, dim> zero_ind;
202 for (unsigned int d = 0; d < dim; ++d)
203 zero_ind[d] = 0;
204
205 std::vector<std::array<std::uint64_t, dim>> ind(points.size(), zero_ind);
206 // manually check effective_dim == 1 and effective_dim == 2
207 if (dim == 3 && effective_dim == 2)
208 {
209 const auto ind2 =
210 inverse_Hilbert_space_filling_curve_effective<dim,
211 Number,
212 2,
213 LongDouble,
214 Integer>(
215 points, bl, extents, valid_extents, min_bits, max_int);
216
217 for (unsigned int i = 0; i < ind.size(); ++i)
218 for (unsigned int d = 0; d < 2; ++d)
219 ind[i][d + 1] = ind2[i][d];
220
221 return ind;
222 }
223 else if (effective_dim == 1)
224 {
225 const auto ind1 =
226 inverse_Hilbert_space_filling_curve_effective<dim,
227 Number,
228 1,
229 LongDouble,
230 Integer>(
231 points, bl, extents, valid_extents, min_bits, max_int);
232
233 for (unsigned int i = 0; i < ind.size(); ++i)
234 ind[i][dim - 1] = ind1[i][0];
235
236 return ind;
237 }
238
239 // we should get here only if effective_dim == 0
240 Assert(effective_dim == 0, ExcInternalError());
241
242 // if the bounding box is degenerate in all dimensions,
243 // can't do much but exit gracefully by setting index according
244 // to the index of each point so that there is no re-ordering
245 for (unsigned int i = 0; i < points.size(); ++i)
246 ind[i][dim - 1] = i;
247
248 return ind;
249 }
250
251
252
253 template <int dim>
254 std::vector<std::array<std::uint64_t, dim>>
256 const std::vector<std::array<std::uint64_t, dim>> &points,
257 const int bits_per_dim)
258 {
259 using Integer = std::uint64_t;
260
261 std::vector<std::array<Integer, dim>> int_points(points);
262
263 std::vector<std::array<Integer, dim>> res(int_points.size());
264
265 // follow
266 // J. Skilling, Programming the Hilbert curve, AIP Conf. Proc. 707, 381
267 // (2004); http://dx.doi.org/10.1063/1.1751381 also see
268 // https://stackoverflow.com/questions/499166/mapping-n-dimensional-value-to-a-point-on-hilbert-curve
269 // https://gitlab.com/octopus-code/octopus/blob/develop/src/grid/hilbert.c
270 // https://github.com/trilinos/Trilinos/blob/master/packages/zoltan/src/hsfc/hsfc_hilbert.c
271 // (Zoltan_HSFC_InvHilbertXd)
272 // https://github.com/aditi137/Hilbert/blob/master/Hilbert/hilbert.cpp
273
274 // now we can map to 1d coordinate stored in Transpose format
275 // adopt AxestoTranspose function from the paper, that
276 // transforms in-place between geometrical axes and Hilbert transpose.
277 // Example: b=5 bits for each of n=3 coordinates.
278 // 15-bit Hilbert integer = A B C D E F G H I J K L M N O is
279 // stored as its Transpose
280 // X[0] = A D G J M X[2]|
281 // X[1] = B E H K N <-------> | /X[1]
282 // X[2] = C F I L O axes |/
283 // high low 0------ X[0]
284
285 // Depth of the Hilbert curve
286 Assert(bits_per_dim <= std::numeric_limits<Integer>::digits,
287 ExcMessage("This integer type can not hold " +
288 std::to_string(bits_per_dim) + " bits."));
289
290 const Integer M = Integer(1) << (bits_per_dim - 1); // largest bit
291
292 for (unsigned int index = 0; index < int_points.size(); ++index)
293 {
294 auto &X = int_points[index];
295 auto &L = res[index];
296
297 // Inverse undo
298 for (Integer q = M; q > 1; q >>= 1)
299 {
300 const Integer p = q - 1;
301 for (unsigned int i = 0; i < dim; ++i)
302 {
303 // invert
304 if (X[i] & q)
305 {
306 X[0] ^= p;
307 }
308 // exchange
309 else
310 {
311 const Integer t = (X[0] ^ X[i]) & p;
312 X[0] ^= t;
313 X[i] ^= t;
314 }
315 }
316 }
317
318 // Gray encode (inverse of decode)
319 for (unsigned int i = 1; i < dim; ++i)
320 X[i] ^= X[i - 1];
321
322 Integer t = 0;
323 for (Integer q = M; q > 1; q >>= 1)
324 if (X[dim - 1] & q)
325 t ^= q - 1;
326 for (unsigned int i = 0; i < dim; ++i)
327 X[i] ^= t;
328
329 // now we need to go from index stored in transpose format to
330 // consecutive format, which is better suited for comparators.
331 // we could interleave into some big unsigned int...
332 // https://www.forceflow.be/2013/10/07/morton-encodingdecoding-through-bit-interleaving-implementations/
333 // https://stackoverflow.com/questions/4431522/given-2-16-bit-ints-can-i-interleave-those-bits-to-form-a-single-32-bit-int
334 // ...but we would loose spatial resolution!
335
336 // interleave using brute force, follow TransposetoLine from
337 // https://github.com/aditi137/Hilbert/blob/master/Hilbert/hilbert.cpp
338 {
339 Integer p = M;
340 unsigned int j = 0;
341 for (unsigned int i = 0; i < dim; ++i)
342 {
343 L[i] = 0;
344 // go through bits using a mask q
345 for (Integer q = M; q > 0; q >>= 1)
346 {
347 if (X[j] & p)
348 L[i] |= q;
349 if (++j == dim)
350 {
351 j = 0;
352 p >>= 1;
353 }
354 }
355 }
356 }
357
358 } // end of the loop over points
359
360 return res;
361 }
362
363
364
365 template <int dim>
366 std::uint64_t
367 pack_integers(const std::array<std::uint64_t, dim> &index,
368 const int bits_per_dim)
369 {
370 using Integer = std::uint64_t;
371
372 AssertIndexRange(bits_per_dim * dim, 65);
373 Assert(bits_per_dim > 0, ExcMessage("bits_per_dim should be positive"));
374
375 const Integer mask = (Integer(1) << bits_per_dim) - 1;
376
377 Integer res = 0;
378 for (unsigned int i = 0; i < dim; ++i)
379 {
380 // take bits_per_dim from each integer and shift them
381 const Integer v = (mask & index[dim - 1 - i]) << (bits_per_dim * i);
382 res |= v;
383 }
384 return res;
385 }
386
387
388
389 std::string
390 compress(const std::string &input)
391 {
392#ifdef DEAL_II_WITH_ZLIB
393 namespace bio = boost::iostreams;
394
395 std::stringstream compressed;
396 std::stringstream origin(input);
397
398 bio::filtering_streambuf<bio::input> out;
399 out.push(bio::gzip_compressor());
400 out.push(origin);
401 bio::copy(out, compressed);
402
403 return compressed.str();
404#else
405 return input;
406#endif
407 }
408
409
410
411 std::string
412 decompress(const std::string &compressed_input)
413 {
414#ifdef DEAL_II_WITH_ZLIB
415 namespace bio = boost::iostreams;
416
417 std::stringstream compressed(compressed_input);
418 std::stringstream decompressed;
419
420 bio::filtering_streambuf<bio::input> out;
421 out.push(bio::gzip_decompressor());
422 out.push(compressed);
423 bio::copy(out, decompressed);
424
425 return decompressed.str();
426#else
427 return compressed_input;
428#endif
429 }
430
431
432
433 std::string
434 encode_base64(const std::vector<unsigned char> &binary_input)
435 {
436 using namespace boost::archive::iterators;
437 using It = base64_from_binary<
438 transform_width<std::vector<unsigned char>::const_iterator, 6, 8>>;
439 auto base64 = std::string(It(binary_input.begin()), It(binary_input.end()));
440 // Add padding.
441 return base64.append((3 - binary_input.size() % 3) % 3, '=');
442 }
443
444
445
446 std::vector<unsigned char>
447 decode_base64(const std::string &base64_input)
448 {
449 using namespace boost::archive::iterators;
450 using It =
451 transform_width<binary_from_base64<std::string::const_iterator>, 8, 6>;
452 auto binary = std::vector<unsigned char>(It(base64_input.begin()),
453 It(base64_input.end()));
454 // Remove padding.
455 auto length = base64_input.size();
456 if (binary.size() > 2 && base64_input[length - 1] == '=' &&
457 base64_input[length - 2] == '=')
458 {
459 binary.erase(binary.end() - 2, binary.end());
460 }
461 else if (binary.size() > 1 && base64_input[length - 1] == '=')
462 {
463 binary.erase(binary.end() - 1, binary.end());
464 }
465 return binary;
466 }
467
468
469
470 std::string
471 int_to_string(const unsigned int value, const unsigned int digits)
472 {
473 return to_string(value, digits);
474 }
475
476
477
478 template <typename number>
479 std::string
480 to_string(const number value, const unsigned int digits)
481 {
482 // For integer data types, use the standard std::to_string()
483 // function. On the other hand, that function is defined in terms
484 // of std::sprintf, which does not use the usual std::iostream
485 // interface and tries to render floating point numbers in awkward
486 // ways (see
487 // https://en.cppreference.com/w/cpp/string/basic_string/to_string). So
488 // resort to boost::lexical_cast for all other types (in
489 // particular for floating point types.
490 std::string lc_string = (std::is_integral<number>::value ?
491 std::to_string(value) :
492 boost::lexical_cast<std::string>(value));
493
494 if ((digits != numbers::invalid_unsigned_int) &&
495 (lc_string.size() < digits))
496 {
497 // We have to add the padding zeroes in front of the number
498 const unsigned int padding_position = (lc_string[0] == '-') ? 1 : 0;
499
500 const std::string padding(digits - lc_string.size(), '0');
501 lc_string.insert(padding_position, padding);
502 }
503
504 return lc_string;
505 }
506
507
508
509 std::string
510 replace_in_string(const std::string &input,
511 const std::string &from,
512 const std::string &to)
513 {
514 if (from.empty())
515 return input;
516
517 std::string out = input;
518 std::string::size_type pos = out.find(from);
519
520 while (pos != std::string::npos)
521 {
522 out.replace(pos, from.size(), to);
523 pos = out.find(from, pos + to.size());
524 }
525 return out;
526 }
527
528 std::string
529 trim(const std::string &input)
530 {
531 std::string::size_type left = 0;
532 std::string::size_type right = input.size() > 0 ? input.size() - 1 : 0;
533
534 for (; left < input.size(); ++left)
535 {
536 if (std::isspace(input[left]) == 0)
537 {
538 break;
539 }
540 }
541
542 for (; right >= left; --right)
543 {
544 if (std::isspace(input[right]) == 0)
545 {
546 break;
547 }
548 }
549
550 return std::string(input, left, right - left + 1);
551 }
552
553
554
555 std::string
556 dim_string(const int dim, const int spacedim)
557 {
558 if (dim == spacedim)
559 return int_to_string(dim);
560 else
561 return int_to_string(dim) + "," + int_to_string(spacedim);
562 }
563
564
565 unsigned int
566 needed_digits(const unsigned int max_number)
567 {
568 if (max_number > 0)
569 return static_cast<int>(
570 std::ceil(std::log10(std::fabs(max_number + 0.1))));
571
572 return 1;
573 }
574
575
576
577 template <typename Number>
578 Number
579 truncate_to_n_digits(const Number number, const unsigned int n_digits)
580 {
581 AssertThrow(n_digits >= 1, ExcMessage("invalid parameter."));
582
583 if (!(std::fabs(number) > std::numeric_limits<Number>::min()))
584 return number;
585
586 const int order =
587 static_cast<int>(std::floor(std::log10(std::fabs(number))));
588
589 const int shift = -order + static_cast<int>(n_digits) - 1;
590
591 Assert(shift <= static_cast<int>(std::floor(
592 std::log10(std::numeric_limits<Number>::max()))),
594 "Overflow. Use a smaller value for n_digits and/or make sure "
595 "that the absolute value of 'number' does not become too small."));
596
597 const Number factor = std::pow(10.0, static_cast<Number>(shift));
598
599 const Number number_cutoff = std::trunc(number * factor) / factor;
600
601 return number_cutoff;
602 }
603
604
605 int
606 string_to_int(const std::string &s_)
607 {
608 // trim whitespace on either side of the text if necessary
609 std::string s = s_;
610 while ((s.size() > 0) && (s[0] == ' '))
611 s.erase(s.begin());
612 while ((s.size() > 0) && (s.back() == ' '))
613 s.erase(s.end() - 1);
614
615 // Now convert and see whether we succeed. Note that strtol only
616 // touches errno if an error occurred, so if we want to check
617 // whether an error happened, we need to make sure that errno==0
618 // before calling strtol since otherwise it may be that the
619 // conversion succeeds and that errno remains at the value it
620 // was before, whatever that was.
621 char *p;
622 errno = 0;
623 const int i = std::strtol(s.c_str(), &p, 10);
624
625 // We have an error if one of the following conditions is true:
626 // - strtol sets errno != 0
627 // - The original string was empty (we could have checked that
628 // earlier already)
629 // - The string has non-zero length and strtol converted the
630 // first part to something useful, but stopped converting short
631 // of the terminating '\0' character. This happens, for example,
632 // if the given string is "1234 abc".
633 AssertThrow(!((errno != 0) || (s.size() == 0) ||
634 ((s.size() > 0) && (*p != '\0'))),
635 ExcMessage("Can't convert <" + s + "> to an integer."));
636
637 return i;
638 }
639
640
641
642 std::vector<int>
643 string_to_int(const std::vector<std::string> &s)
644 {
645 std::vector<int> tmp(s.size());
646 for (unsigned int i = 0; i < s.size(); ++i)
647 tmp[i] = string_to_int(s[i]);
648 return tmp;
649 }
650
651
652
653 double
654 string_to_double(const std::string &s_)
655 {
656 // trim whitespace on either side of the text if necessary
657 std::string s = s_;
658 while ((s.size() > 0) && (s[0] == ' '))
659 s.erase(s.begin());
660 while ((s.size() > 0) && (s.back() == ' '))
661 s.erase(s.end() - 1);
662
663 // Now convert and see whether we succeed. Note that strtol only
664 // touches errno if an error occurred, so if we want to check
665 // whether an error happened, we need to make sure that errno==0
666 // before calling strtol since otherwise it may be that the
667 // conversion succeeds and that errno remains at the value it
668 // was before, whatever that was.
669 char *p;
670 errno = 0;
671 const double d = std::strtod(s.c_str(), &p);
672
673 // We have an error if one of the following conditions is true:
674 // - strtod sets errno != 0
675 // - The original string was empty (we could have checked that
676 // earlier already)
677 // - The string has non-zero length and strtod converted the
678 // first part to something useful, but stopped converting short
679 // of the terminating '\0' character. This happens, for example,
680 // if the given string is "1.234 abc".
681 AssertThrow(!((errno != 0) || (s.size() == 0) ||
682 ((s.size() > 0) && (*p != '\0'))),
683 ExcMessage("Can't convert <" + s + "> to a double."));
684
685 return d;
686 }
687
688
689
690 std::vector<double>
691 string_to_double(const std::vector<std::string> &s)
692 {
693 std::vector<double> tmp(s.size());
694 for (unsigned int i = 0; i < s.size(); ++i)
695 tmp[i] = string_to_double(s[i]);
696 return tmp;
697 }
698
699
700
701 std::vector<std::string>
702 split_string_list(const std::string &s, const std::string &delimiter)
703 {
704 // keep the currently remaining part of the input string in 'tmp' and
705 // keep chopping elements of the list off the front
706 std::string tmp = s;
707
708 // as discussed in the documentation, eat whitespace from the end
709 // of the string
710 while (tmp.size() != 0 && tmp.back() == ' ')
711 tmp.erase(tmp.size() - 1, 1);
712
713 // split the input list until it is empty. since in every iteration
714 // 'tmp' is what's left of the string after the next delimiter,
715 // and since we've stripped trailing space already, 'tmp' will
716 // be empty at one point if 's' ended in a delimiter, even if
717 // there was space after the last delimiter. this matches what's
718 // discussed in the documentation
719 std::vector<std::string> split_list;
720 while (tmp.size() != 0)
721 {
722 std::string name;
723 name = tmp;
724
725 if (name.find(delimiter) != std::string::npos)
726 {
727 name.erase(name.find(delimiter), std::string::npos);
728 tmp.erase(0, tmp.find(delimiter) + delimiter.size());
729 }
730 else
731 tmp = "";
732
733 // strip spaces from this element's front and end
734 while ((name.size() != 0) && (name[0] == ' '))
735 name.erase(0, 1);
736 while (name.size() != 0 && name.back() == ' ')
737 name.erase(name.size() - 1, 1);
738
739 split_list.push_back(name);
740 }
741
742 return split_list;
743 }
744
745
746 std::vector<std::string>
747 split_string_list(const std::string &s, const char delimiter)
748 {
749 std::string d = ",";
750 d[0] = delimiter;
751 return split_string_list(s, d);
752 }
753
754
755 std::vector<std::string>
756 break_text_into_lines(const std::string &original_text,
757 const unsigned int width,
758 const char delimiter)
759 {
760 std::string text = original_text;
761 std::vector<std::string> lines;
762
763 // remove trailing spaces
764 while ((text.size() != 0) && (text.back() == delimiter))
765 text.erase(text.size() - 1, 1);
766
767 // then split the text into lines
768 while (text.size() != 0)
769 {
770 // in each iteration, first remove
771 // leading spaces
772 while ((text.size() != 0) && (text[0] == delimiter))
773 text.erase(0, 1);
774
775 std::size_t pos_newline = text.find_first_of('\n', 0);
776 if (pos_newline != std::string::npos && pos_newline <= width)
777 {
778 std::string line(text, 0, pos_newline);
779 while ((line.size() != 0) && (line.back() == delimiter))
780 line.erase(line.size() - 1, 1);
781 lines.push_back(line);
782 text.erase(0, pos_newline + 1);
783 continue;
784 }
785
786 // if we can fit everything into one
787 // line, then do so. otherwise, we have
788 // to keep breaking
789 if (text.size() < width)
790 {
791 // remove trailing spaces
792 while ((text.size() != 0) && (text.back() == delimiter))
793 text.erase(text.size() - 1, 1);
794 lines.push_back(text);
795 text = "";
796 }
797 else
798 {
799 // starting at position width, find the
800 // location of the previous space, so
801 // that we can break around there
802 int location = std::min<int>(width, text.size() - 1);
803 for (; location > 0; --location)
804 if (text[location] == delimiter)
805 break;
806
807 // if there are no spaces, then try if
808 // there are spaces coming up
809 if (location == 0)
810 for (location = std::min<int>(width, text.size() - 1);
811 location < static_cast<int>(text.size());
812 ++location)
813 if (text[location] == delimiter)
814 break;
815
816 // now take the text up to the found
817 // location and put it into a single
818 // line, and remove it from 'text'
819 std::string line(text, 0, location);
820 while ((line.size() != 0) && (line.back() == delimiter))
821 line.erase(line.size() - 1, 1);
822 lines.push_back(line);
823 text.erase(0, location);
824 }
825 }
826
827 return lines;
828 }
829
830
831
832 bool
833 match_at_string_start(const std::string &name, const std::string &pattern)
834 {
835 if (pattern.size() > name.size())
836 return false;
837
838 for (unsigned int i = 0; i < pattern.size(); ++i)
839 if (pattern[i] != name[i])
840 return false;
841
842 return true;
843 }
844
845
846
847 std::pair<int, unsigned int>
848 get_integer_at_position(const std::string &name, const unsigned int position)
849 {
850 Assert(position < name.size(), ExcInternalError());
851
852 const std::string test_string(name.begin() + position, name.end());
853
854 std::istringstream str(test_string);
855
856 int i;
857 if (str >> i)
858 {
859 // compute the number of
860 // digits of i. assuming it
861 // is less than 8 is likely
862 // ok
863 if (i < 10)
864 return std::make_pair(i, 1U);
865 else if (i < 100)
866 return std::make_pair(i, 2U);
867 else if (i < 1000)
868 return std::make_pair(i, 3U);
869 else if (i < 10000)
870 return std::make_pair(i, 4U);
871 else if (i < 100000)
872 return std::make_pair(i, 5U);
873 else if (i < 1000000)
874 return std::make_pair(i, 6U);
875 else if (i < 10000000)
876 return std::make_pair(i, 7U);
877 else
878 {
879 Assert(false, ExcNotImplemented());
880 return std::make_pair(-1, numbers::invalid_unsigned_int);
881 }
882 }
883 else
884 return std::make_pair(-1, numbers::invalid_unsigned_int);
885 }
886
887
888
889 double
890 generate_normal_random_number(const double a, const double sigma)
891 {
892 // if no noise: return now
893 if (sigma == 0)
894 return a;
895
896 // we would want to use rand(), but that function is not reentrant
897 // in a thread context. one could use rand_r, but this does not
898 // produce reproducible results between threads either (though at
899 // least it is reentrant). these two approaches being
900 // non-workable, use a thread-local random number generator here.
901 // we could use std::mt19937 but doing so results in compiler-dependent
902 // output.
903 static Threads::ThreadLocalStorage<boost::mt19937> random_number_generator;
904 return boost::normal_distribution<>(a,
905 sigma)(random_number_generator.get());
906 }
907
908
909
910 namespace System
911 {
912#ifdef __linux__
913
914 double
916 {
917 std::ifstream cpuinfo;
918 cpuinfo.open("/proc/loadavg");
919
920 AssertThrow(cpuinfo.fail() == false, ExcIO());
921
922 double load;
923 cpuinfo >> load;
924
925 return load;
926 }
927
928#else
929
930 double
932 {
933 return 0.;
934 }
935
936#endif
937
938 const std::string
940 {
942 {
943 case 0:
944 return "disabled";
945 case 128:
946#ifdef __ALTIVEC__
947 return "AltiVec";
948#else
949 return "SSE2";
950#endif
951 case 256:
952 return "AVX";
953 case 512:
954 return "AVX512";
955 default:
956 AssertThrow(false,
958 "Invalid DEAL_II_VECTORIZATION_WIDTH_IN_BITS."));
959 return "ERROR";
960 }
961 }
962
963
964 void
966 {
967 stats.VmPeak = stats.VmSize = stats.VmHWM = stats.VmRSS = 0;
968
969 // parsing /proc/self/stat would be a
970 // lot easier, but it does not contain
971 // VmHWM, so we use /status instead.
972#ifdef __linux__
973 std::ifstream file("/proc/self/status");
974 std::string line;
975 std::string name;
976 while (!file.eof())
977 {
978 file >> name;
979 if (name == "VmPeak:")
980 file >> stats.VmPeak;
981 else if (name == "VmSize:")
982 file >> stats.VmSize;
983 else if (name == "VmHWM:")
984 file >> stats.VmHWM;
985 else if (name == "VmRSS:")
986 {
987 file >> stats.VmRSS;
988 break; // this is always the last entry
989 }
990
991 getline(file, line);
992 }
993#endif
994 }
995
996
997
998 std::string
1000 {
1001#if defined(DEAL_II_HAVE_UNISTD_H) && defined(DEAL_II_HAVE_GETHOSTNAME)
1002 const unsigned int N = 1024;
1003 char hostname[N];
1004 gethostname(&(hostname[0]), N - 1);
1005#else
1006 std::string hostname("unknown");
1007#endif
1008 return hostname;
1009 }
1010
1011
1012
1013 std::string
1015 {
1016 std::time_t time1 = std::time(nullptr);
1017 std::tm * time = std::localtime(&time1);
1018
1019 std::ostringstream o;
1020 o << time->tm_hour << ":" << (time->tm_min < 10 ? "0" : "")
1021 << time->tm_min << ":" << (time->tm_sec < 10 ? "0" : "")
1022 << time->tm_sec;
1023
1024 return o.str();
1025 }
1026
1027
1028
1029 std::string
1031 {
1032 std::time_t time1 = std::time(nullptr);
1033 std::tm * time = std::localtime(&time1);
1034
1035 std::ostringstream o;
1036 o << time->tm_year + 1900 << "/" << time->tm_mon + 1 << "/"
1037 << time->tm_mday;
1038
1039 return o.str();
1040 }
1041
1042
1043
1044 void
1045 posix_memalign(void **memptr, std::size_t alignment, std::size_t size)
1046 {
1047#ifndef DEAL_II_MSVC
1048 const int ierr = ::posix_memalign(memptr, alignment, size);
1049
1050 AssertThrow(ierr == 0, ExcOutOfMemory(size));
1051 AssertThrow(*memptr != nullptr, ExcOutOfMemory(size));
1052#else
1053 // Windows does not appear to have posix_memalign. just use the
1054 // regular malloc in that case
1055 *memptr = malloc(size);
1056 (void)alignment;
1057 AssertThrow(*memptr != 0, ExcOutOfMemory(size));
1058#endif
1059 }
1060
1061
1062
1063 bool
1065 {
1067 }
1068 } // namespace System
1069
1070#ifndef DOXYGEN
1071 template std::string
1072 to_string<int>(int, unsigned int);
1073 template std::string
1074 to_string<long int>(long int, unsigned int);
1075 template std::string
1076 to_string<long long int>(long long int, unsigned int);
1077 template std::string
1078 to_string<unsigned int>(unsigned int, unsigned int);
1079 template std::string
1080 to_string<unsigned long int>(unsigned long int, unsigned int);
1081 template std::string
1082 to_string<unsigned long long int>(unsigned long long int, unsigned int);
1083 template std::string
1084 to_string<float>(float, unsigned int);
1085 template std::string
1086 to_string<double>(double, unsigned int);
1087 template std::string
1088 to_string<long double>(long double, unsigned int);
1089
1090 template double
1091 truncate_to_n_digits(const double, const unsigned int);
1092 template float
1093 truncate_to_n_digits(const float, const unsigned int);
1094
1095 template std::vector<std::array<std::uint64_t, 1>>
1096 inverse_Hilbert_space_filling_curve<1, double>(
1097 const std::vector<Point<1, double>> &,
1098 const int);
1099 template std::vector<std::array<std::uint64_t, 1>>
1100 inverse_Hilbert_space_filling_curve<1>(
1101 const std::vector<std::array<std::uint64_t, 1>> &,
1102 const int);
1103 template std::vector<std::array<std::uint64_t, 2>>
1104 inverse_Hilbert_space_filling_curve<2, double>(
1105 const std::vector<Point<2, double>> &,
1106 const int);
1107 template std::vector<std::array<std::uint64_t, 2>>
1108 inverse_Hilbert_space_filling_curve<2>(
1109 const std::vector<std::array<std::uint64_t, 2>> &,
1110 const int);
1111 template std::vector<std::array<std::uint64_t, 3>>
1112 inverse_Hilbert_space_filling_curve<3, double>(
1113 const std::vector<Point<3, double>> &,
1114 const int);
1115 template std::vector<std::array<std::uint64_t, 3>>
1116 inverse_Hilbert_space_filling_curve<3>(
1117 const std::vector<std::array<std::uint64_t, 3>> &,
1118 const int);
1119
1120 template std::uint64_t
1121 pack_integers<1>(const std::array<std::uint64_t, 1> &, const int);
1122 template std::uint64_t
1123 pack_integers<2>(const std::array<std::uint64_t, 2> &, const int);
1124 template std::uint64_t
1125 pack_integers<3>(const std::array<std::uint64_t, 3> &, const int);
1126
1127#endif
1128
1129} // namespace Utilities
1130
Definition point.h:112
A class that provides a separate storage location on each thread that accesses the object.
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_PACKAGE_VERSION
Definition config.h:26
#define DEAL_II_VECTORIZATION_WIDTH_IN_BITS
Definition config.h:123
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define DEAL_II_PACKAGE_NAME
Definition config.h:24
static ::ExceptionBase & ExcInvalidNumber2StringConversersion(unsigned int arg1, unsigned int arg2)
static ::ExceptionBase & ExcOutOfMemory(std::size_t arg1)
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcInvalidNumber(unsigned int arg1)
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition exceptions.h:533
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcCantConvertString(std::string arg1)
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:510
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
bool job_supports_mpi()
Definition mpi.cc:1057
void get_memory_stats(MemoryStats &stats)
Definition utilities.cc:965
std::string get_hostname()
Definition utilities.cc:999
std::string get_time()
void posix_memalign(void **memptr, std::size_t alignment, std::size_t size)
double get_cpu_load()
Definition utilities.cc:931
const std::string get_current_vectorization_level()
Definition utilities.cc:939
std::string get_date()
std::vector< std::string > split_string_list(const std::string &s, const std::string &delimiter=",")
Definition utilities.cc:702
Number truncate_to_n_digits(const Number number, const unsigned int n_digits)
Definition utilities.cc:579
std::string dim_string(const int dim, const int spacedim)
Definition utilities.cc:556
std::uint64_t pack_integers(const std::array< std::uint64_t, dim > &index, const int bits_per_dim)
Definition utilities.cc:367
std::pair< int, unsigned int > get_integer_at_position(const std::string &name, const unsigned int position)
Definition utilities.cc:848
std::string encode_base64(const std::vector< unsigned char > &binary_input)
Definition utilities.cc:434
std::string replace_in_string(const std::string &input, const std::string &from, const std::string &to)
Definition utilities.cc:510
std::vector< unsigned char > decode_base64(const std::string &base64_input)
Definition utilities.cc:447
std::vector< std::string > break_text_into_lines(const std::string &original_text, const unsigned int width, const char delimiter=' ')
Definition utilities.cc:756
std::string to_string(const number value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition utilities.cc:480
std::vector< std::array< std::uint64_t, dim > > inverse_Hilbert_space_filling_curve(const std::vector< Point< dim, Number > > &points, const int bits_per_dim=64)
Definition utilities.cc:146
std::string compress(const std::string &input)
Definition utilities.cc:390
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition utilities.cc:471
bool match_at_string_start(const std::string &name, const std::string &pattern)
Definition utilities.cc:833
std::string decompress(const std::string &compressed_input)
Definition utilities.cc:412
unsigned int needed_digits(const unsigned int max_number)
Definition utilities.cc:566
double string_to_double(const std::string &s)
Definition utilities.cc:654
std::string dealii_version_string()
Definition utilities.cc:95
std::string trim(const std::string &input)
Definition utilities.cc:529
double generate_normal_random_number(const double a, const double sigma)
Definition utilities.cc:890
int string_to_int(const std::string &s)
Definition utilities.cc:606
static const unsigned int invalid_unsigned_int
Definition types.h:213
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > pow(const ::VectorizedArray< Number, width > &, const Number p)