Loading...
Searching...
No Matches
Foam2Eigen.H
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 License
13 This file is part of ITHACA-FV
14 ITHACA-FV is free software: you can redistribute it and/or modify
15 it under the terms of the GNU Lesser General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18 ITHACA-FV is distributed in the hope that it will be useful,
19 but WITHOUT ANY WARRANTY; without even the implied warranty of
20 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21 GNU Lesser General Public License for more details.
22 You should have received a copy of the GNU Lesser General Public License
23 along with ITHACA-FV. If not, see <http://www.gnu.org/licenses/>.
24 Class
25 foam2eigen
26 Description
27 Converts foam objects into Eigen and viceversa
28 SourceFiles
29 foam2eigen.C
30 \*---------------------------------------------------------------------------*/
31
36
37#ifndef Foam2Eigen_H
38#define Foam2Eigen_H
39
40#include "fvCFD.H"
41#include "IOmanip.H"
42#include "ITHACAassert.H"
43#include "ITHACAutilities.H"
44#include <tuple>
45#include <sys/stat.h>
46#pragma GCC diagnostic push
47#pragma GCC diagnostic ignored "-Wold-style-cast"
48#include <Eigen/Eigen>
49#pragma GCC diagnostic pop
50// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51
52
53/*---------------------------------------------------------------------------*\
54 Class foam2eigen Declaration
55\*---------------------------------------------------------------------------*/
56
59{
60 private:
61
62 public:
63
64 //----------------------------------------------------------------------
86 template <class type_foam_matrix, class type_A, class type_B>
87 static void fvMatrix2Eigen(fvMatrix<type_foam_matrix> foam_matrix, type_A& A,
88 type_B& b);
89
90 //----------------------------------------------------------------------
106 template <class type_foam_matrix, class type_A>
107 static void fvMatrix2EigenM(fvMatrix<type_foam_matrix>& foam_matrix, type_A& A);
108
109 //----------------------------------------------------------------------
123 template <class type_foam_matrix, class type_B>
124 static void fvMatrix2EigenV(fvMatrix<type_foam_matrix>& foam_matrix, type_B& b);
125
126
127 //----------------------------------------------------------------------
141 template<class Type, template<class> class PatchField, class GeoMesh>
142 static Eigen::MatrixXd PtrList2Eigen(
143 PtrList<GeometricField<Type, PatchField, GeoMesh>>& fields,
144 label Nfields = -1);
145
146 //----------------------------------------------------------------------
156 template<template<class> class PatchField, class GeoMesh>
157 static Eigen::VectorXd field2Eigen(GeometricField<tensor, PatchField, GeoMesh>&
158 field);
159
160 //----------------------------------------------------------------------
170 template<template<class> class PatchField, class GeoMesh>
171 static Eigen::VectorXd field2Eigen(GeometricField<vector, PatchField, GeoMesh>&
172 field);
173
174 //----------------------------------------------------------------------
184 template<template<class> class PatchField, class GeoMesh>
185 static Eigen::VectorXd field2Eigen(GeometricField<scalar, PatchField, GeoMesh>&
186 field);
187
188 //----------------------------------------------------------------------
198 template<template<class> class PatchField, class GeoMesh>
199 static Eigen::Map<Eigen::MatrixXd> field2EigenMap(
200 GeometricField<scalar, PatchField, GeoMesh>&
201 field);
202
203 //----------------------------------------------------------------------
214 template<template<class> class PatchField, class GeoMesh>
215 static Eigen::Map<Eigen::MatrixXd> field2EigenMapBC(
216 GeometricField<scalar, PatchField, GeoMesh>&
217 field, int BC_index);
218
219 //----------------------------------------------------------------------
229 template<class Type, class GeoMesh>
230 static Eigen::VectorXd field2Eigen(const
231 DimensionedField<Type, GeoMesh>& field);
232
233 //----------------------------------------------------------------------
243 template<class Type>
244 static Eigen::VectorXd field2Eigen(const Field<Type>& field);
245
246 //----------------------------------------------------------------------
257 template<template<class> class PatchField, class GeoMesh>
258 static List<Eigen::VectorXd> field2EigenBC(
259 GeometricField<scalar, PatchField, GeoMesh>& field);
260
261 //----------------------------------------------------------------------
272 template<template<class> class PatchField, class GeoMesh>
273 static List<Eigen::VectorXd> field2EigenBC(
274 GeometricField<vector, PatchField, GeoMesh>& field);
275
276 //----------------------------------------------------------------------
287 template<template<class> class PatchField, class GeoMesh>
288 static List<Eigen::VectorXd> field2EigenBC(
289 GeometricField<tensor, PatchField, GeoMesh>& field);
290
291
292 //----------------------------------------------------------------------
304 template<template<class> class PatchField, class GeoMesh>
305 static List<Eigen::MatrixXd> PtrList2EigenBC(
306 PtrList<GeometricField<scalar, PatchField, GeoMesh>>& fields,
307 label Nfields = -1);
308
309 //----------------------------------------------------------------------
321 template<template<class> class PatchField, class GeoMesh>
322 static List<Eigen::MatrixXd> PtrList2EigenBC(
323 PtrList<GeometricField<vector, PatchField, GeoMesh>>& fields,
324 label Nfields = -1);
325
326 //----------------------------------------------------------------------
338 template<template<class> class PatchField, class GeoMesh>
339 static List<Eigen::MatrixXd> PtrList2EigenBC(
340 PtrList<GeometricField<tensor, PatchField, GeoMesh>>& fields,
341 label Nfields = -1);
342
343 //----------------------------------------------------------------------
355 template<template<class> class PatchField, class GeoMesh>
356 static GeometricField<scalar, PatchField, GeoMesh> Eigen2field(
357 GeometricField<scalar, PatchField, GeoMesh>& field,
358 Eigen::VectorXd& eigen_vector, bool correctBC = true);
359
360 //----------------------------------------------------------------------
372 template<template<class> class PatchField, class GeoMesh>
373 static GeometricField<vector, PatchField, GeoMesh> Eigen2field(
374 GeometricField<vector, PatchField, GeoMesh>& field,
375 Eigen::VectorXd& eigen_vector, bool correctBC = true);
376
377 //----------------------------------------------------------------------
390 template<template<class> class PatchField, class GeoMesh>
391 static GeometricField<scalar, PatchField, GeoMesh> Eigen2field(
392 GeometricField<scalar, PatchField, GeoMesh>& field,
393 Eigen::VectorXd& eigen_vector, List<Eigen::VectorXd>& eigen_vector_boundary);
394
395 //----------------------------------------------------------------------
408 template<template<class> class PatchField, class GeoMesh>
409 static GeometricField<vector, PatchField, GeoMesh> Eigen2field(
410 GeometricField<vector, PatchField, GeoMesh>& field,
411 Eigen::VectorXd& eigen_vector, List<Eigen::VectorXd>& eigen_vector_boundary);
412
413
414 //----------------------------------------------------------------------
426 template<template<class> class PatchField, class GeoMesh>
427 static GeometricField<tensor, PatchField, GeoMesh> Eigen2field(
428 GeometricField<tensor, PatchField, GeoMesh>& field,
429 Eigen::VectorXd& eigen_vector, bool correctBC = true);
430
431 //----------------------------------------------------------------------
442 template <class Type>
443 static Field<Type> Eigen2field(
444 Field<Type>& field, Eigen::MatrixXd& matrix, bool correctBC = true);
445
446 //----------------------------------------------------------------------
459 template<class Type, template<class> class PatchField, class GeoMesh>
460 static std::tuple<Eigen::MatrixXd, Eigen::VectorXd> projectFvMatrix(
461 fvMatrix<Type>& matrix,
462 PtrList<GeometricField<Type, PatchField, GeoMesh>>& modes, label Nmodes);
463
464
465 //----------------------------------------------------------------------
491 template<class Type, template<class> class PatchField, class GeoMesh>
492 static Eigen::VectorXd projectField(
493 GeometricField<Type, PatchField, GeoMesh>& field,
494 PtrList<GeometricField<Type, PatchField, GeoMesh>>& modes,
495 label Nmodes);
496
497 //----------------------------------------------------------------------
509 template<class Type, template<class> class PatchField, class GeoMesh>
510 static Eigen::MatrixXd MassMatrix(
511 PtrList<GeometricField<Type, PatchField, GeoMesh>>& modes, label Nmodes);
512
513 //----------------------------------------------------------------------
525 template<class Type>
526 static std::tuple<List<Eigen::SparseMatrix<double>>, List<Eigen::VectorXd>>
527 LFvMatrix2LSM(PtrList<fvMatrix<Type>>& MatrixList);
528
529 //--------------------------------------------------------------------------
538 template <class type_matrix>
539 static Eigen::Matrix<type_matrix, Eigen::Dynamic, Eigen::Dynamic>
540 List2EigenMatrix ( List<type_matrix> list );
541
542 //--------------------------------------------------------------------------
551 template <class type_matrix>
552 static List<type_matrix> EigenMatrix2List (
553 Eigen::Matrix<type_matrix, Eigen::Dynamic, Eigen::Dynamic> matrix );
554
555 //----------------------------------------------------------------------
564 template <class type_list>
565 static Eigen::MatrixXd field2Eigen(const List<type_list>& list);
566};
567#endif
Implementation of the assert function for ITHACA-FV.
Header file of the ITHACAutilities namespace.
Class to converts OpenFOAM objects into Eigen and viceversa.
Definition Foam2Eigen.H:59
static GeometricField< scalar, PatchField, GeoMesh > Eigen2field(GeometricField< scalar, PatchField, GeoMesh > &field, Eigen::VectorXd &eigen_vector, bool correctBC=true)
Convert a vector in Eigen format into an OpenFOAM scalar GeometricField.
Definition Foam2Eigen.C:517
static Eigen::VectorXd field2Eigen(GeometricField< tensor, PatchField, GeoMesh > &field)
Convert a vector OpenFOAM field into an Eigen Vector.
Definition Foam2Eigen.C:40
static List< Eigen::VectorXd > field2EigenBC(GeometricField< scalar, PatchField, GeoMesh > &field)
Convert an OpenFOAM scalar field to a List of Eigen Vectors, one for each boundary.
Definition Foam2Eigen.C:244
static Eigen::Map< Eigen::MatrixXd > field2EigenMap(GeometricField< scalar, PatchField, GeoMesh > &field)
Convert a scalar OpenFOAM field into an Eigen Vector.
Definition Foam2Eigen.C:97
static Field< Type > Eigen2field(Field< Type > &field, Eigen::MatrixXd &matrix, bool correctBC=true)
Converts a matrix in Eigen format into an OpenFOAM Field.
static std::tuple< List< Eigen::SparseMatrix< double > >, List< Eigen::VectorXd > > LFvMatrix2LSM(PtrList< fvMatrix< Type > > &MatrixList)
Convert a PtrList of OpenFOAM fvMatrix into a tuple of lists of Eigen Sparse Matrices and source vect...
static Eigen::Matrix< type_matrix, Eigen::Dynamic, Eigen::Dynamic > List2EigenMatrix(List< type_matrix > list)
Convert a Foam List into an Eigen matrix with one column.
static Eigen::VectorXd projectField(GeometricField< Type, PatchField, GeoMesh > &field, PtrList< GeometricField< Type, PatchField, GeoMesh > > &modes, label Nmodes)
Perform the projection of an OpenFOAM field onto a set of modes using the Eigen matrix multiplication...
static Eigen::MatrixXd MassMatrix(PtrList< GeometricField< Type, PatchField, GeoMesh > > &modes, label Nmodes)
Obtain the Mass Matrix from a list of reduced basis.
static List< type_matrix > EigenMatrix2List(Eigen::Matrix< type_matrix, Eigen::Dynamic, Eigen::Dynamic > matrix)
Convert an Eigen matrix with one column into a Foam List.
static void fvMatrix2EigenV(fvMatrix< type_foam_matrix > &foam_matrix, type_B &b)
Convert a ldu OpenFOAM matrix into a source vector b.
static List< Eigen::MatrixXd > PtrList2EigenBC(PtrList< GeometricField< scalar, PatchField, GeoMesh > > &fields, label Nfields=-1)
Convert an OpenFOAM scalar field to a List of Eigen Vectors, one for each boundary.
Definition Foam2Eigen.C:268
static Eigen::MatrixXd PtrList2Eigen(PtrList< GeometricField< Type, PatchField, GeoMesh > > &fields, label Nfields=-1)
Convert a PtrList of snapshots to Eigen matrix (only internal field)
Definition Foam2Eigen.C:649
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.
static std::tuple< Eigen::MatrixXd, Eigen::VectorXd > projectFvMatrix(fvMatrix< Type > &matrix, PtrList< GeometricField< Type, PatchField, GeoMesh > > &modes, label Nmodes)
Obtain the Mass Matrix from a list of reduced basis.
static Eigen::VectorXd field2Eigen(const DimensionedField< Type, GeoMesh > &field)
Convert an OpenFOAM field into an Eigen Vector.
static Eigen::Map< Eigen::MatrixXd > field2EigenMapBC(GeometricField< scalar, PatchField, GeoMesh > &field, int BC_index)
Map a scalar OpenFOAM field boundary into an Eigen Matrix.
Definition Foam2Eigen.C:105
static Eigen::VectorXd field2Eigen(const Field< Type > &field)
Convert an OpenFOAM field into an Eigen Vector.
static Eigen::MatrixXd field2Eigen(const List< type_list > &list)
Function to convert an OpenFOAM list to an Eigen Matrix.
static void fvMatrix2EigenM(fvMatrix< type_foam_matrix > &foam_matrix, type_A &A)
Convert a ldu OpenFOAM matrix into a Eigen Matrix A.
volScalarField & A