Nektar++
Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Member Functions | Static Protected Attributes | Private Attributes | Friends | List of all members
Nektar::SolverUtils::DriverSteadyState Class Reference

#include <DriverSteadyState.h>

Inheritance diagram for Nektar::SolverUtils::DriverSteadyState:
Inheritance graph
[legend]
Collaboration diagram for Nektar::SolverUtils::DriverSteadyState:
Collaboration graph
[legend]

Public Member Functions

void PrintSummarySFD ()
 
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 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 EvalEV_ScalarSFD (const NekDouble &X_input, const NekDouble &Delta_input, const complex< NekDouble > &alpha, NekDouble &MaxEV)
 
void GradientDescentMethod (const complex< NekDouble > &alpha, NekDouble &X_output, NekDouble &Delta_output)
 
void ReadEVfile (int &KrylovSubspaceDim, NekDouble &growthEV, NekDouble &frequencyEV)
 
void ComputeOptimization ()
 
void SetSFDOperator (const NekDouble X_input, const NekDouble Delta_input)
 
- Public Member Functions inherited from Nektar::SolverUtils::DriverArnoldi
SOLVER_UTILS_EXPORT void ArnoldiSummary (std::ostream &out)
 
- Public Member Functions inherited from Nektar::SolverUtils::Driver
virtual ~Driver ()
 Destructor. More...
 
SOLVER_UTILS_EXPORT void InitObject (ostream &out=cout)
 Initialise Object. More...
 
SOLVER_UTILS_EXPORT void Execute (ostream &out=cout)
 Execute driver. More...
 
SOLVER_UTILS_EXPORT Array< OneD, EquationSystemSharedPtrGetEqu ()
 
SOLVER_UTILS_EXPORT Array< OneD, NekDoubleGetRealEvl (void)
 
SOLVER_UTILS_EXPORT Array< OneD, NekDoubleGetImagEvl (void)
 

Static Public Member Functions

static DriverSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession)
 Creates an instance of this class. More...
 
- Static Public Member Functions inherited from Nektar::SolverUtils::DriverModifiedArnoldi
static DriverSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession)
 Creates an instance of this class. More...
 

Static Public Attributes

static std::string className
 Name of the class. More...
 
- Static Public Attributes inherited from Nektar::SolverUtils::DriverModifiedArnoldi
static std::string className
 Name of the class. More...
 

Protected Member Functions

SOLVER_UTILS_EXPORT DriverSteadyState (const LibUtilities::SessionReaderSharedPtr pSession)
 Constructor. More...
 
virtual SOLVER_UTILS_EXPORT ~DriverSteadyState ()
 Destructor. More...
 
virtual SOLVER_UTILS_EXPORT void v_InitObject (ostream &out=cout)
 Second-stage initialisation. More...
 
virtual SOLVER_UTILS_EXPORT void v_Execute (ostream &out=cout)
 Virtual function for solve implementation. More...
 
- Protected Member Functions inherited from Nektar::SolverUtils::DriverModifiedArnoldi
 DriverModifiedArnoldi (const LibUtilities::SessionReaderSharedPtr pSession)
 Constructor. More...
 
virtual ~DriverModifiedArnoldi ()
 Destructor. More...
 
- Protected Member Functions inherited from Nektar::SolverUtils::DriverArnoldi
 DriverArnoldi (const LibUtilities::SessionReaderSharedPtr pSession)
 Constructor. More...
 
virtual ~DriverArnoldi ()
 Destructor. More...
 
void CopyArnoldiArrayToField (Array< OneD, NekDouble > &array)
 Copy Arnoldi storage to fields. More...
 
void CopyFieldToArnoldiArray (Array< OneD, NekDouble > &array)
 Copy fields to Arnoldi storage. More...
 
void CopyFwdToAdj ()
 Copy the forward field to the adjoint system in transient growth calculations. More...
 
void WriteFld (std::string file, std::vector< Array< OneD, NekDouble > > coeffs)
 Write coefficients to file. More...
 
void WriteFld (std::string file, Array< OneD, NekDouble > coeffs)
 
void WriteEvs (ostream &evlout, const int k, const NekDouble real, const NekDouble imag, NekDouble resid=NekConstants::kNekUnsetDouble, bool DumpInverse=true)
 
virtual Array< OneD, NekDoublev_GetRealEvl (void)
 
virtual Array< OneD, NekDoublev_GetImagEvl (void)
 
- Protected Member Functions inherited from Nektar::SolverUtils::Driver
 Driver (const LibUtilities::SessionReaderSharedPtr pSession)
 Initialises EquationSystem class members. More...
 

