OpenLB 1.8.1
Loading...
Searching...
No Matches
olb::VTUsurfaceWriter< T > Class Template Reference

#include <vtuSurfaceWriter.h>

+ Collaboration diagram for olb::VTUsurfaceWriter< T >:

Public Member Functions

 VTUsurfaceWriter (std::string name, CuboidDecomposition3D< T > &cuboidDecomposition, LoadBalancer< T > &loadBalancer, bool binary=false)
 
void createMasterFile ()
 
void write (std::size_t iT)
 writes vtu files and assigns values
 
void addPoint (Vector< T, 3 > &new_position)
 Add single location to output set.
 
void addPoints (std::vector< Vector< T, 3 > > &new_positions)
 Add multiple locations to output set.
 
void addConnectivity (std::vector< std::vector< int > > &connections)
 Define connectivity between previously-added output points, rendering them a surface.
 
void addSTL (STLreader< T > &stlReader)
 Add full STL surface as output surface.
 
void addFunctor (AnalyticalF3D< T, T > &f)
 Schedule analytical functor for output.
 
void addFunctor (SuperF3D< T, T > &f)
 Schedule lattice functor for output (interpolates)
 

Protected Member Functions

void preamblePVD (const std::string &fullNamePVD)
 
void closePVD (const std::string &fullNamePVD)
 
void preambleVTU (const std::string &fullName, int num)
 performes <VTKFile ...>, <ImageData ...> and <PieceExtent ...>
 
void closeVTU (const std::string &fullNamePiece)
 performes </ImageData> and </VTKFile>
 
void dataPVD (int iT, int i, const std::string &fullNamePVD, const std::string &namePiece)
 performes <DataSet timestep= ... file=namePiece >
 
void dataPVDmaster (int iT, int i, const std::string &fullNamePVDMaster, const std::string &namePiece)
 performes <DataSet timestep= ... file=namePiece >
 
void dataArray (const std::string &fullName)
 
void writeFunctor (const std::string &fullName, std::ofstream &fout, SuperF3D< T, T > &f)
 interpolates and writes functors stored at functors
 
void writeAnalyticalFunctor (const std::string &fullName, std::ofstream &fout, AnalyticalF3D< T, T > &f)
 writes functors.
 
void writePosition (std::ofstream &fout)
 writes coordinates of points
 
std::vector< std::vector< T > > removeDuplicates (const std::vector< std::vector< T > > &points)
 removes duplicates from points, needed for STL files
 
std::vector< std::vector< int > > mapTrianglesToIndices (STLreader< T > &stlReader)
 

Static Protected Member Functions

static bool areVectorsEqual (const std::vector< T > &v1, const std::vector< T > &v2, T epsilon=1e-6)
 

Protected Attributes

std::string _name
 
std::vector< SuperF3D< T, T > * > _functors
 
std::vector< AnalyticalF3D< T, T > * > _functorsA
 
std::vector< Vector< T, 3 > > pos
 
std::vector< std::vector< int > > connectivity
 
CuboidDecomposition3D< T > & _cuboidDecomposition
 
LoadBalancer< T > & _loadBalancer
 
bool _binary
 
OstreamManager clout
 

Detailed Description

template<typename T>
class olb::VTUsurfaceWriter< T >

Definition at line 36 of file vtuSurfaceWriter.h.

Constructor & Destructor Documentation

◆ VTUsurfaceWriter()

template<typename T >
olb::VTUsurfaceWriter< T >::VTUsurfaceWriter ( std::string name,
CuboidDecomposition3D< T > & cuboidDecomposition,
LoadBalancer< T > & loadBalancer,
bool binary = false )

Definition at line 36 of file vtuSurfaceWriter.hh.

38 : _name(name)
39 , _cuboidDecomposition(cuboidDecomposition)
40 , _loadBalancer(loadBalancer)
41 , _binary(binary)
42 , clout(std::cout, "VTUsurfaceWriter3D")
43{}
CuboidDecomposition3D< T > & _cuboidDecomposition
LoadBalancer< T > & _loadBalancer

