Commit ecba1490 authored by Dave Moxey's avatar Dave Moxey

Merge remote-tracking branch 'upstream/master' into fix/mc-face-curve

parents 28a36772 5177e572
......@@ -320,10 +320,10 @@ Possible values are:
{Equations} & {\texttt{EQTYPE}} &{Dim.}&{Projections} & Alg.\\
\midrule
Steady Stokes (SS)& \texttt{SteadyStokes} & All & CG &VCS \\
Steady Onseen (SO) & \texttt{SteadyOseen} & All & CG& DS \\
Steady Oseen (SO) & \texttt{SteadyOseen} & All & CG& DS \\
Unsteady Stokes (US) & \texttt{UnsteadyStokes} & All & CG &VCS \\
Steady Linearised NS (SLNS) & \texttt{SteadyLinearisedNS} & All & CG & DS \\
Unsteady Linearised NS (ULNS) & \texttt{UnsteadyLinearisedNS} & All & CG & DS \\
Unsteady Linearised NS (ULNS) & \texttt{UnsteadyLinearisedNS} & All & CG & VCS,DS \\
Unsteady NS (UNS) & \texttt{UnsteadyNavierStokes} & All & CG,CG-DG & VCS \\
\bottomrule
\end{tabular}
......@@ -676,30 +676,64 @@ The following parameters can be specified in the \texttt{PARAMETERS} section of
\item \inltt{Re}: sets the Reynolds number
\item \inltt{Kinvis}: sets the kinematic viscosity $\nu$.
\item \inltt{kdim}: sets the dimension of the Krylov subspace $\kappa$. Can be used in: \inltt{ModifiedArnoldi} and \inltt{Arpack}. Default value: 16.
\item \inltt{evtol}: sets the tolerance of the eigenvalues. Can be used in: \inltt{ModifiedArnoldi} and \inltt{Arpack}. Default value: $10^{-6}$.
\item \inltt{evtol}: sets the tolerance of the eigenvalues. Can be used in: \item \inltt{imagShift}: provide an imaginary shift to the direct sovler eigenvalue problem by the specified value.
ltt{ModifiedArnoldi} and \inltt{Arpack}. Default value: $10^{-6}$.
\item \inltt{nits}: sets the maximum number of iterations. Can be used in: \inltt{ModifiedArnoldi} and \inltt{Arpack}. Default value: 500.
\item \inltt{LZ}: sets the length in the spanswise direction $L_z$. Can be used in \inltt{Homogeneous} set to \inltt{1D}. Default value: 1.
\item \inltt{HomModesZ}: sets the number of planes in the homogeneous directions. Can be used in \inltt{Homogeneous} set to \inltt{1D} and \inltt{ModeType} set to \inltt{MultipleModes}.
\item \inltt{N\_slices}: sets the number of temporal slices for Floquet stability analysis.
\item \inltt{period}: sets the periodicity of the base flow.
\item \inltt{realShift}: provide a real shift to the direct sovler eigenvalue problem by the specified value.
\end{itemize}
\subsection{Functions}
\begin{itemize}
\item To be INserted
\end{itemize}
When using the direct solver for stability analysis it is necessary to specify a Forcing function ``StabilityCoupledLNS'' in the form:
\begin{lstlisting}[style=XMLStyle]
<FORCING>
<FORCE TYPE="StabilityCoupledLNS">
</FORCE>
</FORCING>
\end{lstlisting}
This is required since we need to tell the solver to use the existing
field as a forcing function to the direct matrix inverse as part of
the Arnoldi iteration.
\begin{notebox}
Examples of the set up of the direct solver stability analysis (and
other incompressible Navier-Stokes solvers) can be found in the
regression test directory
\inlsh{NEKTAR/solvers/IncNavierStokesSolver/Tests}. See for example the files
\inlsh{PPF\_R15000\_ModifiedArnoldi\_Shift.tst} and \inlsh{PPF\_R15000\_3D.xml}
noting that some parameters are specified in the .tst files.
\end{notebox}
\section{Steady-state solver Session file configuration}
\label{SectionSFD_XML}
In this section, we detail how to use the steady-state solver (that implements the selective frequency damping method, see Sec. \ref{SectionSFD}).
Two cases are detailed here: the execution of the classical SFD method and the \textit{adaptive} SFD method, where the control coefficient $\chi$ and the filter width $\Delta$ of the SFD method are updated all along the solver execution. For the second case, the parameters of the SFD method do not need to be defined by the user (they will be automatically calculated all along the solver execution) but several session files must be defined in a very specific way.
In this section, we detail how to use the steady-state solver (that
implements the selective frequency damping method, see
Sec. \ref{SectionSFD}). Two cases are detailed here: the execution of
the classical SFD method and the \textit{adaptive} SFD method, where
the control coefficient $\chi$ and the filter width $\Delta$ of the
SFD method are updated all along the solver execution. For the second
case, the parameters of the SFD method do not need to be defined by
the user (they will be automatically calculated all along the solver
execution) but several session files must be defined in a very
specific way.
\subsection{Execution of the classical steady-state solver}
\subsubsection{Solver Info}
The definition of \inltt{Eqtype}, \inltt{TimeIntegrationMethod} and \inltt{Projection} is similar as what is explained in \ref{SectionIncNS_SolverInfo_Stab}. The use of the steady-state solver is enforced through the definition of the \inltt{Driver}: is has to be \inltt{SteadyState}. \inltt{EvolutionOperator} does not need to be defined to run the unadapted SFD method (by default, it is \inltt{Nonlinear}).
The definition of \inltt{Eqtype}, \inltt{TimeIntegrationMethod} and
\inltt{Projection} is similar as what is explained in
\ref{SectionIncNS_SolverInfo_Stab}. The use of the steady-state solver
is enforced through the definition of the \inltt{Driver} which has to be
\inltt{SteadyState}. \inltt{EvolutionOperator} does not need to be
defined to run the unadapted SFD method (by default, it is set to
\inltt{Nonlinear}).
\subsubsection{Parameters}
......
......@@ -291,6 +291,20 @@ the Incompressible Navier-Stokes solver,
</REGION>
\end{lstlisting}
Boundary conditions can also be loaded simultaneously from a file and from an
analytic expression. For example in the scenario where a spatial boundary
condition is read from a file, but needs to be modulated by a time-dependent
analytic expression:
\begin{lstlisting}[style=XMLStyle]
<REGION REF="1">
<D VAR="u" VALUE="USERDEFINEDTYPE="TimeDependent" VALUE="sin(PI*(x-t))"
FILE="bcsfromfiles_u_1.bc" />
</REGION>
\end{lstlisting}
In the case where both \inltt{VALUE} and \inltt{FILE} are specified, the values
are multiplied together to give the final value for the boundary condition.
\subsection{Functions}
Finally, multi-variable functions such as initial conditions and analytic
......@@ -507,4 +521,4 @@ will be the y-axis and the z-axis.
%%% Local Variables:
%%% mode: latex
%%% TeX-master: "../user-guide"
%%% End:
%%% End:
\ No newline at end of file
......@@ -92,6 +92,7 @@ ADD_NEKTAR_TEST(Helmholtz2D_CG_P7_PreconBlock)
ADD_NEKTAR_TEST(Helmholtz2D_CG_P7_PreconDiagonal)
ADD_NEKTAR_TEST(Helmholtz2D_HDG_P7_Modes)
ADD_NEKTAR_TEST(Helmholtz2D_HDG_P7_Modes_AllBCs)
ADD_NEKTAR_TEST(Helmholtz2D_CG_varP_Modes)
ADD_NEKTAR_TEST_LENGTHY(Helmholtz3D_CG_Hex)
......@@ -132,6 +133,7 @@ IF (NEKTAR_USE_MPI)
ADD_NEKTAR_TEST(Helmholtz3D_CG_Prism_iter_ml_par3)
ADD_NEKTAR_TEST_LENGTHY(Helmholtz3D_CG_Hex_AllBCs_xxt_sc_par3)
ADD_NEKTAR_TEST(Helmholtz2D_CG_P14_xxt_per)
ADD_NEKTAR_TEST(Helmholtz2D_CG_varP_Modes_par)
# TODO: This test fails due to a bug with Dirichlet bnd conditions.
# To be resolved in a separate branch.
......
<?xml version="1.0" encoding="utf-8"?>
<test>
<description>Helmholtz 2D CG with P=7</description>
<executable>Helmholtz2D</executable>
<parameters>Helmholtz2D_varP.xml</parameters>
<files>
<file description="Session File">Helmholtz2D_varP.xml</file>
</files>
<metrics>
<metric type="L2" id="1">
<value tolerance="1e-6">0.0045269</value>
</metric>
<metric type="Linf" id="2">
<value tolerance="1e-6">0.00420409</value>
</metric>
</metrics>
</test>
<?xml version="1.0" encoding="utf-8"?>
<test>
<description>Helmholtz 2D CG with P=7, parallel</description>
<executable>Helmholtz2D</executable>
<parameters>Helmholtz2D_varP.xml</parameters>
<processes>2</processes>
<files>
<file description="Session File">Helmholtz2D_varP.xml</file>
</files>
<metrics>
<metric type="L2" id="1">
<value tolerance="1e-6">0.0045269</value>
</metric>
<metric type="Linf" id="2">
<value tolerance="1e-6">0.00420409</value>
</metric>
</metrics>
</test>
This diff is collapsed.
This diff is collapsed.
......@@ -348,7 +348,7 @@ void PtsField::SetPointsPerEdge(const vector< int > nPtsPerEdge)
totPts = totPts * nPtsPerEdge.at(i);
}
ASSERTL0(totPts == m_pts.num_elements(),
ASSERTL0(totPts == m_pts[0].num_elements(),
"nPtsPerEdge does not match total number of points");
m_nPtsPerEdge = nPtsPerEdge;
......
......@@ -55,23 +55,25 @@ namespace Nektar
typedef vector<map<int, int> > DofGraph;
MULTI_REGIONS_EXPORT
pair<int, StdRegions::Orientation> DeterminePeriodicEdgeOrientId(
int meshEdgeId,
StdRegions::Orientation edgeOrient,
const vector<PeriodicEntity> &periodicEdges);
MULTI_REGIONS_EXPORT
StdRegions::Orientation DeterminePeriodicFaceOrient(
StdRegions::Orientation faceOrient1,
StdRegions::Orientation faceOrient2);
/// Constructs mappings for the C0 scalar continuous Galerkin formulation.
class AssemblyMapCG: public AssemblyMap
{
typedef Array<OneD, const ExpListSharedPtr> BndCondExp;
typedef Array<OneD, const SpatialDomains::BoundaryConditionShPtr>
BndCond;
public:
/// Default constructor.
MULTI_REGIONS_EXPORT AssemblyMapCG(
......
......@@ -2435,7 +2435,7 @@
int npoints;
int nbnd = m_bndCondExpansions.num_elements();
MultiRegions::ExpListSharedPtr locExpList;
for (i = 0; i < nbnd; ++i)
{
if (time == 0.0 || m_bndConditions[i]->IsTimeDependent())
......@@ -2446,7 +2446,8 @@
Array<OneD, NekDouble> x0(npoints, 0.0);
Array<OneD, NekDouble> x1(npoints, 0.0);
Array<OneD, NekDouble> x2(npoints, 0.0);
Array<OneD, NekDouble> valuesFile(npoints, 1.0), valuesExp(npoints, 1.0);
locExpList->GetCoords(x0, x1, x2);
if (m_bndConditions[i]->GetBoundaryConditionType()
......@@ -2456,23 +2457,30 @@
SpatialDomains::DirichletBoundaryCondition>(
m_bndConditions[i])->m_filename;
string exprbcs = boost::static_pointer_cast<
SpatialDomains::DirichletBoundaryCondition>(
m_bndConditions[i])->m_expr;
if (filebcs != "")
{
ExtractFileBCs(filebcs, varName, locExpList);
valuesFile = locExpList->GetPhys();
}
else
if (exprbcs != "")
{
LibUtilities::Equation condition = boost::static_pointer_cast<SpatialDomains::
DirichletBoundaryCondition >(
m_bndConditions[i])->m_dirichletCondition;
condition.Evaluate(x0, x1, x2, time,
locExpList->UpdatePhys());
locExpList->FwdTrans_BndConstrained(
locExpList->GetPhys(),
locExpList->UpdateCoeffs());
condition.Evaluate(x0, x1, x2, time, valuesExp);
}
Vmath::Vmul(npoints, valuesExp, 1, valuesFile, 1, locExpList->UpdatePhys(), 1);
locExpList->FwdTrans_BndConstrained(
locExpList->GetPhys(),
locExpList->UpdateCoeffs());
}
else if (m_bndConditions[i]->GetBoundaryConditionType()
== SpatialDomains::eNeumann)
......
......@@ -799,7 +799,7 @@ namespace Nektar
(*m_exp)[eid]->ExtractDataToCoeffs(&fielddata[offset], fielddef->m_numModes,modes_offset,&coeffs[m_coeff_offset[eid] + planes_offset*ncoeffs_per_plane]);
}
}
modes_offset += (*m_exp)[0]->GetNumBases();
modes_offset += (*m_exp)[0]->GetNumBases() + fielddef->m_numHomogeneousDir;
}
}
}
......
This diff is collapsed.
......@@ -40,73 +40,78 @@
namespace Nektar
{
namespace SolverUtils
namespace SolverUtils
{
/// Base class for the development of solvers.
class DriverArnoldi: public Driver
{
public:
friend class MemoryManager<DriverArnoldi>;
SOLVER_UTILS_EXPORT void ArnoldiSummary(std::ostream &out);
protected:
int m_kdim; /// Dimension of Krylov subspace
int m_nvec; /// Number of vectors to test
int m_nits; /// Maxmum number of iterations
NekDouble m_evtol; /// Tolerance of iteratiosn
NekDouble m_period;/// Period of time stepping algorithm
bool m_timeSteppingAlgorithm; /// underlying operator is time stepping
int m_infosteps; /// interval to dump information if required.
int m_nfields;
NekDouble m_realShift;
NekDouble m_imagShift;
int m_negatedOp; /// Operator in solve call is negated
Array<OneD, NekDouble> m_real_evl;
Array<OneD, NekDouble> m_imag_evl;
/// Constructor
DriverArnoldi(const LibUtilities::SessionReaderSharedPtr pSession);
/// Destructor
virtual ~DriverArnoldi();
/// Copy Arnoldi storage to fields.
void CopyArnoldiArrayToField(Array<OneD, NekDouble> &array);
/// Copy fields to Arnoldi storage.
void CopyFieldToArnoldiArray(Array<OneD, NekDouble> &array);
/// Copy the forward field to the adjoint system in transient growth
/// calculations
void CopyFwdToAdj();
/// Write coefficients to file
void WriteFld(std::string file,
std::vector<Array<OneD, NekDouble> > coeffs);
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 void v_InitObject(ostream &out = cout);
virtual Array<OneD, NekDouble> v_GetRealEvl(void)
{
/// Base class for the development of solvers.
class DriverArnoldi: public Driver
{
public:
friend class MemoryManager<DriverArnoldi>;
SOLVER_UTILS_EXPORT void ArnoldiSummary(std::ostream &out);
protected:
int m_kdim; /// Dimension of Krylov subspace
int m_nvec; /// Number of vectors to test
int m_nits; /// Maxmum number of iterations
NekDouble m_evtol;/// Tolerance of iteratiosn
NekDouble m_period;/// Period of time stepping algorithm
bool m_timeSteppingAlgorithm; /// underlying operator is time stepping
int m_infosteps; /// interval to dump information if required.
int m_nfields;
NekDouble m_realShift;
NekDouble m_imagShift;
Array<OneD, NekDouble> m_real_evl;
Array<OneD, NekDouble> m_imag_evl;
/// Constructor
DriverArnoldi(const LibUtilities::SessionReaderSharedPtr pSession);
/// Destructor
virtual ~DriverArnoldi();
/// Copy Arnoldi storage to fields.
void CopyArnoldiArrayToField(Array<OneD, NekDouble> &array);
/// Copy fields to Arnoldi storage.
void CopyFieldToArnoldiArray(Array<OneD, NekDouble> &array);
///Copy the forward field to the adjoint system in transient growth calculations
void CopyFwdToAdj();
// write coefficients to file.
void WriteFld(std::string file, std::vector<Array<OneD, NekDouble> > coeffs);
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);
virtual void v_InitObject(ostream &out = cout);
virtual Array<OneD, NekDouble> v_GetRealEvl(void)
{
return m_real_evl;
}
virtual Array<OneD, NekDouble> v_GetImagEvl(void)
{
return m_imag_evl;
}
};
return m_real_evl;
}
} //end of namespace
#endif //NEKTAR_SOLVERUTILS_DRIVERARNOLDI_H
virtual Array<OneD, NekDouble> v_GetImagEvl(void)
{
return m_imag_evl;
}
};
}
} //end of namespace
#endif //NEKTAR_SOLVERUTILS_DRIVERARNOLDI_H
\ No newline at end of file
This diff is collapsed.
......@@ -42,53 +42,53 @@
namespace Nektar
{
namespace SolverUtils
{
/// Base class for the development of solvers.
class DriverArpack: public DriverArnoldi
{
public:
friend class MemoryManager<DriverArpack>;
/// Creates an instance of this class
static DriverSharedPtr create(const LibUtilities::SessionReaderSharedPtr& pSession) {
DriverSharedPtr p = MemoryManager<DriverArpack>::AllocateSharedPtr(pSession);
p->InitObject();
return p;
}
///Name of the class
static std::string className;
protected:
int m_maxn; //Maximum size of the problem
int m_maxnev; //maximum number of eigenvalues requested
int m_maxncv; //Largest number of basis vector used in Implicitly Restarted Arnoldi
// std::string m_arpackProblemType; //Arpack input for problem type
/// Constructor
DriverArpack( const LibUtilities::SessionReaderSharedPtr pSession);
/// Destructor
virtual ~DriverArpack();
/// Virtual function for initialisation implementation.
virtual void v_InitObject(ostream &out = cout);
/// Virtual function for solve implementation.
virtual void v_Execute(ostream &out = cout);
static std::string arpackProblemTypeLookupIds[];
static std::string arpackProblemTypeDefault;
static std::string driverLookupId;
private:
static std::string ArpackProblemTypeTrans[];
};
namespace SolverUtils
{
/// Base class for the development of solvers.
class DriverArpack: public DriverArnoldi
{
public:
friend class MemoryManager<DriverArpack>;
/// Creates an instance of this class
static DriverSharedPtr create(const LibUtilities::SessionReaderSharedPtr& pSession) {
DriverSharedPtr p = MemoryManager<DriverArpack>::AllocateSharedPtr(pSession);
p->InitObject();
return p;
}
///Name of the class
static std::string className;
protected:
int m_maxn;//Maximum size of the problem
int m_maxnev;//maximum number of eigenvalues requested
int m_maxncv;//Largest number of basis vector used in Implicitly Restarted Arnoldi
/// Constructor
DriverArpack( const LibUtilities::SessionReaderSharedPtr pSession);
/// Destructor
virtual ~DriverArpack();
/// Virtual function for initialisation implementation.
virtual void v_InitObject(ostream &out = cout);
/// Virtual function for solve implementation.
virtual void v_Execute(ostream &out = cout);
static std::string arpackProblemTypeLookupIds[];
static std::string arpackProblemTypeDefault;
static std::string driverLookupId;
private:
static std::string ArpackProblemTypeTrans[];
};
}
} //end of namespace
#endif //NEKTAR_SOLVERUTILS_DRIVERARPACK_H
......
......@@ -41,98 +41,98 @@
namespace Nektar
{
namespace SolverUtils
{
class DriverModifiedArnoldi: public DriverArnoldi
{
public:
friend class MemoryManager<DriverModifiedArnoldi>;
/// Creates an instance of this class
static DriverSharedPtr create(
const LibUtilities::SessionReaderSharedPtr& pSession)
{
DriverSharedPtr p = MemoryManager<DriverModifiedArnoldi>
::AllocateSharedPtr(pSession);
p->InitObject();
return p;
}
///Name of the class
static std::string className;
protected:
/// Constructor
DriverModifiedArnoldi(const LibUtilities::SessionReaderSharedPtr pSession);
/// Destructor
virtual ~DriverModifiedArnoldi();
/// Virtual function for initialisation implementation.
virtual void v_InitObject(ostream &out = cout );
/// Virtual function for solve implementation.
virtual void v_Execute(ostream &out = cout);
private:
/// Generates a new vector in the sequence by applying the linear operator.
void EV_update(Array<OneD, NekDouble> &src,
Array<OneD, NekDouble> &tgt);
/// Generates the upper Hessenberg matrix H and computes its eigenvalues.
void EV_small(Array<OneD, Array<OneD, NekDouble> > &Kseq,
const int ntot,
const Array<OneD, NekDouble> &alpha,
const int kdim,
Array<OneD, NekDouble> &zvec,
Array<OneD, NekDouble> &wr,
Array<OneD, NekDouble> &wi,
NekDouble &resnorm);
/// Tests for convergence of eigenvalues of H.
int EV_test(const int itrn,
const int kdim,
Array<OneD, NekDouble> &zvec,
Array<OneD, NekDouble> &wr,
Array<OneD, NekDouble> &wi,
const NekDouble resnorm,
const int nvec,
ofstream &evlout,
NekDouble &resid0);
/// Sorts a sequence of eigenvectors/eigenvalues by magnitude.
void EV_sort(Array<OneD, NekDouble> &evec,
Array<OneD, NekDouble> &wr,
Array<OneD, NekDouble> &wi,
Array<OneD, NekDouble> &test,
const int dim);
void EV_post(Array<OneD, Array<OneD, NekDouble> > &Tseq,
Array<OneD, Array<OneD, NekDouble> > &Kseq,
const int ntot,
const int kdim,
const int nvec,
Array<OneD, NekDouble> &zvec,
Array<OneD, NekDouble> &wr,
Array<OneD, NekDouble> &wi,
const int icon);
void EV_big(Array<OneD, Array<OneD, NekDouble> > &bvecs,
Array<OneD, Array<OneD, NekDouble> > &evecs,
const int ntot,
const int kdim,
const int nvec,
Array<OneD, NekDouble> &zvec,
Array<OneD, NekDouble> &wr,
Array<OneD, NekDouble> &wi);
static std::string driverLookupId;
};
}
namespace SolverUtils
{
class DriverModifiedArnoldi: public DriverArnoldi
{
public:
friend class MemoryManager<DriverModifiedArnoldi>;
/// Creates an instance of this class
static DriverSharedPtr create(
const LibUtilities::SessionReaderSharedPtr& pSession)
{
DriverSharedPtr p = MemoryManager<DriverModifiedArnoldi>
::AllocateSharedPtr(pSession);
p->InitObject();
return p;
}
///Name of the class
static std::string className;
protected:
/// Constructor
DriverModifiedArnoldi(const LibUtilities::SessionReaderSharedPtr pSession);
/// Destructor
virtual ~DriverModifiedArnoldi();
/// Virtual function for initialisation implementation.
virtual void v_InitObject(ostream &out = cout );
/// Virtual function for solve implementation.
virtual void v_Execute(ostream &out = cout);
private:
/// Generates a new vector in the sequence by applying the linear operator.
void EV_update( Array<OneD, NekDouble> &src,
Array<OneD, NekDouble> &tgt);
/// Generates the upper Hessenberg matrix H and computes its eigenvalues.
void EV_small( Array<OneD, Array<OneD, NekDouble> > &Kseq,
const int ntot,
const Array<OneD, NekDouble> &alpha,
const int kdim,
Array<OneD, NekDouble> &zvec,
Array<OneD, NekDouble> &wr,
Array<OneD, NekDouble> &wi,
NekDouble &resnorm);
/// Tests for convergence of eigenvalues of H.
int EV_test( const int itrn,