Loading...
Searching...
No Matches
ITHACAstream.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ██╗████████╗██╗ ██╗ █████╗ ██████╗ █████╗ ███████╗██╗ ██╗
3 ██║╚══██╔══╝██║ ██║██╔══██╗██╔════╝██╔══██╗ ██╔════╝██║ ██║
4 ██║ ██║ ███████║███████║██║ ███████║█████╗█████╗ ██║ ██║
5 ██║ ██║ ██╔══██║██╔══██║██║ ██╔══██║╚════╝██╔══╝ ╚██╗ ██╔╝
6 ██║ ██║ ██║ ██║██║ ██║╚██████╗██║ ██║ ██║ ╚████╔╝
7 ╚═╝ ╚═╝ ╚═╝ ╚═╝╚═╝ ╚═╝ ╚═════╝╚═╝ ╚═╝ ╚═╝ ╚═══╝
8
9 * In real Time Highly Advanced Computational Applications for Finite Volumes
10 * Copyright (C) 2017 by the ITHACA-FV authors
11-------------------------------------------------------------------------------
12
13License
14 This file is part of ITHACA-FV
15
16 ITHACA-FV is free software: you can redistribute it and/or modify
17 it under the terms of the GNU Lesser General Public License as published by
18 the Free Software Foundation, either version 3 of the License, or
19 (at your option) any later version.
20
21 ITHACA-FV is distributed in the hope that it will be useful,
22 but WITHOUT ANY WARRANTY; without even the implied warranty of
23 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24 GNU Lesser General Public License for more details.
25
26 You should have received a copy of the GNU Lesser General Public License
27 along with ITHACA-FV. If not, see <http://www.gnu.org/licenses/>.
28
29\*---------------------------------------------------------------------------*/
30
31
32#include "ITHACAstream.H"
33
34
37
38// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
39
40namespace ITHACAstream
41{
42
43template<typename Type>
44void exportFvMatrix(fvMatrix<Type>& Matrix, word folder,
45 word MatrixName)
46{
47 Eigen::SparseMatrix<double> A;
48 Eigen::VectorXd b;
49 Foam2Eigen::fvMatrix2Eigen(Matrix, A, b);
50 SaveSparseMatrix(A, folder + "/", "A_" + MatrixName);
51 SaveDenseMatrix(b, folder + "/", "B_" + MatrixName);
52}
53
54template <typename T, int dim>
55void exportMatrix(Eigen::Matrix < T, -1, dim > & matrix,
56 word Name, word type,
57 word folder)
58{
59 std::string message = "The extension \"" + type +
60 "\" was not implemented. Check the list of possible extensions.";
61 M_Assert(type == "python" || type == "matlab"
62 || type == "eigen", message.c_str()
63 );
64 mkDir(folder);
65 word est;
66
67 if (type == "python")
68 {
69 est = ".py";
70 OFstream str(folder + "/" + Name + "_mat" + est);
71 str << Name << "=np.array([";
72
73 for (int i = 0; i < matrix.rows(); i++)
74 {
75 for (int j = 0; j < matrix.cols(); j++)
76 {
77 if (j == 0)
78 {
79 str << "[" << setprecision(10) << matrix(i, j);
80 }
81 else
82 {
83 str << "," << setprecision(10) << matrix(i, j);
84 }
85 }
86
87 if (i != (matrix.rows() - 1))
88 {
89 str << "]," << endl;
90 }
91 }
92
93 str << "]])" << endl;
94 }
95
96 if (type == "matlab")
97 {
98 est = ".m";
99 OFstream str(folder + "/" + Name + "_mat" + est);
100 str << Name << "=[";
101
102 for (int i = 0; i < matrix.rows(); i++)
103 {
104 for (int j = 0; j < matrix.cols(); j++)
105 {
106 str << " " << setprecision(10) << matrix(i, j);
107 }
108
109 if (i != (matrix.rows() - 1))
110 {
111 str << ";" << endl;
112 }
113 }
114
115 str << "];" << endl;
116 }
117
118 if (type == "eigen")
119 {
120 const static Eigen::IOFormat CSVFormat(6, false, ", ", "\n");
121 std::ofstream ofs;
122 ofs.precision(20);
123 ofs.open (folder + "/" + Name + "_mat.txt");
124
125 for (int i = 0; i < matrix.rows(); i++)
126 {
127 for (int j = 0; j < matrix.cols(); j++)
128 {
129 if (j == 0)
130 {
131 ofs << matrix(i, j);
132 }
133 else
134 {
135 ofs << " " << matrix(i, j);
136 }
137 }
138
139 if (i != (matrix.rows() - 1))
140 {
141 ofs << std::endl;
142 }
143 }
144
145 ofs.close();
146 }
147}
148
149template void exportMatrix(Eigen::Matrix < double, -1,
150 -1 > & matrix, word Name, word type,
151 word folder);
152
153template void exportMatrix(Eigen::Matrix < int, -1,
154 -1 > & matrix, word Name, word type,
155 word folder);
156
157template void exportMatrix(Eigen::Matrix < float, -1,
158 -1 > & matrix, word Name, word type,
159 word folder);
160
161template void exportMatrix(Eigen::Matrix < double, -1,
162 1 > & matrix, word Name, word type,
163 word folder);
164
165template void exportMatrix(Eigen::Matrix < int, -1,
166 1 > & matrix, word Name, word type,
167 word folder);
168
169template void exportMatrix(Eigen::Matrix < float, -1,
170 1 > & matrix, word Name, word type,
171 word folder);
172
173void exportMatrix(List <Eigen::MatrixXd>& matrix, word Name,
174 word type, word folder)
175{
176 std::string message = "The extension \"" + type +
177 "\" was not implemented. Check the list of possible extensions.";
178 M_Assert(type == "python" || type == "matlab"
179 || type == "eigen", message.c_str()
180 );
181 mkDir(folder);
182 word est;
183
184 // Python Case
185 if (type == "python")
186 {
187 est = ".py";
188 OFstream str(folder + "/" + Name + "_mat" + est);
189 str << Name << "=np.zeros([" << matrix.size() << "," << matrix[0].rows() <<
190 "," << matrix[0].cols() << "])\n";
191
192 for (int i = 0; i < matrix.size(); i++)
193 {
194 str << Name << "[" << i << ",:,:]=np.array([";
195
196 for (int j = 0; j < matrix[0].rows(); j++)
197 {
198 for (int k = 0; k < matrix[0].cols(); k++)
199 {
200 if ( k == 0)
201 {
202 str << "[" << setprecision(10) << matrix[i](j, k);
203 }
204 else
205 {
206 str << "," << setprecision(10) << matrix[i](j, k);
207 }
208 }
209
210 if (j != (matrix[0].rows() - 1))
211 {
212 str << "]," << endl;
213 }
214 }
215
216 str << "]])\n" << endl;
217 }
218 }
219 // Matlab case
220 else if (type == "matlab")
221 {
222 est = ".m";
223 OFstream str(folder + "/" + Name + "_mat" + est);
224
225 for (int i = 0; i < matrix.size(); i++)
226 {
227 str << Name << "(" << i + 1 << ",:,:)=[";
228
229 for (int j = 0; j < matrix[0].rows(); j++)
230 {
231 for (int k = 0; k < matrix[0].cols(); k++)
232 {
233 str << " " << setprecision(10) << matrix[i](j, k);
234 }
235
236 if (j != (matrix[0].rows() - 1))
237 {
238 str << ";" << endl;
239 }
240 }
241
242 str << "];" << endl;
243 }
244 }
245 else if (type == "eigen")
246 {
247 for (int i = 0; i < matrix.size(); i++)
248 {
249 word Namei = Name + name(i);
250 exportMatrix(matrix[i], Namei, "eigen", folder);
251 }
252 }
253}
254
255void exportVector(Eigen::VectorXd& vector,
256 word Name, word type,
257 word folder)
258{
259 Eigen::MatrixXd matrix = vector;
260 exportMatrix(matrix, Name, type, folder);
261}
262
263template<typename T>
264void exportTensor(Eigen::Tensor<T, 3> tensor, word Name,
265 word type, word folder)
266{
267 std::string message = "The extension \"" + type +
268 "\" was not implemented. Check the list of possible extensions.";
269 M_Assert(type == "python" || type == "matlab"
270 || type == "eigen", message.c_str()
271 );
272 mkDir(folder);
273 word est;
274
275 // Python Case
276 if (type == "python")
277 {
278 est = ".py";
279 OFstream str(folder + "/" + Name + "_mat" + est);
280 str << Name << "=np.zeros([" << tensor.dimension(0) << "," <<
281 Eigen::SliceFromTensor(
282 tensor, 0,
283 0).rows() <<
284 "," << Eigen::SliceFromTensor(tensor, 0,
285 0).cols() << "])\n";
286
287 for (int i = 0; i < tensor.dimension(0); i++)
288 {
289 str << Name << "[" << i << ",:,:]=np.array([";
290
291 for (int j = 0; j < Eigen::SliceFromTensor(tensor, 0, 0).rows(); j++)
292 {
293 for (int k = 0; k < Eigen::SliceFromTensor(tensor, 0, 0).cols(); k++)
294 {
295 if ( k == 0)
296 {
297 str << "[" << setprecision(10) << Eigen::SliceFromTensor(tensor, 0,
298 i)(j, k);
299 }
300 else
301 {
302 str << "," << setprecision(10) << Eigen::SliceFromTensor(tensor, 0,
303 i)(j, k);
304 }
305 }
306
307 if (j != (Eigen::SliceFromTensor(tensor, 0,
308 0).rows() - 1))
309 {
310 str << "]," << endl;
311 }
312 }
313
314 str << "]])\n" << endl;
315 }
316 }
317 // Matlab case
318 else if (type == "matlab")
319 {
320 est = ".m";
321 OFstream str(folder + "/" + Name + "_mat" + est);
322
323 for (int i = 0; i < tensor.dimension(0); i++)
324 {
325 str << Name << "(" << i + 1 << ",:,:)=[";
326
327 for (int j = 0; j < Eigen::SliceFromTensor(tensor, 0,
328 0).rows(); j++)
329 {
330 for (int k = 0; k < Eigen::SliceFromTensor(tensor, 0,
331 0).cols(); k++)
332 {
333 str << " " << setprecision(10) << Eigen::SliceFromTensor(tensor, 0,
334 i)(j, k);
335 }
336
337 if (j != (Eigen::SliceFromTensor(tensor, 0,
338 0).rows() - 1))
339 {
340 str << ";" << endl;
341 }
342 }
343
344 str << "];" << endl;
345 }
346 }
347 else if (type == "eigen")
348 {
349 for (int i = 0; i < tensor.dimension(0); i++)
350 {
351 Eigen::Matrix < T, -1, -1 > matrixAux = Eigen::SliceFromTensor(tensor, 0, i);
352 word Namei = Name + name(i);
353 exportMatrix(matrixAux, Namei, "eigen", folder);
354 }
355 }
356}
357
358
359
360template void exportTensor(Eigen::Tensor<double, 3> tensor,
361 word Name,
362 word type, word folder);
363
364template void exportTensor(Eigen::Tensor<int, 3> tensor,
365 word Name,
366 word type, word folder);
367
368template void exportTensor(Eigen::Tensor<float, 3> tensor,
369 word Name,
370 word type, word folder);
371
372List<Eigen::MatrixXd> readMatrix(word folder, word mat_name)
373{
374 int file_count = 0;
375 DIR* dirp;
376 struct dirent* entry;
377 dirp = opendir(folder.c_str());
378 List <Eigen::MatrixXd> result;
379
380 while ((entry = readdir(dirp)) != NULL)
381 {
382 if (entry->d_type == DT_REG)
383 {
384 file_count++;
385 }
386 }
387
388 closedir (dirp);
389
390 for (int j = 0; j < file_count ; j++)
391 {
392 word matname = folder + "/" + mat_name + name(j) + "_mat.txt";
393 Eigen::MatrixXd temp = readMatrix(matname);
394 result.append(temp);
395 }
396
397 return result;
398}
399
400Eigen::MatrixXd readMatrix(word filename)
401{
402 int cols = 0, rows = 0;
403 double buff[MAXBUFSIZE];
404 // Read numbers from file into buffer.
405 std::ifstream infile;
406 infile.open(filename.c_str());
407 std::string message = "The matrix file \"" + filename +
408 "\" does not exist. Check the existence of the file or the way it is named.";
409 M_Assert(infile.good() != 0, message.c_str()
410 );
411
412 while (! infile.eof())
413 {
414 string line;
415 getline(infile, line);
416 int temp_cols = 0;
417 std::stringstream stream(line);
418
419 while (! stream.eof())
420 {
421 stream >> buff[cols * rows + temp_cols++];
422 }
423
424 if (temp_cols == 0)
425 {
426 continue;
427 }
428
429 if (cols == 0)
430 {
431 cols = temp_cols;
432 }
433
434 rows++;
435 }
436
437 infile.close();
438 rows--;
439 // Populate matrix with numbers.
440 Eigen::MatrixXd result(rows, cols);
441
442 for (int i = 0; i < rows; i++)
443 {
444 for (int j = 0; j < cols; j++)
445 {
446 result(i, j) = buff[ cols * i + j ];
447 }
448 }
449
450 return result;
451}
452
453template<class Type, template<class> class PatchField, class GeoMesh>
454void read_fields(
455 PtrList<GeometricField<Type, PatchField, GeoMesh>>& Lfield, word Name,
456 fileName casename, int first_snap, int n_snap)
457{
458 ITHACAparameters* para(ITHACAparameters::getInstance());
459 fvMesh& mesh = para->mesh;
460
461 if (!Pstream::parRun())
462 {
463 Info << "######### Reading the Data for " << Name << " #########" << endl;
464 fileName rootpath(".");
465 int last_s;
466 Foam::Time runTime2(Foam::Time::controlDictName, rootpath, casename);
467
468 if (first_snap >= runTime2.times().size())
469 {
470 Info << "Error the index of the first snapshot must be smaller than the number of snapshots"
471 << endl;
472 exit(0);
473 }
474
475 if (n_snap == 0)
476 {
477 last_s = runTime2.times().size();
478 }
479 else
480 {
481 last_s = min(runTime2.times().size(), n_snap + 2);
482 }
483
484 for (int i = 2 + first_snap; i < last_s + first_snap; i++)
485 {
486 GeometricField<Type, PatchField, GeoMesh> tmp_field(
487 IOobject
488 (
489 Name,
490 casename + runTime2.times()[i].name(),
491 mesh,
492 IOobject::MUST_READ
493 ),
494 mesh
495 );
496 Lfield.append(tmp_field.clone());
497 printProgress(double(i + 1) / (last_s + first_snap));
498 }
499
500 std::cout << std::endl;
501 }
502 else
503 {
504 Info << "######### Reading the Data for " << Name << " #########" <<
505 endl;
506 word timename(mesh.time().rootPath() + "/" +
507 mesh.time().caseName() );
508 timename = timename.substr(0, timename.find_last_of("\\/"));
509 timename = timename + "/" + casename + "processor" + name(Pstream::myProcNo());
510 int last_s = numberOfFiles(casename,
511 "processor" + name(Pstream::myProcNo()) + "/");
512
513 if (first_snap > last_s)
514 {
515 Info << "Error the index of the first snapshot must be smaller than the number of snapshots"
516 << endl;
517 exit(0);
518 }
519
520 if (n_snap == 0)
521 {
522 }
523 else
524 {
525 last_s = min(last_s, n_snap + 2);
526 }
527
528 for (int i = 1 + first_snap; i < last_s + first_snap; i++)
529 {
530 GeometricField<Type, PatchField, GeoMesh> tmp_field(
531 IOobject
532 (
533 Name,
534 timename + "/" + name(i),
535 mesh,
536 IOobject::MUST_READ
537 ),
538 mesh
539 );
540 Lfield.append(tmp_field.clone());
541 printProgress(double(i + 1) / (last_s + first_snap));
542 }
543
544 Info << endl;
545 }
546}
547
548template<class Type, template<class> class PatchField, class GeoMesh>
549void read_fields(
550 PtrList<GeometricField<Type, PatchField, GeoMesh>>& Lfield,
551 GeometricField<Type, PatchField, GeoMesh>& field,
552 fileName casename, int first_snap, int n_snap)
553{
554 if (!Pstream::parRun())
555 {
556 Info << "######### Reading the Data for " << field.name() << " #########" <<
557 endl;
558 fileName rootpath(".");
559 Foam::Time runTime2(Foam::Time::controlDictName, rootpath, casename);
560 int last_s;
561
562 if (first_snap >= runTime2.times().size())
563 {
564 Info << "Error the index of the first snapshot must be smaller than the number of snapshots"
565 << endl;
566 exit(0);
567 }
568
569 if (n_snap == 0)
570 {
571 last_s = runTime2.times().size();
572 }
573 else
574 {
575 last_s = min(runTime2.times().size(), n_snap + 2);
576 }
577
578 for (int i = 2 + first_snap; i < last_s + first_snap; i++)
579 {
580 GeometricField<Type, PatchField, GeoMesh> tmp_field(
581 IOobject
582 (
583 field.name(),
584 casename + runTime2.times()[i].name(),
585 field.mesh(),
586 IOobject::MUST_READ
587 ),
588 field.mesh()
589 );
590 Lfield.append(tmp_field.clone());
591 printProgress(double(i + 1) / (last_s + first_snap));
592 }
593
594 std::cout << std::endl;
595 }
596 else
597 {
598 Info << "######### Reading the Data for " << field.name() << " #########" <<
599 endl;
600 word timename(field.mesh().time().rootPath() + "/" +
601 field.mesh().time().caseName() );
602 timename = timename.substr(0, timename.find_last_of("\\/"));
603 timename = timename + "/" + casename + "processor" + name(Pstream::myProcNo());
604 int last_s = numberOfFiles(casename,
605 "processor" + name(Pstream::myProcNo()) + "/");
606
607 if (first_snap > last_s)
608 {
609 Info << "Error the index of the first snapshot must be smaller than the number of snapshots"
610 << endl;
611 exit(0);
612 }
613
614 if (n_snap == 0)
615 {
616 }
617 else
618 {
619 last_s = min(last_s, n_snap + 2);
620 }
621
622 for (int i = 1 + first_snap; i < last_s + first_snap; i++)
623 {
624 GeometricField<Type, PatchField, GeoMesh> tmp_field(
625 IOobject
626 (
627 field.name(),
628 timename + "/" + name(i),
629 field.mesh(),
630 IOobject::MUST_READ
631 ),
632 field.mesh()
633 );
634 Lfield.append(tmp_field.clone());
635 printProgress(double(i + 1) / (last_s + first_snap));
636 }
637
638 Info << endl;
639 }
640}
641
642template<class Type, template<class> class PatchField, class GeoMesh>
643void readMiddleFields(
644 PtrList<GeometricField<Type, PatchField, GeoMesh>>& Lfield,
645 GeometricField<Type, PatchField, GeoMesh>& field, fileName casename)
646{
647 int par = 1;
648 M_Assert(ITHACAutilities::check_folder(casename + name(par)) != 0,
649 "No parameter dependent solutions stored into Offline folder");
650
651 while (ITHACAutilities::check_folder(casename + name(par)))
652 {
653 read_fields(Lfield, field, casename + name(par) + "/");
654 par++;
655 }
656}
657
658template<class Type, template<class> class PatchField, class GeoMesh>
659void readConvergedFields(
660 PtrList<GeometricField<Type, PatchField, GeoMesh>>& Lfield,
661 GeometricField<Type, PatchField, GeoMesh>& field,
662 fileName casename)
663{
664 int par = 1;
665 M_Assert(ITHACAutilities::check_folder(casename + name(par)) != 0,
666 "No parameter dependent solutions stored into Offline folder");
667 std::cout << "######### Reading the Data for " << field.name() << " #########"
668 << std::endl;
669
670 while (ITHACAutilities::check_folder(casename + name(par)))
671 {
672 int last = 1;
673
674 while (ITHACAutilities::check_folder(casename + name(par) + "/" + name(last)))
675 {
676 last++;
677 }
678
679 GeometricField<Type, PatchField, GeoMesh> tmpField(
680 IOobject
681 (
682 field.name(),
683 casename + name(par) + "/" + name(last - 1),
684 field.mesh(),
685 IOobject::MUST_READ
686 ),
687 field.mesh()
688 );
689 Lfield.append(tmpField.clone());
690 par++;
691 }
692}
693
694template<class Type, template<class> class PatchField, class GeoMesh>
695void read_last_fields(
696 PtrList<GeometricField<Type, PatchField, GeoMesh>>& Lfield,
697 const GeometricField<Type, PatchField, GeoMesh>& field,
698 const fileName casename)
699{
700 if (!Pstream::parRun())
701 {
702 Info << "######### Reading the Data for " << field.name() << " #########" <<
703 endl;
704 fileName rootpath(".");
705 Foam::Time runTime2(Foam::Time::controlDictName, rootpath, casename);
706 int last_s (runTime2.times().size());
707
708 #if defined(OFVER) && (OFVER >= 2212)
709 Lfield.emplace_back
710 (
711 IOobject
712 (
713 field.name(),
714 casename + runTime2.times()[last_s-1].name(),
715 field.mesh(),
716 IOobject::MUST_READ
717 ),
718 field.mesh()
719 );
720 #else
721 auto tfld =
722 autoPtr<GeometricField<Type, PatchField, GeoMesh>>::New
723 (
724 IOobject
725 (
726 field.name(),
727 casename + runTime2.times()[last_s-1].name(),
728 field.mesh(),
729 IOobject::MUST_READ
730 ),
731 field.mesh()
732 );
733 Lfield.append(std::move(tfld));
734 #endif
735
736 std::cout << std::endl;
737 }
738 else
739 {
740 Info << "######### Reading the Data for " << field.name() << " #########" <<
741 endl;
742 word timename(field.mesh().time().rootPath() + "/" +
743 field.mesh().time().caseName() );
744 timename = timename.substr(0, timename.find_last_of("\\/"));
745 timename = timename + "/" + casename + "processor" + name(Pstream::myProcNo());
746 int last_s = numberOfFiles(casename,
747 "processor" + name(Pstream::myProcNo()) + "/");
748
749 #if defined(OFVER) && (OFVER >= 2212)
750 Lfield.emplace_back
751 (
752 IOobject
753 (
754 field.name(),
755 timename + "/" + name(last_s-1),
756 field.mesh(),
757 IOobject::MUST_READ
758 ),
759 field.mesh()
760 );
761 #else
762 auto tfld =
763 autoPtr<GeometricField<Type, PatchField, GeoMesh>>::New
764 (
765 IOobject
766 (
767 field.name(),
768 timename + "/" + name(last_s-1),
769 field.mesh(),
770 IOobject::MUST_READ
771 ),
772 field.mesh()
773 );
774 Lfield.append(std::move(tfld));
775 #endif
776
777 Info << endl;
778 }
779}
780
781template<class Type, template<class> class PatchField, class GeoMesh>
782void readLastFields(
783 PtrList<GeometricField<Type, PatchField, GeoMesh>>& Lfield,
784 const GeometricField<Type, PatchField, GeoMesh>& field, const fileName casename)
785{
786 int par = 1;
787 M_Assert(ITHACAutilities::check_folder(casename + name(par)) != 0,
788 "No parameter dependent solutions stored into Offline folder");
789
790 while (ITHACAutilities::check_folder(casename + name(par)))
791 {
792 read_last_fields(Lfield, field, casename + name(par) + "/");
793 par++;
794 }
795}
796
797int numberOfFiles(word folder, word MatrixName, word ext)
798{
799 int number_of_files = 0;
800 std::ifstream in;
801 in.open(folder + "/" + MatrixName + name(0) + ext,
802 std::ios::in | std::ios::binary);
803 Info << folder + "/" + MatrixName + name(0) + ext << endl;
804
805 while (in.good())
806 {
807 in.close();
808 number_of_files++;
809 in.open(folder + "/" + MatrixName + name(number_of_files) + ext,
810 std::ios::in | std::ios::binary);
811 }
812
813 in.close();
814 return number_of_files;
815}
816
817template<class Type, template<class> class PatchField, class GeoMesh>
818void exportFields(
819 PtrList<GeometricField<Type, PatchField, GeoMesh>>& field,
820 word folder, word fieldname)
821{
822 ITHACAutilities::createSymLink(folder);
823 Info << "######### Exporting the Data for " << fieldname << " #########" <<
824 endl;
825
826 for (int j = 0; j < field.size() ; j++)
827 {
828 exportSolution(field[j], name(j + 1), folder, fieldname);
829 printProgress(double(j + 1) / field.size());
830 }
831
832 std::cout << std::endl;
833}
834
835template void exportFields(
836 PtrList<GeometricField<scalar, fvPatchField, volMesh>>& field,
837 word folder, word fieldname);
838template void exportFields(
839 PtrList<GeometricField<scalar, fvsPatchField, surfaceMesh>>& field,
840 word folder, word fieldname);
841template void exportFields(
842 PtrList<GeometricField<vector, fvPatchField, volMesh>>& field,
843 word folder, word fieldname);
844template void exportFields(
845 PtrList<GeometricField<tensor, fvPatchField, volMesh>>& field,
846 word folder, word fieldname);
847
848template<class Type, template<class> class PatchField, class GeoMesh>
849void exportSolution(GeometricField<Type, PatchField, GeoMesh>& s,
850 fileName subfolder, fileName folder,
851 word fieldName)
852{
853 if (!Pstream::parRun())
854 {
855 mkDir(folder + "/" + subfolder);
856 ITHACAutilities::createSymLink(folder);
857 GeometricField<Type, PatchField, GeoMesh> act(fieldName, s);
858 fileName fieldname = folder + "/" + subfolder + "/" + fieldName;
859 OFstream os(fieldname);
860 act.writeHeader(os);
861 os << act << endl;
862 }
863 else
864 {
865 mkDir(folder + "/processor" + name(Pstream::myProcNo()) + "/" + subfolder);
866 ITHACAutilities::createSymLink(folder);
867 GeometricField<Type, PatchField, GeoMesh> act(fieldName, s);
868 fileName fieldname = folder + "/processor" + name(Pstream::myProcNo()) + "/" +
869 subfolder + "/" + fieldName;
870 std::cout << fieldname << std::endl;
871 OFstream os(fieldname);
872 act.writeHeader(os);
873 os << act << endl;
874 }
875}
876
877template void exportSolution(
878 GeometricField<scalar, fvPatchField, volMesh>& s,
879 fileName subfolder, fileName folder,
880 word fieldName);
881template void exportSolution(
882 GeometricField<vector, fvPatchField, volMesh>& s,
883 fileName subfolder, fileName folder,
884 word fieldName);
885template void exportSolution(
886 GeometricField<tensor, fvPatchField, volMesh>& s,
887 fileName subfolder, fileName folder,
888 word fieldName);
889template void exportSolution(
890 GeometricField<scalar, fvsPatchField, surfaceMesh>& s,
891 fileName subfolder, fileName folder,
892 word fieldName);
893
894template<class Type, template<class> class PatchField, class GeoMesh>
895void exportSolution(GeometricField<Type, PatchField, GeoMesh>& s,
896 fileName subfolder, fileName folder)
897{
898 if (!Pstream::parRun())
899 {
900 mkDir(folder + "/" + subfolder);
901 ITHACAutilities::createSymLink(folder);
902 fileName fieldname = folder + "/" + subfolder + "/" + s.name();
903 OFstream os(fieldname);
904 s.writeHeader(os);
905 os << s << endl;
906 }
907 else
908 {
909 mkDir(folder + "/processor" + name(Pstream::myProcNo()) + "/" + subfolder);
910 ITHACAutilities::createSymLink(folder);
911 fileName fieldname = folder + "/processor" + name(Pstream::myProcNo()) + "/" +
912 subfolder + "/" + s.name();
913 OFstream os(fieldname);
914 s.writeHeader(os);
915 os << s << endl;
916 }
917}
918
919template void exportSolution(
920 GeometricField<scalar, fvPatchField, volMesh>& s,
921 fileName subfolder, fileName folder);
922template void exportSolution(
923 GeometricField<vector, fvPatchField, volMesh>& s,
924 fileName subfolder, fileName folder);
925template void exportSolution(
926 GeometricField<tensor, fvPatchField, volMesh>& s,
927 fileName subfolder, fileName folder);
928template void exportSolution(
929 GeometricField<scalar, fvsPatchField, surfaceMesh>& s,
930 fileName subfolder, fileName folder);
931
932template void exportSolution(
933 GeometricField<scalar, pointPatchField, pointMesh>& s,
934 fileName subfolder, fileName folder);
935template void exportSolution(
936 GeometricField<vector, pointPatchField, pointMesh>& s,
937 fileName subfolder, fileName folder);
938template void exportSolution(
939 GeometricField<tensor, pointPatchField, pointMesh>& s,
940 fileName subfolder, fileName folder);
941
942void writePoints(pointField points, fileName folder,
943 fileName subfolder)
944{
945 if (!Pstream::parRun())
946 {
947 mkDir(folder + "/" + subfolder);
948 ITHACAutilities::createSymLink(folder);
949 fileName fieldname = folder + "/" + subfolder + "/" + "points";
950 OFstream os(fieldname);
951 os << "FoamFile \n { \n version 2.0; \n format ascii; \n class vectorField; \n location ""1 / polyMesh""; \n object points; \n }"
952 << endl;
953 os << points << endl;
954 }
955 else
956 {
957 mkDir(folder + "/processor" + name(Pstream::myProcNo()) + "/" + subfolder);
958 ITHACAutilities::createSymLink(folder);
959 fileName fieldname = folder + "/processor" + name(Pstream::myProcNo()) + "/" +
960 subfolder + "/" + "points";
961 OFstream os(fieldname);
962 os << "FoamFile \n { \n version 2.0; \n format ascii; \n class vectorField; \n location ""1 / polyMesh""; \n object points; \n }"
963 << endl;
964 os << points << endl;
965 }
966}
967
968void printProgress(double percentage)
969{
970 int val = static_cast<int>(percentage * 100);
971 int lpad = static_cast<int> (percentage * PBWIDTH);
972 int rpad = PBWIDTH - lpad;
973
974 if (Pstream::master())
975 {
976 printf ("\r%3d%% [%.*s%*s]", val, lpad, PBSTR, rpad, "");
977 fflush (stdout);
978 }
979}
980
981template<typename T>
982void save(const List<Eigen::SparseMatrix<T>>& MatrixList, word folder,
983 word MatrixName)
984{
985 mkDir(folder);
986
987 for (int i = 0; i < MatrixList.size(); i++)
988 {
989 word fileName = folder + "/" + MatrixName + name(i) + ".npz";
990 cnpy::save(MatrixList[i], fileName);
991 }
992}
993
994template<typename T>
995void load(List<Eigen::SparseMatrix<T>>& MatrixList, word folder,
996 word MatrixName)
997{
998 int number_of_files = numberOfFiles(folder, MatrixName, ".npz");
999 std::cout << "Reading the Matrix " + folder + "/" + MatrixName << std::endl;
1000 M_Assert(number_of_files != 0,
1001 "Check if the file you are trying to read exists" );
1002 MatrixList.resize(0);
1003 Eigen::SparseMatrix<T> A;
1004
1005 for (int i = 0; i < number_of_files; i++)
1006 {
1007 cnpy::load(A, folder + "/" + MatrixName + name(i) + ".npz");
1008 MatrixList.append(A);
1009 }
1010}
1011
1012template void read_fields(PtrList<volScalarField>& Lfield,
1013 word Name,
1014 fileName casename, int first_snap, int n_snap);
1015template void read_fields(PtrList<volVectorField>& Lfield,
1016 word Name,
1017 fileName casename, int first_snap, int n_snap);
1018template void read_fields(PtrList<volTensorField>& Lfield,
1019 word Name,
1020 fileName casename, int first_snap, int n_snap);
1021template void read_fields(PtrList<surfaceScalarField>& Lfield,
1022 word Name,
1023 fileName casename, int first_snap, int n_snap);
1024template void read_fields(PtrList<surfaceVectorField>& Lfield,
1025 word Name,
1026 fileName casename, int first_snap, int n_snap);
1027template void read_fields(PtrList<volScalarField>& Lfield,
1028 volScalarField& field, fileName casename, int first_snap, int n_snap);
1029template void read_fields(PtrList<volVectorField>& Lfield,
1030 volVectorField& field, fileName casename, int first_snap, int n_snap);
1031template void read_fields(PtrList<volTensorField>& Lfield,
1032 volTensorField& field, fileName casename, int first_snap, int n_snap);
1033template void read_fields(PtrList<surfaceScalarField>& Lfield,
1034 surfaceScalarField& field, fileName casename, int first_snap, int n_snap);
1035template void read_fields(PtrList<surfaceVectorField>& Lfield,
1036 surfaceVectorField& field, fileName casename, int first_snap, int n_snap);
1037template void readMiddleFields(PtrList<volScalarField>& Lfield,
1038 volScalarField& field, fileName casename);
1039template void readMiddleFields(PtrList<volVectorField>& Lfield,
1040 volVectorField& field, fileName casename);
1041template void readMiddleFields(PtrList<volTensorField>& Lfield,
1042 volTensorField& field, fileName casename);
1043template void readMiddleFields(PtrList<surfaceScalarField>&
1044 Lfield, surfaceScalarField& field, fileName casename);
1045template void readMiddleFields(PtrList<surfaceVectorField>&
1046 Lfield, surfaceVectorField& field, fileName casename);
1047template void readConvergedFields(PtrList<volScalarField>& Lfield,
1048 volScalarField& field, fileName casename);
1049template void readConvergedFields(PtrList<volVectorField>& Lfield,
1050 volVectorField& field, fileName casename);
1051template void readConvergedFields(PtrList<volTensorField>& Lfield,
1052 volTensorField& field, fileName casename);
1053template void readConvergedFields(PtrList<surfaceScalarField>&
1054 Lfield, surfaceScalarField& field, fileName casename);
1055template void readConvergedFields(PtrList<surfaceVectorField>&
1056 Lfield, surfaceVectorField& field, fileName casename);
1057
1058template void read_last_fields(PtrList<volScalarField>& Lfield,
1059 const volScalarField& field, const fileName casename);
1060template void read_last_fields(PtrList<volVectorField>& Lfield,
1061 const volVectorField& field, const fileName casename);
1062template void read_last_fields(PtrList<volTensorField>& Lfield,
1063 const volTensorField& field, const fileName casename);
1064template void read_last_fields(PtrList<surfaceScalarField>& Lfield,
1065 const surfaceScalarField& field, const fileName casename);
1066template void read_last_fields(PtrList<surfaceVectorField>& Lfield,
1067 const surfaceVectorField& field, const fileName casename);
1068template void readLastFields(PtrList<volScalarField>& Lfield,
1069 const volScalarField& field, const fileName casename);
1070template void readLastFields(PtrList<volVectorField>& Lfield,
1071 const volVectorField& field, const fileName casename);
1072template void readLastFields(PtrList<volTensorField>& Lfield,
1073 const volTensorField& field, const fileName casename);
1074template void readLastFields(PtrList<surfaceScalarField>&
1075 Lfield, const surfaceScalarField& field, const fileName casename);
1076template void readLastFields(PtrList<surfaceVectorField>&
1077 Lfield, const surfaceVectorField& field, const fileName casename);
1078
1079template<typename T>
1080void exportList(T& list, word folder, word filename)
1081{
1082 mkDir(folder);
1083 word fieldname = folder + filename;
1084 OFstream os(fieldname);
1085
1086 for (int i = 0; i < list.size(); i++)
1087 {
1088 os << list[i] << endl;
1089 }
1090}
1091
1092template void exportList(Field<scalar>& list, word folder,
1093 word filename);
1094template void exportList(Field<vector>& list, word folder,
1095 word filename);
1096template void exportList(Field<tensor>& list, word folder,
1097 word filename);
1098
1099template void save(const List<Eigen::SparseMatrix<double>>& MatrixList,
1100 word folder, word MatrixName);
1101
1102template void load(List<Eigen::SparseMatrix<double>>& MatrixList, word folder,
1103 word MatrixName);
1104
1105}
mesh
Foam::fvMesh & mesh
Definition createMesh.H:47
M_Assert
#define M_Assert(Expr, Msg)
Definition ITHACAassert.H:46
PBWIDTH
#define PBWIDTH
Definition ITHACAregularization.H:71
MAXBUFSIZE
#define MAXBUFSIZE
Definition ITHACAregularization.H:69
PBSTR
#define PBSTR
Definition ITHACAregularization.H:70
ITHACAstream.H
Header file of the ITHACAstream class, it contains the implementation of several methods for input ou...
Foam2Eigen::fvMatrix2Eigen
static void fvMatrix2Eigen(fvMatrix< type_foam_matrix > foam_matrix, type_A &A, type_B &b)
Convert a FvMatrix OpenFOAM matrix (Linear System) into a Eigen Matrix A and a source vector b.
ITHACAparameters
Class for the definition of some general parameters, the parameters must be defined from the file ITH...
Definition ITHACAparameters.H:41
ITHACAparameters::mesh
fvMesh & mesh
Mesh object defined locally.
Definition ITHACAparameters.H:70
ITHACAparameters::getInstance
static ITHACAparameters * getInstance()
Gets an instance of ITHACAparameters, to be used if the instance is already existing.
Definition ITHACAparameters.C:55
T
volScalarField & T
Definition createT.H:46
A
volScalarField & A
Definition createConstants.H:144
Eigen::SliceFromTensor
Matrix< VectorType, Dynamic, Dynamic > SliceFromTensor(Eigen::Tensor< VectorType, 3 > &tensor, label dim, label index1)
Definition EigenFunctions.H:477
ITHACAstream
Namespace for input-output manipulation.
Definition ITHACAstream.C:41
ITHACAstream::writePoints
void writePoints(pointField points, fileName folder, fileName subfolder)
Write points of a mesh to a file.
Definition ITHACAstream.C:942
ITHACAstream::SaveSparseMatrix
void SaveSparseMatrix(Eigen::SparseMatrix< T, Nrows, IND > &m, word folder, word MatrixName)
Export an Eigen sparse matrix into bynary format file.
ITHACAstream::numberOfFiles
int numberOfFiles(word folder, word MatrixName, word ext)
Count the number of files with a certain name prefix, used by reading functions of list of matrices.
Definition ITHACAstream.C:797
ITHACAstream::readConvergedFields
void readConvergedFields(PtrList< GeometricField< Type, PatchField, GeoMesh > > &Lfield, GeometricField< Type, PatchField, GeoMesh > &field, fileName casename)
Function to read a list of volVectorField from name of the field including only converged snapshots.
Definition ITHACAstream.C:659
ITHACAstream::save
void save(const List< Eigen::SparseMatrix< T > > &MatrixList, word folder, word MatrixName)
Definition ITHACAstream.C:982
ITHACAstream::exportFields
void exportFields(PtrList< GeometricField< Type, PatchField, GeoMesh > > &field, word folder, word fieldname)
Function to export a scalar of vector field.
Definition ITHACAstream.C:818
ITHACAstream::load
void load(List< Eigen::SparseMatrix< T > > &MatrixList, word folder, word MatrixName)
Definition ITHACAstream.C:995
ITHACAstream::exportSolution
void exportSolution(GeometricField< Type, PatchField, GeoMesh > &s, fileName subfolder, fileName folder, word fieldName)
Export a field to file in a certain folder and subfolder.
Definition ITHACAstream.C:849
ITHACAstream::read_last_fields
void read_last_fields(PtrList< GeometricField< Type, PatchField, GeoMesh > > &Lfield, const GeometricField< Type, PatchField, GeoMesh > &field, const fileName casename)
Function to read a list of fields from the name of the field including only the last snapshot.
Definition ITHACAstream.C:695
ITHACAstream::readMiddleFields
void readMiddleFields(PtrList< GeometricField< Type, PatchField, GeoMesh > > &Lfield, GeometricField< Type, PatchField, GeoMesh > &field, fileName casename)
Funtion to read a list of volVectorField from name of the field including all the intermediate snapsh...
Definition ITHACAstream.C:643
ITHACAstream::exportMatrix
void exportMatrix(Eigen::Matrix< T, -1, dim > &matrix, word Name, word type, word folder)
Export the reduced matrices in numpy (type=python), matlab (type=matlab) and txt (type=eigen) format ...
Definition ITHACAstream.C:55
ITHACAstream::exportList
void exportList(T &list, word folder, word filename)
Export a list to file.
Definition ITHACAstream.C:1080
ITHACAstream::SaveDenseMatrix
void SaveDenseMatrix(MatrixType &Matrix, word folder, word MatrixName)
Save a dense matrix to a binary format file.
Definition ITHACAstream.H:655
ITHACAstream::exportFvMatrix
void exportFvMatrix(fvMatrix< Type > &Matrix, word folder, word MatrixName)
Export an Fv Matrix to folder together with its source term.
Definition ITHACAstream.C:44
ITHACAstream::exportVector
void exportVector(Eigen::VectorXd &vector, word Name, word type, word folder)
Definition ITHACAstream.C:255
ITHACAstream::readMatrix
List< Eigen::MatrixXd > readMatrix(word folder, word mat_name)
Read a three dimensional matrix from a txt file in Eigen format.
Definition ITHACAstream.C:372
ITHACAstream::read_fields
void read_fields(PtrList< GeometricField< Type, PatchField, GeoMesh > > &Lfield, word Name, fileName casename, int first_snap, int n_snap)
Function to read a list of fields from the name of the field and casename.
Definition ITHACAstream.C:454
ITHACAstream::readLastFields
void readLastFields(PtrList< GeometricField< Type, PatchField, GeoMesh > > &Lfield, const GeometricField< Type, PatchField, GeoMesh > &field, const fileName casename)
Funtion to read a list of volVectorField from name of the field including only the last snapshots.
Definition ITHACAstream.C:782
ITHACAstream::exportTensor
void exportTensor(Eigen::Tensor< T, 3 > tensor, word Name, word type, word folder)
Export the reduced tensor in numpy (tipo=python), matlab (tipo=matlab) and txt (tipo=eigen) format.
Definition ITHACAstream.C:264
ITHACAstream::printProgress
void printProgress(double percentage)
Print progress bar given the percentage.
Definition ITHACAstream.C:968
ITHACAutilities::check_folder
bool check_folder(word folder)
Checks if a folder exists.
Definition ITHACAsystem.C:89
ITHACAutilities::createSymLink
void createSymLink(word folder)
Creates symbolic links to 0, system and constant.
Definition ITHACAsystem.C:30
cnpy::save
void save(const Eigen::Matrix< T, -1, dim > &mat, const std::string fname)
cnpy::load
Eigen::Matrix< T, -1, dim > load(Eigen::Matrix< T, -1, dim > &mat, const std::string fname, std::string order="rowMajor")
i
label i
Definition pEqn.H:46
Generated by   doxygen 1.11.0