Commit 8a902627 authored by Chris Cantwell's avatar Chris Cantwell
Browse files

Merge branch 'feature/NewHOBCs' into 'master'

Feature/newhobcs

Added in a new convective like outflow boundary condition as well as a velocity correction scheme using a weak pressure Poisson equations. 

To complement this I have also tried to restructure the description in the user guide to use the Guermond and Shen approach to motivate both the standard velocity correction scheme and the weak pressure form.  

See merge request !713
parents f57dfa9d 3e801554
......@@ -48,6 +48,8 @@ v4.4.0
- Add ability to simulate additional scalar fields (!624)
- Improve performance when using homogeneous dealiasing (!622)
- Fix linearised advection for full 3D cases (!708)
- Added a weak pressure formulation following Guermond & Shen (!713)
- Added a convective like outflow boundary condition from Dong (!713)
**FieldConvert:**
- Allow equi-spaced output for 1D and 2DH1D fields (!613)
......
......@@ -345,6 +345,16 @@ year={2011}
year={2014}
}
@article{Dong15,
title={A convective-like energy-stable open boundary condition for simulation of incompressible flows},
author={S. Dong},
journal={Journal of Computational Physics},
volume={302},
pages={300-328},
year={2015}
}
@article{Ko07,
title = {Vectorized Matlab codes for linear two-dimensional elasticity},
author = {Koko, Jonas},
......@@ -441,3 +451,12 @@ year={2011}
publisher = {Springer London}
}
@Article{BaPlGrSh16,
author = {Bao, Y. and Palacios, R. and Graham, M. and Sherwin, S.J. },
title = {Generalized “thick” strip modelling for vortex-induced vibration of long flexible cylinders},
journal = {J. Comp. Phys},
year = {2016},
volume = {321},
pages = {1079-1097},
}
......@@ -75,7 +75,7 @@ namespace ErrorUtil
NekError(const std::string& message) : std::runtime_error(message) {}
};
inline static void Error(ErrType type, const char *routine, int lineNumber, const char *msg, unsigned int level)
inline static void Error(ErrType type, const char *routine, int lineNumber, const char *msg, unsigned int level, bool DoComm = false)
{
// The user of outStream is primarily for the unit tests.
// The unit tests often generate errors on purpose to make sure
......@@ -94,11 +94,14 @@ namespace ErrorUtil
// the correct rank. Messages are only printed on rank zero.
int rank = 0;
#if defined(NEKTAR_USE_MPI)
int flag;
MPI_Initialized(&flag);
if(flag)
int flag = 0;
if(DoComm)
{
MPI_Comm_rank(MPI_COMM_WORLD,&rank);
MPI_Initialized(&flag);
if(flag)
{
MPI_Comm_rank(MPI_COMM_WORLD,&rank);
}
}
#endif
......@@ -138,9 +141,12 @@ namespace ErrorUtil
}
}
#if defined(NEKTAR_USE_MPI)
if (flag)
if(DoComm)
{
MPI_Barrier(MPI_COMM_WORLD);
if (flag)
{
MPI_Barrier(MPI_COMM_WORLD);
}
}
#endif
throw NekError(baseMsg);
......@@ -185,6 +191,10 @@ namespace ErrorUtil
#define NEKERROR(type, msg) \
ErrorUtil::Error(type, __FILE__, __LINE__, msg, 0);
#define ROOTONLY_NEKERROR(type, msg) \
ErrorUtil::Error(type, __FILE__, __LINE__, msg, 0,true);
#define ASSERTL0(condition,msg) \
if(!(condition)) \
{ \
......@@ -211,7 +221,7 @@ namespace ErrorUtil
#define WARNINGL1(condition,msg) \
if(!(condition)) \
{ \
ErrorUtil::Error(ErrorUtil::ewarning, __FILE__, __LINE__, msg, 0); \
ErrorUtil::Error(ErrorUtil::ewarning, __FILE__, __LINE__, msg, 1); \
}
#else //defined(NEKTAR_DEBUG) || defined(NEKTAR_FULLDEBUG)
......@@ -233,7 +243,7 @@ namespace ErrorUtil
#define WARNINGL2(condition,msg) \
if(!(condition)) \
{ \
ErrorUtil::Error(ErrorUtil::ewarning, __FILE__, __LINE__, msg, 0); \
ErrorUtil::Error(ErrorUtil::ewarning, __FILE__, __LINE__, msg, 2); \
}
#else //NEKTAR_FULLDEBUG
......
......@@ -1103,6 +1103,8 @@ namespace Nektar
"Not set up for non boundary-interior expansions");
ASSERTL1(inoutmat->GetRows() == inoutmat->GetColumns(),
"Assuming that input matrix was square");
ASSERTL1(GetBasisType(0) == LibUtilities::eModified_A,
"Method set up only for modified modal bases curretly");
int i,j;
int id1,id2;
......@@ -1132,6 +1134,9 @@ namespace Nektar
// - if inoutmat.m_rows() == v_NCoeffs() it is a full
// matrix system
// - if inoutmat.m_rows() == v_GetNverts() it is a vertex space
// preconditioner.
// - if inoutmat.m_rows() == v_NumBndCoeffs() it is a
// boundary CG system
......@@ -1144,6 +1149,32 @@ namespace Nektar
{
GetFaceToElementMap(face,GetForient(face),map,sign);
}
else if (rows == GetNverts())
{
int nfvert = faceExp->GetNverts();
// Need to find where linear vertices are in facemat
Array<OneD, unsigned int> linmap;
Array<OneD, int> linsign;
// Use a linear expansion to get correct mapping
GetLinStdExp()->GetFaceToElementMap(face,GetForient(face),linmap, linsign);
// zero out sign map to remove all other modes
sign = Array<OneD, int> (order_f,0);
map = Array<OneD, unsigned int>(order_f,(unsigned int)0);
int fmap;
// Reset sign map to only have contribution from vertices
for(i = 0; i < nfvert; ++i)
{
fmap = faceExp->GetVertexMap(i,true);
sign[fmap] = 1;
// need to reset map
map[fmap] = linmap[i];
}
}
else if(rows == NumBndryCoeffs())
{
int nbndry = NumBndryCoeffs();
......
......@@ -567,6 +567,19 @@ namespace Nektar
}
StdRegions::StdExpansionSharedPtr HexExp::v_GetLinStdExp(void) const
{
LibUtilities::BasisKey bkey0(m_base[0]->GetBasisType(),
2, m_base[0]->GetPointsKey());
LibUtilities::BasisKey bkey1(m_base[1]->GetBasisType(),
2, m_base[1]->GetPointsKey());
LibUtilities::BasisKey bkey2(m_base[2]->GetBasisType(),
2, m_base[2]->GetPointsKey());
return MemoryManager<StdRegions::StdHexExp>
::AllocateSharedPtr( bkey0, bkey1, bkey2);
}
/**
* \brief Retrieves the physical coordinates of a given set of
* reference coordinates.
......
......@@ -160,6 +160,9 @@ namespace Nektar
LOCAL_REGIONS_EXPORT virtual
StdRegions::StdExpansionSharedPtr v_GetStdExp(void) const;
LOCAL_REGIONS_EXPORT virtual
StdRegions::StdExpansionSharedPtr v_GetLinStdExp(void) const;
LOCAL_REGIONS_EXPORT virtual void v_ExtractDataToCoeffs(
const NekDouble *data,
......
......@@ -669,6 +669,17 @@ namespace Nektar
m_nodalPointsKey.GetPointsType());
}
StdRegions::StdExpansionSharedPtr NodalTriExp::v_GetLinStdExp(void) const
{
LibUtilities::BasisKey bkey0(m_base[0]->GetBasisType(),
2, m_base[0]->GetPointsKey());
LibUtilities::BasisKey bkey1(m_base[1]->GetBasisType(),
2, m_base[1]->GetPointsKey());
return MemoryManager<StdRegions::StdNodalTriExp>
::AllocateSharedPtr( bkey0, bkey1, m_nodalPointsKey.GetPointsType());
}
DNekMatSharedPtr NodalTriExp::v_GenMatrix(const StdRegions::StdMatrixKey &mkey)
{
DNekMatSharedPtr returnval;
......
......@@ -178,6 +178,8 @@ namespace Nektar
const StdRegions::StdMatrixKey &mkey);
virtual StdRegions::StdExpansionSharedPtr v_GetStdExp(void) const;
virtual StdRegions::StdExpansionSharedPtr v_GetLinStdExp(void) const;
virtual DNekMatSharedPtr v_GenMatrix(const StdRegions::StdMatrixKey &mkey);
......
......@@ -456,6 +456,19 @@ namespace Nektar
m_base[2]->GetBasisKey());
}
StdRegions::StdExpansionSharedPtr PrismExp::v_GetLinStdExp(void) const
{
LibUtilities::BasisKey bkey0(m_base[0]->GetBasisType(),
2, m_base[0]->GetPointsKey());
LibUtilities::BasisKey bkey1(m_base[1]->GetBasisType(),
2, m_base[1]->GetPointsKey());
LibUtilities::BasisKey bkey2(m_base[2]->GetBasisType(),
2, m_base[2]->GetPointsKey());
return MemoryManager<StdRegions::StdPrismExp>
::AllocateSharedPtr( bkey0, bkey1, bkey2);
}
/**
* @brief Get the coordinates #coords at the local coordinates
* #Lcoords.
......
......@@ -135,6 +135,10 @@ namespace Nektar
LOCAL_REGIONS_EXPORT virtual
StdRegions::StdExpansionSharedPtr v_GetStdExp(void) const;
LOCAL_REGIONS_EXPORT virtual
StdRegions::StdExpansionSharedPtr v_GetLinStdExp(void) const;
LOCAL_REGIONS_EXPORT virtual int v_GetCoordim();
LOCAL_REGIONS_EXPORT virtual void v_ExtractDataToCoeffs(
const NekDouble *data,
......
......@@ -309,6 +309,19 @@ namespace Nektar
m_base[2]->GetBasisKey());
}
StdRegions::StdExpansionSharedPtr PyrExp::v_GetLinStdExp(void) const
{
LibUtilities::BasisKey bkey0(m_base[0]->GetBasisType(),
2, m_base[0]->GetPointsKey());
LibUtilities::BasisKey bkey1(m_base[1]->GetBasisType(),
2, m_base[1]->GetPointsKey());
LibUtilities::BasisKey bkey2(m_base[2]->GetBasisType(),
2, m_base[2]->GetPointsKey());
return MemoryManager<StdRegions::StdPyrExp>
::AllocateSharedPtr( bkey0, bkey1, bkey2);
}
/*
* @brief Get the coordinates #coords at the local coordinates
* #Lcoords
......
......@@ -104,6 +104,9 @@ namespace Nektar
LOCAL_REGIONS_EXPORT virtual
StdRegions::StdExpansionSharedPtr v_GetStdExp(void) const;
LOCAL_REGIONS_EXPORT virtual
StdRegions::StdExpansionSharedPtr v_GetLinStdExp(void) const;
LOCAL_REGIONS_EXPORT virtual void v_GetCoord(
const Array<OneD, const NekDouble> &Lcoords,
Array<OneD, NekDouble> &coords);
......
......@@ -616,6 +616,18 @@ namespace Nektar
m_base[1]->GetBasisKey());
}
StdRegions::StdExpansionSharedPtr QuadExp::v_GetLinStdExp(void) const
{
LibUtilities::BasisKey bkey0(m_base[0]->GetBasisType(),
2, m_base[0]->GetPointsKey());
LibUtilities::BasisKey bkey1(m_base[1]->GetBasisType(),
2, m_base[1]->GetPointsKey());
return MemoryManager<StdRegions::StdQuadExp>
::AllocateSharedPtr( bkey0, bkey1);
}
void QuadExp::v_GetCoords(
Array<OneD, NekDouble> &coords_0,
Array<OneD, NekDouble> &coords_1,
......
......@@ -146,6 +146,9 @@ namespace Nektar
LOCAL_REGIONS_EXPORT virtual
StdRegions::StdExpansionSharedPtr v_GetStdExp(void) const;
LOCAL_REGIONS_EXPORT virtual
StdRegions::StdExpansionSharedPtr v_GetLinStdExp(void) const;
LOCAL_REGIONS_EXPORT virtual void v_GetCoord(
const Array<OneD, const NekDouble> &Lcoords,
Array<OneD, NekDouble> &coords);
......
......@@ -792,6 +792,15 @@ cout<<"deps/dx ="<<inarray_d0[i]<<" deps/dy="<<inarray_d1[i]<<endl;
::AllocateSharedPtr(m_base[0]->GetBasisKey());
}
StdRegions::StdExpansionSharedPtr SegExp::v_GetLinStdExp(void) const
{
LibUtilities::BasisKey bkey0(m_base[0]->GetBasisType(),
2, m_base[0]->GetPointsKey());
return MemoryManager<StdRegions::StdSegExp>
::AllocateSharedPtr( bkey0);
}
int SegExp::v_GetCoordim()
{
return m_geom->GetCoordim();
......
......@@ -166,6 +166,9 @@ namespace Nektar
LOCAL_REGIONS_EXPORT virtual
StdRegions::StdExpansionSharedPtr v_GetStdExp(void) const;
LOCAL_REGIONS_EXPORT virtual
StdRegions::StdExpansionSharedPtr v_GetLinStdExp(void) const;
LOCAL_REGIONS_EXPORT virtual int v_GetCoordim();
LOCAL_REGIONS_EXPORT virtual void v_SetCoeffsToOrientation(
......
......@@ -558,6 +558,21 @@ namespace Nektar
m_base[2]->GetBasisKey());
}
StdRegions::StdExpansionSharedPtr TetExp::v_GetLinStdExp(void) const
{
LibUtilities::BasisKey bkey0(m_base[0]->GetBasisType(),
2, m_base[0]->GetPointsKey());
LibUtilities::BasisKey bkey1(m_base[1]->GetBasisType(),
2, m_base[1]->GetPointsKey());
LibUtilities::BasisKey bkey2(m_base[2]->GetBasisType(),
2, m_base[2]->GetPointsKey());
return MemoryManager<StdRegions::StdTetExp>
::AllocateSharedPtr( bkey0, bkey1, bkey2);
}
int TetExp::v_GetCoordim()
{
return m_geom->GetCoordim();
......
......@@ -131,6 +131,9 @@ namespace Nektar
LOCAL_REGIONS_EXPORT virtual
StdRegions::StdExpansionSharedPtr v_GetStdExp(void) const;
LOCAL_REGIONS_EXPORT virtual
StdRegions::StdExpansionSharedPtr v_GetLinStdExp(void) const;
LOCAL_REGIONS_EXPORT virtual int v_GetCoordim();
LOCAL_REGIONS_EXPORT virtual void v_ExtractDataToCoeffs(
......
......@@ -583,6 +583,17 @@ namespace Nektar
m_base[1]->GetBasisKey());
}
StdRegions::StdExpansionSharedPtr TriExp::v_GetLinStdExp(void) const
{
LibUtilities::BasisKey bkey0(m_base[0]->GetBasisType(),
2, m_base[0]->GetPointsKey());
LibUtilities::BasisKey bkey1(m_base[1]->GetBasisType(),
2, m_base[1]->GetPointsKey());
return MemoryManager<StdRegions::StdTriExp>
::AllocateSharedPtr( bkey0, bkey1);
}
void TriExp::v_GetCoord(const Array<OneD, const NekDouble> &Lcoords,
Array<OneD,NekDouble> &coords)
{
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment