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 "PETScInterface/romPETSC.hpp" |
21 |
|
|
|
22 |
|
|
namespace felisce { |
23 |
|
|
|
24 |
|
✗ |
RomPETSC::RomPETSC(): |
25 |
|
✗ |
m_fstransient(nullptr), |
26 |
|
✗ |
m_buildBasis(false), |
27 |
|
✗ |
m_buildMass(false), |
28 |
|
✗ |
m_buildInverse(false), |
29 |
|
✗ |
m_isFirstStep(false) |
30 |
|
|
{ |
31 |
|
✗ |
m_setPreconditionerOption = SAME_NONZERO_PATTERN; |
32 |
|
|
} |
33 |
|
|
|
34 |
|
✗ |
RomPETSC::~RomPETSC() { |
35 |
|
✗ |
if(m_buildBasis) { |
36 |
|
✗ |
m_romBasis.destroy(); |
37 |
|
✗ |
m_romMatrix.destroy(); |
38 |
|
|
} |
39 |
|
✗ |
if(m_buildMass) { |
40 |
|
✗ |
m_mass.destroy(); |
41 |
|
|
} |
42 |
|
✗ |
if(m_buildInverse) { |
43 |
|
✗ |
m_inverseMatrix.destroy(); |
44 |
|
|
} |
45 |
|
✗ |
if (m_isFirstStep) { |
46 |
|
✗ |
m_romRHS.destroy(); |
47 |
|
✗ |
m_romSolution.destroy(); |
48 |
|
|
} |
49 |
|
|
} |
50 |
|
|
|
51 |
|
|
// initialize the pod basis parameters (passing AO& ao too - useful to find petsc dof indeces mapping) |
52 |
|
✗ |
void RomPETSC::initializeRom(FelisceTransient::Pointer fstransient,MPI_Comm& comm, AO& ao) { |
53 |
|
✗ |
m_fstransient = fstransient; |
54 |
|
✗ |
m_verbose = FelisceParam::verbose(); |
55 |
|
✗ |
m_MPIcomm = comm; |
56 |
|
✗ |
m_ao = ao; |
57 |
|
✗ |
m_filenameBasis = FelisceParam::instance().podBasisFile; |
58 |
|
✗ |
m_dimRomBasis = FelisceParam::instance().dimRomBasis; |
59 |
|
✗ |
m_dimFullRomBasis = m_dimRomBasis; |
60 |
|
✗ |
m_useInverse = true; |
61 |
|
|
} |
62 |
|
|
|
63 |
|
|
// initialize the pod basis parameters |
64 |
|
✗ |
void RomPETSC::initializeRom(FelisceTransient::Pointer fstransient,MPI_Comm& comm) { |
65 |
|
✗ |
m_fstransient = fstransient; |
66 |
|
✗ |
m_verbose = FelisceParam::verbose(); |
67 |
|
✗ |
m_MPIcomm = comm; |
68 |
|
✗ |
m_filenameBasis = FelisceParam::instance().podBasisFile; |
69 |
|
✗ |
m_dimRomBasis = FelisceParam::instance().dimRomBasis; |
70 |
|
✗ |
m_dimFullRomBasis = m_dimRomBasis; |
71 |
|
✗ |
m_useInverse = true; |
72 |
|
|
} |
73 |
|
|
|
74 |
|
|
// Intialize basis matrix from matrix |
75 |
|
✗ |
void RomPETSC::initializeRomBasis(const PetscMatrix& basis) { |
76 |
|
|
felInt lengthBasis; |
77 |
|
✗ |
basis.getSize(&lengthBasis,&m_dimRomBasis); |
78 |
|
|
|
79 |
|
✗ |
m_romBasis.duplicateFrom(basis,MAT_COPY_VALUES); |
80 |
|
✗ |
m_romBasis.assembly(MAT_FINAL_ASSEMBLY); |
81 |
|
|
|
82 |
|
✗ |
m_buildBasis = true; |
83 |
|
|
|
84 |
|
|
} |
85 |
|
|
|
86 |
|
|
// Intialize basis matrix from matrix |
87 |
|
✗ |
void RomPETSC::initializeRomBasis(const std::vector<PetscVector*>& basis) { |
88 |
|
✗ |
m_dimRomBasis = basis.size(); |
89 |
|
|
felInt lengthBasis; |
90 |
|
✗ |
(*basis[0]).getSize(&lengthBasis); |
91 |
|
✗ |
m_romBasis.createDense(m_MPIcomm,PETSC_DECIDE,PETSC_DECIDE,lengthBasis,m_dimRomBasis,FELISCE_PETSC_NULLPTR); |
92 |
|
✗ |
m_romBasis.setFromOptions(); |
93 |
|
|
|
94 |
|
✗ |
double tmpValue[lengthBasis]; |
95 |
|
✗ |
felInt ia[lengthBasis]; |
96 |
|
✗ |
for (felInt i=0; i< lengthBasis; i++) { |
97 |
|
✗ |
ia[i] = i; |
98 |
|
|
} |
99 |
|
✗ |
AOApplicationToPetsc(m_ao,lengthBasis,ia); |
100 |
|
|
int ja; |
101 |
|
✗ |
for (int j=0; j<m_dimRomBasis; j++) { |
102 |
|
✗ |
(*basis[j]).getValues(lengthBasis,ia,tmpValue); |
103 |
|
✗ |
ja = j; |
104 |
|
✗ |
AOApplicationToPetsc(m_ao,1,&ja); |
105 |
|
✗ |
m_romBasis.setValues(lengthBasis,ia,1,&ja,tmpValue,INSERT_VALUES); |
106 |
|
|
} |
107 |
|
✗ |
m_romBasis.assembly(MAT_FINAL_ASSEMBLY); |
108 |
|
|
|
109 |
|
✗ |
m_buildBasis = true; |
110 |
|
|
|
111 |
|
|
} |
112 |
|
|
|
113 |
|
|
// read POD basis from file |
114 |
|
|
// format: |
115 |
|
|
// # of elements |
116 |
|
|
// # of dof |
117 |
|
|
// element1 || ...|| elementK |
118 |
|
|
// element1 || ...|| elementK |
119 |
|
|
// ... |
120 |
|
|
// Intialize basis matrix from file |
121 |
|
✗ |
void RomPETSC::initializeRomBasis() { |
122 |
|
✗ |
m_filenameBasis = FelisceParam::instance().inputDirectory + "/" + FelisceParam::instance().podBasisFile; |
123 |
|
✗ |
if (m_verbose) { |
124 |
|
✗ |
std::cout << " Rom::initBasis: Reading from file "<< m_filenameBasis << std::endl; |
125 |
|
|
} |
126 |
|
✗ |
std::ifstream basisFile(m_filenameBasis.c_str(),std::ios_base::in); |
127 |
|
✗ |
if (!basisFile) { |
128 |
|
✗ |
FEL_ERROR(" Rom::initBasis(): Impossible to read file: "+m_filenameBasis+" (POD basis)"); |
129 |
|
|
} |
130 |
|
|
|
131 |
|
|
felInt lengthBasis; |
132 |
|
✗ |
basisFile >> lengthBasis; // read n. of rows |
133 |
|
✗ |
basisFile >> m_dimFullRomBasis; // read n. of columns (total) |
134 |
|
✗ |
if (m_verbose) { |
135 |
|
✗ |
std::cout << "Read matrix of " << lengthBasis << " rows and " << m_dimFullRomBasis << " columns.\n"; |
136 |
|
|
} |
137 |
|
|
|
138 |
|
✗ |
if (m_dimFullRomBasis<m_dimRomBasis) { |
139 |
|
✗ |
std::cout << " requested dimBasis = " << m_dimRomBasis << ", available elements = "<< m_dimFullRomBasis << std::endl; |
140 |
|
✗ |
FEL_ERROR(" Rom::initBasis: not enough basis element available !"); |
141 |
|
|
} |
142 |
|
|
|
143 |
|
✗ |
m_romBasis.createDense(m_MPIcomm, PETSC_DECIDE, PETSC_DECIDE, lengthBasis, m_dimRomBasis, FELISCE_PETSC_NULLPTR); |
144 |
|
✗ |
m_romBasis.setFromOptions(); |
145 |
|
|
|
146 |
|
✗ |
double tmpValue[m_dimRomBasis]; |
147 |
|
✗ |
felInt ia[m_dimRomBasis]; |
148 |
|
|
felInt ja; |
149 |
|
|
char line[1024]; |
150 |
|
✗ |
double k = 1.0; |
151 |
|
✗ |
for (felInt j=0; j<lengthBasis; j++) { //rows |
152 |
|
✗ |
ja = j; |
153 |
|
✗ |
if ((j+1.0)>=(k*0.1*lengthBasis)) { |
154 |
|
✗ |
std::cout << "... " << 10*k << "% "; |
155 |
|
✗ |
k=k+1.0; |
156 |
|
|
} |
157 |
|
✗ |
for (felInt i=0; i< m_dimRomBasis; i++) { // columns |
158 |
|
✗ |
ia[i] = i; |
159 |
|
✗ |
basisFile >> tmpValue[i]; |
160 |
|
|
} |
161 |
|
✗ |
basisFile.getline(line, 1024, '\n'); |
162 |
|
✗ |
AOApplicationToPetsc(m_ao,1,&ja); |
163 |
|
✗ |
m_romBasis.setValues(1,&ja,m_dimRomBasis,ia,tmpValue,INSERT_VALUES); |
164 |
|
|
} |
165 |
|
✗ |
m_romBasis.assembly(MAT_FINAL_ASSEMBLY); |
166 |
|
|
|
167 |
|
✗ |
m_buildBasis = true; |
168 |
|
|
} |
169 |
|
|
|
170 |
|
✗ |
void RomPETSC::gatherVector(PetscVector& parVector, PetscVector& seqVector) { |
171 |
|
✗ |
PetscVector vout; |
172 |
|
✗ |
parVector.scatterToAll(vout,INSERT_VALUES,SCATTER_FORWARD); |
173 |
|
✗ |
seqVector.duplicateFrom(vout); |
174 |
|
✗ |
seqVector.copyFrom(vout); |
175 |
|
✗ |
vout.destroy(); |
176 |
|
|
} |
177 |
|
|
|
178 |
|
|
// given a full std::vector f, compute \hat f = Phi^t f |
179 |
|
✗ |
void RomPETSC::fullToReducedRHS(PetscVector& fullVector) { |
180 |
|
✗ |
if (m_fstransient->iteration == 1) { |
181 |
|
✗ |
m_romRHS.create(m_MPIcomm); |
182 |
|
✗ |
m_romRHS.setSizes(PETSC_DECIDE, m_dimRomBasis); |
183 |
|
✗ |
m_romRHS.setFromOptions(); |
184 |
|
✗ |
m_romRHS.zeroEntries(); |
185 |
|
|
} |
186 |
|
|
|
187 |
|
✗ |
multTranspose(m_romBasis,fullVector,m_romRHS); |
188 |
|
|
|
189 |
|
✗ |
m_isFirstStep = true; |
190 |
|
|
} |
191 |
|
|
|
192 |
|
|
// given a full matrix A, compute \hat A = Phi^t A \Phi |
193 |
|
✗ |
void RomPETSC::fullToReducedMatrix(PetscMatrix& fullMatrix) { |
194 |
|
✗ |
PetscMatrix tmpMat; |
195 |
|
|
// m_tmpMat = fullMatrix * m_romBasis |
196 |
|
✗ |
matMult(fullMatrix, m_romBasis, MAT_INITIAL_MATRIX, PETSC_DEFAULT, tmpMat); |
197 |
|
|
// m_romMatrix = m_romBasis^T * tmpMat |
198 |
|
✗ |
transposeMatMult(m_romBasis, tmpMat, MAT_INITIAL_MATRIX, PETSC_DEFAULT, m_romMatrix); |
199 |
|
|
|
200 |
|
✗ |
tmpMat.destroy(); |
201 |
|
|
|
202 |
|
✗ |
if (m_useInverse) |
203 |
|
✗ |
calculateInverse(); |
204 |
|
|
|
205 |
|
|
} |
206 |
|
|
|
207 |
|
|
// solve the reduced linear system |
208 |
|
✗ |
void RomPETSC::solveALP() { |
209 |
|
✗ |
if (m_fstransient->iteration == 1) { |
210 |
|
✗ |
m_romSolution.duplicateFrom(m_romRHS); |
211 |
|
|
} |
212 |
|
✗ |
m_romSolution.zeroEntries(); |
213 |
|
|
|
214 |
|
✗ |
if (m_verbose) { |
215 |
|
✗ |
std::cout << " Rom:: solveROM() - solve reduced model "<< std::endl; |
216 |
|
|
} |
217 |
|
|
|
218 |
|
✗ |
mult(m_inverseMatrix, m_romRHS, m_romSolution); |
219 |
|
|
|
220 |
|
|
} |
221 |
|
|
|
222 |
|
|
|
223 |
|
|
// solve the reduced linear system |
224 |
|
✗ |
void RomPETSC::solveROM(PetscVector& solution) { |
225 |
|
✗ |
if (m_fstransient->iteration == 1) { |
226 |
|
✗ |
m_romSolution.duplicateFrom(m_romRHS); |
227 |
|
|
} |
228 |
|
✗ |
m_romSolution.zeroEntries(); |
229 |
|
|
|
230 |
|
✗ |
if (m_verbose) { |
231 |
|
✗ |
std::cout << " Rom:: solveROM() - solve reduced model "<< std::endl; |
232 |
|
|
} |
233 |
|
|
|
234 |
|
✗ |
if (!m_useInverse) { |
235 |
|
|
|
236 |
|
✗ |
PCFactorSetReuseFill(m_pc, PETSC_TRUE); |
237 |
|
|
|
238 |
|
✗ |
KSPSetOperators(m_ksp,m_romMatrix.toPetsc(),m_romMatrix.toPetsc()); |
239 |
|
✗ |
if(strcmp(FelisceParam::instance().setPreconditionerOption[0].c_str(), "SAME_PRECONDITIONER") == 0) |
240 |
|
✗ |
KSPSetReusePreconditioner(m_ksp, PETSC_TRUE); |
241 |
|
|
else |
242 |
|
✗ |
KSPSetReusePreconditioner(m_ksp, PETSC_FALSE); |
243 |
|
|
|
244 |
|
✗ |
KSPSetFromOptions(m_ksp); |
245 |
|
✗ |
KSPSolve(m_ksp,m_romRHS.toPetsc(),m_romSolution.toPetsc()); |
246 |
|
|
|
247 |
|
|
PetscInt its; |
248 |
|
✗ |
KSPGetIterationNumber(m_ksp,&its); |
249 |
|
✗ |
if (m_verbose) { |
250 |
|
✗ |
PetscPrintf(m_MPIcomm, "\n/================Information about solver==============/\n"); |
251 |
|
✗ |
PetscPrintf(m_MPIcomm, "Solver use: <%s>.\n",FelisceParam::instance().solver[0].c_str()); |
252 |
|
✗ |
PetscPrintf(m_MPIcomm, "Preconditioner use: <%s>.\n",FelisceParam::instance().preconditioner[0].c_str()); |
253 |
|
✗ |
PetscPrintf(m_MPIcomm, "Preconditioner re-use option: <%s>.\n",FelisceParam::instance().setPreconditionerOption[0].c_str()); |
254 |
|
✗ |
if (m_verbose > 1 ) { |
255 |
|
✗ |
PetscPrintf(m_MPIcomm, "Relative tolerance: %lf.\n",FelisceParam::instance().relativeTolerance[0]); |
256 |
|
✗ |
PetscPrintf(m_MPIcomm, "Absolute tolerance: %lf.\n",FelisceParam::instance().absoluteTolerance[0]); |
257 |
|
✗ |
PetscPrintf(m_MPIcomm, "Maximum of iterations: %d.\n",FelisceParam::instance().maxIteration[0]); |
258 |
|
|
} |
259 |
|
✗ |
PetscPrintf(m_MPIcomm, "Number of iterations to inverse matrix: %d\n", its); |
260 |
|
|
} |
261 |
|
|
} else { |
262 |
|
✗ |
mult(m_inverseMatrix, m_romRHS, m_romSolution); |
263 |
|
|
} |
264 |
|
|
|
265 |
|
✗ |
if (m_verbose) { |
266 |
|
✗ |
std::cout << " Rom:: solveROM() - get full solution "<< std::endl; |
267 |
|
|
} |
268 |
|
✗ |
mult(m_romBasis,m_romSolution,solution); |
269 |
|
|
} |
270 |
|
|
|
271 |
|
✗ |
void RomPETSC::setSolver(KSP ksp, PC pc, MatStructure setPreconditionerOption) |
272 |
|
|
{ |
273 |
|
|
IGNORE_UNUSED_ARGUMENT(ksp); |
274 |
|
|
IGNORE_UNUSED_ARGUMENT(pc); |
275 |
|
|
IGNORE_UNUSED_ARGUMENT(setPreconditionerOption); |
276 |
|
|
|
277 |
|
✗ |
const auto& r_instance = FelisceParam::instance(); |
278 |
|
✗ |
KSPCreate(m_MPIcomm,&m_ksp); |
279 |
|
✗ |
PCCreate(m_MPIcomm,&m_pc); |
280 |
|
✗ |
const auto& solver = r_instance.solver[0]; |
281 |
|
✗ |
if ( solver == "preonly" || solver == "mumps" || solver == "superlu_dist" ) |
282 |
|
✗ |
KSPSetType(m_ksp, KSPPREONLY); |
283 |
|
|
else |
284 |
|
✗ |
KSPSetType(m_ksp, solver.c_str()); |
285 |
|
✗ |
PCSetFromOptions(m_pc); |
286 |
|
✗ |
KSPGetPC(m_ksp,&m_pc); |
287 |
|
✗ |
PCSetType(m_pc,r_instance.preconditioner[0].c_str()); |
288 |
|
✗ |
KSPSetTolerances(m_ksp,r_instance.relativeTolerance[0],r_instance.absoluteTolerance[0],PETSC_DEFAULT,r_instance.maxIteration[0]); |
289 |
|
|
|
290 |
|
|
// TODO: MATSOLVERPASTIX once it works |
291 |
|
✗ |
if ( solver == "preonly" || solver == "mumps" ) { |
292 |
|
✗ |
PCFactorSetMatSolverType(m_pc, MATSOLVERMUMPS); |
293 |
|
|
} else |
294 |
|
|
#if SUPERLU_DIST_COMPILED |
295 |
|
✗ |
if ( solver == "superlu_dist" ) { |
296 |
|
✗ |
PCFactorSetMatSolverType(m_pc, MATSOLVERSUPERLU_DIST); |
297 |
|
|
} else |
298 |
|
|
#else |
299 |
|
|
if ( solver == "superlu_dist" ) { |
300 |
|
|
FEL_WARNING("SuperLUDist not available, we proceed using MUMPS") |
301 |
|
|
PCFactorSetMatSolverType(m_pc, MATSOLVERMUMPS); |
302 |
|
|
} else |
303 |
|
|
#endif |
304 |
|
|
{ |
305 |
|
✗ |
if ( r_instance.initSolverWithPreviousSolution[0] ) { |
306 |
|
✗ |
KSPSetInitialGuessNonzero(m_ksp, PETSC_TRUE); |
307 |
|
|
} |
308 |
|
|
} |
309 |
|
|
|
310 |
|
|
{ |
311 |
|
|
// TODO: to be removed |
312 |
|
✗ |
PetscOptionsSetValue(NULL, "-mat_mumps_icntl_20", "0"); |
313 |
|
|
} |
314 |
|
|
} |
315 |
|
|
|
316 |
|
✗ |
void RomPETSC::calculateInverse() { |
317 |
|
|
//Create the identityMatrix |
318 |
|
✗ |
PetscMatrix identityMatrix; |
319 |
|
✗ |
identityMatrix.createDense(m_MPIcomm, PETSC_DECIDE, PETSC_DECIDE, m_dimRomBasis, m_dimRomBasis, FELISCE_PETSC_NULLPTR); |
320 |
|
✗ |
identityMatrix.setFromOptions(); |
321 |
|
✗ |
identityMatrix.zeroEntries(); |
322 |
|
✗ |
m_inverseMatrix.duplicateFrom(identityMatrix,MAT_COPY_VALUES); |
323 |
|
✗ |
identityMatrix.shift(1.0); |
324 |
|
|
|
325 |
|
|
// ILU factorization of_romMatrix |
326 |
|
|
IS perm; |
327 |
|
✗ |
ISCreateStride(m_MPIcomm, m_dimRomBasis, 0, 1, &perm); |
328 |
|
|
MatFactorInfo iluinfo; |
329 |
|
✗ |
iluinfo.fill = 1.0; |
330 |
|
✗ |
m_romMatrix.luFactor(perm, perm, &iluinfo); |
331 |
|
|
|
332 |
|
|
// Solve m_romMatrix * m_inverseMatrix = identityMatrix |
333 |
|
✗ |
MatMatSolve(m_romMatrix.toPetsc(), identityMatrix.toPetsc(), m_inverseMatrix.toPetsc()); |
334 |
|
|
|
335 |
|
✗ |
identityMatrix.destroy(); |
336 |
|
|
|
337 |
|
✗ |
m_buildInverse = true; |
338 |
|
|
|
339 |
|
|
} |
340 |
|
|
|
341 |
|
|
|
342 |
|
✗ |
void RomPETSC::readCase(int& precisionSol,std::vector<std::string>& variable, felInt& snapSize) { |
343 |
|
|
// Solution files. |
344 |
|
✗ |
std::string caseFileName = FelisceParam::instance().resultDir + "/" + FelisceParam::instance().prefixName[0] + ".case"; |
345 |
|
|
|
346 |
|
✗ |
std::ifstream casefile(caseFileName.c_str(), std::ios_base::in); |
347 |
|
|
char line[1024]; |
348 |
|
✗ |
std::string aux; |
349 |
|
✗ |
if ( !casefile ) { |
350 |
|
✗ |
std::cout << " ERROR: Can not open file "+caseFileName+"."<< std::endl; |
351 |
|
✗ |
exit(1); |
352 |
|
|
} else { |
353 |
|
✗ |
if (FelisceParam::verbose() > 0) { |
354 |
|
✗ |
std::cout << "Reading " << caseFileName << std::endl; |
355 |
|
|
} |
356 |
|
|
} |
357 |
|
|
|
358 |
|
✗ |
casefile.getline(line, 1024, ':'); |
359 |
|
✗ |
casefile.getline(line, 1024, ':'); |
360 |
|
✗ |
casefile.getline(line, 1024, ' '); |
361 |
|
✗ |
casefile.getline(line, 1024, ' '); |
362 |
|
✗ |
casefile.getline(line, 1024, '\n'); |
363 |
|
✗ |
std::string geoFileName = FelisceParam::instance().resultDir + "/" + line; |
364 |
|
✗ |
std::ifstream geoFile(geoFileName.c_str(),std::ios_base::in); |
365 |
|
✗ |
for (int i=0; i<5; i++) { |
366 |
|
✗ |
geoFile.getline(line,1024,'\n'); |
367 |
|
|
} |
368 |
|
✗ |
geoFile >> snapSize; |
369 |
|
✗ |
geoFile.close(); |
370 |
|
|
|
371 |
|
✗ |
bool endCase = false; |
372 |
|
✗ |
while (!endCase) { |
373 |
|
✗ |
casefile.getline(line, 1024, ':'); |
374 |
|
✗ |
casefile.getline(line, 1024, ' '); |
375 |
|
✗ |
casefile.getline(line, 1024, ' '); |
376 |
|
✗ |
casefile.getline(line, 1024, ' '); |
377 |
|
✗ |
casefile.getline(line, 1024, '.'); |
378 |
|
✗ |
variable.emplace_back(line); |
379 |
|
✗ |
casefile.getline(line, 1024, '.'); |
380 |
|
✗ |
aux = line; |
381 |
|
✗ |
precisionSol = aux.length(); |
382 |
|
✗ |
casefile.getline(line, 1024, '\n'); |
383 |
|
✗ |
casefile.getline(line, 1024, ' '); |
384 |
|
✗ |
aux = line; |
385 |
|
✗ |
if (aux == "TIME\ntime") { |
386 |
|
✗ |
endCase = true; |
387 |
|
|
} |
388 |
|
|
} |
389 |
|
✗ |
casefile.close(); |
390 |
|
|
|
391 |
|
✗ |
if (m_verbose) { |
392 |
|
✗ |
std::cout << "Variables:\n"; |
393 |
|
✗ |
for (std::size_t k=0; k < variable.size(); k++) { |
394 |
|
✗ |
std::cout << variable[k] << std::endl; |
395 |
|
|
} |
396 |
|
✗ |
std::cout << "Numbers of nodes per variable : " << snapSize << std::endl; |
397 |
|
|
} |
398 |
|
|
|
399 |
|
|
} |
400 |
|
|
|
401 |
|
✗ |
void RomPETSC::readEnsight(double* scalarValue,std::string sclFile, felInt snapSize) { |
402 |
|
✗ |
std::string snapFileName = FelisceParam::instance().resultDir + sclFile + ".scl"; |
403 |
|
✗ |
std::ifstream snap(snapFileName.c_str(), std::ios_base::in); |
404 |
|
|
|
405 |
|
✗ |
if ( !snap ) { |
406 |
|
✗ |
std::cout << " ERROR: Can not open file "+snapFileName+"."<< std::endl; |
407 |
|
|
} else { |
408 |
|
✗ |
if (FelisceParam::verbose() > 1) { |
409 |
|
✗ |
std::cout << "Reading " << snapFileName << std::endl; |
410 |
|
|
} |
411 |
|
|
} |
412 |
|
|
char newline[1024]; |
413 |
|
✗ |
snap.getline(newline, 1024, '\n'); |
414 |
|
✗ |
for (felInt j=0; j < snapSize; j++) { |
415 |
|
✗ |
snap >> scalarValue[j]; |
416 |
|
|
} |
417 |
|
✗ |
snap.close(); |
418 |
|
|
} |
419 |
|
|
|
420 |
|
|
|
421 |
|
|
} |
422 |
|
|
|