67 UnsteadySystem::UnsteadySystem(
85 m_session->MatchSolverInfo(
"DIFFUSIONADVANCEMENT",
"Explicit",
87 m_session->MatchSolverInfo(
"ADVECTIONADVANCEMENT",
"Explicit",
89 m_session->MatchSolverInfo(
"REACTIONADVANCEMENT",
"Explicit",
93 if (
m_session->DefinesSolverInfo(
"TIMEINTEGRATIONMETHOD"))
97 "TIMEINTEGRATIONMETHOD"));
105 boost::lexical_cast<std::string>(
m_time);
112 LibUtilities::FilterMap::const_iterator x;
114 for (x = f.begin(); x != f.end(); ++x)
142 TimeStability = 2.784;
164 "No CFL control implementation for this time"
165 "integration scheme");
168 return TimeStability;
181 int nfields =
m_fields.num_elements();
185 for (i = 0; i < nfields; ++i)
189 nvariables = nfields;
199 for(i = 0; i < nfields; ++i)
213 for(i = 0; i < nvariables; ++i)
235 "Final condition not unique: "
236 "fintime > 0.0 and Nsteps > 0");
240 "Timestep not unique: timestep > 0.0 & CFL > 0.0");
247 "Only one of IO_CheckTime and IO_CheckSteps "
251 bool doCheckTime =
false;
297 if (
m_session->GetComm()->GetRank() == 0 &&
300 cout <<
"Steps: " << setw(8) << left << step+1 <<
" "
301 <<
"Time: " << setw(12) << left <<
m_time;
305 cout <<
" Time-step: " << setw(12)
310 ss << cpuTime <<
"s";
311 cout <<
" CPU Time: " << setw(8) << left
317 for (i = 0; i < nvariables; ++i)
334 for (i = 0; i < nvariables; ++i)
336 if (
Vmath::Nnan(fields[i].num_elements(), fields[i], 1) > 0)
341 m_session->GetComm()->AllReduce(nanFound,
344 "NaN found during time integration.");
359 vector<bool> transformed(nfields,
false);
360 for(i = 0; i < nfields; i++)
368 transformed[i] =
true;
372 for(i = 0; i < nfields; i++)
396 if (
m_session->GetComm()->GetRank() == 0)
404 if (
m_session->GetSolverInfo(
"Driver") !=
"SteadyState")
406 cout <<
"Time-integration : " << intTime <<
"s" << endl;
413 for(i = 0; i < nfields; i++)
426 for(i = 0; i < nvariables; ++i)
466 if(
m_session->DefinesSolverInfo(
"AdvectionType"))
469 m_session->GetSolverInfo(
"AdvectionType"));
476 ==
"SteadyAdvectionDiffusionReaction")
505 outfile.open(
"solution1D.txt");
508 outfile << scientific << setw (17) << setprecision(16) << x[i]
509 <<
" " << solution1D[0][i] << endl;
511 outfile << endl << endl;
520 "This function is not implemented for this equation.");
529 "This function is not implemented for this equation.");
538 int nvariables =
m_fields.num_elements();
539 int nqvar = uflux.num_elements();
550 for (j = 0; j < nqvar; ++j)
552 for (i = 0; i < nvariables ; ++i)
555 m_fields[i]->GetFwdBwdTracePhys(ufield[i], Fwd, Bwd);
577 if(
m_fields[0]->GetBndCondExpansions().num_elements())
608 int nvariables =
m_fields.num_elements();
609 int nqvar = qfield.num_elements();
624 for (
int i = 0; i < nvariables; ++i)
627 for (
int j = 0; j < nqvar; ++j)
630 m_fields[i]->GetFwdBwdTracePhys(qfield[j][i],qFwd,qBwd);
654 m_fields[i]->GetFwdBwdTracePhys(ufield[i], Fwd, Bwd);
661 -1.0 * C11, uterm, 1,
671 if (
m_fields[0]->GetBndCondExpansions().num_elements())
692 if (
m_session->DefinesFunction(
"InitialConditions"))
694 for (
int i = 0; i <
m_fields.num_elements(); ++i)
699 "InitialConditions",
m_session->GetVariable(i));
705 "InitialConditions",
m_session->GetVariable(i));
707 fs::path pfilename(filename);
710 if(fs::is_directory(pfilename))
712 fs::path metafile(
"Info.xml");
713 fs::path fullpath = pfilename / metafile;
744 int i, e, npoints, id1, id2;
747 int nbnd =
m_fields[var]->GetBndCondExpansions().num_elements();
754 m_fields[var]->ExtractTracePhys(physfield, uplus);
755 for (i = 0; i < nbnd; ++i)
759 GetBndCondExpansions()[i]->GetExpSize();
763 m_session->GetFunction(
"InitialConditions", 0);
766 GetBndCondExpansions()[i]->GetNpoints();
773 m_fields[var]->GetBndCondExpansions()[i]->GetCoords(x0,x1,x2);
774 ifunc->Evaluate(x0,x1,x2,time,BDphysics);
777 for (e = 0; e < numBDEdge ; ++e)
781 GetBndCondExpansions()[i]->GetExp(e)->GetNumPoints(0);
784 GetBndCondExpansions()[i]->GetPhys_Offset(e);
788 GetBndCondTraceToGlobalTraceMap(cnt++));
791 if (
m_fields[var]->GetBndConditions()[i]->
796 &penaltyflux[id2], 1);
799 else if ((
m_fields[var]->GetBndConditions()[i])->
804 &penaltyflux[id2], 1);
823 int i, e, npoints, id1, id2;
824 int nbnd =
m_fields[var]->GetBndCondExpansions().num_elements();
831 m_fields[var]->ExtractTracePhys(physfield,qtemp);
833 for (i = 0; i < nbnd; ++i)
836 GetBndCondExpansions()[i]->GetExpSize();
840 m_session->GetFunction(
"InitialConditions", 0);
843 GetBndCondExpansions()[i]->GetNpoints();
850 m_fields[var]->GetBndCondExpansions()[i]->GetCoords(x0,x1,x2);
851 ifunc->Evaluate(x0,x1,x2,time,BDphysics);
854 for (e = 0; e < numBDEdge ; ++e)
857 GetBndCondExpansions()[i]->GetExp(e)->GetNumPoints(0);
860 GetBndCondExpansions()[i]->GetPhys_Offset(e);
864 GetBndCondTraceToGlobalTraceMap(cnt++));
868 if(
m_fields[var]->GetBndConditions()[i]->
874 &penaltyflux[id2], 1);
877 else if((
m_fields[var]->GetBndConditions()[i])->
883 &penaltyflux[id2], 1);
912 ASSERTL0(
false,
"Not defined for this class");
935 int phystot =
m_fields[0]->GetTotPoints();
936 int nvel = vel.num_elements();
942 for(
int n = 1; n < nvel; ++n)
944 Vmath::Vvtvp(phystot,vel[n],1,vel[n],1,varcoeff,1,varcoeff,1);
948 for(
int i = 0; i <
m_fields[0]->GetNumElmts(); ++i)
950 int offset =
m_fields[0]->GetPhys_Offset(i);
951 int nq =
m_fields[0]->GetExp(i)->GetTotPoints();
956 for(
int n = 0; n <
m_fields[0]->GetExp(i)->GetNumBases(); ++n)
959 m_fields[0]->GetExp(i)->GetBasisNumModes(n));
965 Vmath::Smul(nq,h,varcoeff+offset,1,tmp = varcoeff+offset,1);
virtual SOLVER_UTILS_EXPORT bool v_SteadyStateCheck(int step)
virtual SOLVER_UTILS_EXPORT void v_DoSolve()
Solves an unsteady problem.
virtual SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &l)
Virtual function for generating summary information.
#define ASSERTL0(condition, msg)
bool m_explicitReaction
Indicates if explicit or implicit treatment of reaction is used.
virtual SOLVER_UTILS_EXPORT void v_NumFluxforVector(const Array< OneD, Array< OneD, NekDouble > > &ufield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &qfield, Array< OneD, Array< OneD, NekDouble > > &qflux)
SOLVER_UTILS_EXPORT void CheckForRestartTime(NekDouble &time)
virtual SOLVER_UTILS_EXPORT bool v_PreIntegrate(int step)
A base class for describing how to solve specific equations.
bool m_homoInitialFwd
Flag to determine if simulation should start in homogeneous forward transformed state.
Adams-Bashforth Forward multi-step scheme of order 2.
NekDouble m_time
Current time of simulation.
bool m_explicitDiffusion
Indicates if explicit or implicit treatment of diffusion is used.
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
Runge-Kutta multi-stage scheme 4th order explicit (old name)
LibUtilities::TimeIntegrationSchemeOperators m_ode
The time integration scheme operators to use.
NekDouble m_timestep
Time step size.
std::vector< std::pair< std::string, std::string > > SummaryList
void WeakPenaltyforVector(const int var, const int dir, const Array< OneD, const NekDouble > &physfield, Array< OneD, NekDouble > &penaltyflux, NekDouble C11, NekDouble time=0.0)
Classical RungeKutta4 method (new name for eClassicalRungeKutta4)
SOLVER_UTILS_EXPORT void Checkpoint_Output(const int n)
Write checkpoint file of m_fields.
virtual SOLVER_UTILS_EXPORT void v_AppendOutput1D(Array< OneD, Array< OneD, NekDouble > > &solution1D)
Print the solution at each solution point in a txt file.
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
virtual SOLVER_UTILS_EXPORT NekDouble v_GetTimeStep(const Array< OneD, const Array< OneD, NekDouble > > &inarray)
Return the timestep to be used for the next step in the time-marching loop.
NekDouble m_checktime
Time between checkpoints.
const char *const TimeIntegrationMethodMap[]
SOLVER_UTILS_EXPORT NekDouble MaxTimeStepEstimator()
Get the maximum timestep estimator for cfl control.
Nonlinear SSP RungeKutta3 explicit.
int m_nchk
Number of checkpoints written so far.
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
int m_checksteps
Number of steps between checkpoints.
NekDouble m_fintime
Finish time of the simulation.
virtual SOLVER_UTILS_EXPORT void v_InitObject()
Initialisation object for EquationSystem.
SOLVER_UTILS_EXPORT void SVVVarDiffCoeff(const Array< OneD, Array< OneD, NekDouble > > vel, StdRegions::VarCoeffMap &varCoeffMap)
Evaluate the SVV diffusion coefficient according to Moura's paper where it should proportional to h t...
int m_steps
Number of steps to take.
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Array holding trace normals for DG simulations in the forwards direction.
virtual SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &s)
Print a summary of time stepping parameters.
Classical RungeKutta2 method (new name for eMidpoint)
static const NekDouble kNekZeroTol
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
bool m_explicitAdvection
Indicates if explicit or implicit treatment of advection is used.
virtual SOLVER_UTILS_EXPORT ~UnsteadySystem()
Destructor.
virtual SOLVER_UTILS_EXPORT bool v_PostIntegrate(int step)
Adams-Bashforth Forward multi-step scheme of order 1.
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
void AddSummaryItem(SummaryList &l, const std::string &name, const std::string &value)
Adds a summary item to the summary info list.
int Nnan(int n, const T *x, const int incx)
Return number of NaN elements of x.
std::vector< std::pair< std::string, FilterParams > > FilterMap
Nonlinear SSP RungeKutta2 explicit (surrogate for eRungeKutta2_ImprovedEuler)
Improved RungeKutta2 explicit (old name meaning Heun's method)
int m_spacedim
Spatial dimension (>= expansion dim).
LibUtilities::FieldMetaDataMap m_fieldMetaDataMap
Map to identify relevant solver info to dump in output fields.
std::string PortablePath(const boost::filesystem::path &path)
create portable path on different platforms for boost::filesystem path
TimeIntegrationWrapperFactory & GetTimeIntegrationWrapperFactory()
virtual SOLVER_UTILS_EXPORT void v_InitObject()
Init object for UnsteadySystem class.
SOLVER_UTILS_EXPORT NekDouble GetTimeStep()
NekDouble m_cflSafetyFactor
CFL safety factor (comprise between 0 to 1).
virtual SOLVER_UTILS_EXPORT void v_NumFluxforScalar(const Array< OneD, Array< OneD, NekDouble > > &ufield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &uflux)
boost::shared_ptr< Equation > EquationSharedPtr
SOLVER_UTILS_EXPORT void SetBoundaryConditions(NekDouble time)
Evaluates the boundary conditions at the given time.
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.
virtual SOLVER_UTILS_EXPORT void v_NumericalFlux(Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &numflux)
SOLVER_UTILS_EXPORT int GetPhys_Offset(int n)
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
SOLVER_UTILS_EXPORT int GetNpoints()
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
SOLVER_UTILS_EXPORT int GetTraceNpoints()
LibUtilities::TimeIntegrationWrapperSharedPtr m_intScheme
Wrapper to the time integration scheme.
LibUtilities::FieldIOSharedPtr m_fld
Field input/output.
int m_initialStep
Number of the step where the simulation should begin.
FilterFactory & GetFilterFactory()
std::vector< FilterSharedPtr > m_filters
int m_infosteps
Number of time steps between outputting status information.
SOLVER_UTILS_EXPORT void SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
Initialise the data in the dependent fields.
void WeakPenaltyforScalar(const int var, const Array< OneD, const NekDouble > &physfield, Array< OneD, NekDouble > &penaltyflux, NekDouble time=0.0)
midpoint method (old name)
NekDouble TimePerTest(unsigned int n)
Returns amount of seconds per iteration in a test with n iterations.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
static FieldMetaDataMap NullFieldMetaDataMap
LibUtilities::TimeIntegrationSolutionSharedPtr m_intSoln
virtual SOLVER_UTILS_EXPORT void v_DoInitialise()
Sets up initial conditions.
enum HomogeneousType m_HomogeneousType
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
std::vector< int > m_intVariables