ryujin 2.1.1 revision aeaed7a1ebbcf627d0fd44135960de4dec297929
checkpointing.h
Go to the documentation of this file.
1//
2// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
3// Copyright (C) 2020 - 2022 by the ryujin authors
4//
5
6#pragma once
7
9#include "offline_data.h"
10
11#include <deal.II/base/utilities.h>
12#include <deal.II/distributed/solution_transfer.h>
13
14#include <boost/archive/binary_iarchive.hpp>
15#include <boost/archive/binary_oarchive.hpp>
16#include <boost/core/demangle.hpp>
17
18#include <filesystem>
19#include <fstream>
20#include <string>
21
22namespace ryujin
23{
24 namespace Checkpointing
25 {
26 template <int dim>
28 std::is_same<typename Discretization<dim>::Triangulation,
29 dealii::parallel::distributed::Triangulation<dim>>::value;
30
38 template <int dim>
39 void load_mesh(Discretization<dim> &discretization,
40 const std::string &base_name)
41 {
42 if constexpr (have_distributed_triangulation<dim>) {
43 discretization.refinement() = 0; /* do not refine */
44 discretization.prepare();
45 discretization.triangulation().load(base_name + "-checkpoint.mesh");
46 } else {
47 AssertThrow(false, dealii::ExcNotImplemented());
48 __builtin_trap();
49 }
50 }
51
52
60 template <int dim, typename Number, int n_comp, int simd_length>
62 const std::string &base_name,
64 Number &t,
65 unsigned int &output_cycle,
66 const MPI_Comm &mpi_communicator)
67 {
68 if constexpr (have_distributed_triangulation<dim>) {
69 const auto &dof_handler = offline_data.dof_handler();
70
71 /* Create temporary scalar component vectors: */
72
73 const auto &scalar_partitioner = offline_data.scalar_partitioner();
74
75 using scalar_type = typename OfflineData<dim, Number>::scalar_type;
76 std::array<scalar_type, n_comp> state_vector;
77 for (auto &it : state_vector) {
78 it.reinit(scalar_partitioner);
79 }
80
81 /* Create SolutionTransfer object, attach state vector and deserialize:
82 */
83
84 dealii::parallel::distributed::SolutionTransfer<dim, scalar_type>
85 solution_transfer(dof_handler);
86
87 std::vector<scalar_type *> ptr_state;
88 std::transform(state_vector.begin(),
89 state_vector.end(),
90 std::back_inserter(ptr_state),
91 [](auto &it) { return &it; });
92
93 solution_transfer.deserialize(ptr_state);
94
95 unsigned int d = 0;
96 for (auto &it : state_vector) {
97 U.insert_component(it, d++);
98 }
99 U.update_ghost_values();
100
101 /* Read in and broadcast metadata: */
102
103 std::string name = base_name + "-checkpoint";
104
105 if (dealii::Utilities::MPI::this_mpi_process(mpi_communicator) == 0) {
106 std::string meta = name + ".metadata";
107
108 std::ifstream file(meta, std::ios::binary);
109 boost::archive::binary_iarchive ia(file);
110 ia >> t >> output_cycle;
111 }
112
113 int ierr;
114 if constexpr (std::is_same_v<Number, double>)
115 ierr = MPI_Bcast(&t, 1, MPI_DOUBLE, 0, mpi_communicator);
116 else
117 ierr = MPI_Bcast(&t, 1, MPI_FLOAT, 0, mpi_communicator);
118 AssertThrowMPI(ierr);
119
120 ierr = MPI_Bcast(&output_cycle, 1, MPI_UNSIGNED, 0, mpi_communicator);
121 AssertThrowMPI(ierr);
122
123 ierr = MPI_Barrier(mpi_communicator);
124 AssertThrowMPI(ierr);
125
126 } else {
127 AssertThrow(false, dealii::ExcNotImplemented());
128 __builtin_trap();
129 }
130 }
131
132
143 template <int dim, typename Number, int n_comp, int simd_length>
144 void
146 const std::string &base_name,
148 const Number t,
149 const unsigned int output_cycle,
150 const MPI_Comm &mpi_communicator)
151 {
152 if constexpr (have_distributed_triangulation<dim>) {
153 const auto &triangulation =
154 offline_data.discretization().triangulation();
155 const auto &dof_handler = offline_data.dof_handler();
156
157 /* Copy state into scalar component vectors: */
158
159 const auto &scalar_partitioner = offline_data.scalar_partitioner();
160
161 using scalar_type = typename OfflineData<dim, Number>::scalar_type;
162 std::array<scalar_type, n_comp> state_vector;
163 unsigned int d = 0;
164 for (auto &it : state_vector) {
165 it.reinit(scalar_partitioner);
166 U.extract_component(it, d++);
167 }
168
169 /* Create SolutionTransfer object, attach state vector and write out: */
170
171 dealii::parallel::distributed::SolutionTransfer<dim, scalar_type>
172 solution_transfer(dof_handler);
173
174 std::vector<const scalar_type *> ptr_state;
175 std::transform(state_vector.begin(),
176 state_vector.end(),
177 std::back_inserter(ptr_state),
178 [](auto &it) { return &it; });
179 solution_transfer.prepare_for_serialization(ptr_state);
180
181 std::string name = base_name + "-checkpoint";
182
183 if (dealii::Utilities::MPI::this_mpi_process(mpi_communicator) == 0) {
184 for (const std::string suffix :
185 {".mesh", ".mesh_fixed.data", ".mesh.info", ".metadata"})
186 if (std::filesystem::exists(name + suffix))
187 std::filesystem::rename(name + suffix, name + suffix + "~");
188 }
189
190 triangulation.save(name + ".mesh");
191
192 /* Metadata: */
193
194 if (dealii::Utilities::MPI::this_mpi_process(mpi_communicator) == 0) {
195 std::string meta = name + ".metadata";
196 std::ofstream file(meta, std::ios::binary | std::ios::trunc);
197 boost::archive::binary_oarchive oa(file);
198 oa << t << output_cycle;
199 }
200
201 const int ierr = MPI_Barrier(mpi_communicator);
202 AssertThrowMPI(ierr);
203
204 } else {
205 AssertThrow(false, dealii::ExcNotImplemented());
206 __builtin_trap();
207 }
208 }
209 } // namespace Checkpointing
210} // namespace ryujin
auto & triangulation() const
void insert_component(const scalar_type &scalar_vector, unsigned int component)
void extract_component(scalar_type &scalar_vector, unsigned int component) const
auto & dof_handler() const
Definition: offline_data.h:93
const auto & scalar_partitioner() const
Definition: offline_data.h:105
auto & discretization() const
Definition: offline_data.h:230
void write_checkpoint(const OfflineData< dim, Number > &offline_data, const std::string &base_name, const MultiComponentVector< Number, n_comp, simd_length > &U, const Number t, const unsigned int output_cycle, const MPI_Comm &mpi_communicator)
void load_mesh(Discretization< dim > &discretization, const std::string &base_name)
Definition: checkpointing.h:39
void load_state_vector(const OfflineData< dim, Number > &offline_data, const std::string &base_name, MultiComponentVector< Number, n_comp, simd_length > &U, Number &t, unsigned int &output_cycle, const MPI_Comm &mpi_communicator)
Definition: checkpointing.h:61
constexpr bool have_distributed_triangulation
Definition: checkpointing.h:27