Static Protected Attributes

static std::string driverLookupId
 
- Static Protected Attributes inherited from Nektar::SolverUtils::Driver
static std::string evolutionOperatorLookupIds []
 
static std::string evolutionOperatorDef
 
static std::string driverDefault
 

Private Attributes

int m_stepCounter
 
int m_Check
 
int m_Check_BaseFlow
 
NekDouble TOL
 
int m_infosteps
 
int m_checksteps
 
int NumVar_SFD
 
Timer timer
 
NekDouble cpuTime
 
NekDouble totalTime
 
NekDouble elapsed
 
std::ofstream m_file
 
NekDouble m_Delta
 Definition of the SFD parameters: More...
 
NekDouble m_X
 
NekDouble m_dt
 
NekDouble M11
 Definition of the SFD operator. More...
 
NekDouble M12
 
NekDouble M21
 
NekDouble M22
 
NekDouble Diff_q_qBar
 
NekDouble Diff_q1_q0
 
bool FlowPartiallyConverged
 For adaptive SFD method. More...
 
NekDouble AdaptiveTOL
 
NekDouble AdaptiveTime
 
int m_NonConvergingStepsCounter
 
NekDouble GrowthRateEV
 
NekDouble FrequencyEV
 

Friends

class MemoryManager< DriverSteadyState >
 

Additional Inherited Members

- Protected Attributes inherited from Nektar::SolverUtils::DriverArnoldi
int m_kdim
 
int m_nvec
 Dimension of Krylov subspace. More...
 
int m_nits
 Number of vectors to test. More...
 
NekDouble m_evtol
 Maxmum number of iterations. More...
 
NekDouble m_period
 Tolerance of iteratiosn. More...
 
bool m_timeSteppingAlgorithm
 Period of time stepping algorithm. More...
 
int m_infosteps
 underlying operator is time stepping More...
 
int m_nfields
 interval to dump information if required. More...
 
NekDouble m_realShift
 
NekDouble m_imagShift
 
int m_negatedOp
 
Array< OneD, NekDoublem_real_evl
 Operator in solve call is negated. More...
 
Array< OneD, NekDoublem_imag_evl
 
- Protected Attributes inherited from Nektar::SolverUtils::Driver
LibUtilities::CommSharedPtr m_comm
 Communication object. More...
 
LibUtilities::SessionReaderSharedPtr m_session
 Session reader object. More...
 
LibUtilities::SessionReaderSharedPtr session_LinNS
 I the Coupling between SFD and arnoldi. More...
 
Array< OneD, EquationSystemSharedPtrm_equ
 Equation system to solve. More...
 
int m_nequ
 number of equations More...
 
enum EvolutionOperatorType m_EvolutionOperator
 Evolution Operator. More...
 

Detailed Description

Definition at line 49 of file DriverSteadyState.h.

Constructor & Destructor Documentation

Nektar::SolverUtils::DriverSteadyState::DriverSteadyState ( const LibUtilities::SessionReaderSharedPtr  pSession)
protected

Constructor.

Definition at line 54 of file DriverSteadyState.cpp.

56  : DriverModifiedArnoldi(pSession)
57 {
58 }
DriverModifiedArnoldi(const LibUtilities::SessionReaderSharedPtr pSession)
Constructor.
Nektar::SolverUtils::DriverSteadyState::~DriverSteadyState ( )
protectedvirtual

Destructor.

Definition at line 64 of file DriverSteadyState.cpp.

65 {
66 }

Member Function Documentation

void Nektar::SolverUtils::DriverSteadyState::ComputeOptimization ( )

Definition at line 356 of file DriverSteadyState.cpp.

References GradientDescentMethod(), m_Delta, Nektar::SolverUtils::DriverArnoldi::m_kdim, m_X, ReadEVfile(), and SetSFDOperator().

Referenced by v_Execute().

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 }
void ReadEVfile(int &KrylovSubspaceDim, NekDouble &growthEV, NekDouble &frequencyEV)
double NekDouble
void GradientDescentMethod(const complex< NekDouble > &alpha, NekDouble &X_output, NekDouble &Delta_output)
void SetSFDOperator(const NekDouble X_input, const NekDouble Delta_input)
NekDouble m_Delta
Definition of the SFD parameters:
void Nektar::SolverUtils::DriverSteadyState::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 
)

Encapsulated SFD method

Definition at line 336 of file DriverSteadyState.cpp.

References M11, M12, M21, M22, Nektar::SolverUtils::Driver::m_equ, Nektar::SolverUtils::Driver::m_nequ, and Vmath::Svtvp().

