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