PeriDEM 0.3.0
PeriDEM -- Peridynamics-based high-fidelity model for granular media
Loading...
Searching...
No Matches
Namespaces | Data Structures | Functions
rw Namespace Reference

Collection of methods and database related to reading and writing. More...

Namespaces

namespace  reader
 Collection of methods and database related to reading.
 
namespace  writer
 Collection of methods and database related to writing.
 

Data Structures

struct  PvtuPointArray
 

Functions

void writeOutput (data::ModelData &data)
 Write VTU/PVD (and optional strain VTU / particle-location CSV).
 
void writePvdCollectionFile (const std::string &pvd_path, const std::vector< std::pair< double, std::string > > &time_and_vtu_relative_path)
 Write a ParaView VTK collection (.pvd) that lists VTU files with timesteps.
 
void writePvtuCollectionFile (const std::string &pvtu_path, const std::vector< std::string > &piece_vtu_relative_paths, const std::vector< PvtuPointArray > &point_arrays={})
 Write a ParaView parallel VTU (.pvtu) that lists per-rank .vtu pieces.
 

Detailed Description

Collection of methods and database related to reading and writing.

This namespace provides methods and data members specific to reading and writing of the mesh data and simulation data.

Function Documentation

◆ writeOutput()

void rw::writeOutput ( data::ModelData &  data)

Write VTU/PVD (and optional strain VTU / particle-location CSV).

Definition at line 32 of file particleOutput.cpp.