Referenced by v_Execute().

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 }
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
Array< OneD, EquationSystemSharedPtr > m_equ
Equation system to solve.
Definition: Driver.h:94
NekDouble M11
Definition of the SFD operator.
int m_nequ
number of equations
Definition: Driver.h:97
void Nektar::SolverUtils::DriverSteadyState::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 
)

This routine evaluates |q-qBar|_inf (and |q1-q0|_inf) and writes the values in "ConvergenceHistory.txt"

Definition at line 580 of file DriverSteadyState.cpp.

References cpuTime, elapsed, Nektar::SolverUtils::Driver::m_comm, m_Delta, Nektar::SolverUtils::Driver::m_equ, m_file, Nektar::SolverUtils::Driver::m_nequ, m_stepCounter, m_X, NumVar_SFD, Nektar::Timer::Start(), Nektar::Timer::Stop(), Nektar::Timer::TimePerTest(), timer, and totalTime.

Referenced by v_Execute().

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 }
void Stop()
Definition: Timer.cpp:62
Array< OneD, EquationSystemSharedPtr > m_equ
Equation system to solve.
Definition: Driver.h:94
LibUtilities::CommSharedPtr m_comm
Communication object.
Definition: Driver.h:85
void Start()
Definition: Timer.cpp:51
NekDouble m_Delta
Definition of the SFD parameters:
NekDouble TimePerTest(unsigned int n)
Returns amount of seconds per iteration in a test with n iterations.
Definition: Timer.cpp:108
int m_nequ
number of equations
Definition: Driver.h:97
static DriverSharedPtr Nektar::SolverUtils::DriverSteadyState::create ( const LibUtilities::SessionReaderSharedPtr pSession)
inlinestatic

Creates an instance of this class.

Definition at line 55 of file DriverSteadyState.h.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr().

56  {
58  ::AllocateSharedPtr(pSession);
59  p->InitObject();
60  return p;
61  }
boost::shared_ptr< Driver > DriverSharedPtr
A shared pointer to a Driver object.
Definition: Driver.h:52
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
void Nektar::SolverUtils::DriverSteadyState::EvalEV_ScalarSFD ( const NekDouble X_input,
const NekDouble Delta_input,
const complex< NekDouble > &  alpha,
NekDouble MaxEV 
)

This routine evaluates the maximum eigenvalue of the SFD system when applied to the 1D model u(n+1) = alpha*u(n)

Definition at line 472 of file DriverSteadyState.cpp.

Referenced by GradientDescentMethod().

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 }
double NekDouble
void Nektar::SolverUtils::DriverSteadyState::GradientDescentMethod ( const complex< NekDouble > &  alpha,
NekDouble X_output,
NekDouble Delta_output 
)

This routine implements a gradient descent method to find the parameters X end Delta which give the minimum eigenlavue of the SFD problem applied to the scalar case u(n+1) = *u(n).

Definition at line 386 of file DriverSteadyState.cpp.

References EvalEV_ScalarSFD().

Referenced by ComputeOptimization(), and v_Execute().

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 }
void EvalEV_ScalarSFD(const NekDouble &X_input, const NekDouble &Delta_input, const complex< NekDouble > &alpha, NekDouble &MaxEV)
double NekDouble
void Nektar::SolverUtils::DriverSteadyState::PrintSummarySFD ( )

Definition at line 634 of file DriverSteadyState.cpp.

References AdaptiveTime, AdaptiveTOL, Nektar::SolverUtils::eAdaptiveSFD, m_Delta, Nektar::SolverUtils::Driver::m_EvolutionOperator, m_X, and TOL.

Referenced by v_Execute().

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 }
enum EvolutionOperatorType m_EvolutionOperator
Evolution Operator.
Definition: Driver.h:100
NekDouble m_Delta
Definition of the SFD parameters:
void Nektar::SolverUtils::DriverSteadyState::ReadEVfile ( int &  KrylovSubspaceDim,
NekDouble growthEV,
NekDouble frequencyEV 
)

Definition at line 508 of file DriverSteadyState.cpp.

References ASSERTL0, and Nektar::SolverUtils::Driver::m_session.

Referenced by ComputeOptimization().

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 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
double NekDouble
LibUtilities::SessionReaderSharedPtr m_session
Session reader object.
Definition: Driver.h:88
void Nektar::SolverUtils::DriverSteadyState::SetSFDOperator ( const NekDouble  X_input,
const NekDouble  Delta_input 
)

