Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
DriverSteadyState.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File DriverSteadyState.cpp
4 //
5 // For more information, please see: http://www.nektar.info
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: Incompressible Navier Stokes solver
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
38 
39 using namespace std;
40 
41 namespace Nektar
42 {
43 namespace SolverUtils
44 {
45 
46 string DriverSteadyState::className
48  "SteadyState", DriverSteadyState::create);
49 string DriverSteadyState::driverLookupId
50  = LibUtilities::SessionReader::RegisterEnumValue(
51  "Driver","SteadyState",0);
52 
53 /**
54  *
55  */
56 DriverSteadyState::DriverSteadyState(
58  : DriverModifiedArnoldi(pSession)
59 {
60 }
61 
62 
63 /**
64  *
65  */
67 {
68 }
69 
70 
71 /**
72  *
73  */
75 {
77 }
78 
79 
80 void DriverSteadyState::v_Execute(ostream &out)
81 
82 {
83  // With a loop over "DoSolve", this Driver implements the
84  // "encaplulated" Selective Frequency Damping method
85  // to find the steady state of a flow above the critical Reynolds
86  // number.
87  m_equ[m_nequ - 1]->PrintSummary(out);
88  m_equ[m_nequ - 1]->DoInitialise();
89 
90  m_session->LoadParameter("IO_InfoSteps", m_infosteps, 1000);
91  m_session->LoadParameter("IO_CheckSteps", m_checksteps, 100000);
92  m_session->LoadParameter("ControlCoeff",m_X, 1);
93  m_session->LoadParameter("FilterWidth", m_Delta, 2);
94 
95  // To evaluate optimum SFD parameters if growth rate provided in the
96  // xml file
97  m_session->LoadParameter("GrowthRateEV", GrowthRateEV, 0.0);
98 
99  // To evaluate optimum SFD parameters if frequency provided in the xml
100  // file
101  m_session->LoadParameter("FrequencyEV", FrequencyEV, 0.0);
102 
103  // Determine when SFD method is converged
104  m_session->LoadParameter("TOL", TOL, 1.0e-08);
105 
106  // Used only for the Adaptive SFD method
107  m_session->LoadParameter("AdaptiveTOL", AdaptiveTOL, 1.0e-02);
108 
109  // Used only for the Adaptive SFD method
110  m_session->LoadParameter("AdaptiveTime", AdaptiveTime, 50.0*m_Delta);
111 
112  if (m_comm->GetRank() == 0)
113  {
114  PrintSummarySFD();
115  }
116 
117  timer.Start();
118 
119  // Definition of shared pointer (used only for the Adaptive SFD method)
121 
122  // Condition necessary to run SFD for the compressible case
123  NumVar_SFD = m_equ[m_nequ - 1]->UpdateFields()[0]->GetCoordim(0);
124  if (m_session->GetSolverInfo("EqType") == "EulerCFE" ||
125  m_session->GetSolverInfo("EqType") == "NavierStokesCFE")
126  {
127  // Number of variables for the compressible equations
128  NumVar_SFD += 2;
129  }
130  if(m_session->DefinesSolverInfo("HOMOGENEOUS"))
131  {
132  if (m_session->GetSolverInfo("HOMOGENEOUS") == "1D")
133  {
134  NumVar_SFD += 1;
135  }
136  }
137 
138  // We store the time step
139  m_dt = m_equ[m_nequ - 1]->GetTimeStep();
140 
141  // Evaluate optimum SFD parameters if dominent EV given by xml file
142  if (GrowthRateEV != 0.0 && FrequencyEV != 0.0)
143  {
144  cout << "Besed on the dominant EV given in the xml file,"
145  << "a 1D model is used to evaluate the optumum parameters"
146  << "of the SFD method:" << endl;
147  complex<NekDouble> EV = polar(exp(GrowthRateEV), FrequencyEV);
149  }
150 
151 
152  // We set up the elements of the operator of the encapsulated
153  // formulation of the selective frequencive damping method
155 
156  // m_steps is set to 1. Then "m_equ[m_nequ - 1]->DoSolve()" will run
157  // for only one time step
158  m_equ[m_nequ - 1]->SetSteps(1);
159  ofstream m_file("ConvergenceHistory.txt", ios::out | ios::trunc);
160 
161  Array<OneD, Array<OneD, NekDouble> > q0(NumVar_SFD);
162  Array<OneD, Array<OneD, NekDouble> > q1(NumVar_SFD);
163  Array<OneD, Array<OneD, NekDouble> > qBar0(NumVar_SFD);
164  Array<OneD, Array<OneD, NekDouble> > qBar1(NumVar_SFD);
165 
166  for(int i = 0; i < NumVar_SFD; ++i)
167  {
168  q0[i] = Array<OneD, NekDouble> (m_equ[m_nequ-1]->GetTotPoints(),
169  0.0); //q0 is initialised
170  qBar0[i] = Array<OneD, NekDouble> (m_equ[m_nequ-1]->GetTotPoints(),
171  0.0);
172  m_equ[m_nequ - 1]->CopyFromPhysField(i, qBar0[i]);
173  }
174 
175  ///Definition of variables used in this algorithm
176  m_stepCounter = 0;
177  m_Check = 0;
178  m_Check_BaseFlow = 1;
180  Diff_q_qBar = 1.0;
181  Diff_q1_q0 = 1.0;
182  cpuTime = 0.0;
183  elapsed = 0.0;
184  totalTime = 0.0;
185  FlowPartiallyConverged = false;
186 
187  while (max(Diff_q_qBar, Diff_q1_q0) > TOL)
188  {
189  ///Call the Navier-Stokes solver for one time step
190  m_equ[m_nequ - 1]->DoSolve();
191 
192  for(int i = 0; i < NumVar_SFD; ++i)
193  {
194  ///Copy the current flow field into q0
195  m_equ[m_nequ - 1]->CopyFromPhysField(i, q0[i]);
196 
197  ///Apply the linear operator to the outcome of the solver
198  ComputeSFD(i, q0, qBar0, q1, qBar1);
199 
200  ///Update qBar
201  qBar0[i] = qBar1[i];
202 
203  ///Copy the output of the SFD method into the current flow field
204  m_equ[m_nequ - 1]->CopyToPhysField(i, q1[i]);
205  }
206 
208  {
210 
211  ///Loop for the adaptive SFD method
213  FlowPartiallyConverged == false)
214  {
215 
217 
218  if (Diff_q_qBar < AdaptiveTOL)
219  {
220  if (m_comm->GetRank() == 0)
221  {
222  cout << "\n\t The SFD method is converging: we compute "
223  << "stability analysis using the 'partially "
224  << "converged' steady state as base flow:\n" << endl;
225  }
226 
227  m_equ[m_nequ - 1]->Checkpoint_BaseFlow(m_Check_BaseFlow);
229 
230  A->GetAdvObject()->SetBaseFlow(q0);
232 
233  if (m_comm->GetRank() == 0)
234  {
235  // Compute the update of the SFD parameters only on
236  // one processor
238  }
239  else
240  {
241  // On all the other processors, the parameters are set
242  // to 0
243  m_X = 0;
244  m_Delta = 0;
245  }
246  // The we give to all the processors the value of X and
247  // Delta of the first processor
250 
251  FlowPartiallyConverged = true;
252  }
254  >= AdaptiveTime)
255  {
256  if (m_comm->GetRank() == 0)
257  {
258  cout << "\n\t We compute stability analysis using"
259  << " the current flow field as base flow:\n"
260  << endl;
261  }
262 
263  m_equ[m_nequ - 1]->Checkpoint_BaseFlow(m_Check_BaseFlow);
265 
266  A->GetAdvObject()->SetBaseFlow(q0);
268 
269  if (m_comm->GetRank() == 0)
270  {
271  // Compute the update of the SFD parameters only on
272  // one processor
274  }
275  else
276  {
277  // On all the other processors, the parameters are set
278  // to 0
279  m_X = 0;
280  m_Delta = 0;
281  }
282  // The we give to all the processors the value of X and
283  // Delta of the first processor
286 
288  }
289  }
290  }
291 
293  {
294  m_Check++;
295  m_equ[m_nequ - 1]->Checkpoint_Output(m_Check);
296  }
297  m_stepCounter++;
298  }
299 
300  m_file.close();
301 
302  ///We save the final solution into a .fld file
303  m_equ[m_nequ - 1]->Output();
304 
305  for(int j = 0; j < m_equ[m_nequ - 1]->GetNvariables(); ++j)
306  {
307  NekDouble vL2Error = m_equ[m_nequ - 1]->L2Error(j,false);
308  NekDouble vLinfError = m_equ[m_nequ - 1]->LinfError(j);
309  if (m_comm->GetRank() == 0)
310  {
311  out << "L 2 error (variable " << m_equ[m_nequ - 1]->GetVariable(j)
312  << ") : " << vL2Error << endl;
313  out << "L inf error (variable " << m_equ[m_nequ - 1]->GetVariable(j)
314  << ") : " << vLinfError << endl;
315  }
316  }
317 }
318 
319 
320 /**
321  * This routine defines the encapsulated SFD operator with first-order
322  * splitting and exact resolution of the second subproblem.
323  *(See http://scitation.aip.org/content/aip/journal/pof2/26/3/10.1063/1.4867482 for details)
324  */
326  const NekDouble Delta_input)
327 {
328  NekDouble X = X_input*m_dt;
329  NekDouble Delta = Delta_input/m_dt;
330  NekDouble coeff = 1.0/(1.0 + X*Delta);
331  M11 = coeff*(1.0 + X*Delta*exp(-(X + 1.0/Delta)));
332  M12 = coeff*(X*Delta*(1.0 - exp(-(X + 1.0/Delta))));
333  M21 = coeff*(1.0 - exp(-(X + 1.0/Delta)));
334  M22 = coeff*(X*Delta + exp(-(X + 1.0/Delta)));
335 }
336 
337 
339  const Array<OneD, const Array<OneD, NekDouble> > &q0,
340  const Array<OneD, const Array<OneD, NekDouble> > &qBar0,
343 {
344  q1[i] = Array<OneD, NekDouble> (m_equ[m_nequ - 1]->GetTotPoints(),0.0);
345  qBar1[i] = Array<OneD, NekDouble> (m_equ[m_nequ - 1]->GetTotPoints(),0.0);
346 
347  ///Encapsulated SFD method
348  Vmath::Svtvp(q1[i].num_elements(), M11, q0[i], 1, q1[i], 1, q1[i], 1 );
349  Vmath::Svtvp(q1[i].num_elements(), M12, qBar0[i], 1, q1[i], 1, q1[i], 1 );
350 
351  Vmath::Svtvp(qBar1[i].num_elements(), M21, q0[i], 1, qBar1[i], 1,
352  qBar1[i], 1 );
353  Vmath::Svtvp(qBar1[i].num_elements(), M22, qBar0[i], 1, qBar1[i], 1,
354  qBar1[i], 1 );
355 }
356 
357 
359 {
360  NekDouble growthEV(0.0);
361  NekDouble frequencyEV(0.0);
362 
363  // m_kdim is the dimension of Krylov subspace (defined in the xml file and
364  // used in DriverArnoldi.cpp)
365  ReadEVfile(m_kdim, growthEV, frequencyEV);
366 
367  cout << "\n\tgrowthEV = " << growthEV << endl;
368  cout << "\tfrequencyEV = " << frequencyEV << endl;
369 
370  complex<NekDouble> ApproxEV = polar(exp(growthEV), frequencyEV);
371 
372  NekDouble X_new = m_X;
373  NekDouble Delta_new = m_Delta;
374 
375  GradientDescentMethod(ApproxEV, X_new, Delta_new);
376 
377  m_X = X_new;
378  m_Delta = Delta_new;
379 
381 }
382 
383 /**
384  * This routine implements a gradient descent method to find the parameters X
385  * end Delta which give the minimum eigenlavue of the SFD problem applied to
386  * the scalar case u(n+1) = \alpha*u(n).
387  */
389  const complex<NekDouble> &alpha,
390  NekDouble &X_output,
391  NekDouble &Delta_output)
392 {
393  cout << "\n\tWe enter the Gradient Descent Method [...]" << endl;
394  bool OptParmFound = false;
395  bool Descending = true;
396  NekDouble X_input = X_output;
397  NekDouble Delta_input = Delta_output;
398 
399  NekDouble X_init = X_output;
400  NekDouble Delta_init = Delta_output;
401  int stepCounter(0);
402 
403  NekDouble F0(0.0);
404  NekDouble F0x(0.0);
405  NekDouble F0y(0.0);
406  NekDouble F1(0.0);
407  NekDouble dx = 0.00000001;
408  NekDouble dirx(0.0);
409  NekDouble diry(0.0);
410  NekDouble s(0.0);
411  NekDouble CurrentMin = 1.0;
412 
413  while (OptParmFound == false)
414  {
415  Descending = true;
416  EvalEV_ScalarSFD(X_input, Delta_input, alpha, F0);
417  EvalEV_ScalarSFD(X_input + dx, Delta_input, alpha, F0x);
418  EvalEV_ScalarSFD(X_input, Delta_input + dx, alpha, F0y);
419  dirx = (F0 - F0x);
420  diry = (F0 - F0y);
421  s = abs(0.000001/dirx);
422  X_output = X_input + s*dirx;
423  Delta_output = Delta_input + s*diry;
424  F1 = F0;
425 
426  while (Descending == true)
427  {
428  CurrentMin = F1;
429  X_input = X_output;
430  Delta_input = Delta_output;
431  EvalEV_ScalarSFD(X_output, Delta_output, alpha, F1);
432 
433  if (F1 > CurrentMin)
434  {
435  Descending = false;
436  }
437  else
438  {
439  s = s + s*0.01;
440  X_output = X_input + s*dirx;
441  Delta_output = Delta_input + s*diry;
442  }
443 
444  if(stepCounter > 9999999)
445  {
446  //We are stuck in this loop..
447  //Then we restart it with different initail conditions
448  Descending = false;
449  X_input = X_init;
450  Delta_init = Delta_init + Delta_init*0.1;
451  Delta_input = Delta_init;
452  stepCounter = 0;
453  }
454  stepCounter++;
455  }
456 
457  if (abs(F0-F1) < dx)
458  {
459  cout << "\tThe Gradient Descent Method has converged!" << endl;
460  EvalEV_ScalarSFD(X_output, Delta_output, alpha, F1);
461  cout << "\n\tThe updated parameters are: X_tilde = " << X_output
462  << " and Delta_tilde = " << Delta_output << endl;
463  OptParmFound = true;
464  }
465 
466  }
467 }
468 
469 
470 /**
471  * This routine evaluates the maximum eigenvalue of the SFD system when applied
472  * to the 1D model u(n+1) = alpha*u(n)
473  */
475  const NekDouble &X_input,
476  const NekDouble &Delta_input,
477  const complex<NekDouble> &alpha,
478  NekDouble &MaxEV)
479 {
480  NekDouble A11 = ( 1.0 + X_input * Delta_input *
481  exp(-(X_input + 1.0/Delta_input)) )/(1.0 + X_input*Delta_input);
482  NekDouble A12 = ( X_input*Delta_input - X_input * Delta_input *
483  exp(-(X_input + 1.0/Delta_input)) )/(1.0 + X_input*Delta_input);
484  NekDouble A21 = ( 1.0 - 1.0 *
485  exp(-(X_input + 1.0/Delta_input)) )/(1 + X_input*Delta_input);
486  NekDouble A22 = ( X_input*Delta_input + 1.0 *
487  exp(-(X_input + 1.0/Delta_input)) )/(1.0 + X_input*Delta_input);
488 
489  complex<NekDouble> B11 = alpha;
490  NekDouble B12 = 0.0;
491  NekDouble B21 = 0.0;
492  NekDouble B22 = 1.0;
493 
494  complex<NekDouble> a = A11*B11 + A12*B21;
495  NekDouble b = A11*B12 + A12*B22;
496  complex<NekDouble> c = A21*B11 + A22*B21;
497  NekDouble d = A21*B12 + A22*B22;
498 
499  complex<NekDouble> delt = sqrt((a-d)*(a-d) + 4.0*b*c);
500  complex<NekDouble> lambda_1 = 0.5*(a+d + delt);
501  complex<NekDouble> lambda_2 = 0.5*(a+d - delt);
502 
503  NekDouble NORM_1 = abs(lambda_1);
504  NekDouble NORM_2 = abs(lambda_2);
505 
506  MaxEV = max(NORM_1, NORM_2);
507 }
508 
509 
511  int &KrylovSubspaceDim,
512  NekDouble &growthEV,
513  NekDouble &frequencyEV)
514 {
515  // This routine reads the .evl file written by the Arnoldi algorithm
516  // (written in September 2014)
517  std::string line;
518  int NumLinesInFile = 0;
519  std::string EVfileName = m_session->GetSessionName() + ".evl";
520  std::ifstream EVfile(EVfileName.c_str());
521  ASSERTL0(EVfile.good(), "Cannot open .evl file.");
522 
523  if(EVfile)
524  {
525  // This block counts the total number of lines of the .evl file
526  // We keep going util we reach the end of the file
527  while(getline(EVfile, line))
528  {
529  ++NumLinesInFile;
530  }
531 
532  // It may happen that the Stability method that have produced the .elv
533  // file converges in less than m_kdim iterations. In this case,
534  // KrylovSubspaceDim has to be changed here
535  if(NumLinesInFile < KrylovSubspaceDim*2.0
536  + KrylovSubspaceDim*(KrylovSubspaceDim+1.0)/2.0)
537  {
538  for(int i = 1; i <= KrylovSubspaceDim; ++i)
539  {
540  if(NumLinesInFile == i*2.0 + i*(i+1.0)/2.0)
541  {
542  KrylovSubspaceDim = i;
543  }
544  }
545  }
546 
547  // go back to the beginning of the file
548  EVfile.clear();
549  EVfile.seekg(0, ios::beg);
550 
551  // We now want to go to the line where the most unstable eigenlavue was
552  // written
553  for(int i = 0; i < (NumLinesInFile - KrylovSubspaceDim); ++i)
554  {
555  std::getline(EVfile, line);
556  cout << "Discard line: " << line << endl;
557  }
558 
559  std::vector<std::string> tokens;
560  std::getline(EVfile, line);
561  boost::algorithm::split(tokens, line,
562  boost::is_any_of("\t "),boost::token_compress_on);
563 
564  ASSERTL0(tokens.size() >= 5,
565  "Unexpected formatting of .evl file while reading line:\n"
566  + line);
567  growthEV = boost::lexical_cast<NekDouble>(tokens[4]);
568  frequencyEV = boost::lexical_cast<NekDouble>(tokens[5]);
569  }
570  else
571  {
572  cout << "An error occurred when opening the .evl file" << endl;
573  }
574  EVfile.close();
575 }
576 
577 
578 /**
579  * This routine evaluates |q-qBar|_inf (and |q1-q0|_inf) and writes the values
580  * in "ConvergenceHistory.txt"
581  */
583  const Array<OneD, const Array<OneD, NekDouble> > &qBar1,
584  const Array<OneD, const Array<OneD, NekDouble> > &q0,
585  NekDouble &MaxNormDiff_q_qBar,
586  NekDouble &MaxNormDiff_q1_q0)
587 {
588  Array<OneD, NekDouble > NormDiff_q_qBar(NumVar_SFD, 1.0);
589  Array<OneD, NekDouble > NormDiff_q1_q0(NumVar_SFD, 1.0);
590  MaxNormDiff_q_qBar=0.0;
591  MaxNormDiff_q1_q0=0.0;
592 
593  for(int i = 0; i < NumVar_SFD; ++i)
594  {
595  NormDiff_q_qBar[i] = m_equ[m_nequ - 1]->LinfError(i, qBar1[i]);
596  NormDiff_q1_q0[i] = m_equ[m_nequ - 1]->LinfError(i, q0[i]);
597 
598  if (MaxNormDiff_q_qBar < NormDiff_q_qBar[i])
599  {
600  MaxNormDiff_q_qBar = NormDiff_q_qBar[i];
601  }
602  if (MaxNormDiff_q1_q0 < NormDiff_q1_q0[i])
603  {
604  MaxNormDiff_q1_q0 = NormDiff_q1_q0[i];
605  }
606  }
607 
608  timer.Stop();
610  cpuTime += elapsed;
611  totalTime += elapsed;
612 
613  if (m_comm->GetRank() == 0)
614  {
615  cout << "SFD - Step: " << left << m_stepCounter+1
616  << ";\tTime: " << left << m_equ[m_nequ - 1]->GetFinalTime()
617  << ";\tCPU time = " << left << cpuTime << " s"
618  << ";\tTot time = " << left << totalTime << " s"
619  << ";\tX = " << left << m_X
620  << ";\tDelta = " << left << m_Delta
621  << ";\t|q-qBar|inf = " << left << MaxNormDiff_q_qBar << endl;
622  std::ofstream m_file( "ConvergenceHistory.txt", std::ios_base::app);
623  m_file << m_stepCounter+1 << "\t"
624  << m_equ[m_nequ - 1]->GetFinalTime() << "\t"
625  << totalTime << "\t"
626  << MaxNormDiff_q_qBar << "\t"
627  << MaxNormDiff_q1_q0 << endl;
628  m_file.close();
629  }
630 
631  cpuTime = 0.0;
632  timer.Start();
633 }
634 
635 
637 {
638  cout << "\n====================================="
639  "=================================="
640  << endl;
641  cout << "Parameters for the SFD method:" << endl;
642  cout << "\tControl Coefficient: X = " << m_X << endl;
643  cout << "\tFilter Width: Delta = " << m_Delta << endl;
644  cout << "The simulation is stopped when |q-qBar|inf < " << TOL << endl;
646  {
647  cout << "\nWe use the adaptive SFD method:" << endl;
648  cout << " The parameters are updated every " << AdaptiveTime
649  << " time units;" << endl;
650  cout << " until |q-qBar|inf becomes smaller than " << AdaptiveTOL
651  << endl;
652  }
653  cout << "====================================="
654  "==================================\n" << endl;
655 }
656 
657 }
658 }
659 
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
boost::shared_ptr< AdvectionSystem > AdvectionSystemSharedPtr
Shared pointer to an AdvectionSystem class.
void ConvergenceHistory(const Array< OneD, const Array< OneD, NekDouble > > &qBar1, const Array< OneD, const Array< OneD, NekDouble > > &q0, NekDouble &MaxNormDiff_q_qBar, NekDouble &MaxNormDiff_q1_q0)
void Svtvp(int n, const T alpha, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
svtvp (scalar times vector plus vector): z = alpha*x + y
Definition: Vmath.cpp:471
STL namespace.
virtual SOLVER_UTILS_EXPORT void v_InitObject(std::ostream &out=std::cout)
Second-stage initialisation.
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
virtual void v_InitObject(std::ostream &out=std::cout)
Virtual function for initialisation implementation.
void ComputeSFD(const int i, const Array< OneD, const Array< OneD, NekDouble > > &q0, const Array< OneD, const Array< OneD, NekDouble > > &qBar0, Array< OneD, Array< OneD, NekDouble > > &q1, Array< OneD, Array< OneD, NekDouble > > &qBar1)
void Stop()
Definition: Timer.cpp:62
void ReadEVfile(int &KrylovSubspaceDim, NekDouble &growthEV, NekDouble &frequencyEV)
enum EvolutionOperatorType m_EvolutionOperator
Evolution Operator.
Definition: Driver.h:100
void GradientDescentMethod(const std::complex< NekDouble > &alpha, NekDouble &X_output, NekDouble &Delta_output)
double NekDouble
void EvalEV_ScalarSFD(const NekDouble &X_input, const NekDouble &Delta_input, const std::complex< NekDouble > &alpha, NekDouble &MaxEV)
virtual SOLVER_UTILS_EXPORT ~DriverSteadyState()
Destructor.
void SetSFDOperator(const NekDouble X_input, const NekDouble Delta_input)
Array< OneD, EquationSystemSharedPtr > m_equ
Equation system to solve.
Definition: Driver.h:94
NekDouble M11
Definition of the SFD operator.
DriverFactory & GetDriverFactory()
Definition: Driver.cpp:66
LibUtilities::CommSharedPtr m_comm
Communication object.
Definition: Driver.h:85
void Start()
Definition: Timer.cpp:51
NekDouble m_Delta
Definition of the SFD parameters:
virtual void v_Execute(std::ostream &out=std::cout)
Virtual function for solve implementation.
virtual SOLVER_UTILS_EXPORT void v_Execute(std::ostream &out=std::cout)
Virtual function for solve implementation.
A base class for PDEs which include an advection component.
NekDouble TimePerTest(unsigned int n)
Returns amount of seconds per iteration in a test with n iterations.
Definition: Timer.cpp:108
LibUtilities::SessionReaderSharedPtr m_session
Session reader object.
Definition: Driver.h:88
int m_nequ
number of equations
Definition: Driver.h:97
bool FlowPartiallyConverged
For adaptive SFD method.
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, tDescription pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:215