32 {
33
34
35 // write out % completion of simulation at 10% interval
36 {
37 float p = float(data.d_n) * 100. / data.d_modelDeck_p->d_Nt;
38 int m = std::max(1, int(data.d_modelDeck_p->d_Nt / 10));
39 if (data.d_n % m == 0 && int(p) > 0)
40 util::io::log(0, std::format("{}: Simulation {}% complete\n",
41 data.d_name, int(p)));
42 ;
43 }
44
45 util::io::log(2, std::format("{}: Output step = {}, time = {:.6f} \n",
46 data.d_name, data.d_n, data.d_time));
47
48 if (data.d_outputDeck_p->d_debug > 0 and data.getKeyData("debug_once") < 0) {
49
50 data.setKeyData("debug_once", 1);
51
52 size_t nt = 1;
53 auto tabS = util::io::getTabS(nt);
54 std::ostringstream oss;
55 oss << tabS << "*******************************************\n";
56 oss << tabS << "Debug various input decks\n\n\n";
57 oss << data.d_modelDeck_p->printStr(nt + 1);
58 oss << data.d_particleDeck_p->printStr(nt + 1);
59 oss << data.d_outputDeck_p->printStr(nt + 1);
60 oss << data.d_restartDeck_p->printStr(nt + 1);
61 oss << data.d_testDeck_p->printStr(nt + 1);
62 oss << data.d_bcDeck_p->printStr(nt + 1);
63 oss << tabS << "\n\n*******************************************\n";
64 oss << tabS << "Debug particle data\n\n\n";
65 oss << tabS << "Number of particles = " << data.d_particlesListTypeAll.size() << std::endl;
66 oss << tabS << "Number of particle zones = " << data.d_zInfo.size() << std::endl;
67 for (auto zone : data.d_zInfo) {
68 oss << tabS << "zone of d_zInfo: " << util::io::printStr(zone)
69 << std::endl;
70 }
71
72 // wall info
73 oss << tabS << "Number of walls = " << data.d_particlesListTypeWall.size() << std::endl;
74 for (auto &d_wall : data.d_particlesListTypeWall)
75 oss << tabS << "Number of nodes in wall " << d_wall->d_id
76 << " is " << d_wall->getNumNodes() << std::endl;
77
78 oss << tabS << "h_min = " << data.d_hMin << ", h_max = " << data.d_hMax << std::endl;
79
80 util::io::log(2, oss);
81 } // end of debug
82
83 size_t dt_out = data.d_outputDeck_p->d_dtOutCriteria;
84 const size_t frame = data.d_n / dt_out;
85 const std::string tag =
86 data.d_outputDeck_p->d_tagPPFile.empty()
87 ? std::to_string(frame)
88 : data.d_outputDeck_p->d_tagPPFile + "_" + std::to_string(frame);
89 const std::string path = data.d_outputDeck_p->d_path;
90 const int mpi_size = util::parallel::mpiSize();
91 const int mpi_rank = util::parallel::mpiRank();
92 const bool parallel_pieces = mpi_size > 1;
93
94 // Per-rank piece (or single file). No solution gather for I/O.
95 std::string piece_stem = path + "output_" + tag;
96 if (parallel_pieces)
97 piece_stem += "_r" + std::to_string(mpi_rank);
98
99 {
100 auto writer = rw::writer::VtkParticleWriter(piece_stem);
101 writer.appendMeshParallelPiece(&data, data.d_outputDeck_p->d_outTags);
102 writer.addTimeStep(data.d_time);
103 writer.close();
104 }
105
106 if (data.d_outputDeck_p->d_outFormat == "vtu" &&
107 data.d_outputDeck_p->d_pvdCollection) {
108 if (parallel_pieces) {
109 // Ensure all piece files are on disk before rank 0 writes the .pvtu.
110 MPI_Barrier(util::parallel::mpiComm());
111 if (mpi_rank == 0) {
112 std::vector<std::string> pieces;
113 pieces.reserve(static_cast<size_t>(mpi_size));
114 for (int r = 0; r < mpi_size; ++r)
115 pieces.push_back("output_" + tag + "_r" + std::to_string(r) + ".vtu");
116 // Must match appendPointArraysForNodes: ParaView only lists arrays
117 // declared in PPointData of the .pvtu (not discovered from pieces).
118 const auto &tags = data.d_outputDeck_p->d_outTags;
119 std::vector<rw::PvtuPointArray> point_arrays;
120 auto add_arr = [&](const char *name, int ncomp) {
121 point_arrays.push_back({name, "Float64", ncomp});
122 };
123 if (util::methods::isTagInList("Displacement", tags))
124 add_arr("Displacement", 3);
125 if (util::methods::isTagInList("Velocity", tags))
126 add_arr("Velocity", 3);
127 if (util::methods::isTagInList("Force_Density", tags))
128 add_arr("Force_Density", 3);
129 if (util::methods::isTagInList("Force", tags))
130 add_arr("Force", 3);
131 if (util::methods::isTagInList("Damage_Z", tags) && !data.d_Z.empty())
132 add_arr("Damage_Z", 1);
133 if (util::methods::isTagInList("Damage", tags) && !data.d_phi.empty())
134 add_arr("Damage", 1);
135 if (util::methods::isTagInList("Damage_Bond", tags) && !data.d_phiBond.empty())
136 add_arr("Damage_Bond", 1);
137 if (util::methods::isTagInList("Particle_ID", tags))
138 add_arr("Particle_ID", 1);
139 const std::string pvtu_name = "output_" + tag + ".pvtu";
140 rw::writePvtuCollectionFile(path + pvtu_name, pieces, point_arrays);
141 data.d_pvdParticleEntries.push_back({data.d_time, pvtu_name});
142 rw::writePvdCollectionFile(path + "output.pvd",
143 data.d_pvdParticleEntries);
144 }
145 } else {
146 data.d_pvdParticleEntries.push_back(
147 {data.d_time, "output_" + tag + ".vtu"});
148 rw::writePvdCollectionFile(path + "output.pvd",
149 data.d_pvdParticleEntries);
150 }
151 }
152
153 if (util::methods::isTagInList("Strain_Stress", data.d_outputDeck_p->d_outTags)) {
154
155 // compute current position of quadrature points and strain/stress data
156 {
157 // if particle mat data is not computed, compute them
158 if (data.d_particlesMatDataList.empty()) {
159 for (auto &p: data.d_particlesListTypeAll) {
160 data.d_particlesMatDataList.push_back(p->getMaterial()->computeMaterialProperties(
161 p->getMeshP()->getDimension()));
162 }
163 }
164
165 for (auto &p: data.d_particlesListTypeAll) {
166
167 const auto particle_mesh_p = p->getMeshP();
168
169 mesh::getCurrentQuadPoints(particle_mesh_p.get(), data.d_xRef, data.d_u, data.d_xQuadCur,
170 p->d_globStart,
171 p->d_globQuadStart,
172 data.d_modelDeck_p->d_quadOrder);
173
174 auto isPlaneStrain = p->d_material_p->isPlaneStrain();
175 mesh::getStrainStress(particle_mesh_p.get(), data.d_xRef, data.d_u,
176 isPlaneStrain,
177 data.d_strain, data.d_stress,
178 p->d_globStart,
179 p->d_globQuadStart,
180 data.d_particlesMatDataList[p->getId()].d_nu,
181 data.d_particlesMatDataList[p->getId()].d_lambda,
182 data.d_particlesMatDataList[p->getId()].d_mu,
183 true,
184 data.d_modelDeck_p->d_quadOrder);
185 } // for loop over particles
186 } // compute strain/stress block
187
188 // Strain field is global; only rank 0 writes (legacy path). Prefer primary
189 // particle VTU/PVTU for parallel visualization.
190 if (mpi_rank == 0) {
191 std::string out_filename = path + "output_strain_" + tag;
192 auto writer1 = rw::writer::VtkParticleWriter(out_filename);
193 writer1.appendStrainStress(&data);
194 writer1.addTimeStep(data.d_time);
195 writer1.close();
196
197 if (data.d_outputDeck_p->d_outFormat == "vtu" &&
198 data.d_outputDeck_p->d_pvdCollection) {
199 data.d_pvdStrainEntries.push_back(
200 {data.d_time, "output_strain_" + tag + ".vtu"});
201 rw::writePvdCollectionFile(path + "output_strain.pvd",
202 data.d_pvdStrainEntries);
203 }
204 }
205 }
206
207 // output particle locations to csv file (rank 0)
208 if (util::methods::isTagInList("Particle_Locations",
209 data.d_outputDeck_p->d_outTags) &&
210 mpi_rank == 0) {
211
212 std::string out_filename = path + "particle_locations_" + tag + ".csv";
213 std::ofstream oss(out_filename);
214 oss << "i, x, y, z, r\n";
215 for (const auto &p : data.d_particlesListTypeAll) {
216 auto xc = p->getXCenter();
217 oss << p->d_id << ", " << xc.d_x << ", " << xc.d_y << ", " << xc.d_z
218 << ", " << p->d_geom_p->boundingRadius() << "\n";
219 }
220 oss.close();
221 }
222
223}
rw::writer::VtkParticleWriter
A vtk writer for simple point data and complex fem mesh data.
Definition vtkParticleWriter.h:28
data
Definition contact.h:20
mesh::getCurrentQuadPoints
void getCurrentQuadPoints(const mesh::Mesh *mesh_p, const std::vector< util::Point > &xRef, const std::vector< util::Point > &u, std::vector< util::Point > &xQuadCur, size_t iNodeStart, size_t iQuadStart, size_t quadOrder)
Get current location of quadrature points of elements in the mesh. This function expects mesh has ele...
Definition meshUtil.cpp:305
mesh::getStrainStress
void getStrainStress(const mesh::Mesh *mesh_p, const std::vector< util::Point > &xRef, const std::vector< util::Point > &u, bool isPlaneStrain, std::vector< util::SymMatrix3 > &strain, std::vector< util::SymMatrix3 > &stress, size_t iNodeStart, size_t iStrainStart, double nu, double lambda, double mu, bool computeStress, size_t quadOrder)
Strain and stress at quadrature points in the mesh.
Definition meshUtil.cpp:373
rw::writePvdCollectionFile
void writePvdCollectionFile(const std::string &pvd_path, const std::vector< std::pair< double, std::string > > &time_and_vtu_relative_path)
Write a ParaView VTK collection (.pvd) that lists VTU files with timesteps.
Definition pvdCollectionWriter.cpp:18
rw::writePvtuCollectionFile
void writePvtuCollectionFile(const std::string &pvtu_path, const std::vector< std::string > &piece_vtu_relative_paths, const std::vector< PvtuPointArray > &point_arrays)
Write a ParaView parallel VTU (.pvtu) that lists per-rank .vtu pieces.
Definition pvdCollectionWriter.cpp:37
util::io::getTabS
std::string getTabS(int nt)
Returns tab spaces of a given size.
Definition io.h:82
util::io::printStr
std::string printStr(const T &msg, int nt=print_default_tab)
Returns formatted string for output.
Definition io.h:96
util::io::log
void log(std::ostringstream &oss, bool screen_out=false, int printMpiRank=print_default_mpi_rank)
Global method to log the message.
Definition io.cpp:41
util::methods::isTagInList
bool isTagInList(const std::string &tag, const std::vector< std::string > &tags)
Returns true if tag is found in the list of tags.
Definition vecMethods.h:284
util::parallel::mpiSize
int mpiSize()
Get size (number) of processors.
Definition parallelUtil.cpp:88
util::parallel::mpiRank
int mpiRank()
get rank (id) of this processor
Definition parallelUtil.cpp:96
util::parallel::mpiComm
MPI_Comm mpiComm()
Get MPI comm.
Definition parallelUtil.cpp:104

