Loading [MathJax]/extensions/tex2jax.js
PolyFEM
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
StateOutput.cpp
Go to the documentation of this file.
1#include <polyfem/State.hpp>
2
5
7
8#include <filesystem>
9
10namespace polyfem
11{
12 void State::compute_errors(const Eigen::MatrixXd &sol)
13 {
14 if (!args["output"]["advanced"]["compute_error"])
15 return;
16
17 double tend = 0;
18
19 if (!args["time"].is_null())
20 {
21 tend = args["time"]["tend"];
22 }
23
25 }
26
27 std::string State::root_path() const
28 {
29 if (utils::is_param_valid(args, "root_path"))
30 return args["root_path"].get<std::string>();
31 return "";
32 }
33
34 std::string State::resolve_input_path(const std::string &path, const bool only_if_exists) const
35 {
36 return utils::resolve_path(path, root_path(), only_if_exists);
37 }
38
39 std::string State::resolve_output_path(const std::string &path) const
40 {
41 if (output_dir.empty() || path.empty() || std::filesystem::path(path).is_absolute())
42 {
43 return path;
44 }
45 return std::filesystem::weakly_canonical(std::filesystem::path(output_dir) / path).string();
46 }
47
48 void State::save_timestep(const double time, const int t, const double t0, const double dt, const Eigen::MatrixXd &sol, const Eigen::MatrixXd &pressure)
49 {
50 if (args["output"]["advanced"]["save_time_sequence"] && !(t % args["output"]["paraview"]["skip_frame"].get<int>()))
51 {
52 logger().trace("Saving VTU...");
53 POLYFEM_SCOPED_TIMER("Saving VTU");
54 const std::string step_name = args["output"]["advanced"]["timestep_prefix"];
55
57 resolve_output_path(fmt::format(step_name + "{:d}.vtu", t)),
58 *this, sol, pressure, time, dt,
60
62 resolve_output_path(args["output"]["paraview"]["file_name"]),
63 [step_name](int i) { return fmt::format(step_name + "{:d}.vtm", i); },
64 t, t0, dt, args["output"]["paraview"]["skip_frame"].get<int>());
65 }
66 }
67
68 void State::save_json(const Eigen::MatrixXd &sol)
69 {
70 const std::string out_path = resolve_output_path(args["output"]["json"]);
71 if (!out_path.empty())
72 {
73 std::ofstream out(out_path);
74 if (!out.is_open())
75 {
76 logger().error("Unable to save simulation JSON to {}", out_path);
77 return;
78 }
79 save_json(sol, out);
80 out.close();
81 }
82 }
83
84 void State::save_json(const Eigen::MatrixXd &sol, std::ostream &out)
85 {
86 if (!mesh)
87 {
88 logger().error("Load the mesh first!");
89 return;
90 }
91 if (sol.size() <= 0)
92 {
93 logger().error("Solve the problem first!");
94 return;
95 }
96
97 logger().info("Saving json...");
98
99 using json = nlohmann::json;
100 json j;
103 assembler->name(), iso_parametric(), args["output"]["advanced"]["sol_at_node"],
104 j);
105 out << j.dump(4) << std::endl;
106 }
107
108 void State::save_subsolve(const int i, const int t, const Eigen::MatrixXd &sol, const Eigen::MatrixXd &pressure)
109 {
110 if (!args["output"]["advanced"]["save_solve_sequence_debug"].get<bool>())
111 return;
112
113 double dt = 1;
114 if (!args["time"].is_null())
115 dt = args["time"]["dt"];
116
118 resolve_output_path(fmt::format("solve_{:d}.vtu", i)),
119 *this, sol, pressure, t, dt,
121 }
122
123 void State::export_data(const Eigen::MatrixXd &sol, const Eigen::MatrixXd &pressure)
124 {
125 if (!mesh)
126 {
127 logger().error("Load the mesh first!");
128 return;
129 }
130 if (n_bases <= 0)
131 {
132 logger().error("Build the bases first!");
133 return;
134 }
135 // if (rhs.size() <= 0)
136 // {
137 // logger().error("Assemble the rhs first!");
138 // return;
139 // }
140 if (sol.size() <= 0)
141 {
142 logger().error("Solve the problem first!");
143 return;
144 }
145
146 // Export vtu mesh of solution + wire mesh of deformed input
147 // + mesh colored with the bases
148 const std::string vis_mesh_path = resolve_output_path(args["output"]["paraview"]["file_name"]);
149 const std::string nodes_path = resolve_output_path(args["output"]["data"]["nodes"]);
150 const std::string solution_path = resolve_output_path(args["output"]["data"]["solution"]);
151 const std::string stress_path = resolve_output_path(args["output"]["data"]["stress_mat"]);
152 const std::string mises_path = resolve_output_path(args["output"]["data"]["mises"]);
153
154 double tend = args.value("tend", 1.0);
155 double dt = 1;
156 if (!args["time"].is_null())
157 dt = args["time"]["dt"];
158
160 *this, sol, pressure,
161 !args["time"].is_null(),
162 tend, dt,
163 io::OutGeometryData::ExportOptions(args, mesh->is_linear(), problem->is_scalar()),
164 vis_mesh_path,
165 nodes_path,
166 solution_path,
167 stress_path,
168 mises_path,
170
171 if (assembler->name() == "Electrostatics")
172 {
173 std::shared_ptr<assembler::Electrostatics> electrostatics_assembler = std::dynamic_pointer_cast<assembler::Electrostatics>(assembler);
174 double energy = electrostatics_assembler->compute_stored_energy(mesh->is_volume(), n_bases, bases, geom_bases(), ass_vals_cache, 0, sol);
175 double capacitance = 2 * energy;
176 logger().info("Capacitance computation: {}", capacitance);
177 }
178 }
179
180 void State::save_restart_json(const double t0, const double dt, const int t) const
181 {
182 const std::string restart_json_path = args["output"]["restart_json"];
183 if (restart_json_path.empty())
184 return;
185
186 json restart_json;
187 restart_json["root_path"] = root_path();
188 restart_json["common"] = root_path();
189 restart_json["time"] = {{"t0", t0 + dt * t}};
190
191 restart_json["space"] = R"({
192 "remesh": {
193 "collapse": {
194 "abs_max_edge_length": -1,
195 "rel_max_edge_length": -1
196 }
197 }
198 })"_json;
199 restart_json["space"]["remesh"]["collapse"]["abs_max_edge_length"] = std::min(
200 args["space"]["remesh"]["collapse"]["abs_max_edge_length"].get<double>(),
201 starting_min_edge_length * args["space"]["remesh"]["collapse"]["rel_max_edge_length"].get<double>());
202 restart_json["space"]["remesh"]["collapse"]["rel_max_edge_length"] = std::numeric_limits<float>::max();
203
204 std::string rest_mesh_path = args["output"]["data"]["rest_mesh"].get<std::string>();
205 if (!rest_mesh_path.empty())
206 {
207 rest_mesh_path = resolve_output_path(fmt::format(args["output"]["data"]["rest_mesh"], t));
208
209 std::vector<json> patch;
210 if (args["geometry"].is_array())
211 {
212 const std::vector<json> in_geometry = args["geometry"];
213 for (int i = 0; i < in_geometry.size(); ++i)
214 {
215 if (!in_geometry[i]["is_obstacle"].get<bool>())
216 {
217 patch.push_back({
218 {"op", "remove"},
219 {"path", fmt::format("/geometry/{}", i)},
220 });
221 }
222 }
223
224 const int remaining_geometry = in_geometry.size() - patch.size();
225 assert(remaining_geometry >= 0);
226
227 patch.push_back({
228 {"op", "add"},
229 {"path", fmt::format("/geometry/{}", remaining_geometry > 0 ? "0" : "-")},
230 {"value",
231 {
232 // TODO: this does not set the surface selections
233 {"mesh", rest_mesh_path},
234 }},
235 });
236 }
237 else
238 {
239 assert(args["geometry"].is_object());
240 patch.push_back({
241 {"op", "remove"},
242 {"path", "/geometry"},
243 });
244 patch.push_back({
245 {"op", "replace"},
246 {"path", "/geometry"},
247 {"value",
248 {
249 // TODO: this does not set the surface selections
250 {"mesh", rest_mesh_path},
251 }},
252 });
253 }
254
255 restart_json["patch"] = patch;
256 }
257
258 restart_json["input"] = {{
259 "data",
260 {
261 {"state", resolve_output_path(fmt::format(args["output"]["data"]["state"], t))},
262 },
263 }};
264
265 std::ofstream file(resolve_output_path(fmt::format(restart_json_path, t)));
266 file << restart_json;
267 }
268} // namespace polyfem
#define POLYFEM_SCOPED_TIMER(...)
Definition Timer.hpp:10
io::OutRuntimeData timings
runtime statistics
Definition State.hpp:590
int n_bases
number of bases
Definition State.hpp:178
int n_pressure_bases
number of pressure bases
Definition State.hpp:180
assembler::AssemblyValsCache ass_vals_cache
used to store assembly values for small problems
Definition State.hpp:196
const std::vector< basis::ElementBases > & geom_bases() const
Get a constant reference to the geometry mapping bases.
Definition State.hpp:223
std::string root_path() const
Get the root path for the state (e.g., args["root_path"] or ".")
std::shared_ptr< assembler::Assembler > assembler
assemblers
Definition State.hpp:155
void compute_errors(const Eigen::MatrixXd &sol)
computes all errors
void save_subsolve(const int i, const int t, const Eigen::MatrixXd &sol, const Eigen::MatrixXd &pressure)
saves a subsolve when save_solve_sequence_debug is true
io::OutGeometryData out_geom
visualization stuff
Definition State.hpp:588
std::string resolve_output_path(const std::string &path) const
Resolve output path relative to output_dir if the path is not absolute.
std::string output_dir
Directory for output files.
Definition State.hpp:586
void save_timestep(const double time, const int t, const double t0, const double dt, const Eigen::MatrixXd &sol, const Eigen::MatrixXd &pressure)
saves a timestep
void save_json(const Eigen::MatrixXd &sol, std::ostream &out)
saves the output statistic to a stream
double starting_min_edge_length
Definition State.hpp:593
std::unique_ptr< mesh::Mesh > mesh
current mesh, it can be a Mesh2D or Mesh3D
Definition State.hpp:471
std::shared_ptr< assembler::Problem > problem
current problem, it contains rhs and bc
Definition State.hpp:168
json args
main input arguments containing all defaults
Definition State.hpp:101
void save_restart_json(const double t0, const double dt, const int t) const
Save a JSON sim file for restarting the simulation at time t.
io::OutStatsData stats
Other statistics.
Definition State.hpp:592
std::vector< basis::ElementBases > bases
FE bases, the size is #elements.
Definition State.hpp:171
bool iso_parametric() const
check if using iso parametric bases
Definition State.cpp:462
void export_data(const Eigen::MatrixXd &sol, const Eigen::MatrixXd &pressure)
saves all data on the disk according to the input params
std::string resolve_input_path(const std::string &path, const bool only_if_exists=false) const
Resolve input path relative to root_path() if the path is not absolute.
bool is_contact_enabled() const
does the simulation have contact
Definition State.hpp:556
Eigen::VectorXi disc_orders
vector of discretization orders, used when not all elements have the same degree, one per element
Definition State.hpp:190
void save_pvd(const std::string &name, const std::function< std::string(int)> &vtu_names, int time_steps, double t0, double dt, int skip_frame=1) const
save a PVD of a time dependent simulation
Definition OutData.cpp:2339
void save_vtu(const std::string &path, const State &state, const Eigen::MatrixXd &sol, const Eigen::MatrixXd &pressure, const double t, const double dt, const ExportOptions &opts, const bool is_contact_enabled) const
saves the vtu file for time t
Definition OutData.cpp:1123
void export_data(const State &state, const Eigen::MatrixXd &sol, const Eigen::MatrixXd &pressure, const bool is_time_dependent, const double tend_in, const double dt, const ExportOptions &opts, const std::string &vis_mesh_path, const std::string &nodes_path, const std::string &solution_path, const std::string &stress_path, const std::string &mises_path, const bool is_contact_enabled) const
exports everytihng, txt, vtu, etc
Definition OutData.cpp:940
void save_json(const nlohmann::json &args, const int n_bases, const int n_pressure_bases, const Eigen::MatrixXd &sol, const mesh::Mesh &mesh, const Eigen::VectorXi &disc_orders, const assembler::Problem &problem, const OutRuntimeData &runtime, const std::string &formulation, const bool isoparametric, const int sol_at_node_id, nlohmann::json &j)
saves the output statistic to a json object
Definition OutData.cpp:2815
void compute_errors(const int n_bases, const std::vector< polyfem::basis::ElementBases > &bases, const std::vector< polyfem::basis::ElementBases > &gbases, const polyfem::mesh::Mesh &mesh, const assembler::Problem &problem, const double tend, const Eigen::MatrixXd &sol)
compute errors
Definition OutData.cpp:2582
std::string resolve_path(const std::string &path, const std::string &input_file_path, const bool only_if_exists=false)
bool is_param_valid(const json &params, const std::string &key)
Determine if a key exists and is non-null in a json object.
spdlog::logger & logger()
Retrieves the current logger.
Definition Logger.cpp:44
nlohmann::json json
Definition Common.hpp:9