Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
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 std::complex< NekDouble > &alpha, NekDouble &MaxEV)
 
void GradientDescentMethod (const std::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 (std::ostream &out=std::cout)
 Initialise Object. More...
 
SOLVER_UTILS_EXPORT void Execute (std::ostream &out=std::cout)
 Execute driver. More...
 
SOLVER_UTILS_EXPORT Array
< OneD,
EquationSystemSharedPtr
GetEqu ()
 
SOLVER_UTILS_EXPORT Array
< OneD, NekDouble
GetRealEvl (void)
 
SOLVER_UTILS_EXPORT Array
< OneD, NekDouble
GetImagEvl (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 (std::ostream &out=std::cout)
 Second-stage initialisation. More...
 
virtual SOLVER_UTILS_EXPORT void v_Execute (std::ostream &out=std::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 (std::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,
EquationSystemSharedPtr
m_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 56 of file DriverSteadyState.cpp.

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

Destructor.

Definition at line 66 of file DriverSteadyState.cpp.

67 {
68 }

Member Function Documentation

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

Definition at line 358 of file DriverSteadyState.cpp.

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

Referenced by v_Execute().

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 }
void ReadEVfile(int &KrylovSubspaceDim, NekDouble &growthEV, NekDouble &frequencyEV)
void GradientDescentMethod(const std::complex< NekDouble > &alpha, NekDouble &X_output, NekDouble &Delta_output)
double NekDouble
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 338 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().

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 }
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 582 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().

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 }
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 std::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 474 of file DriverSteadyState.cpp.

Referenced by GradientDescentMethod().

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 }
double NekDouble
void Nektar::SolverUtils::DriverSteadyState::GradientDescentMethod ( const std::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 388 of file DriverSteadyState.cpp.

References EvalEV_ScalarSFD().

Referenced by ComputeOptimization(), and v_Execute().

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

Definition at line 636 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().

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 }
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 510 of file DriverSteadyState.cpp.

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

Referenced by ComputeOptimization().

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 }
#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 325 of file DriverSteadyState.cpp.

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

Referenced by ComputeOptimization(), and v_Execute().

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

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)
120  AdvectionSystemSharedPtr A = m_equ[0]->as<AdvectionSystem>();
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 }
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)
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 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:
virtual void v_Execute(std::ostream &out=std::cout)
Virtual function for solve implementation.
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 ( std::ostream &  out = std::cout)
protectedvirtual

Second-stage initialisation.

Reimplemented from Nektar::SolverUtils::DriverModifiedArnoldi.

Definition at line 74 of file DriverSteadyState.cpp.

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

75 {
77 }
virtual void v_InitObject(std::ostream &out=std::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().