References mesh::getCurrentQuadPoints(), mesh::getStrainStress(), util::io::getTabS(), util::methods::isTagInList(), util::io::log(), util::parallel::mpiComm(), util::parallel::mpiRank(), util::parallel::mpiSize(), util::io::printStr(), writePvdCollectionFile(), and writePvtuCollectionFile().

Referenced by PeriDEMModel::output().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ writePvdCollectionFile()

void rw::writePvdCollectionFile ( const std::string &  pvd_path,
const std::vector< std::pair< double, std::string > > &  time_and_vtu_relative_path 
)

Write a ParaView VTK collection (.pvd) that lists VTU files with timesteps.

Open the .pvd in ParaView to animate all listed .vtu snapshots. Each entry's second string must be the VTU filename only (or a path relative to the .pvd file).

Definition at line 18 of file pvdCollectionWriter.cpp.

20 {
21
22 std::ofstream os(pvd_path);
23 if (!os)
24 throw std::runtime_error("writePvdCollectionFile: could not open " + pvd_path);
25
26 os << std::setprecision(17);
27 os << "<?xml version=\"1.0\"?>\n";
28 os << "<VTKFile type=\"Collection\" version=\"0.1\" byte_order=\"LittleEndian\">\n";
29 os << " <Collection>\n";
30 for (const auto &e : time_and_vtu_relative_path) {
31 os << " <DataSet timestep=\"" << e.first << "\" file=\"" << e.second << "\"/>\n";
32 }
33 os << " </Collection>\n";
34 os << "</VTKFile>\n";
35}