Member Function Documentation

◆ addConnectivity()

template<typename T >
void olb::VTUsurfaceWriter< T >::addConnectivity ( std::vector< std::vector< int > > & connections)

Define connectivity between previously-added output points, rendering them a surface.

Definition at line 72 of file vtuSurfaceWriter.hh.

73{
74 for (unsigned long i = 0; i < connections.size(); i++) {
75 connectivity.push_back(connections[i]);
76 }
77}
std::vector< std::vector< int > > connectivity

◆ addFunctor() [1/2]

template<typename T >
void olb::VTUsurfaceWriter< T >::addFunctor ( AnalyticalF3D< T, T > & f)

Schedule analytical functor for output.

Definition at line 46 of file vtuSurfaceWriter.hh.

47{
48 _functorsA.push_back(&f);
49}
std::vector< AnalyticalF3D< T, T > * > _functorsA

◆ addFunctor() [2/2]

template<typename T >
void olb::VTUsurfaceWriter< T >::addFunctor ( SuperF3D< T, T > & f)

Schedule lattice functor for output (interpolates)

Definition at line 52 of file vtuSurfaceWriter.hh.

53{
54 _functors.push_back(&f);
55}
std::vector< SuperF3D< T, T > * > _functors

◆ addPoint()

template<typename T >
void olb::VTUsurfaceWriter< T >::addPoint ( Vector< T, 3 > & new_position)

Add single location to output set.

Definition at line 66 of file vtuSurfaceWriter.hh.

67{
68 pos.push_back(new_position);
69}
std::vector< Vector< T, 3 > > pos

◆ addPoints()

template<typename T >
void olb::VTUsurfaceWriter< T >::addPoints ( std::vector< Vector< T, 3 > > & new_positions)

Add multiple locations to output set.

Definition at line 58 of file vtuSurfaceWriter.hh.

59{
60 for (unsigned long i = 0; i < new_positions.size(); i++) {
61 pos.push_back(new_positions[i]);
62 }
63}

◆ addSTL()

template<typename T >
void olb::VTUsurfaceWriter< T >::addSTL ( STLreader< T > & stlReader)

Add full STL surface as output surface.

Definition at line 509 of file vtuSurfaceWriter.hh.

510{
511 STLmesh<T> mesh = stlReader.getMesh();
512 std::vector<std::vector<T>> points;
513
514 // Collect points from triangles
515 for (const STLtriangle<T>& triangle : mesh.getTriangles()) {
516 std::vector<T> pt = {triangle.point[0][0], triangle.point[0][1], triangle.point[0][2]};
517 points.push_back(pt);
518 pt = {triangle.point[1][0], triangle.point[1][1], triangle.point[1][2]};
519 points.push_back(pt);
520 pt = {triangle.point[2][0], triangle.point[2][1], triangle.point[2][2]};
521 points.push_back(pt);
522 }
523 // Remove duplicate points
524 points = removeDuplicates(points);
525
526 // Add points to pos
527 for (const std::vector<T>& point : points) {
528 pos.push_back(Vector<T, 3>(point[0], point[1], point[2]));
529 }
530
531 std::vector<std::vector<int>> indexedTriangles = mapTrianglesToIndices(stlReader);
532 for (const auto& tri : indexedTriangles) {
533 connectivity.push_back(tri);
534 }
535}
std::vector< std::vector< T > > removeDuplicates(const std::vector< std::vector< T > > &points)
removes duplicates from points, needed for STL files
std::vector< std::vector< int > > mapTrianglesToIndices(STLreader< T > &stlReader)
T triangle(Vector< T, 2 > p, Vector< T, 2 > a, Vector< T, 2 > b, Vector< T, 2 > c) any_platform
Exact signed distance to the surface of two-dimensional triangle.
Definition sdf.h:152

