preCICE v3.1.2
Loading...
Searching...
No Matches
ExportXML.cpp
Go to the documentation of this file.
1#include "io/ExportXML.hpp"
2#include <Eigen/Core>
3#include <algorithm>
4#include <filesystem>
5#include <fstream>
6#include <memory>
7#include <string>
8#include "io/Export.hpp"
10#include "mesh/Data.hpp"
11#include "mesh/Edge.hpp"
12#include "mesh/Mesh.hpp"
14#include "mesh/Tetrahedron.hpp"
15#include "mesh/Triangle.hpp"
16#include "mesh/Vertex.hpp"
17#include "utils/Helpers.hpp"
18#include "utils/IntraComm.hpp"
19#include "utils/assertion.hpp"
20
21namespace precice::io {
22
24 const std::string &name,
25 const std::string &location,
26 const mesh::Mesh & mesh)
27{
28 PRECICE_TRACE(name, location, mesh.getName());
30 if (not location.empty())
33 writeParallelFile(name, location, mesh);
34 }
35 if (mesh.nVertices() > 0) { // only procs at the coupling interface should write output (for performance reasons)
36 writeSubFile(name, location, mesh);
37 }
38}
39
41{
44 const bool isThreeDim = (mesh.getDimensions() == 3);
45 for (const mesh::PtrData &data : mesh.data()) {
46 int dataDimensions = data->getDimensions();
47 const bool hasGradient = data->hasGradient();
48 PRECICE_ASSERT(dataDimensions >= 1);
49 std::string dataName = data->getName();
50 if (dataDimensions == 1) {
52 if (hasGradient) {
53 _vectorDataNames.push_back(dataName + "_gradient");
54 }
55 } else {
57 if (hasGradient) {
58 _vectorDataNames.push_back(dataName + "_dx");
59 _vectorDataNames.push_back(dataName + "_dy");
60 if (isThreeDim) {
61 _vectorDataNames.push_back(dataName + "_dz");
62 }
63 }
64 }
65 }
66}
67
69 const std::string &name,
70 const std::string &location,
71 const mesh::Mesh & mesh) const
72{
73 namespace fs = std::filesystem;
74 fs::path outfile(location);
75 outfile = outfile / fs::path(name + getParallelExtension());
76 std::ofstream outParallelFile(outfile.string(), std::ios::trunc);
77
78 PRECICE_CHECK(outParallelFile, "{} export failed to open primary file \"{}\"", getVTKFormat(), outfile.generic_string());
79
80 const auto formatType = getVTKFormat();
81 outParallelFile << "<?xml version=\"1.0\"?>\n";
82 outParallelFile << "<VTKFile type=\"P" << formatType << "\" version=\"0.1\" byte_order=\"";
83 outParallelFile << (utils::isMachineBigEndian() ? "BigEndian\">" : "LittleEndian\">") << '\n';
84 outParallelFile << " <P" << formatType << " GhostLevel=\"0\">\n";
85
86 outParallelFile << " <PPoints>\n";
87 outParallelFile << " <PDataArray type=\"Float64\" Name=\"Position\" NumberOfComponents=\"" << 3 << "\"/>\n";
88 outParallelFile << " </PPoints>\n";
89
90 writeParallelCells(outParallelFile);
91
92 writeParallelData(outParallelFile);
93
94 const auto &offsets = mesh.getVertexOffsets();
95 PRECICE_ASSERT(offsets.size() > 0);
96 if (offsets[0] > 0) {
97 outParallelFile << " <Piece Source=\"" << name << "_" << 0 << getPieceExtension() << "\"/>\n";
98 }
99 for (size_t rank : utils::IntraComm::allSecondaryRanks()) {
100 PRECICE_ASSERT(rank < offsets.size());
101 if (offsets[rank] - offsets[rank - 1] > 0) {
102 // only non-empty subfiles
103 outParallelFile << " <Piece Source=\"" << name << "_" << rank << getPieceExtension() << "\"/>\n";
104 }
105 }
106
107 outParallelFile << " </P" << formatType << ">\n";
108 outParallelFile << "</VTKFile>\n";
109
110 outParallelFile.close();
111}
112
113namespace {
114std::string getPieceSuffix()
115{
116 if (!utils::IntraComm::isParallel()) {
117 return "";
118 }
119 return "_" + std::to_string(utils::IntraComm::getRank());
120}
121} // namespace
122
124 const std::string &name,
125 const std::string &location,
126 const mesh::Mesh & mesh) const
127{
128 namespace fs = std::filesystem;
129 fs::path outfile(location);
130 outfile /= fs::path(name + getPieceSuffix() + getPieceExtension());
131 std::ofstream outSubFile(outfile.string(), std::ios::trunc);
132
133 PRECICE_CHECK(outSubFile, "{} export failed to open secondary file \"{}\"", getVTKFormat(), outfile.generic_string());
134
135 const auto formatType = getVTKFormat();
136 outSubFile << "<?xml version=\"1.0\"?>\n";
137 outSubFile << "<VTKFile type=\"" << formatType << "\" version=\"0.1\" byte_order=\"";
138 outSubFile << (utils::isMachineBigEndian() ? "BigEndian\">" : "LittleEndian\">") << '\n';
139
140 outSubFile << " <" << formatType << ">\n";
141 outSubFile << " <Piece " << getPieceAttributes(mesh) << "> \n";
142 exportPoints(outSubFile, mesh);
143
144 // Write Mesh
145 exportConnectivity(outSubFile, mesh);
146
147 // Write data
148 exportData(outSubFile, mesh);
149
150 outSubFile << " </Piece>\n";
151 outSubFile << " </" << formatType << "> \n";
152 outSubFile << "</VTKFile>\n";
153
154 outSubFile.close();
155}
156
157void ExportXML::exportGradient(const mesh::PtrData data, const int spaceDim, std::ostream &outFile) const
158{
159 const auto & gradients = data->gradients();
160 const int dataDimensions = data->getDimensions();
162 if (dataDimensions == 1) {
163 suffices = {"_gradient"};
164 } else if (spaceDim == 2) {
165 suffices = {"_dx", "_dy"};
166 } else if (spaceDim == 3) {
167 suffices = {"_dx", "_dy", "_dz"};
168 }
169 int counter = 0; // Counter for multicomponent
170 for (const auto &suffix : suffices) {
171 const std::string dataName(data->getName());
172 outFile << " <DataArray type=\"Float64\" Name=\"" << dataName << suffix << "\" NumberOfComponents=\"" << 3;
173 outFile << "\" format=\"ascii\">\n";
174 outFile << " ";
175 for (int i = counter; i < gradients.cols(); i += spaceDim) { // Loop over vertices
176 int j = 0;
177 for (; j < gradients.rows(); j++) { // Loop over components
178 outFile << gradients.coeff(j, i) << " ";
179 }
180 if (j < 3) { // If 2D data add additional zero as third component
181 outFile << "0.0"
182 << " ";
183 }
184 }
185 outFile << '\n'
186 << " </DataArray>\n";
187 counter++; // Increment counter for next component
188 }
189}
190
192 std::ostream & outFile,
193 const mesh::Mesh &mesh) const
194{
195 outFile << " <PointData Scalars=\"Rank ";
196 for (const auto &scalarDataName : _scalarDataNames) {
197 outFile << scalarDataName << ' ';
198 }
199 outFile << "\" Vectors=\"";
200 for (const auto &vectorDataName : _vectorDataNames) {
201 outFile << vectorDataName << ' ';
202 }
203 outFile << "\">\n";
204
205 // Export the current rank
206 outFile << " <DataArray type=\"UInt32\" Name=\"Rank\" NumberOfComponents=\"1\" format=\"ascii\">\n";
207 outFile << " ";
208 const auto rank = utils::IntraComm::getRank();
209 for (size_t count = 0; count < mesh.nVertices(); ++count) {
210 outFile << rank << ' ';
211 }
212 outFile << "\n </DataArray>\n";
213
214 for (const mesh::PtrData &data : mesh.data()) { // Plot vertex data
215 Eigen::VectorXd &values = data->values();
216 int dataDimensions = data->getDimensions();
217 std::string dataName(data->getName());
218 int numberOfComponents = (dataDimensions == 2) ? 3 : dataDimensions;
219 const bool hasGradient = data->hasGradient();
220 outFile << " <DataArray type=\"Float64\" Name=\"" << dataName << "\" NumberOfComponents=\"" << numberOfComponents;
221 outFile << "\" format=\"ascii\">\n";
222 outFile << " ";
223 if (dataDimensions > 1) {
224 Eigen::VectorXd viewTemp(dataDimensions);
225 for (size_t count = 0; count < mesh.nVertices(); count++) {
226 size_t offset = count * dataDimensions;
227 for (int i = 0; i < dataDimensions; i++) {
228 viewTemp[i] = values(offset + i);
229 }
230 for (int i = 0; i < dataDimensions; i++) {
231 outFile << viewTemp[i] << ' ';
232 }
233 if (dataDimensions == 2) {
234 outFile << "0.0" << ' '; // 2D data needs to be 3D for vtk
235 }
236 outFile << ' ';
237 }
238 } else if (dataDimensions == 1) {
239 for (size_t count = 0; count < mesh.nVertices(); count++) {
240 outFile << values(count) << ' ';
241 }
242 }
243 outFile << '\n'
244 << " </DataArray>\n";
245 if (hasGradient) {
246 exportGradient(data, dataDimensions, outFile);
247 }
248 }
249 outFile << " </PointData> \n";
250}
251
253 const Eigen::VectorXd &position,
254 std::ostream & outFile)
255{
256 outFile << " ";
257 for (int i = 0; i < position.size(); i++) {
258 outFile << position(i) << " ";
259 }
260 if (position.size() == 2) {
261 outFile << 0.0 << " "; // also for 2D scenario, vtk needs 3D data
262 }
263 outFile << '\n';
264}
265
267 const mesh::Triangle &triangle,
268 std::ostream & outFile)
269{
270 outFile << triangle.vertex(0).getID() << " ";
271 outFile << triangle.vertex(1).getID() << " ";
272 outFile << triangle.vertex(2).getID() << " ";
273}
274
276 const mesh::Tetrahedron &tetra,
277 std::ostream & outFile)
278{
279 outFile << tetra.vertex(0).getID() << " ";
280 outFile << tetra.vertex(1).getID() << " ";
281 outFile << tetra.vertex(2).getID() << " ";
282 outFile << tetra.vertex(3).getID() << " ";
283}
284
286 const mesh::Edge &edge,
287 std::ostream & outFile)
288{
289 outFile << edge.vertex(0).getID() << " ";
290 outFile << edge.vertex(1).getID() << " ";
291}
292
294 std::ostream & outFile,
295 const mesh::Mesh &mesh) const
296{
297 outFile << " <Points> \n";
298 outFile << " <DataArray type=\"Float64\" Name=\"Position\" NumberOfComponents=\"" << 3 << "\" format=\"ascii\"> \n";
299 for (const mesh::Vertex &vertex : mesh.vertices()) {
300 writeVertex(vertex.getCoords(), outFile);
301 }
302 outFile << " </DataArray>\n";
303 outFile << " </Points> \n\n";
304}
305
307{
308 // write scalar data names
309 out << " <PPointData Scalars=\"Rank ";
310 for (const auto &scalarDataName : _scalarDataNames) {
311 out << scalarDataName << ' ';
312 }
313 // write vector data names
314 out << "\" Vectors=\"";
315 for (const auto &vectorDataName : _vectorDataNames) {
316 out << vectorDataName << ' ';
317 }
318 out << "\">\n";
319
320 out << " <PDataArray type=\"Int32\" Name=\"Rank\" NumberOfComponents=\"1\"/>\n";
321
322 for (const auto &scalarDataName : _scalarDataNames) {
323 out << " <PDataArray type=\"Float64\" Name=\"" << scalarDataName << "\" NumberOfComponents=\"" << 1 << "\"/>\n";
324 }
325
326 for (const auto &vectorDataName : _vectorDataNames) {
327 out << " <PDataArray type=\"Float64\" Name=\"" << vectorDataName << "\" NumberOfComponents=\"" << 3 << "\"/>\n";
328 }
329 out << " </PPointData>\n";
330}
331
332} // namespace precice::io
std::ostream & out
#define PRECICE_TRACE(...)
Definition LogMacros.hpp:95
#define PRECICE_CHECK(check,...)
Definition LogMacros.hpp:35
std::string name
#define PRECICE_ASSERT(...)
Definition assertion.hpp:87
virtual std::string getPieceExtension() const =0
virtual std::string getPieceAttributes(const mesh::Mesh &mesh) const =0
void writeParallelData(std::ostream &out) const
static void writeLine(const mesh::Edge &edge, std::ostream &outFile)
virtual std::string getVTKFormat() const =0
void exportData(std::ostream &outFile, const mesh::Mesh &mesh) const
virtual void exportConnectivity(std::ostream &outFile, const mesh::Mesh &mesh) const =0
void exportPoints(std::ostream &outFile, const mesh::Mesh &mesh) const
void writeSubFile(const std::string &name, const std::string &location, const mesh::Mesh &mesh) const
Writes the sub file for each rank.
void doExport(const std::string &name, const std::string &location, const mesh::Mesh &mesh) override
Does export. Has to be implemented in subclass.
Definition ExportXML.cpp:23
static void writeTriangle(const mesh::Triangle &triangle, std::ostream &outFile)
std::vector< std::string > _scalarDataNames
List of names of all scalar data on mesh.
Definition ExportXML.hpp:51
virtual std::string getParallelExtension() const =0
static void writeTetrahedron(const mesh::Tetrahedron &tetra, std::ostream &outFile)
std::vector< std::string > _vectorDataNames
List of names of all vector data on mesh.
Definition ExportXML.hpp:54
void exportGradient(const mesh::PtrData data, const int dataDim, std::ostream &outFile) const
void writeParallelFile(const std::string &name, const std::string &location, const mesh::Mesh &mesh) const
Writes the primary file (called only by the primary rank)
Definition ExportXML.cpp:68
void processDataNamesAndDimensions(const mesh::Mesh &mesh)
Stores scalar and vector data names in string vectors Needed for writing primary file and sub files.
Definition ExportXML.cpp:40
static void writeVertex(const Eigen::VectorXd &position, std::ostream &outFile)
virtual void writeParallelCells(std::ostream &out) const =0
Linear edge of a mesh, defined by two Vertex objects.
Definition Edge.hpp:16
Vertex & vertex(int i)
Returns the edge's vertex with index 0 or 1.
Definition Edge.hpp:77
Container and creator for meshes.
Definition Mesh.hpp:39
int getDimensions() const
Definition Mesh.cpp:98
VertexContainer & vertices()
Returns modifieable container holding all vertices.
Definition Mesh.cpp:53
const std::string & getName() const
Returns the name of the mesh, as set in the config file.
Definition Mesh.cpp:218
std::size_t nVertices() const
Returns the number of vertices.
Definition Mesh.cpp:63
const VertexOffsets & getVertexOffsets() const
Definition Mesh.hpp:248
const DataContainer & data() const
Allows access to all data.
Definition Mesh.cpp:170
Tetrahedron of a mesh, defined by 4 vertices.
Vertex & vertex(int i)
Returns tetrahedron vertex with index 0, 1, 2 or 3.
Triangle of a mesh, defined by three vertices.
Definition Triangle.hpp:27
Vertex & vertex(int i)
Returns triangle vertex with index 0, 1 or 2.
Definition Triangle.hpp:143
Vertex of a mesh.
Definition Vertex.hpp:16
VertexID getID() const
Returns the unique (among vertices of one mesh on one processor) ID of the vertex.
Definition Vertex.hpp:111
static Rank getRank()
Current rank.
Definition IntraComm.cpp:42
static bool isPrimary()
True if this process is running the primary rank.
Definition IntraComm.cpp:52
static auto allSecondaryRanks()
Returns an iterable range over salve ranks [1, _size)
Definition IntraComm.hpp:37
T clear(T... args)
T close(T... args)
T create_directories(T... args)
T empty(T... args)
T generic_string(T... args)
provides Import and Export of the coupling mesh and data.
bool isMachineBigEndian()
Returns true if machine is big-endian needed for parallel vtk output.
Definition Helpers.cpp:5
T push_back(T... args)
T to_string(T... args)