Loading [MathJax]/extensions/tex2jax.js
PolyFEM
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
OutData.hpp
Go to the documentation of this file.
1#pragma once
2
3#include <polyfem/Common.hpp>
4
6
8
10
12
13#include <paraviewo/ParaviewWriter.hpp>
14#include <paraviewo/VTUWriter.hpp>
15#include <paraviewo/HDF5VTUWriter.hpp>
16
18
19#include <Eigen/Dense>
20
21#include <ipc/potentials/normal_adhesion_potential.hpp>
22#include <ipc/potentials/tangential_adhesion_potential.hpp>
23
24namespace polyfem
25{
26 class State;
27}
28
29namespace polyfem::io
30{
33 {
34 public:
37 {
38 std::vector<std::string> fields; // fields to export, empty means all
39
40 bool volume;
41 bool surface;
42 bool wire;
43 bool points;
48 bool forces;
50
61 bool nodes;
62
65
67
72 ExportOptions(const json &args,
73 const bool is_mesh_linear,
74 const bool is_problem_scalar);
75
78 inline std::string file_extension() const { return use_hdf5 ? ".hdf" : ".vtu"; }
79
80 bool export_field(const std::string &field) const;
81 };
82
92 static void extract_boundary_mesh(
93 const mesh::Mesh &mesh,
94 const int n_bases,
95 const std::vector<basis::ElementBases> &bases,
96 const std::vector<mesh::LocalBoundary> &total_local_boundary,
97 Eigen::MatrixXd &node_positions,
98 Eigen::MatrixXi &boundary_edges,
99 Eigen::MatrixXi &boundary_triangles,
100 std::vector<Eigen::Triplet<double>> &displacement_map_entries);
101
105 void init_sampler(const polyfem::mesh::Mesh &mesh, const double vismesh_rel_area);
106
110 void build_grid(const polyfem::mesh::Mesh &mesh, const double spacing);
111
126 void export_data(
127 const State &state,
128 const Eigen::MatrixXd &sol,
129 const Eigen::MatrixXd &pressure,
130 const bool is_time_dependent,
131 const double tend_in,
132 const double dt,
133 const ExportOptions &opts,
134 const std::string &vis_mesh_path,
135 const std::string &nodes_path,
136 const std::string &solution_path,
137 const std::string &stress_path,
138 const std::string &mises_path,
139 const bool is_contact_enabled) const;
140
150 void save_vtu(const std::string &path,
151 const State &state,
152 const Eigen::MatrixXd &sol,
153 const Eigen::MatrixXd &pressure,
154 const double t,
155 const double dt,
156 const ExportOptions &opts,
157 const bool is_contact_enabled) const;
158
167 void save_volume(const std::string &path,
168 const State &state,
169 const Eigen::MatrixXd &sol,
170 const Eigen::MatrixXd &pressure,
171 const double t,
172 const double dt,
173 const ExportOptions &opts) const;
174
184 void save_surface(const std::string &export_surface,
185 const State &state,
186 const Eigen::MatrixXd &sol,
187 const Eigen::MatrixXd &pressure,
188 const double t,
189 const double dt_in,
190 const ExportOptions &opts,
191 const bool is_contact_enabled) const;
192
203 const std::string &export_surface,
204 const State &state,
205 const Eigen::MatrixXd &sol,
206 const Eigen::MatrixXd &pressure,
207 const double t,
208 const double dt_in,
209 const ExportOptions &opts,
210 const bool is_contact_enabled) const;
211
218 void save_wire(const std::string &name,
219 const State &state,
220 const Eigen::MatrixXd &sol,
221 const double t,
222 const ExportOptions &opts) const;
223
229 void save_points(
230 const std::string &path,
231 const State &state,
232 const Eigen::MatrixXd &sol,
233 const ExportOptions &opts) const;
234
242 void save_pvd(const std::string &name, const std::function<std::string(int)> &vtu_names,
243 int time_steps, double t0, double dt, int skip_frame = 1) const;
244
245 private:
248
250 Eigen::MatrixXd grid_points;
252 Eigen::MatrixXi grid_points_to_elements;
254 Eigen::MatrixXd grid_points_bc;
255
271 const mesh::Mesh &mesh,
272 const std::vector<basis::ElementBases> &bases,
273 const std::vector<basis::ElementBases> &gbases,
274 const std::vector<mesh::LocalBoundary> &total_local_boundary,
275 const Eigen::MatrixXd &solution,
276 const int problem_dim,
277 Eigen::MatrixXd &boundary_vis_vertices,
278 Eigen::MatrixXd &boundary_vis_local_vertices,
279 Eigen::MatrixXi &boundary_vis_elements,
280 Eigen::MatrixXi &boundary_vis_elements_ids,
281 Eigen::MatrixXi &boundary_vis_primitive_ids,
282 Eigen::MatrixXd &boundary_vis_normals,
283 Eigen::MatrixXd &displaced_boundary_vis_normals) const;
284
299 void build_vis_mesh(
300 const mesh::Mesh &mesh,
301 const Eigen::VectorXi &disc_orders,
302 const std::vector<basis::ElementBases> &gbases,
303 const std::map<int, Eigen::MatrixXd> &polys,
304 const std::map<int, std::pair<Eigen::MatrixXd, Eigen::MatrixXi>> &polys_3d,
305 const bool boundary_only,
306 Eigen::MatrixXd &points,
307 Eigen::MatrixXi &tets,
308 Eigen::MatrixXi &el_id,
309 Eigen::MatrixXd &discr) const;
310
322 const mesh::Mesh &mesh,
323 const Eigen::VectorXi &disc_orders,
324 const std::vector<basis::ElementBases> &bases,
325 Eigen::MatrixXd &points,
326 std::vector<std::vector<int>> &elements,
327 Eigen::MatrixXi &el_id,
328 Eigen::MatrixXd &discr) const;
329
331 const State &state,
332 const Eigen::MatrixXd &points,
333 const ExportOptions &opts,
334 const std::string &name,
335 const Eigen::VectorXd &field,
336 paraviewo::ParaviewWriter &writer) const;
337 };
338
365
368 {
369 public:
371 Eigen::Vector4d spectrum;
372
376
378 double mesh_size;
383
386
390
392 double max_angle;
395
398
419
428 void compute_errors(const int n_bases,
429 const std::vector<polyfem::basis::ElementBases> &bases,
430 const std::vector<polyfem::basis::ElementBases> &gbases,
431 const polyfem::mesh::Mesh &mesh,
432 const assembler::Problem &problem,
433 const double tend,
434 const Eigen::MatrixXd &sol);
435
438 void compute_mesh_stats(const polyfem::mesh::Mesh &mesh);
439
447 void compute_mesh_size(const polyfem::mesh::Mesh &mesh_in, const std::vector<polyfem::basis::ElementBases> &bases_in, const int n_samples, const bool use_curved_mesh_size);
448
450 void reset();
451
455 void count_flipped_elements(const polyfem::mesh::Mesh &mesh, const std::vector<polyfem::basis::ElementBases> &gbases);
456
459
473 void save_json(const nlohmann::json &args,
474 const int n_bases, const int n_pressure_bases,
475 const Eigen::MatrixXd &sol,
476 const mesh::Mesh &mesh,
477 const Eigen::VectorXi &disc_orders,
478 const assembler::Problem &problem,
479 const OutRuntimeData &runtime,
480 const std::string &formulation,
481 const bool isoparametric,
482 const int sol_at_node_id,
483 nlohmann::json &j);
484 };
485
487 {
488 public:
489 EnergyCSVWriter(const std::string &path, const solver::SolveData &solve_data);
491
492 void write(const int i, const Eigen::MatrixXd &sol);
493
494 protected:
495 std::ofstream file;
497 };
498
500 {
501 public:
502 RuntimeStatsCSVWriter(const std::string &path, const State &state, const double t0, const double dt);
504
505 void write(const int t, const double forward, const double remeshing, const double global_relaxation, const Eigen::MatrixXd &sol);
506
507 protected:
508 std::ofstream file;
509 const State &state;
510 const double t0;
511 const double dt;
515 };
516} // namespace polyfem::io
main class that contains the polyfem solver and all its state
Definition State.hpp:79
const solver::SolveData & solve_data
Definition OutData.hpp:496
void write(const int i, const Eigen::MatrixXd &sol)
Definition OutData.cpp:2965
Utilies related to export of geometry.
Definition OutData.hpp:33
void build_vis_mesh(const mesh::Mesh &mesh, const Eigen::VectorXi &disc_orders, const std::vector< basis::ElementBases > &gbases, const std::map< int, Eigen::MatrixXd > &polys, const std::map< int, std::pair< Eigen::MatrixXd, Eigen::MatrixXi > > &polys_3d, const bool boundary_only, Eigen::MatrixXd &points, Eigen::MatrixXi &tets, Eigen::MatrixXi &el_id, Eigen::MatrixXd &discr) const
builds visualzation mesh, upsampled mesh used for visualization the visualization mesh is a dense mes...
Definition OutData.cpp:645
void build_high_order_vis_mesh(const mesh::Mesh &mesh, const Eigen::VectorXi &disc_orders, const std::vector< basis::ElementBases > &bases, Eigen::MatrixXd &points, std::vector< std::vector< int > > &elements, Eigen::MatrixXi &el_id, Eigen::MatrixXd &discr) const
builds high-der visualzation mesh per element all disconnected it also retuns the mapping to element ...
Definition OutData.cpp:788
void save_volume_vector_field(const State &state, const Eigen::MatrixXd &points, const ExportOptions &opts, const std::string &name, const Eigen::VectorXd &field, paraviewo::ParaviewWriter &writer) const
Definition OutData.cpp:1734
Eigen::MatrixXd grid_points_bc
grid mesh boundaries
Definition OutData.hpp:254
Eigen::MatrixXd grid_points
grid mesh points to export solution sampled on a grid
Definition OutData.hpp:250
void build_vis_boundary_mesh(const mesh::Mesh &mesh, const std::vector< basis::ElementBases > &bases, const std::vector< basis::ElementBases > &gbases, const std::vector< mesh::LocalBoundary > &total_local_boundary, const Eigen::MatrixXd &solution, const int problem_dim, Eigen::MatrixXd &boundary_vis_vertices, Eigen::MatrixXd &boundary_vis_local_vertices, Eigen::MatrixXi &boundary_vis_elements, Eigen::MatrixXi &boundary_vis_elements_ids, Eigen::MatrixXi &boundary_vis_primitive_ids, Eigen::MatrixXd &boundary_vis_normals, Eigen::MatrixXd &displaced_boundary_vis_normals) const
builds the boundary mesh for visualization
Definition OutData.cpp:397
void save_points(const std::string &path, const State &state, const Eigen::MatrixXd &sol, const ExportOptions &opts) const
saves the nodal values
Definition OutData.cpp:2286
void build_grid(const polyfem::mesh::Mesh &mesh, const double spacing)
builds the grid to export the solution
Definition OutData.cpp:2352
void save_contact_surface(const std::string &export_surface, const State &state, const Eigen::MatrixXd &sol, const Eigen::MatrixXd &pressure, const double t, const double dt_in, const ExportOptions &opts, const bool is_contact_enabled) const
saves the surface vtu file for for constact quantites, eg contact or friction forces
Definition OutData.cpp:1945
static void extract_boundary_mesh(const mesh::Mesh &mesh, const int n_bases, const std::vector< basis::ElementBases > &bases, const std::vector< mesh::LocalBoundary > &total_local_boundary, Eigen::MatrixXd &node_positions, Eigen::MatrixXi &boundary_edges, Eigen::MatrixXi &boundary_triangles, std::vector< Eigen::Triplet< double > > &displacement_map_entries)
extracts the boundary mesh
Definition OutData.cpp:145
void save_volume(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
saves the volume vtu file
Definition OutData.cpp:1202
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_wire(const std::string &name, const State &state, const Eigen::MatrixXd &sol, const double t, const ExportOptions &opts) const
saves the wireframe
Definition OutData.cpp:2096
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 init_sampler(const polyfem::mesh::Mesh &mesh, const double vismesh_rel_area)
unitalize the ref element sampler
Definition OutData.cpp:2347
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_surface(const std::string &export_surface, const State &state, const Eigen::MatrixXd &sol, const Eigen::MatrixXd &pressure, const double t, const double dt_in, const ExportOptions &opts, const bool is_contact_enabled) const
saves the surface vtu file for for surface quantites, eg traction forces
Definition OutData.cpp:1759
Eigen::MatrixXi grid_points_to_elements
grid mesh mapping to fe elements
Definition OutData.hpp:252
utils::RefElementSampler ref_element_sampler
used to sample the solution
Definition OutData.hpp:247
stores all runtime data
Definition OutData.hpp:341
double loading_mesh_time
time to load the mesh
Definition OutData.hpp:346
double total_time()
computes total time
Definition OutData.hpp:360
double assembling_stiffness_mat_time
time to assembly
Definition OutData.hpp:350
double assigning_rhs_time
time to computing the rhs
Definition OutData.hpp:354
double assembling_mass_mat_time
time to assembly mass
Definition OutData.hpp:352
double building_basis_time
time to construct the basis
Definition OutData.hpp:344
double solving_time
time to solve
Definition OutData.hpp:356
double computing_poly_basis_time
time to build the polygonal/polyhedral bases
Definition OutData.hpp:348
all stats from polyfem
Definition OutData.hpp:368
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
int n_flipped
number of flipped elements, compute only when using count_flipped_els (false by default)
Definition OutData.hpp:397
int simplex_count
statiscs on the mesh (simplices)
Definition OutData.hpp:400
json solver_info
information of the solver, eg num iteration, time, errors, etc the informations varies depending on t...
Definition OutData.hpp:375
int simple_singular_count
statiscs on the mesh (irregular quad/hex part of the mesh), see Polyspline paper for desciption
Definition OutData.hpp:406
Eigen::Vector4d spectrum
spectrum of the stiffness matrix, enable only if POLYSOLVE_WITH_SPECTRA is ON (off by default)
Definition OutData.hpp:371
void count_flipped_elements(const polyfem::mesh::Mesh &mesh, const std::vector< polyfem::basis::ElementBases > &gbases)
counts the number of flipped elements
Definition OutData.cpp:2537
double l2_err
errors, lp_err is in fact an L8 error
Definition OutData.hpp:385
int undefined_count
statiscs on the mesh (not quad/hex simplex), see Polyspline paper for desciption
Definition OutData.hpp:416
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
void compute_mesh_size(const polyfem::mesh::Mesh &mesh_in, const std::vector< polyfem::basis::ElementBases > &bases_in, const int n_samples, const bool use_curved_mesh_size)
computes the mesh size, it samples every edges n_samples times uses curved_mesh_size (false by defaul...
Definition OutData.cpp:2450
double average_edge_length
avg edge lenght
Definition OutData.hpp:382
double max_angle
statiscs on angle, compute only when using p_ref (false by default)
Definition OutData.hpp:392
int non_regular_boundary_count
statiscs on the mesh (irregular boundary quad/hex part of the mesh), see Polyspline paper for descipt...
Definition OutData.hpp:412
int non_regular_count
statiscs on the mesh (irregular quad/hex part of the mesh), see Polyspline paper for desciption
Definition OutData.hpp:414
int regular_boundary_count
statiscs on the mesh (regular quad/hex boundary part of the mesh), see Polyspline paper for desciptio...
Definition OutData.hpp:404
long long nn_zero
non zeros and sytem matrix size num dof is the total dof in the system
Definition OutData.hpp:389
double mesh_size
max edge lenght
Definition OutData.hpp:378
void reset()
clears all stats
Definition OutData.cpp:2528
int multi_singular_boundary_count
statiscs on the mesh (irregular boundary quad/hex part of the mesh), see Polyspline paper for descipt...
Definition OutData.hpp:418
void compute_mesh_stats(const polyfem::mesh::Mesh &mesh)
compute stats (counts els type, mesh lenght, etc), step 1 of solve
Definition OutData.cpp:2743
double sigma_max
statiscs on tri/tet quality, compute only when using p_ref (false by default)
Definition OutData.hpp:394
int boundary_count
statiscs on the mesh (boundary quads/hexs), see Polyspline paper for desciption
Definition OutData.hpp:410
int regular_count
statiscs on the mesh (regular quad/hex part of the mesh), see Polyspline paper for desciption
Definition OutData.hpp:402
double min_edge_length
min edge lenght
Definition OutData.hpp:380
int multi_singular_count
statiscs on the mesh (irregular quad/hex part of the mesh), see Polyspline paper for desciption
Definition OutData.hpp:408
void write(const int t, const double forward, const double remeshing, const double global_relaxation, const Eigen::MatrixXd &sol)
Definition OutData.cpp:2991
Abstract mesh class to capture 2d/3d conforming and non-conforming meshes.
Definition Mesh.hpp:39
class to store time stepping data
Definition SolveData.hpp:59
nlohmann::json json
Definition Common.hpp:9
std::string file_extension() const
return the extension of the output paraview files depending on use_hdf5
Definition OutData.hpp:78
std::vector< std::string > fields
Definition OutData.hpp:38
bool export_field(const std::string &field) const
Definition OutData.cpp:1080