Referenced by writeOutput().

Here is the caller graph for this function:

◆ writePvtuCollectionFile()

void rw::writePvtuCollectionFile ( const std::string &  pvtu_path,
const std::vector< std::string > &  piece_vtu_relative_paths,
const std::vector< PvtuPointArray > &  point_arrays = {} 
)

Write a ParaView parallel VTU (.pvtu) that lists per-rank .vtu pieces.

Each string in piece_vtu_relative_paths is a path relative to the .pvtu (typically output_N_r0.vtu, …). point_arrays must list every PointData array present in the pieces; ParaView will not expose fields omitted here.

Definition at line 37 of file pvdCollectionWriter.cpp.

40 {
41
42 std::ofstream os(pvtu_path);
43 if (!os)
44 throw std::runtime_error("writePvtuCollectionFile: could not open " +
45 pvtu_path);
46
47 os << "<?xml version=\"1.0\"?>\n";
48 os << "<VTKFile type=\"PUnstructuredGrid\" version=\"0.1\" "
49 "byte_order=\"LittleEndian\">\n";
50 os << " <PUnstructuredGrid GhostLevel=\"0\">\n";
51 os << " <PPoints>\n";
52 os << " <PDataArray type=\"Float32\" NumberOfComponents=\"3\"/>\n";
53 os << " </PPoints>\n";
54 if (!point_arrays.empty()) {
55 os << " <PPointData>\n";
56 for (const auto &a : point_arrays) {
57 os << " <PDataArray type=\"" << a.type << "\" Name=\"" << a.name
58 << "\"";
59 if (a.number_of_components > 1)
60 os << " NumberOfComponents=\"" << a.number_of_components << "\"";
61 os << "/>\n";
62 }
63 os << " </PPointData>\n";
64 }
65 for (const auto &piece : piece_vtu_relative_paths) {
66 os << " <Piece Source=\"" << piece << "\"/>\n";
67 }
68 os << " </PUnstructuredGrid>\n";
69 os << "</VTKFile>\n";
70}

Referenced by writeOutput().

Here is the caller graph for this function: