Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties 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 
89  m_session->LoadParameter("IO_InfoSteps", m_infosteps, 1000);
90  m_session->LoadParameter("IO_CheckSteps", m_checksteps, 100000);
91  m_session->LoadParameter("ControlCoeff",m_X, 1);
92  m_session->LoadParameter("FilterWidth", m_Delta, 2);
93 
94  // To evaluate optimum SFD parameters if growth rate provided in the
95  // xml file
96  m_session->LoadParameter("GrowthRateEV", GrowthRateEV, 0.0);
97 
98  // To evaluate optimum SFD parameters if frequency provided in the xml
99  // file
100  m_session->LoadParameter("FrequencyEV", FrequencyEV, 0.0);
101 
102  // Determine when SFD method is converged
103  m_session->LoadParameter("TOL", TOL, 1.0e-08);
104 
105  // Used only for the Adaptive SFD method
106  m_session->LoadParameter("AdaptiveTOL", AdaptiveTOL, 1.0e-02);
107 
108  // Used only for the Adaptive SFD method
109  m_session->LoadParameter("AdaptiveTime", AdaptiveTime, 50.0*m_Delta);
110 
111  if (m_comm->GetRank() == 0)
112  {
113  PrintSummarySFD();
114  }
115 
116  timer.Start();
117 
118  // Definition of shared pointer (used only for the Adaptive SFD method)
120 
121  // Condition necessary to run SFD for the compressible case
122  NumVar_SFD = m_equ[m_nequ - 1]->UpdateFields()[0]->GetCoordim(0);
123  if (m_session->GetSolverInfo("EqType") == "EulerCFE" ||
124  m_session->GetSolverInfo("EqType") == "NavierStokesCFE")
125  {
126  // Number of variables for the compressible equations
127  NumVar_SFD += 2;
128  }
129  if(m_session->DefinesSolverInfo("HOMOGENEOUS"))
130  {
131  if (m_session->GetSolverInfo("HOMOGENEOUS") == "1D")
132  {
133  NumVar_SFD += 1;
134  }
135  }
136 
137  // We store the time step
138  m_dt = m_equ[m_nequ - 1]->GetTimeStep();
139 
140  // Evaluate optimum SFD parameters if dominent EV given by xml file
141  if (GrowthRateEV != 0.0 && FrequencyEV != 0.0)
142  {
143  cout << "Besed on the dominant EV given in the xml file,"
144  << "a 1D model is used to evaluate the optumum parameters"
145  << "of the SFD method:" << endl;
146  complex<NekDouble> EV = polar(exp(GrowthRateEV), FrequencyEV);
148  }
149 
150 
151  // We set up the elements of the operator of the encapsulated
152  // formulation of the selective frequencive damping method
154 
155  // m_steps is set to 1. Then "m_equ[m_nequ - 1]->DoSolve()" will run
156  // for only one time step
157  m_equ[m_nequ - 1]->SetSteps(1);
158  ofstream m_file("ConvergenceHistory.txt", ios::out | ios::trunc);
159 
160  Array<OneD, Array<OneD, NekDouble> > q0(NumVar_SFD);
161  Array<OneD, Array<OneD, NekDouble> > q1(NumVar_SFD);
162  Array<OneD, Array<OneD, NekDouble> > qBar0(NumVar_SFD);
163  Array<OneD, Array<OneD, NekDouble> > qBar1(NumVar_SFD);
164 
165  for(int i = 0; i < NumVar_SFD; ++i)
166  {
167  q0[i] = Array<OneD, NekDouble> (m_equ[m_nequ-1]->GetTotPoints(),
168  0.0); //q0 is initialised
169  qBar0[i] = Array<OneD, NekDouble> (m_equ[m_nequ-1]->GetTotPoints(),
170  0.0);
171  m_equ[m_nequ - 1]->CopyFromPhysField(i, qBar0[i]);
172  }
173 
174  ///Definition of variables used in this algorithm
175  m_stepCounter = 0;
176  m_Check = 0;
177  m_Check_BaseFlow = 1;
179  Diff_q_qBar = 1.0;
180  Diff_q1_q0 = 1.0;
181  cpuTime = 0.0;
182  elapsed = 0.0;
183  totalTime = 0.0;
184  FlowPartiallyConverged = false;
185 
186  while (max(Diff_q_qBar, Diff_q1_q0) > TOL)
187  {
188  ///Call the Navier-Stokes solver for one time step
189  m_equ[m_nequ - 1]->DoSolve();
190 
191  for(int i = 0; i < NumVar_SFD; ++i)
192  {
193  ///Copy the current flow field into q0
194  m_equ[m_nequ - 1]->CopyFromPhysField(i, q0[i]);
195 
196  ///Apply the linear operator to the outcome of the solver
197  ComputeSFD(i, q0, qBar0, q1, qBar1);
198 
199  ///Update qBar
200  qBar0[i] = qBar1[i];
201 
202  ///Copy the output of the SFD method into the current flow field
203  m_equ[m_nequ - 1]->CopyToPhysField(i, q1[i]);
204  }
205 
207  {
209 
210  ///Loop for the adaptive SFD method
212  FlowPartiallyConverged == false)
213  {
214 
216 
217  if (Diff_q_qBar < AdaptiveTOL)
218  {
219  if (m_comm->GetRank() == 0)
220  {
221  cout << "\n\t The SFD method is converging: we compute "
222  << "stability analysis using the 'partially "
223  << "converged' steady state as base flow:\n" << endl;
224  }
225 
226  m_equ[m_nequ - 1]->Checkpoint_BaseFlow(m_Check_BaseFlow);
228 
229  A->GetAdvObject()->SetBaseFlow(q0,m_equ[0]->UpdateFields());
231 
232  if (m_comm->GetRank() == 0)
233  {
234  // Compute the update of the SFD parameters only on
235  // one processor
237  }
238  else
239  {
240  // On all the other processors, the parameters are set
241  // to 0
242  m_X = 0;
243  m_Delta = 0;
244  }
245  // The we give to all the processors the value of X and
246  // Delta of the first processor
249 
250  FlowPartiallyConverged = true;
251  }
253  >= AdaptiveTime)
254  {
255  if (m_comm->GetRank() == 0)
256  {
257  cout << "\n\t We compute stability analysis using"
258  << " the current flow field as base flow:\n"
259  << endl;
260  }
261 
262  m_equ[m_nequ - 1]->Checkpoint_BaseFlow(m_Check_BaseFlow);
264 
265  A->GetAdvObject()->SetBaseFlow(q0,m_equ[0]->UpdateFields());
267 
268  if (m_comm->GetRank() == 0)
269  {
270  // Compute the update of the SFD parameters only on
271  // one processor
273  }
274  else
275  {
276  // On all the other processors, the parameters are set
277  // to 0
278  m_X = 0;
279  m_Delta = 0;
280  }
281  // The we give to all the processors the value of X and
282  // Delta of the first processor
285 
287  }
288  }
289  }
290 
292  {
293  m_Check++;
294  m_equ[m_nequ - 1]->Checkpoint_Output(m_Check);
295  }
296  m_stepCounter++;
297  }
298 
299  m_file.close();
300 
301  ///We save the final solution into a .fld file
302  m_equ[m_nequ - 1]->Output();
303 
304  for(int j = 0; j < m_equ[m_nequ - 1]->GetNvariables(); ++j)
305  {
306  NekDouble vL2Error = m_equ[m_nequ - 1]->L2Error(j,false);
307  NekDouble vLinfError = m_equ[m_nequ - 1]->LinfError(j);
308  if (m_comm->GetRank() == 0)
309  {
310  out << "L 2 error (variable " << m_equ[m_nequ - 1]->GetVariable(j)
311  << ") : " << vL2Error << endl;
312  out << "L inf error (variable " << m_equ[m_nequ - 1]->GetVariable(j)
313  << ") : " << vLinfError << endl;
314  }
315  }
316 }
317 
318 
319 /**
320  * This routine defines the encapsulated SFD operator with first-order
321  * splitting and exact resolution of the second subproblem.
322  *(See http://scitation.aip.org/content/aip/journal/pof2/26/3/10.1063/1.4867482 for details)
323  */
325  const NekDouble Delta_input)
326 {
327  NekDouble X = X_input*m_dt;
328  NekDouble Delta = Delta_input/m_dt;
329  NekDouble coeff = 1.0/(1.0 + X*Delta);
330  M11 = coeff*(1.0 + X*Delta*exp(-(X + 1.0/Delta)));
331  M12 = coeff*(X*Delta*(1.0 - exp(-(X + 1.0/Delta))));
332  M21 = coeff*(1.0 - exp(-(X + 1.0/Delta)));
333  M22 = coeff*(X*Delta + exp(-(X + 1.0/Delta)));
334 }
335 
336 
338  const Array<OneD, const Array<OneD, NekDouble> > &q0,
339  const Array<OneD, const Array<OneD, NekDouble> > &qBar0,
342 {
343  q1[i] = Array<OneD, NekDouble> (m_equ[m_nequ - 1]->GetTotPoints(),0.0);
344  qBar1[i] = Array<OneD, NekDouble> (m_equ[m_nequ - 1]->GetTotPoints(),0.0);
345 
346  ///Encapsulated SFD method
347  Vmath::Svtvp(q1[i].num_elements(), M11, q0[i], 1, q1[i], 1, q1[i], 1 );
348  Vmath::Svtvp(q1[i].num_elements(), M12, qBar0[i], 1, q1[i], 1, q1[i], 1 );
349 
350  Vmath::Svtvp(qBar1[i].num_elements(), M21, q0[i], 1, qBar1[i], 1,
351  qBar1[i], 1 );
352  Vmath::Svtvp(qBar1[i].num_elements(), M22, qBar0[i], 1, qBar1[i], 1,
353  qBar1[i], 1 );
354 }
355 
356 
358 {
359  NekDouble growthEV(0.0);
360  NekDouble frequencyEV(0.0);
361 
362  // m_kdim is the dimension of Krylov subspace (defined in the xml file and
363  // used in DriverArnoldi.cpp)
364  ReadEVfile(m_kdim, growthEV, frequencyEV);
365 
366  cout << "\n\tgrowthEV = " << growthEV << endl;
367  cout << "\tfrequencyEV = " << frequencyEV << endl;
368 
369  complex<NekDouble> ApproxEV = polar(exp(growthEV), frequencyEV);
370 
371  NekDouble X_new = m_X;
372  NekDouble Delta_new = m_Delta;
373 
374  GradientDescentMethod(ApproxEV, X_new, Delta_new);
375 
376  m_X = X_new;
377  m_Delta = Delta_new;
378 
380 }
381 
382 /**
383  * This routine implements a gradient descent method to find the parameters X
384  * end Delta which give the minimum eigenlavue of the SFD problem applied to
385  * the scalar case u(n+1) = \alpha*u(n).
386  */
388  const complex<NekDouble> &alpha,
389  NekDouble &X_output,
390  NekDouble &Delta_output)
391 {
392  cout << "\n\tWe enter the Gradient Descent Method [...]" << endl;
393  bool OptParmFound = false;
394  bool Descending = true;
395  NekDouble X_input = X_output;
396  NekDouble Delta_input = Delta_output;
397 
398  NekDouble X_init = X_output;
399  NekDouble Delta_init = Delta_output;
400  int stepCounter(0);
401 
402  NekDouble F0(0.0);
403  NekDouble F0x(0.0);
404  NekDouble F0y(0.0);
405  NekDouble F1(0.0);
406  NekDouble dx = 0.00000001;
407  NekDouble dirx(0.0);
408  NekDouble diry(0.0);
409  NekDouble s(0.0);
410  NekDouble CurrentMin = 1.0;
411 
412  while (OptParmFound == false)
413  {
414  Descending = true;
415  EvalEV_ScalarSFD(X_input, Delta_input, alpha, F0);
416  EvalEV_ScalarSFD(X_input + dx, Delta_input, alpha, F0x);
417  EvalEV_ScalarSFD(X_input, Delta_input + dx, alpha, F0y);
418  dirx = (F0 - F0x);
419  diry = (F0 - F0y);
420  s = abs(0.000001/dirx);
421  X_output = X_input + s*dirx;
422  Delta_output = Delta_input + s*diry;
423  F1 = F0;
424 
425  while (Descending == true)
426  {
427  CurrentMin = F1;
428  X_input = X_output;
429  Delta_input = Delta_output;
430  EvalEV_ScalarSFD(X_output, Delta_output, alpha, F1);
431 
432  if (F1 > CurrentMin)
433  {
434  Descending = false;
435  }
436  else
437  {
438  s = s + s*0.01;
439  X_output = X_input + s*dirx;
440  Delta_output = Delta_input + s*diry;
441  }
442 
443  if(stepCounter > 9999999)
444  {
445  //We are stuck in this loop..
446  //Then we restart it with different initail conditions
447  Descending = false;
448  X_input = X_init;
449  Delta_init = Delta_init + Delta_init*0.1;
450  Delta_input = Delta_init;
451  stepCounter = 0;
452  }
453  stepCounter++;
454  }
455 
456  if (abs(F0-F1) < dx)
457  {
458  cout << "\tThe Gradient Descent Method has converged!" << endl;
459  EvalEV_ScalarSFD(X_output, Delta_output, alpha, F1);
460  cout << "\n\tThe updated parameters are: X_tilde = " << X_output
461  << " and Delta_tilde = " << Delta_output << endl;
462  OptParmFound = true;
463  }
464 
465  }
466 }
467 
468 
469 /**
470  * This routine evaluates the maximum eigenvalue of the SFD system when applied
471  * to the 1D model u(n+1) = alpha*u(n)
472  */
474  const NekDouble &X_input,
475  const NekDouble &Delta_input,
476  const complex<NekDouble> &alpha,
477  NekDouble &MaxEV)
478 {
479  NekDouble A11 = ( 1.0 + X_input * Delta_input *
480  exp(-(X_input + 1.0/Delta_input)) )/(1.0 + X_input*Delta_input);
481  NekDouble A12 = ( X_input*Delta_input - X_input * Delta_input *
482  exp(-(X_input + 1.0/Delta_input)) )/(1.0 + X_input*Delta_input);
483  NekDouble A21 = ( 1.0 - 1.0 *
484  exp(-(X_input + 1.0/Delta_input)) )/(1 + X_input*Delta_input);
485  NekDouble A22 = ( X_input*Delta_input + 1.0 *
486  exp(-(X_input + 1.0/Delta_input)) )/(1.0 + X_input*Delta_input);
487 
488  complex<NekDouble> B11 = alpha;
489  NekDouble B12 = 0.0;
490  NekDouble B21 = 0.0;
491  NekDouble B22 = 1.0;
492 
493  complex<NekDouble> a = A11*B11 + A12*B21;
494  NekDouble b = A11*B12 + A12*B22;
495  complex<NekDouble> c = A21*B11 + A22*B21;
496  NekDouble d = A21*B12 + A22*B22;
497 
498  complex<NekDouble> delt = sqrt((a-d)*(a-d) + 4.0*b*c);
499  complex<NekDouble> lambda_1 = 0.5*(a+d + delt);
500  complex<NekDouble> lambda_2 = 0.5*(a+d - delt);
501 
502  NekDouble NORM_1 = abs(lambda_1);
503  NekDouble NORM_2 = abs(lambda_2);
504 
505  MaxEV = max(NORM_1, NORM_2);
506 }
507 
508 
510  int &KrylovSubspaceDim,
511  NekDouble &growthEV,
512  NekDouble &frequencyEV)
513 {
514  // This routine reads the .evl file written by the Arnoldi algorithm
515  // (written in September 2014)
516  std::string line;
517  int NumLinesInFile = 0;
518  std::string EVfileName = m_session->GetSessionName() + ".evl";
519  std::ifstream EVfile(EVfileName.c_str());
520  ASSERTL0(EVfile.good(), "Cannot open .evl file.");
521 
522  if(EVfile)
523  {
524  // This block counts the total number of lines of the .evl file
525  // We keep going util we reach the end of the file
526  while(getline(EVfile, line))
527  {
528  ++NumLinesInFile;
529  }
530 
531  // It may happen that the Stability method that have produced the .elv
532  // file converges in less than m_kdim iterations. In this case,
533  // KrylovSubspaceDim has to be changed here
534  if(NumLinesInFile < KrylovSubspaceDim*2.0
535  + KrylovSubspaceDim*(KrylovSubspaceDim+1.0)/2.0)
536  {
537  for(int i = 1; i <= KrylovSubspaceDim; ++i)
538  {
539  if(NumLinesInFile == i*2.0 + i*(i+1.0)/2.0)
540  {
541  KrylovSubspaceDim = i;
542  }
543  }
544  }
545 
546  // go back to the beginning of the file
547  EVfile.clear();
548  EVfile.seekg(0, ios::beg);
549 
550  // We now want to go to the line where the most unstable eigenlavue was
551  // written
552  for(int i = 0; i < (NumLinesInFile - KrylovSubspaceDim); ++i)
553  {
554  std::getline(EVfile, line);
555  cout << "Discard line: " << line << endl;
556  }
557 
558  std::vector<std::string> tokens;
559  std::getline(EVfile, line);
560  boost::algorithm::split(tokens, line,
561  boost::is_any_of("\t "),boost::token_compress_on);
562 
563  ASSERTL0(tokens.size() >= 5,
564  "Unexpected formatting of .evl file while reading line:\n"
565  + line);
566  growthEV = boost::lexical_cast<NekDouble>(tokens[4]);
567  frequencyEV = boost::lexical_cast<NekDouble>(tokens[5]);
568  }
569  else
570  {
571  cout << "An error occurred when opening the .evl file" << endl;
572  }
573  EVfile.close();
574 }
575 
576 
577 /**
578  * This routine evaluates |q-qBar|_inf (and |q1-q0|_inf) and writes the values
579  * in "ConvergenceHistory.txt"
580  */
582  const Array<OneD, const Array<OneD, NekDouble> > &qBar1,
583  const Array<OneD, const Array<OneD, NekDouble> > &q0,
584  NekDouble &MaxNormDiff_q_qBar,
585  NekDouble &MaxNormDiff_q1_q0)
586 {
587  Array<OneD, NekDouble > NormDiff_q_qBar(NumVar_SFD, 1.0);
588  Array<OneD, NekDouble > NormDiff_q1_q0(NumVar_SFD, 1.0);
589  MaxNormDiff_q_qBar=0.0;
590  MaxNormDiff_q1_q0=0.0;
591 
592  for(int i = 0; i < NumVar_SFD; ++i)
593  {
594  NormDiff_q_qBar[i] = m_equ[m_nequ - 1]->LinfError(i, qBar1[i]);
595  NormDiff_q1_q0[i] = m_equ[m_nequ - 1]->LinfError(i, q0[i]);
596 
597  if (MaxNormDiff_q_qBar < NormDiff_q_qBar[i])
598  {
599  MaxNormDiff_q_qBar = NormDiff_q_qBar[i];
600  }
601  if (MaxNormDiff_q1_q0 < NormDiff_q1_q0[i])
602  {
603  MaxNormDiff_q1_q0 = NormDiff_q1_q0[i];
604  }
605  }
606 
607  timer.Stop();
609  cpuTime += elapsed;
610  totalTime += elapsed;
611 
612  if (m_comm->GetRank() == 0)
613  {
614  cout << "SFD - Step: " << left << m_stepCounter+1
615  << ";\tTime: " << left << m_equ[m_nequ - 1]->GetFinalTime()
616  << ";\tCPU time = " << left << cpuTime << " s"
617  << ";\tTot time = " << left << totalTime << " s"
618  << ";\tX = " << left << m_X
619  << ";\tDelta = " << left << m_Delta
620  << ";\t|q-qBar|inf = " << left << MaxNormDiff_q_qBar << endl;
621  std::ofstream m_file( "ConvergenceHistory.txt", std::ios_base::app);
622  m_file << m_stepCounter+1 << "\t"
623  << m_equ[m_nequ - 1]->GetFinalTime() << "\t"
624  << totalTime << "\t"
625  << MaxNormDiff_q_qBar << "\t"
626  << MaxNormDiff_q1_q0 << endl;
627  m_file.close();
628  }
629 
630  cpuTime = 0.0;
631  timer.Start();
632 }
633 
634 
636 {
637  cout << "\n====================================="
638  "=================================="
639  << endl;
640  cout << "Parameters for the SFD method:" << endl;
641  cout << "\tControl Coefficient: X = " << m_X << endl;
642  cout << "\tFilter Width: Delta = " << m_Delta << endl;
643  cout << "The simulation is stopped when |q-qBar|inf < " << TOL << endl;
645  {
646  cout << "\nWe use the adaptive SFD method:" << endl;
647  cout << " The parameters are updated every " << AdaptiveTime
648  << " time units;" << endl;
649  cout << " until |q-qBar|inf becomes smaller than " << AdaptiveTOL
650  << endl;
651  }
652  cout << "====================================="
653  "==================================\n" << endl;
654 }
655 
656 }
657 }
658 
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
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:485
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