References olb::STLreader< T >::getMesh(), and olb::STLmesh< T >::getTriangles().

+ Here is the call graph for this function:

◆ areVectorsEqual()

template<typename T >
bool olb::VTUsurfaceWriter< T >::areVectorsEqual ( const std::vector< T > & v1,
const std::vector< T > & v2,
T epsilon = 1e-6 )
staticprotected

Definition at line 452 of file vtuSurfaceWriter.hh.

453{
454 for (std::size_t i = 0; i < v1.size(); ++i) {
455 if (std::fabs(v1[i] - v2[i]) > epsilon) {
456 return false;
457 }
458 }
459 return true;
460}

◆ closePVD()

template<typename T >
void olb::VTUsurfaceWriter< T >::closePVD ( const std::string & fullNamePVD)
protected

Definition at line 199 of file vtuSurfaceWriter.hh.

200{
201 int rank = 0;
202#ifdef PARALLEL_MODE_MPI
203 rank = singleton::mpi().getRank();
204#endif
205 if (rank == 0) {
206
207 std::ofstream fout(fullNamePVD.c_str(), std::ios::app);
208 if (!fout) {
209 clout << "Error: could not open " << fullNamePVD << std::endl;
210 }
211 fout << "</Collection>\n";
212 fout << "</VTKFile>\n";
213 fout.close();
214 }
215}
int getRank() const
Returns the process ID.
MpiManager & mpi()

References olb::singleton::MpiManager::getRank(), and olb::singleton::mpi().

+ Here is the call graph for this function:

◆ closeVTU()

template<typename T >
void olb::VTUsurfaceWriter< T >::closeVTU ( const std::string & fullNamePiece)
protected

performes </ImageData> and </VTKFile>

Definition at line 130 of file vtuSurfaceWriter.hh.

131{
132 std::ofstream fout(fullNamePiece.c_str(), std::ios::app);
133 if (!fout) {
134 clout << "Error: could not open " << fullNamePiece << std::endl;
135 }
136 fout << "</UnstructuredGrid>\n";
137 fout << "</VTKFile>\n";
138 fout.close();
139}

◆ createMasterFile()

template<typename T >
void olb::VTUsurfaceWriter< T >::createMasterFile ( )

Definition at line 171 of file vtuSurfaceWriter.hh.

172{
173 int rank = 0;
174#ifdef PARALLEL_MODE_MPI
175 rank = singleton::mpi().getRank();
176#endif
177 if (rank == 0) {
178 std::string fullNamePVDmaster = singleton::directories().getVtkOutDir() + createFileName(_name) + ".pvd";
179 preamblePVD(fullNamePVDmaster);
180 closePVD(fullNamePVDmaster);
181 }
182}
void closePVD(const std::string &fullNamePVD)
void preamblePVD(const std::string &fullNamePVD)
std::string getVtkOutDir() const
Definition singleton.h:103
Directories & directories()
Definition singleton.h:162
std::string createFileName(std::string name)
for .pvd masterFile
Definition fileName.hh:34

References olb::createFileName(), olb::singleton::directories(), olb::singleton::MpiManager::getRank(), olb::singleton::Directories::getVtkOutDir(), and olb::singleton::mpi().

+ Here is the call graph for this function:

◆ dataArray()

template<typename T >
void olb::VTUsurfaceWriter< T >::dataArray ( const std::string & fullName)
protected

Definition at line 218 of file vtuSurfaceWriter.hh.

