GCC Code Coverage Report

Directory: ./
File: Solver/linearProblemReducedSteklov.hpp
Date: 2024-04-14 07:32:34
Exec Total Coverage
Lines: 6 16 37.5%
Branches: 1 2 50.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:
13 //
14
15 #ifndef _LINEARPROBLEMREDUCEDSTEKLOV_HPP
16 #define _LINEARPROBLEMREDUCEDSTEKLOV_HPP
17
18 // System includes
19
20 // External includes
21
22 // Project includes
23 #include "Solver/linearProblem.hpp"
24 #include "FiniteElement/elementField.hpp"
25 #include "Solver/reducedSteklov.hpp"
26 #include "Geometry/curvatures.hpp"
27
28 namespace felisce {
29 /*!
30 \class LinearProblemReducedSteklov
31 \authors M.Aletti
32 \date 2015
33 \brief Commmon interface for different linear problems to be used for steklov reduction
34 */
35 template <class volumeProblem>
36 class ReducedSteklov;
37
38 template<class volumeProblem> class LinearProblemReducedSteklov :
39 public volumeProblem {
40
41 friend class ReducedSteklov<volumeProblem>;
42
43 public:
44 using LinearProblem::sequential;
45 using LinearProblem::parallel;
46 enum OptionForSteklov {
47 LOAD_IT= 0,
48 COMPUTE_IT= 1,
49 COMPUTE_AND_SAVE_IT= 2
50 };
51 enum imgType { dirichlet, neumann };
52 // ====================================================================
53 // Simple constructor/destructor =
54 // ====================================================================
55 LinearProblemReducedSteklov();
56 ~LinearProblemReducedSteklov() override;
57 // ====================================================================
58 // These functions let the model build the steklov oeprator =
59 // and use it. =
60 // Since the internal mapping towards the interface is complicated =
61 // the idea is to pass a petscVector defined in the volume =
62 // and do all the crazy mapping into the function =
63 // ====================================================================
64 void assembleFullRankSteklov();
65 void assembleLowRankSteklov(felInt imesh);
66 // ====================================================================
67 // function that read a given petscVector coming =
68 // from the NS problem =
69 // and save it in the corresponding dof of =
70 // the petscVector <whereToSave> =
71 // ====================================================================
72 // scalar case
73 void readerInterface( const std::map<int, std::map<int,int> > & externalMap,
74 PetscVector& seqExternalVec,
75 const std::map<int,int> &mapBetweenLabels,
76 std::string whereToSave,
77 felInt iUnkWTW=0, /* where to write */
78 felInt iCompWTW=0); /* where to write */
79 void readerInterface( const std::map<int, std::map<int,std::map<int,int> > > & externalMap,
80 PetscVector& seqExternalVec,
81 const std::map<int,int> &mapBetweenLabels,
82 std::string whereToSave,
83 felInt iUnkWTW, /*where to write*/
84 felInt iUCR);/*where to read (iUnknownComp of the other linearProblem)*/
85 // ===================================================================
86 // export functions
87 // ===================================================================
88 void exportOutputSteklov(FelisceTransient::Pointer fstransient,int iFixedPointIteration);
89 void exportSteklovInputOutput( PetscVector& Input, PetscVector& Output, FelisceTransient::Pointer fstransient, int iFixedPointIteration);
90 ✗ inline void exportAllEig(felInt imesh) {m_reducedSteklov->exportAllEig(imesh); }
91
92 virtual void initializeDofBoundaryAndBD2VolMaps()=0;
93 328 virtual void derivedProblemAssemble( FlagMatrixRHS /*flag*/) {};
94 ✗ virtual double userInletPressure() {return 1.;};
95
96 2560 inline DofBoundary& dofBD(std::size_t iBD) { return m_dofBD[iBD]; };
97
98 protected:
99 ✗ virtual void computeTheConstantResponse() {};
100
101 void readerInterfaceComponentWise( const std::map<int,std::map<int,int> > & externalMap,
102 PetscVector& seqExternalVec,
103 const std::map<int,int> & mapBetweenLabels,
104 std::string whereToSave,
105 felInt iUnknownWTW,
106 felInt iCompWTW);
107
108 virtual void userSetNullSpace(PetscVector& v, int k);
109 // ====================================================================
110 // This function is used to apply the boundary condition =
111 // at the interface =
112 // ====================================================================
113 ✗ void userElementComputeNaturalBoundaryCondition( const std::vector<Point*>& /*elemPoint*/, const std::vector<felInt>& /*elemIdPoint*/,felInt& /*iel*/, int /*label*/) override{}
114 // ====================================================================
115 // This function allocate the petsc matrix to store =
116 // the steklov operator and apply steklov allow you to use it =
117 // ====================================================================
118 void applySteklov( PetscVector& in, PetscVector& out, FelisceTransient::Pointer fstransient, int nfp);
119 // ====================================================================
120 // This flag is std::set to false in the constructor =
121 // it is used in the boundary condition loop =
122 // to read data from the correct std::vector in case of the computation =
123 // of the Steklov operator =
124 // it is std::set to true in the meanwhile that we are computing S-Operator=
125 // ====================================================================
126 bool m_computingSteklov;
127 ReducedSteklov<volumeProblem>* m_reducedSteklov;
128 PetscMatrix m_fullSteklov; // - the matrix (parallel and dense)
129
130 //! This method assembles the stiffness matrix defined with the laplacian operator on the boundary
131 PetscMatrix assembleLaplacianBoundary();
132 //! This method assembles the mass matrix defined on the boundary
133 PetscMatrix assembleMassBoundary();
134
135 // This std::set of virtual functions should be reimplemented by the derived linear problem since they differ from scalar
136 // to std::vector, for instance.
137 ✗ virtual void initPerETMASS(){};
138 ✗ virtual void initPerETLAP(){};
139
140 ✗ virtual void updateFE(const std::vector<Point*>& /*elemPoint*/, const std::vector<felInt>& /*elemIdPoint*/){};
141
142 ✗ virtual void massMatrixComputer(felInt /*ielSupportDof*/){};
143 ✗ virtual void laplacianMatrixComputer(felInt /*ielSupportDof*/){};
144
145 void setValueMatrixBD(felInt ielSupportDof);
146 std::vector<felInt> m_globPosColumn;
147 std::vector<felInt> m_globPosRow;
148 std::vector<double> m_matrixValues;
149
150 Curvatures m_curv;
151 ElementField m_normalField;
152
153 void assembleVolumeSystem( FlagMatrixRHS flagMatrixRHS);
154 virtual void useSteklovDataBegin();
155 virtual void useSteklovDataEnd();
156 std::vector<double> m_tmpValues;
157 void setSteklovData( PetscVector& bdInput );
158 PetscVector getSteklovImg();
159
160 std::vector<int> m_interfaceLabels;
161
162 PetscMatrix m_aux;
163
164 imgType m_imgType;
165
166 bool m_useSteklovData = false;
167
168 ChronoInstance::Pointer m_chronoRS;
169
170 // Dof boundary utilities
171 std::vector<DofBoundary> m_dofBD = std::vector<DofBoundary>(1, DofBoundary());
172
173 public:
174 double comulatedSolveTime(){
175 return m_chronoRS->diff_cumul();
176 }
177
178 ✗ ChronoInstance::Pointer chronoRS(){return m_chronoRS;}
179
180 4 void initChronoRS(){
181
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
 4 if (m_chronoRS == nullptr) {
182 4 m_chronoRS = felisce::make_shared<ChronoInstance>();
183 }
184 4 }
185 };
186 }
187
188 #include "linearProblemReducedSteklov.hxx"
189
190 #endif
191