This routine defines the encapsulated SFD operator with first-order splitting and exact resolution of the second subproblem. (See http://scitation.aip.org/content/aip/journal/pof2/26/3/10.1063/1.4867482 for details)

Definition at line 323 of file DriverSteadyState.cpp.

References M11, M12, M21, M22, and m_dt.

Referenced by ComputeOptimization(), and v_Execute().

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 }
double NekDouble
NekDouble M11
Definition of the SFD operator.
void Nektar::SolverUtils::DriverSteadyState::v_Execute ( ostream &  out = cout)
protectedvirtual

Virtual function for solve implementation.

Definition of variables used in this algorithm

Call the Navier-Stokes solver for one time step

Copy the current flow field into q0

Apply the linear operator to the outcome of the solver

Update qBar

Copy the output of the SFD method into the current flow field

Loop for the adaptive SFD method

We save the final solution into a .fld file

Reimplemented from Nektar::SolverUtils::DriverModifiedArnoldi.

Definition at line 78 of file DriverSteadyState.cpp.

References AdaptiveTime, AdaptiveTOL, ComputeOptimization(), ComputeSFD(), ConvergenceHistory(), cpuTime, Diff_q1_q0, Diff_q_qBar, Nektar::SolverUtils::eAdaptiveSFD, elapsed, FlowPartiallyConverged, FrequencyEV, GradientDescentMethod(), GrowthRateEV, m_Check, m_Check_BaseFlow, m_checksteps, Nektar::SolverUtils::Driver::m_comm, m_Delta, m_dt, Nektar::SolverUtils::Driver::m_equ, Nektar::SolverUtils::Driver::m_EvolutionOperator, m_file, m_infosteps, Nektar::SolverUtils::Driver::m_nequ, m_NonConvergingStepsCounter, Nektar::SolverUtils::Driver::m_session, m_stepCounter, m_X, NumVar_SFD, PrintSummarySFD(), Nektar::LibUtilities::ReduceSum, SetSFDOperator(), Nektar::Timer::Start(), timer, TOL, totalTime, and Nektar::SolverUtils::DriverModifiedArnoldi::v_Execute().

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)
118  AdvectionSystemSharedPtr A = m_equ[0]->as<AdvectionSystem>();
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 }
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 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)
virtual void v_Execute(ostream &out=cout)
Virtual function for solve implementation.
enum EvolutionOperatorType m_EvolutionOperator
Evolution Operator.
Definition: Driver.h:100
double NekDouble
void GradientDescentMethod(const complex< NekDouble > &alpha, NekDouble &X_output, NekDouble &Delta_output)
void SetSFDOperator(const NekDouble X_input, const NekDouble Delta_input)
Array< OneD, EquationSystemSharedPtr > m_equ
Equation system to solve.
Definition: Driver.h:94
LibUtilities::CommSharedPtr m_comm
Communication object.
Definition: Driver.h:85
void Start()
Definition: Timer.cpp:51
NekDouble m_Delta
Definition of the SFD parameters:
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.
void Nektar::SolverUtils::DriverSteadyState::v_InitObject ( ostream &  out = cout)
protectedvirtual

Second-stage initialisation.

Reimplemented from Nektar::SolverUtils::DriverModifiedArnoldi.

Definition at line 72 of file DriverSteadyState.cpp.

References Nektar::SolverUtils::DriverModifiedArnoldi::v_InitObject().

73 {
75 }
virtual void v_InitObject(ostream &out=cout)
Virtual function for initialisation implementation.

Friends And Related Function Documentation

friend class MemoryManager< DriverSteadyState >
friend

Definition at line 52 of file DriverSteadyState.h.

Member Data Documentation

NekDouble Nektar::SolverUtils::DriverSteadyState::AdaptiveTime
private

Definition at line 152 of file DriverSteadyState.h.

Referenced by PrintSummarySFD(), and v_Execute().

NekDouble Nektar::SolverUtils::DriverSteadyState::AdaptiveTOL
private

Definition at line 151 of file DriverSteadyState.h.

Referenced by PrintSummarySFD(), and v_Execute().

string Nektar::SolverUtils::DriverSteadyState::className
static
Initial value:

Name of the class.

Definition at line 64 of file DriverSteadyState.h.

NekDouble Nektar::SolverUtils::DriverSteadyState::cpuTime
private

Definition at line 130 of file DriverSteadyState.h.

Referenced by ConvergenceHistory(), and v_Execute().

NekDouble Nektar::SolverUtils::DriverSteadyState::Diff_q1_q0
private

Definition at line 147 of file DriverSteadyState.h.

