GCC Code Coverage Report


Directory: ./
File: Solver/eigenProblemPOD.hpp
Date: 2024-04-14 07:32:34
Exec Total Coverage
Lines: 0 10 0.0%
Branches: 0 10 0.0%

Line Branch Exec Source
1 // ______ ______ _ _ _____ ______
2 // | ____| ____| | (_)/ ____| | ____|
3 // | |__ | |__ | | _| (___ ___| |__
4 // | __| | __| | | | |\___ \ / __| __|
5 // | | | |____| |____| |____) | (__| |____
6 // |_| |______|______|_|_____/ \___|______|
7 // Finite Elements for Life Sciences and Engineering
8 //
9 // License: LGL2.1 License
10 // FELiScE default license: LICENSE in root folder
11 //
12 // Main authors: E.Schenone D. Lombardi
13 //
14
15 #ifndef _EIGENPROBLEMPOD_HPP
16 #define _EIGENPROBLEMPOD_HPP
17
18 // System includes
19
20 // External includes
21
22 // Project includes
23 #include "Core/felisceParam.hpp"
24 #include "Solver/eigenProblem.hpp"
25
26 namespace felisce {
27 /*!
28 \class EigenProblemPOD
29 \authors E.Schenone D.Lombardi
30 \date 20/09/2013
31 \brief Manage functions of a POD eigenvalues problem.
32 */
33
34 class EigenProblemPOD:
35 public EigenProblem {
36 public:
37 // Constructor
38 //============
39 EigenProblemPOD();
40 ~EigenProblemPOD() override;
41
42 void initialize(const GeometricMeshRegion::Pointer& mesh, FelisceTransient::Pointer fstransient, MPI_Comm& comm);
43
44 //Define Physical Variable associate to the problem
45 //=================================================
46 void initPerElementType() override;
47
48 // Assemble Matrix
49 //================
50 void computeElementArray(const std::vector<Point*>& elemPoint, const std::vector<felInt>& elemIdPoint, felInt& iel, FlagMatrixRHS flagMatrixRHS = FlagMatrixRHS::matrix_and_rhs) override;
51 void readSnapshot(IO& io);
52 void computeAutocorr();
53
54 //! Visualize modes as a variable
55 void writeMode(const int iter=0);
56 void writeEnsightSolution(const int iter);
57
58 //! Solve eigenproblem to build rom basis
59 void buildSolver() override;
60 virtual void solve();
61
62 // Reduced model functions
63 //========================
64 // Initialization
65 void readData(IO& io) override;
66 void getFiberDirection(felInt iel, int iUnknown, std::vector<double>& elemFiber);
67 void initializeROM();
68
69 void readEnsightFile(double* vec, std::string varName, int idIter, felInt size);
70 void readEigenValueFromFile();
71 void computeGamma();
72 void computeMatrixZeta();
73
74 void computeTensor();
75 int secondOrderGlobalIndex(int& i, int& j);
76 int thirdOrderGlobalIndex(int& i, int& j, int& k);
77 int fourthOrderGlobalIndex(int& i, int& j, int& k, int& h);
78
79 void initializeSolution();
80 void setFhNf0(std::vector<double>& valuef0);
81 void setIapp(std::vector<double>& iApp, int& idS);
82
83 // Compute solution
84 void checkSolution(PetscVector& solution, PetscVector& refSol, double & error);
85 void checkBasis();
86 void checkOrthonormality();
87
88
89 //Update functions
90 void MGS(std::vector<PetscVector>&);
91 void updateBeta();
92
93
94 // Projection functions from FE space to RO one and viceversa
95 void projectOnBasis(PetscVector& vectorDof, double* vectorBasis, bool flagMass);
96 void projectOnDof(double* vectorBasis, PetscVector& vectorDof, int size);
97 void scalarProduct(PetscVector& firstVec, PetscVector& secondVec, double& prod,bool flagMass);
98
99
100 // Print&Auxiliary functions
101 //==========================
102 void viewPOD(double* vToPrint, int vSize, std::ofstream outStream);
103 void viewPOD(double* vToPrint, int vSize, std::string fName);
104
105
106 // Acces functions
107 //================
108 double eigenValueFirst() {
109 return m_eigenValue[0];
110 }
111 double* eigenValue() {
112 return m_eigenValue;
113 }
114 double* beta() {
115 return m_beta;
116 }
117 double* mu() {
118 return m_mu;
119 }
120 double* gamma() {
121 return m_gamma;
122 }
123 double* eta() {
124 return m_eta;
125 }
126 int dimRomBasis() {
127 return m_dimRomBasis;
128 }
129
130 template <class Templ_functor> void evalFunctionOnDof(const Templ_functor& func, std::vector<double>& array);
131
132 // writing ECG
133 virtual void initializeECG() {}
134 virtual void writeECG(int /*iteration*/) {}
135 virtual void updateEcgOperator() {}
136
137 protected:
138 CurrentFiniteElement* m_fePotTransMemb;
139 felInt m_ipotTransMemb;
140 felInt m_ipotExtraCell;
141
142 std::vector<PetscVector> m_snapshot;
143 PetscMatrix m_autocorr;
144
145 PetscVector m_U_0;
146 PetscVector m_U_0_seq;
147 PetscVector m_W_0;
148 PetscVector m_FhNf0;
149 ElementField m_elemFieldU0;
150 double* m_eigenValue;
151 std::vector<PetscVector> m_basis;
152 std::vector<PetscVector> m_basisPOD;
153 std::vector<PetscVector> m_eigenvector;
154 std::vector<PetscVector> m_zeta;
155
156 felInt m_dimRomBasis;
157 felInt m_numSnapshot;
158 felInt m_samplingFreq;
159
160 double* m_beta;
161 double* m_mu;
162 double* m_gamma;
163 double* m_eta;
164 double* m_modeMean;
165 double* m_fiber;
166 double** m_source;
167 bool m_solutionInitialized;
168
169 int m_idK;
170 int m_idG;
171 int m_idGs;
172
173
174 double* m_tensorCs;
175 double* m_tensorE;
176 double* m_tensorT;
177 double* m_tensorTs;
178 double* m_tensorY;
179 int m_size2;
180 int m_size3;
181 int m_size4;
182 int m_numSource;
183 int m_problem;
184
185
186 ElementField m_elemFieldFhNf0;
187
188 };
189
190
191 template <class Templ_functor> void EigenProblemPOD::evalFunctionOnDof(const Templ_functor& func, std::vector<double>& array) {
192 Point pt;
193 array.resize(m_numDof);
194 felInt position;
195 for( felInt i = 0; i < m_numDof; i++) {
196 findCoordinateWithIdDof(i,pt);
197 position = i;
198 AOApplicationToPetsc(m_ao,1,&position);
199 array[position] = func(pt.x(),pt.y(),pt.z());
200 }
201 }
202
203
204 }
205
206 #endif
207