219{
220
221 std::ofstream fout(fullName.c_str(), std::ios::app);
222 if (!fout) {
223 clout << "Error: could not open " << fullName << std::endl;
224 }
225
226 for (unsigned long i = 0; i < _functors.size(); i++) {
227 writeFunctor(fullName, fout, *_functors[i]);
228 }
229
230 for (unsigned long i = 0; i < _functorsA.size(); i++) {
231 writeAnalyticalFunctor(fullName, fout, *_functorsA[i]);
232 }
233 // master writes points and connectivity
234 int rank = 0;
235#ifdef PARALLEL_MODE_MPI
236 rank = singleton::mpi().getRank();
237#endif
238 if (rank == 0) {
239 fout << "</PointData>" << std::endl;
240 fout << "<CellData /> " << std::endl;
241 fout << "<Cells>" << std::endl;
242 fout << "<DataArray type=\"Int32\" Name=\"connectivity\" format=\"ascii\">" << std::endl;
243 // if there are no connections between the points, connectivity simply counts up
244 if (connectivity.size() == 0) {
245 for (int i = 0; i < 1; ++i) {
246 for (unsigned long j = 0; j < pos.size(); ++j) {
247 fout << j << " ";
248 }
249 }
250 }
251 // if there are connections, connectivity is written out
252 else {
253 for (int i = 0; i < 1; ++i) {
254 for (unsigned long j = 0; j < connectivity.size(); ++j) {
255 for (unsigned long k = 0; k < connectivity[j].size(); ++k) {
256 fout << connectivity[j][k] << " ";
257 }
258 if (j % 2 == 0) {
259 fout << std::endl;
260 }
261 }
262 }
263 }
264
265 fout << "</DataArray>" << std::endl;
266 fout << "<DataArray type=\"Int32\" Name=\"offsets\" format=\"ascii\">" << std::endl;
267 // if there are no connections, offsets just count up
268 if (connectivity.size() == 0) {
269 for (std::size_t i = 0; i < 1; ++i) {
270 for (unsigned long j = 1; j <= connectivity.size(); ++j) {
271 fout << j << " ";
272 if (j % 6 == 0) {
273 fout << std::endl;
274 }
275 }
276 }
277 // if connections are present, offset is multiplied by 3
278 }
279 else { //connectivity.size() != 0
280 for (std::size_t i = 0; i < 1; ++i) {
281 for (unsigned long j = 1; j <= connectivity.size(); ++j) {
282 fout << j * 3 << " ";
283 if (j % 6 == 0) {
284 fout << std::endl;
285 }
286 }
287 }
288 }
289 // Connected points have the type 5, unconnected points have the type 3
290 fout << "</DataArray>" << std::endl;
291 fout << "<DataArray type=\"UInt8\" Name=\"types\" format=\"ascii\">" << std::endl;
292 if (connectivity.size() != 0) {
293 for (unsigned long j = 0; j < connectivity.size(); ++j) {
294 fout << 5 << " ";
295 if (j % 6 == 0) {
296 fout << std::endl;
297 }
298 }
299 }
300 else {
301 for (unsigned long j = 0; j < pos.size(); ++j) {
302 fout << 3 << " ";
303 if (j % 6 == 0) {
304 fout << std::endl;
305 }
306 }
307 }
308 // write out all coordinates of points
309 fout << "</DataArray>" << std::endl;
310 fout << "</Cells>" << std::endl;
311 fout << "<Points>" << std::endl;
312 writePosition(fout);
313 fout << "</Points>" << std::endl;
314 fout << "</Piece>" << std::endl;
315
316 fout.close();
317 }
318}
void writeAnalyticalFunctor(const std::string &fullName, std::ofstream &fout, AnalyticalF3D< T, T > &f)
writes functors.
void writePosition(std::ofstream &fout)
writes coordinates of points
void writeFunctor(const std::string &fullName, std::ofstream &fout, SuperF3D< T, T > &f)
interpolates and writes functors stored at functors

References olb::singleton::MpiManager::getRank(), and olb::singleton::mpi().

+ Here is the call graph for this function:

◆ dataPVD()

template<typename T >
void olb::VTUsurfaceWriter< T >::dataPVD ( int iT,
int i,
const std::string & fullNamePVD,
const std::string & namePiece )
protected