Referenced by v_Execute().

NekDouble Nektar::SolverUtils::DriverSteadyState::Diff_q_qBar
private

Definition at line 146 of file DriverSteadyState.h.

Referenced by v_Execute().

string Nektar::SolverUtils::DriverSteadyState::driverLookupId
staticprotected
Initial value:

Definition at line 118 of file DriverSteadyState.h.

NekDouble Nektar::SolverUtils::DriverSteadyState::elapsed
private

Definition at line 132 of file DriverSteadyState.h.

Referenced by ConvergenceHistory(), and v_Execute().

bool Nektar::SolverUtils::DriverSteadyState::FlowPartiallyConverged
private

For adaptive SFD method.

Definition at line 150 of file DriverSteadyState.h.

Referenced by v_Execute().

NekDouble Nektar::SolverUtils::DriverSteadyState::FrequencyEV
private

Definition at line 155 of file DriverSteadyState.h.

Referenced by v_Execute().

NekDouble Nektar::SolverUtils::DriverSteadyState::GrowthRateEV
private

Definition at line 154 of file DriverSteadyState.h.

Referenced by v_Execute().

NekDouble Nektar::SolverUtils::DriverSteadyState::M11
private

Definition of the SFD operator.

Definition at line 141 of file DriverSteadyState.h.

Referenced by ComputeSFD(), and SetSFDOperator().

NekDouble Nektar::SolverUtils::DriverSteadyState::M12
private

Definition at line 142 of file DriverSteadyState.h.

Referenced by ComputeSFD(), and SetSFDOperator().

NekDouble Nektar::SolverUtils::DriverSteadyState::M21
private

Definition at line 143 of file DriverSteadyState.h.

Referenced by ComputeSFD(), and SetSFDOperator().

NekDouble Nektar::SolverUtils::DriverSteadyState::M22
private

Definition at line 144 of file DriverSteadyState.h.

Referenced by ComputeSFD(), and SetSFDOperator().

int Nektar::SolverUtils::DriverSteadyState::m_Check
private

Definition at line 122 of file DriverSteadyState.h.

Referenced by v_Execute().

int Nektar::SolverUtils::DriverSteadyState::m_Check_BaseFlow
private

Definition at line 123 of file DriverSteadyState.h.

Referenced by v_Execute().

int Nektar::SolverUtils::DriverSteadyState::m_checksteps
private

Definition at line 126 of file DriverSteadyState.h.

Referenced by v_Execute().

NekDouble Nektar::SolverUtils::DriverSteadyState::m_Delta
private

Definition of the SFD parameters:

Definition at line 136 of file DriverSteadyState.h.

Referenced by ComputeOptimization(), ConvergenceHistory(), PrintSummarySFD(), and v_Execute().

NekDouble Nektar::SolverUtils::DriverSteadyState::m_dt
private

Definition at line 138 of file DriverSteadyState.h.

Referenced by SetSFDOperator(), and v_Execute().

std::ofstream Nektar::SolverUtils::DriverSteadyState::m_file
private

Definition at line 133 of file DriverSteadyState.h.

Referenced by ConvergenceHistory(), and v_Execute().

int Nektar::SolverUtils::DriverSteadyState::m_infosteps
private

Definition at line 125 of file DriverSteadyState.h.

Referenced by v_Execute().

int Nektar::SolverUtils::DriverSteadyState::m_NonConvergingStepsCounter
private

Definition at line 153 of file DriverSteadyState.h.

Referenced by v_Execute().

int Nektar::SolverUtils::DriverSteadyState::m_stepCounter
private

Definition at line 121 of file DriverSteadyState.h.

Referenced by ConvergenceHistory(), and v_Execute().

NekDouble Nektar::SolverUtils::DriverSteadyState::m_X
private
int Nektar::SolverUtils::DriverSteadyState::NumVar_SFD
private

Definition at line 127 of file DriverSteadyState.h.

Referenced by ConvergenceHistory(), and v_Execute().

Timer Nektar::SolverUtils::DriverSteadyState::timer
private

Definition at line 129 of file DriverSteadyState.h.

Referenced by ConvergenceHistory(), and v_Execute().

NekDouble Nektar::SolverUtils::DriverSteadyState::TOL
private

Definition at line 124 of file DriverSteadyState.h.

Referenced by PrintSummarySFD(), and v_Execute().

NekDouble Nektar::SolverUtils::DriverSteadyState::totalTime
private

Definition at line 131 of file DriverSteadyState.h.

Referenced by ConvergenceHistory(), and v_Execute().