GCC Code Coverage Report

Directory: ./
File: Model/SSTModel.cpp
Date: 2024-04-14 07:32:34
Exec Total Coverage
Lines: 0 88 0.0%
Branches: 0 42 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
13 //
14
15 // System includes
16
17 // External includes
18
19 // Project includes
20 #include "Model/SSTModel.hpp"
21
22 namespace felisce {
23 ✗ SSTModel::SSTModel():Model() {
24 ✗ m_name = "SST";
25 }
26
27 ✗ SSTModel::~SSTModel() {
28 ✗ m_solExtrapolate.destroy();
29 }
30
31 ✗ void SSTModel::initializeDerivedModel() {
32 ✗ m_bdfNS.defineOrder(FelisceParam::instance().orderBdfNS);
33 ✗ m_linearProblem[0]->initializeTimeScheme(&m_bdfNS);
34 }
35
36 ✗ void SSTModel::preAssemblingMatrixRHS(std::size_t iProblem) {
37 ✗ if(iProblem == 0) {
38 ✗ m_U_0.duplicateFrom(m_linearProblem[iProblem]->vector());
39 ✗ m_U_0.set(0.0);
40
41 ✗ m_solExtrapolate.duplicateFrom(m_linearProblem[iProblem]->vector());
42 ✗ m_solExtrapolate.set( 0.);
43
44 ✗ const auto num_dofs = m_linearProblem[iProblem]->numDofPerUnknown(0);
45 ✗ felInt idGlobalDof[num_dofs];
46 ✗ double valueByDofU_0[num_dofs];
47
48 ✗ for ( felInt i = 0; i < num_dofs; i++) {
49 ✗ idGlobalDof[i] = i;
50 ✗ valueByDofU_0[i] = 0.;
51 }
52
53 //Use global mapping for parallel computation (in sequential: (before mapping idGlobalDof) == (after mapping idGlobalDof)).
54 ✗ AOApplicationToPetsc(m_linearProblem[iProblem]->ao(),num_dofs,idGlobalDof);
55
56 // "add" values in Petsc vectors.
57 ✗ m_U_0.setValues(num_dofs,idGlobalDof,valueByDofU_0,INSERT_VALUES);
58
59 //assemble _U_0.
60 ✗ m_U_0.assembly();
61
62 //Initialize solution for solver.
63 ✗ m_linearProblem[iProblem]->solution().copyFrom(m_U_0);
64
65 //First assembly loop in iteration 0 to build static matrix.
66 ✗ m_linearProblem[iProblem]->assembleMatrixRHS(MpiInfo::rankProc());
67 }
68 }
69
70 ✗ void SSTModel::postAssemblingMatrixRHS(std::size_t iProblem) {
71 ✗ if(iProblem == 0) {
72 // _Matrix = _A + _Matrix (add static matrix to the dynamic matrix to build
73 // complete matrix of the system
74 ✗ m_linearProblem[iProblem]->addMatrixRHS();
75 }
76 }
77
78
79
80 ✗ void SSTModel::forward() {
81 // Write solution for postprocessing (if required)
82 ✗ writeSolution();
83
84 // Advance time step.
85 ✗ updateTime();
86
87 // Print time information
88 ✗ printNewTimeIterationBanner();
89
90
91 /*----------------------------
92 Navier Stokes equation
93 ----------------------------*/
94
95 // update the bdf scheme
96 ✗ if ( m_fstransient->iteration == 1) {
97 // I use m_sol because from Model::setInitialCondition, I have a different _U_0
98 ✗ m_bdfNS.initialize( m_linearProblem[0]->solution() );
99 ✗ m_bdfNS.extrapolate(m_solExtrapolate);
100 ✗ m_linearProblem[0]->initExtrapol(m_solExtrapolate);
101 } else {
102 ✗ m_bdfNS.update(m_linearProblem[0]->solution());
103 ✗ m_bdfNS.extrapolate(m_solExtrapolate);
104 ✗ m_linearProblem[0]->updateExtrapol(m_solExtrapolate);
105 }
106
107 ✗ m_bdfNS.computeRHSTime(m_fstransient->timeStep);
108 ✗ m_linearProblem[0]->gatherVectorBeforeAssembleMatrixRHS();
109
110 //Assembly loop on elements.
111 ✗ m_linearProblem[0]->assembleMatrixRHS(MpiInfo::rankProc());
112
113 //Specific operations before solve the system.
114 ✗ postAssemblingMatrixRHS(0);
115
116 //Apply boundary conditions.
117 ✗ m_linearProblem[0]->finalizeEssBCTransient();
118 ✗ m_linearProblem[0]->applyBC(FelisceParam::instance().essentialBoundaryConditionsMethod, MpiInfo::rankProc());
119
120 //Solve linear system.
121 ✗ PetscPrintf(PETSC_COMM_WORLD, "Solving the Navier-Stokes equation\n");
122 ✗ m_linearProblem[0]->solve(MpiInfo::rankProc(), MpiInfo::numProc());
123 ✗ m_linearProblem[0]->gatherSolution();
124
125 /*----------------------------
126 First ARD equation (omega)
127 ----------------------------*/
128 //Assembly loop on elements.
129 ✗ m_linearProblem[1]->assembleMatrixRHS(MpiInfo::rankProc());
130
131 //Apply boundary conditions.
132 ✗ m_linearProblem[1]->finalizeEssBCTransient();
133 ✗ m_linearProblem[1]->applyBC(FelisceParam::instance().essentialBoundaryConditionsMethod, MpiInfo::rankProc());
134
135 /*----------------------------
136 Second ARD equation (k)
137 ----------------------------*/
138 //Assembly loop on elements.
139 ✗ m_linearProblem[2]->assembleMatrixRHS(MpiInfo::rankProc());
140
141 //Apply boundary conditions.
142 ✗ m_linearProblem[2]->finalizeEssBCTransient();
143 ✗ m_linearProblem[2]->applyBC(FelisceParam::instance().essentialBoundaryConditionsMethod, MpiInfo::rankProc());
144
145 /*----------------------------
146 Solve the two ARD equations
147 ----------------------------*/
148 ✗ PetscPrintf(PETSC_COMM_WORLD, "Solving the Omega equation\n");
149 ✗ m_linearProblem[1]->solve(MpiInfo::rankProc(), MpiInfo::numProc());
150 ✗ m_linearProblem[1]->gatherSolution();
151
152 ✗ PetscPrintf(PETSC_COMM_WORLD, "Solving the K equation\n");
153 ✗ m_linearProblem[2]->solve(MpiInfo::rankProc(), MpiInfo::numProc());
154 ✗ m_linearProblem[2]->gatherSolution();
155
156 // k should be positive, just to be sure.
157 //VecAbs(m_linearProblem[2]->seqSolution());
158 //VecAbs(m_linearProblem[2]->solution());
159 }
160
161
162 ✗ void SSTModel::setExternalVec() {
163 // std::set external vectors for linear problems
164 // (to be used in assembling routines)
165
166 // Navier-Stokes
167 // k
168 ✗ m_linearProblem[0]->pushBackExternalVec(m_linearProblem[2]->sequentialSolution());
169 ✗ m_linearProblem[0]->pushBackExternalAO(m_linearProblem[2]->ao());
170 ✗ m_linearProblem[0]->pushBackExternalDof(m_linearProblem[2]->dof());
171
172 // omega
173 ✗ m_linearProblem[0]->pushBackExternalVec(m_linearProblem[1]->sequentialSolution());
174 ✗ m_linearProblem[0]->pushBackExternalAO(m_linearProblem[1]->ao());
175 ✗ m_linearProblem[0]->pushBackExternalDof(m_linearProblem[1]->dof());
176
177 // y
178 //m_linearProblem[0]->pushBackExternalVec(m_distance);
179 //m_linearProblem[0]->pushBackExternalAO(m_linearProblem[1]->ao());
180 //m_linearProblem[0]->pushBackExternalDof(m_linearProblem[1]->dof());
181
182
183 // Advection Reaction Diffusion for omega
184 // u
185 ✗ m_linearProblem[1]->pushBackExternalVec(m_linearProblem[0]->sequentialSolution());
186 ✗ m_linearProblem[1]->pushBackExternalAO(m_linearProblem[0]->ao());
187 ✗ m_linearProblem[1]->pushBackExternalDof(m_linearProblem[0]->dof());
188
189 // k
190 ✗ m_linearProblem[1]->pushBackExternalVec(m_linearProblem[2]->sequentialSolution());
191 ✗ m_linearProblem[1]->pushBackExternalAO(m_linearProblem[2]->ao());
192 ✗ m_linearProblem[1]->pushBackExternalDof(m_linearProblem[2]->dof());
193
194 // y
195 //m_linearProblem[1]->pushBackExternalVec(m_distance);
196 //m_linearProblem[1]->pushBackExternalAO(m_linearProblem[1]->ao());
197 //m_linearProblem[1]->pushBackExternalDof(m_linearProblem[1]->dof());
198
199
200 // Advection Reaction Diffusion for k
201 // u
202 ✗ m_linearProblem[2]->pushBackExternalVec(m_linearProblem[0]->sequentialSolution());
203 ✗ m_linearProblem[2]->pushBackExternalAO(m_linearProblem[0]->ao());
204 ✗ m_linearProblem[2]->pushBackExternalDof(m_linearProblem[0]->dof());
205
206 // omega
207 ✗ m_linearProblem[2]->pushBackExternalVec(m_linearProblem[1]->sequentialSolution());
208 ✗ m_linearProblem[2]->pushBackExternalAO(m_linearProblem[1]->ao());
209 ✗ m_linearProblem[2]->pushBackExternalDof(m_linearProblem[1]->dof());
210
211 // y
212 //m_linearProblem[2]->pushBackExternalVec(m_distance);
213 //m_linearProblem[2]->pushBackExternalAO(m_linearProblem[1]->ao());
214 //m_linearProblem[2]->pushBackExternalDof(m_linearProblem[1]->dof());
215 }
216
217
218 ✗ int SSTModel::getNstate() const {
219 ✗ return m_linearProblem[0]->numDof();
220 }
221
222 ✗ void SSTModel::getState(double* & state) {
223 ✗ m_linearProblem[0]->getSolution(state, MpiInfo::numProc(), MpiInfo::rankProc());
224 }
225
226 ✗ void SSTModel::getState_swig(double* data, felInt numDof) {
227 double * state;
228 ✗ m_linearProblem[0]->getSolution(state, MpiInfo::numProc(), MpiInfo::rankProc());
229 ✗ for (felInt i=0 ; i<numDof ; i++) {
230 ✗ data[i] = state[i];
231 }
232 }
233
234 ✗ void SSTModel::setState(double* & state) {
235 ✗ m_linearProblem[0]->setSolution(state, MpiInfo::numProc(), MpiInfo::rankProc());
236 }
237 }
238