performes <DataSet timestep= ... file=namePiece >

Definition at line 142 of file vtuSurfaceWriter.hh.

143{
144 std::ofstream fout(fullNamePVD.c_str(), std::ios::app);
145 if (!fout) {
146 clout << "Error: could not open " << fullNamePVD << std::endl;
147 }
148 fout << "<DataSet timestep=\"" << iT << "\" " << "group=\"\" part=\" " << i << "\" " << "file=\"" << namePiece
149 << "\"/>\n";
150 fout.close();
151}

◆ dataPVDmaster()

template<typename T >
void olb::VTUsurfaceWriter< T >::dataPVDmaster ( int iT,
int i,
const std::string & fullNamePVDMaster,
const std::string & namePiece )
protected

performes <DataSet timestep= ... file=namePiece >

Definition at line 154 of file vtuSurfaceWriter.hh.

156{
157 std::ofstream fout(fullNamePVDMaster.c_str(), std::ios::in | std::ios::out | std::ios::ate);
158 if (fout) {
159 fout.seekp(-25, std::ios::end); // jump -25 from the end of file to overwrite closePVD
160 fout << "<DataSet timestep=\"" << iT << "\" " << "group=\"\" part=\" " << i << "\" " << "file=\"" << namePiece
161 << "\"/>\n";
162 fout.close();
163 closePVD(fullNamePVDMaster);
164 }
165 else {
166 clout << "Error: could not open " << fullNamePVDMaster << std::endl;
167 }
168}

◆ mapTrianglesToIndices()

template<typename T >
std::vector< std::vector< int > > olb::VTUsurfaceWriter< T >::mapTrianglesToIndices ( STLreader< T > & stlReader)
protected

Definition at line 463 of file vtuSurfaceWriter.hh.

464{
465 STLmesh<T> mesh = stlReader.getMesh();
466
467 struct VectorHash {
468 std::size_t operator()(const std::vector<T>& v) const
469 {
470 return std::hash<T>()(v[0]) ^ std::hash<T>()(v[1]) ^ std::hash<T>()(v[2]);
471 }
472 };
473
474 struct VectorEqual {
475 bool operator()(const std::vector<T>& a, const std::vector<T>& b) const { return areVectorsEqual(a, b); }
476 };
477
478 std::unordered_map<std::vector<T>, int, VectorHash, VectorEqual> pointIndexMap;
479 std::vector<std::vector<int>> indexedTriangles;
480
481 // Create a map from points to their indices
482 for (std::size_t i = 0; i < pos.size(); ++i) {
483 pointIndexMap[{pos[i][0], pos[i][1], pos[i][2]}] = i;
484 }
485
486 // Map each triangle's points to their indices
487 for (const STLtriangle<T>& triangle : mesh.getTriangles()) {
488 std::vector<int> indices(3);
489 for (int j = 0; j < 3; ++j) {
490 std::vector<T> pt = {triangle.point[j][0], triangle.point[j][1], triangle.point[j][2]};
491 auto it = pointIndexMap.find(pt);
492 if (it != pointIndexMap.end()) {
493 indices[j] = it->second;
494 }
495 else {
496 // Handle the case where the point is not found (should not happen if removeDuplicates worked correctly)
497 indices[j] = -1; // Invalid index
498 }
499 }
500 if (indices[0] != -1 && indices[1] != -1 && indices[2] != -1) {
501 indexedTriangles.push_back(indices);
502 }
503 }
504
505 return indexedTriangles;
506}

References olb::STLreader< T >::getMesh(), and olb::STLmesh< T >::getTriangles().

+ Here is the call graph for this function:

◆ preamblePVD()

template<typename T >
void olb::VTUsurfaceWriter< T >::preamblePVD ( const std::string & fullNamePVD)
protected

Definition at line 185 of file vtuSurfaceWriter.hh.

186{
187 std::ofstream fout(fullNamePVD.c_str(), std::ios::trunc);
188 if (!fout) {
189 clout << "Error: could not open " << fullNamePVD << std::endl;
190 }
191 fout << "<?xml version=\"1.0\"?>\n";
192 fout << "<VTKFile type=\"Collection\" version=\"1.0\" "
193 << "byte_order=\"LittleEndian\">\n"
194 << "<Collection>\n";
195 fout.close();
196}

◆ preambleVTU()

template<typename T >
void olb::VTUsurfaceWriter< T >::preambleVTU ( const std::string & fullName,
int num )
protected

performes <VTKFile ...>, <ImageData ...> and <PieceExtent ...>

Definition at line 115 of file vtuSurfaceWriter.hh.

116{
117 std::ofstream fout(fullName.c_str(), std::ios::trunc);
118 if (!fout) {
119 clout << "Error: could not open " << fullName << std::endl;
120 }
121 fout << "<?xml version=\"1.0\"?>" << std::endl << std::flush;
122 fout << "<VTKFile type=\"UnstructuredGrid\" version=\"1.0\" byte_order=\"LittleEndian\">" << std::endl;
123 fout << "<UnstructuredGrid>" << std::endl;
124 fout << "<Piece NumberOfPoints=\"" << num << "\" NumberOfCells=\"" << connectivity.size() << "\">" << std::endl;
125 fout << "<PointData Vectors=\"Particles\">" << std::endl;
126 fout.close();
127}

◆ removeDuplicates()

template<typename T >
std::vector< std::vector< T > > olb::VTUsurfaceWriter< T >::removeDuplicates ( const std::vector< std::vector< T > > & points)
protected

removes duplicates from points, needed for STL files

Definition at line 432 of file vtuSurfaceWriter.hh.

433{
434 struct VectorHash {
435 std::size_t operator()(const std::vector<T>& v) const
436 {
437 return std::hash<T>()(v[0]) ^ std::hash<T>()(v[1]) ^ std::hash<T>()(v[2]);
438 }
439 };
440
441 struct VectorEqual {
442 bool operator()(const std::vector<T>& a, const std::vector<T>& b) const
443 {
444 return a[0] == b[0] && a[1] == b[1] && a[2] == b[2];
445 }
446 };
447 std::unordered_set<std::vector<T>, VectorHash, VectorEqual> uniquePoints(points.begin(), points.end());
448 return std::vector<std::vector<T>>(uniquePoints.begin(), uniquePoints.end());
449}

◆ write()

template<typename T >
void olb::VTUsurfaceWriter< T >::write ( std::size_t iT)

writes vtu files and assigns values

Definition at line 81 of file vtuSurfaceWriter.hh.

82{
83 int rank = 0;
84#ifdef PARALLEL_MODE_MPI
85 rank = singleton::mpi().getRank();
86#endif
87
88 for (SuperF3D<T, T>* f : _functors) {
89 f->getSuperStructure().communicate();
90 }
91 std::string fullNamePVDmaster = singleton::directories().getVtkOutDir() + createFileName(_name) + ".pvd";
92 std::string namePiece = "data/" + createFileName(_name, iT, 0) + ".vtu";
93
94 if (rank == 0) {
95 // puts name of .vti piece to a .pvd file [fullNamePVD]
96 // adds a namePiece to master.pvd file.
97 // To do so we overwrite closePVD() and add new entry.
98 dataPVDmaster(iT, 1, fullNamePVDmaster, namePiece);
99 }
100
101 //VTU data files
102
103 std::string fullNameVTU = singleton::directories().getVtkOutDir() + "data/" + createFileName(_name, iT, 0) + ".vtu";
104 if (rank == 0) {
105 preambleVTU(fullNameVTU, pos.size());
106 }
107 //writes data arrays into vtu file. All ranks are needed for this
108 this->dataArray(fullNameVTU);
109 if (rank == 0) {
110 closeVTU(fullNameVTU);
111 }
112}
void closeVTU(const std::string &fullNamePiece)
performes </ImageData> and </VTKFile>
void dataPVDmaster(int iT, int i, const std::string &fullNamePVDMaster, const std::string &namePiece)
performes <DataSet timestep= ... file=namePiece >
void preambleVTU(const std::string &fullName, int num)
performes <VTKFile ...>, <ImageData ...> and <PieceExtent ...>
void dataArray(const std::string &fullName)

References olb::createFileName(), olb::singleton::directories(), olb::singleton::MpiManager::getRank(), olb::singleton::Directories::getVtkOutDir(), and olb::singleton::mpi().

+ Here is the call graph for this function:

◆ writeAnalyticalFunctor()

template<typename T >
void olb::VTUsurfaceWriter< T >::writeAnalyticalFunctor ( const std::string & fullName,
std::ofstream & fout,
AnalyticalF3D< T, T > & f )
protected

writes functors.

Definition at line 340 of file vtuSurfaceWriter.hh.

342{
343 std::size_t dimF = f.getTargetDim();
344
345 //Header DataArray, ascii or Base64 binary
346 int rank = 0;
347#ifdef PARALLEL_MODE_MPI
348 rank = singleton::mpi().getRank();
349#endif
350
351 if (rank == 0) {
352 if (!_binary) {
353 fout << "<DataArray type=\"Float32\" Name=\"" << f.getName() << "\" NumberOfComponents=\"" << dimF << "\">"
354 << std::endl;
355 }
356 else if (_binary) {
357 fout << "<DataArray type=\"Float32\" Name=\"" << f.getName()
358 << "\" format=\"binary\" encoding=\"base64\" NumberOfComponents=\"" << dimF << "\">" << std::endl;
359 }
360 }
361
362 if (_binary) {
363 fout.close();
364 }
365 std::size_t numPoints = pos.size();
366 std::ofstream ofstr(fullName.c_str(), std::ios::out | std::ios::app | std::ios::binary);
367 size_t binarySize = size_t(numPoints * dimF * sizeof(float));
368 Base64Encoder<unsigned int> sizeEncoder(ofstr, 1);
369 unsigned int uintBinarySize = (unsigned int)binarySize;
370 if (_binary) {
371 sizeEncoder.encode(&uintBinarySize, 1);
372 }
373 Base64Encoder<float> dataEncoder(ofstr, numPoints);
374
375 std::vector<T> vals(dimF * numPoints, 0);
376
377 T* val = new T[dimF];
378
379 for (std::size_t j = 0; j < numPoints; ++j) {
380 if (_loadBalancer.isLocal(_cuboidDecomposition.getC(pos[j]).value())) {
381 f(val, pos[j].data());
382
383 for (std::size_t i = 0; i < dimF; ++i) {
384 vals[i * numPoints + j] = val[i];
385 }
386 }
387 }
388#ifdef PARALLEL_MODE_MPI
389 // reduce the values to the master, then master writes them to the file
390 std::vector<T> valsCummulated(dimF * numPoints, 0);
391 singleton::mpi().reduceVect(vals, valsCummulated, MPI_SUM);
392
393 if (rank == 0) {
394 for (std::size_t j = 0; j < numPoints; ++j) {
395 for (std::size_t i = 0; i < dimF; ++i) {
396 if (!_binary) {
397 fout << valsCummulated[i * numPoints + j] << " ";
398 }
399 else {
400 const float helper = float(valsCummulated[j]);
401 dataEncoder.encode(&helper, 1);
402 }
403 }
404 }
405 }
406#endif
407#ifndef PARALLEL_MODE_MPI
408 // write the values to the file
409 for (std::size_t j = 0; j < numPoints; ++j) {
410 for (std::size_t i = 0; i < dimF; ++i) {
411 if (!_binary) {
412 fout << vals[i * numPoints + j] << " ";
413 }
414 else {
415 const float helper = float(vals[i * numPoints + j]);
416 dataEncoder.encode(&helper, 1);
417 }
418 }
419 }
420#endif
421
422 if (_binary) {
423 ofstr.close();
424 fout.open(fullName.c_str(), std::ios::out | std::ios::app);
425 }
426 if (rank == 0) {
427 fout << "</DataArray>" << std::endl;
428 }
429}
void reduceVect(std::vector< T > &sendVal, std::vector< T > &recvVal, MPI_Op op, int root=0, MPI_Comm comm=MPI_COMM_WORLD)
Element-per-element reduction of a vector of data.

References olb::Base64Encoder< T >::encode(), olb::GenericF< T, S >::getName(), olb::singleton::MpiManager::getRank(), olb::GenericF< T, S >::getTargetDim(), olb::singleton::mpi(), and olb::singleton::MpiManager::reduceVect().

+ Here is the call graph for this function:

◆ writeFunctor()

template<typename T >
void olb::VTUsurfaceWriter< T >::writeFunctor ( const std::string & fullName,
std::ofstream & fout,
SuperF3D< T, T > & f )
protected

interpolates and writes functors stored at functors

Definition at line 333 of file vtuSurfaceWriter.hh.

334{
335 AnalyticalFfromSuperF3D<T> interpolateF(f, false);
336 writeAnalyticalFunctor(fullName, fout, interpolateF);
337}

◆ writePosition()

template<typename T >
void olb::VTUsurfaceWriter< T >::writePosition ( std::ofstream & fout)
protected

writes coordinates of points

Definition at line 321 of file vtuSurfaceWriter.hh.

322{
323 fout << "<DataArray type=\"Float32\" Name=\"Point\" NumberOfComponents=\"" << 3 << "\">" << std::endl;
324
325 int num = pos.size();
326 for (int j = 0; j < num; ++j) {
327 fout << pos[j][0] << " " << pos[j][1] << " " << pos[j][2] << " " << std::endl;
328 }
329 fout << "</DataArray>" << std::endl;
330}

Member Data Documentation

◆ _binary

template<typename T >
bool olb::VTUsurfaceWriter< T >::_binary
protected

Definition at line 45 of file vtuSurfaceWriter.h.

◆ _cuboidDecomposition

template<typename T >
CuboidDecomposition3D<T>& olb::VTUsurfaceWriter< T >::_cuboidDecomposition
protected

Definition at line 43 of file vtuSurfaceWriter.h.

◆ _functors

template<typename T >
std::vector<SuperF3D<T, T>*> olb::VTUsurfaceWriter< T >::_functors
protected

Definition at line 39 of file vtuSurfaceWriter.h.

◆ _functorsA

template<typename T >
std::vector<AnalyticalF3D<T, T>*> olb::VTUsurfaceWriter< T >::_functorsA
protected

Definition at line 40 of file vtuSurfaceWriter.h.

◆ _loadBalancer

template<typename T >
LoadBalancer<T>& olb::VTUsurfaceWriter< T >::_loadBalancer
protected

Definition at line 44 of file vtuSurfaceWriter.h.

◆ _name

template<typename T >
std::string olb::VTUsurfaceWriter< T >::_name
protected

Definition at line 38 of file vtuSurfaceWriter.h.

◆ clout

template<typename T >
OstreamManager olb::VTUsurfaceWriter< T >::clout
protected

Definition at line 47 of file vtuSurfaceWriter.h.

◆ connectivity

template<typename T >
std::vector<std::vector<int> > olb::VTUsurfaceWriter< T >::connectivity
protected

Definition at line 42 of file vtuSurfaceWriter.h.

◆ pos

template<typename T >
std::vector<Vector<T, 3> > olb::VTUsurfaceWriter< T >::pos
protected

Definition at line 41 of file vtuSurfaceWriter.h.


The documentation for this class was generated from the following files: