Skip to content

GitLab

Commit 00e500bf authored by Chris Cantwell's avatar Chris Cantwell
Browse files

Merge branch 'feature/fieldconvert' of /opt/gitlab/repositories/nektar

parents 5f4c0946 54f0ecbd
Expand all Hide whitespace changes
Inline Side-by-side
library/LibUtilities/BasicUtils/FieldIO.cpp
View file @ 00e500bf
......@@ -76,8 +76,21 @@ namespace Nektar
{
#ifdef NEKTAR_USE_MPI
int size;
MPI_Comm_size( MPI_COMM_WORLD, &size );
ASSERTL0(size == 1, "This function is not available in parallel.");
int init;
MPI_Initialized(&init);
// If MPI has been initialised we can check the number of processes
// and, if > 1, tell the user he should not be running this
// function in parallel. If it is not initialised, we do not
// initialise it here, and assume the user knows what they are
// doing.
if (init)
{
MPI_Comm_size( MPI_COMM_WORLD, &size );
ASSERTL0(size == 1,
"This static function is not available in parallel. Please"
"instantiate a FieldIO object for parallel use.");
}
#endif
CommSharedPtr c = GetCommFactory().CreateInstance("Serial", 0, 0);
FieldIO f(c);
......@@ -99,8 +112,21 @@ namespace Nektar
{
#ifdef NEKTAR_USE_MPI
int size;
MPI_Comm_size( MPI_COMM_WORLD, &size );
ASSERTL0(size == 1, "This function is not available in parallel.");
int init;
MPI_Initialized(&init);
// If MPI has been initialised we can check the number of processes
// and, if > 1, tell the user he should not be running this
// function in parallel. If it is not initialised, we do not
// initialise it here, and assume the user knows what they are
// doing.
if (init)
{
MPI_Comm_size( MPI_COMM_WORLD, &size );
ASSERTL0(size == 1,
"This static function is not available in parallel. Please"
"instantiate a FieldIO object for parallel use.");
}
#endif
CommSharedPtr c = GetCommFactory().CreateInstance("Serial", 0, 0);
FieldIO f(c);
......@@ -396,7 +422,10 @@ namespace Nektar
ASSERTL0(loadOkay1, errstr.str());
ImportFieldDefs(doc1, fielddefs, false);
ImportFieldData(doc1, fielddefs, fielddata);
if(fielddata != NullVectorNekDoubleVector)
{
ImportFieldData(doc1, fielddefs, fielddata);
}
}
}
......@@ -438,7 +467,10 @@ namespace Nektar
ASSERTL0(loadOkay1, errstr.str());
ImportFieldDefs(doc1, fielddefs, false);
ImportFieldData(doc1, fielddefs, fielddata);
if(fielddata != NullVectorNekDoubleVector)
{
ImportFieldData(doc1, fielddefs, fielddata);
}
}
}
}
......@@ -457,7 +489,10 @@ namespace Nektar
ImportFieldMetaData(doc,fieldmetadatamap);
ImportFieldDefs(doc, fielddefs, false);
ImportFieldData(doc, fielddefs, fielddata);
if(fielddata != NullVectorNekDoubleVector)
{
ImportFieldData(doc, fielddefs, fielddata);
}
}
}
......
library/LibUtilities/BasicUtils/FieldIO.h
View file @ 00e500bf
......@@ -51,6 +51,7 @@
#include <boost/iostreams/copy.hpp>
#include <boost/iostreams/filter/zlib.hpp>
#include <boost/iostreams/filtering_stream.hpp>
#include <boost/assign/list_of.hpp>
namespace Nektar
{
......@@ -63,7 +64,7 @@ namespace Nektar
typedef std::map<std::string, std::string> FieldMetaDataMap;
static FieldMetaDataMap NullFieldMetaDataMap;
static std::vector<std::vector< NekDouble> > NullVectorNekDoubleVector = boost::assign::list_of(NullNekDoubleVector);
struct FieldDefinitions
{
......@@ -141,7 +142,7 @@ namespace Nektar
LIB_UTILITIES_EXPORT void Import(
const std::string& infilename,
std::vector<FieldDefinitionsSharedPtr> &fielddefs,
std::vector<std::vector<NekDouble> > &fielddata,
std::vector<std::vector<NekDouble> > &fielddata = NullVectorNekDoubleVector,
FieldMetaDataMap &fieldinfomap = NullFieldMetaDataMap,
const Array<OneD, int> ElementiDs = NullInt1DArray);
......@@ -165,7 +166,7 @@ namespace Nektar
LIB_UTILITIES_EXPORT void Import(
const std::string& infilename,
std::vector<FieldDefinitionsSharedPtr> &fielddefs,
std::vector<std::vector<NekDouble> > &fielddata,
std::vector<std::vector<NekDouble> > &fielddata = NullVectorNekDoubleVector,
FieldMetaDataMap &fieldinfomap = NullFieldMetaDataMap,
const Array<OneD, int> ElementiDs = NullInt1DArray);
......
library/LibUtilities/BasicUtils/ParseUtils.hpp
View file @ 00e500bf
......@@ -124,7 +124,7 @@ namespace Nektar
space_p).full;
}
static bool GenerateUnOrderedVector(const char *const str, std::vector<NekDouble> &vec)
static bool GenerateUnOrderedVector(const char *const str, std::vector<NekDouble> &vec)
{
// Functors used to parse the sequence.
fctor5 functor5(&vec);
......
library/LibUtilities/BasicUtils/SessionReader.cpp
View file @ 00e500bf
......@@ -179,6 +179,14 @@ namespace Nektar
// Create communicator
CreateComm(argc, argv);
// If running in parallel change the default global sys solution
// type.
if (m_comm->GetSize() > 1)
{
m_solverInfoDefaults["GLOBALSYSSOLN"] =
"IterativeStaticCond";
}
}
......@@ -215,6 +223,14 @@ namespace Nektar
{
m_comm = pComm;
}
// If running in parallel change the default global sys solution
// type.
if (m_comm->GetSize() > 1)
{
m_solverInfoDefaults["GLOBALSYSSOLN"] =
"IterativeStaticCond";
}
}
......@@ -1321,14 +1337,6 @@ namespace Nektar
}
m_comm = GetCommFactory().CreateInstance(vCommModule,argc,argv);
// If running in parallel change the default global sys solution
// type.
if (m_comm->GetSize() > 1)
{
m_solverInfoDefaults["GLOBALSYSSOLN"] =
"IterativeStaticCond";
}
}
}
......
library/LibUtilities/Communication/CommMpi.cpp
View file @ 00e500bf
......@@ -52,11 +52,16 @@ namespace Nektar
CommMpi::CommMpi(int narg, char* arg[])
: Comm(narg,arg)
{
int init = 0;
MPI_Initialized(&init);
ASSERTL0(!init, "MPI has already been initialised.");
int retval = MPI_Init(&narg, &arg);
if (retval != MPI_SUCCESS)
{
ASSERTL0(false, "Failed to initialise MPI");
}
m_comm = MPI_COMM_WORLD;
MPI_Comm_size( m_comm, &m_size );
MPI_Comm_rank( m_comm, &m_rank );
......
library/MultiRegions/ContField2D.cpp
View file @ 00e500bf
......@@ -682,6 +682,28 @@ namespace Nektar
Vmath::Vcopy(nDir, tmp, 1, outarray, 1);
}
void ContField2D::v_FillBndCondFromField(void)
{
NekDouble sign;
int bndcnt = 0;
const Array<OneD,const int> &bndMap =
m_locToGloMap->GetBndCondCoeffsToGlobalCoeffsMap();
Array<OneD, NekDouble> tmp(m_locToGloMap->GetNumGlobalCoeffs());
LocalToGlobal(m_coeffs,tmp);
// Now fill in all other Dirichlet coefficients.
for(int i = 0; i < m_bndCondExpansions.num_elements(); ++i)
{
Array<OneD, NekDouble>& coeffs = m_bndCondExpansions[i]->UpdateCoeffs();
for(int j = 0; j < (m_bndCondExpansions[i])->GetNcoeffs(); ++j)
{
sign = m_locToGloMap->GetBndCondCoeffsToGlobalCoeffsSign(bndcnt);
coeffs[j] = sign * tmp[bndMap[bndcnt++]];
}
}
}
/**
* This operation is evaluated as:
......
library/MultiRegions/ContField2D.h
View file @ 00e500bf
......@@ -209,6 +209,7 @@ namespace Nektar
/// Impose the Dirichlet Boundary Conditions on outarray
MULTI_REGIONS_EXPORT virtual void v_ImposeDirichletConditions(Array<OneD,NekDouble>& outarray);
MULTI_REGIONS_EXPORT virtual void v_FillBndCondFromField();
/// Gathers the global coefficients \f$\boldsymbol{\hat{u}}_g\f$
/// from the local coefficients \f$\boldsymbol{\hat{u}}_l\f$.
......
library/MultiRegions/ContField3D.cpp
View file @ 00e500bf
......@@ -317,7 +317,7 @@ namespace Nektar
Array<OneD, NekDouble> &outarray)
{
int bndcnt=0;
const Array<OneD,const int>& map = m_locToGloMap->GetBndCondCoeffsToGlobalCoeffsMap();
const Array<OneD,const int>& map = m_locToGloMap->GetBndCondCoeffsToGlobalCoeffsMap();
NekDouble sign;
for(int i = 0; i < m_bndCondExpansions.num_elements(); ++i)
......@@ -465,6 +465,29 @@ namespace Nektar
Vmath::Vcopy(nDir, tmp, 1, outarray, 1);
}
void ContField3D::v_FillBndCondFromField(void)
{
NekDouble sign;
int bndcnt = 0;
const Array<OneD,const int> &bndMap =
m_locToGloMap->GetBndCondCoeffsToGlobalCoeffsMap();
Array<OneD, NekDouble> tmp(m_locToGloMap->GetNumGlobalCoeffs());
LocalToGlobal(m_coeffs,tmp);
// Now fill in all other Dirichlet coefficients.
for(int i = 0; i < m_bndCondExpansions.num_elements(); ++i)
{
Array<OneD, NekDouble>& coeffs = m_bndCondExpansions[i]->UpdateCoeffs();
for(int j = 0; j < (m_bndCondExpansions[i])->GetNcoeffs(); ++j)
{
sign = m_locToGloMap->GetBndCondCoeffsToGlobalCoeffsSign(bndcnt);
coeffs[j] = sign * tmp[bndMap[bndcnt++]];
}
}
}
void ContField3D::v_LocalToGlobal(void)
{
......
library/MultiRegions/ContField3D.h
View file @ 00e500bf
......@@ -158,6 +158,8 @@ namespace Nektar
/// Impose the Dirichlet Boundary Conditions on outarray
virtual void v_ImposeDirichletConditions(Array<OneD,NekDouble>& outarray);
virtual void v_FillBndCondFromField();
virtual void v_LocalToGlobal(void);
......
library/MultiRegions/DisContField2D.cpp
View file @ 00e500bf
......@@ -100,7 +100,7 @@ namespace Nektar
GenerateBoundaryConditionExpansion(graph2D,bcs,variable,
DeclareCoeffPhysArrays);
if(DeclareCoeffPhysArrays)
if (DeclareCoeffPhysArrays)
{
EvaluateBoundaryConditions();
}
......@@ -1429,6 +1429,11 @@ namespace Nektar
m_traceMap->UniversalTraceAssemble(Fwd);
m_traceMap->UniversalTraceAssemble(Bwd);
}
void DisContField2D::v_FillBndCondFromField(void)
{
}
void DisContField2D::v_ExtractTracePhys(
Array<OneD, NekDouble> &outarray)
......@@ -1456,7 +1461,7 @@ namespace Nektar
{
// Loop over elemente and collect forward expansion
int nexp = GetExpSize();
int n,e,offset,phys_offset;
int n, e, offset, phys_offset;
Array<OneD,NekDouble> e_tmp;
Array<OneD, Array<OneD, StdRegions::StdExpansionSharedPtr> >
&elmtToTrace = m_traceMap->GetElmtToTrace();
......@@ -1485,7 +1490,7 @@ namespace Nektar
const Array<OneD, const NekDouble> &Fy,
Array<OneD, NekDouble> &outarray)
{
int e,n,offset, t_offset;
int e, n, offset, t_offset;
Array<OneD, NekDouble> e_outarray;
Array<OneD, Array<OneD, StdRegions::StdExpansionSharedPtr> >
&elmtToTrace = m_traceMap->GetElmtToTrace();
......@@ -1528,7 +1533,7 @@ namespace Nektar
const Array<OneD, const NekDouble> &Fn,
Array<OneD, NekDouble> &outarray)
{
int e,n,offset, t_offset;
int e, n, offset, t_offset;
Array<OneD, NekDouble> e_outarray;
Array<OneD, Array<OneD, StdRegions::StdExpansionSharedPtr> >
&elmtToTrace = m_traceMap->GetElmtToTrace();
......@@ -1546,7 +1551,7 @@ namespace Nektar
}
}
}
/**
* @brief Add trace contributions into elemental coefficient spaces.
......@@ -1580,16 +1585,16 @@ namespace Nektar
Array<OneD, Array<OneD, StdRegions::StdExpansionSharedPtr> >
&elmtToTrace = m_traceMap->GetElmtToTrace();
for(n = 0; n < GetExpSize(); ++n)
for (n = 0; n < GetExpSize(); ++n)
{
offset = GetCoeff_Offset(n);
for(e = 0; e < (*m_exp)[n]->GetNedges(); ++e)
for (e = 0; e < (*m_exp)[n]->GetNedges(); ++e)
{
t_offset = GetTrace()->GetPhys_Offset(
elmtToTrace[n][e]->GetElmtId());
// Evaluate upwind flux less local edge
if(IsLeftAdjacentEdge(n,e))
if (IsLeftAdjacentEdge(n, e))
{
(*m_exp)[n]->AddEdgeNormBoundaryInt(
e, elmtToTrace[n][e], Fwd+t_offset,
......@@ -1648,17 +1653,19 @@ namespace Nektar
}
LocalRegions::Expansion1DSharedPtr exp1d;
for(cnt = n = 0; n < m_bndCondExpansions.num_elements(); ++n)
for (cnt = n = 0; n < m_bndCondExpansions.num_elements(); ++n)
{
for(i = 0; i < m_bndCondExpansions[n]->GetExpSize(); ++i, ++cnt)
for (i = 0; i < m_bndCondExpansions[n]->GetExpSize();
++i, ++cnt)
{
exp1d = LocalRegions::Expansion1D::FromStdExp(m_bndCondExpansions[n]->GetExp(i));
exp1d = LocalRegions::Expansion1D::FromStdExp(
m_bndCondExpansions[n]->GetExp(i));
// Use edge to element map from MeshGraph2D.
SpatialDomains::ElementEdgeVectorSharedPtr tmp =
graph2D->GetElementsFromEdge(exp1d->GetGeom1D());
ElmtID[cnt] = globalIdMap[(*tmp)[0]->
m_Element->GetGlobalID()];
m_Element->GetGlobalID()];
EdgeID[cnt] = (*tmp)[0]->m_EdgeIndx;
}
}
......@@ -2119,60 +2126,62 @@ namespace Nektar
MultiRegions::ExpListSharedPtr locExpList;
for(i = 0; i < nbnd; ++i)
for (i = 0; i < nbnd; ++i)
{
if(time == 0.0 || m_bndConditions[i]->GetUserDefined() ==
SpatialDomains::eTimeDependent)
if (time == 0.0 ||
m_bndConditions[i]->GetUserDefined() ==
SpatialDomains::eTimeDependent)
{
locExpList = m_bndCondExpansions[i];
npoints = locExpList->GetNpoints();
Array<OneD,NekDouble> x0(npoints,0.0);
Array<OneD,NekDouble> x1(npoints,0.0);
Array<OneD,NekDouble> x2(npoints,0.0);
Array<OneD, NekDouble> x0(npoints, 0.0);
Array<OneD, NekDouble> x1(npoints, 0.0);
Array<OneD, NekDouble> x2(npoints, 0.0);
// Homogeneous input case for x2.
if(x2_in == NekConstants::kNekUnsetDouble)
if (x2_in == NekConstants::kNekUnsetDouble)
{
locExpList->GetCoords(x0,x1,x2);
}
else
{
locExpList->GetCoords(x0,x1,x2);
Vmath::Fill(npoints,x2_in,x2,1);
locExpList->GetCoords(x0, x1, x2);
Vmath::Fill(npoints, x2_in, x2, 1);
}
if(m_bndConditions[i]->GetBoundaryConditionType()
== SpatialDomains::eDirichlet)
if (m_bndConditions[i]->GetBoundaryConditionType()
== SpatialDomains::eDirichlet)
{
string filebcs = boost::static_pointer_cast<
SpatialDomains::DirichletBoundaryCondition>(
m_bndConditions[i])->m_filename;
if(filebcs != "")
if (filebcs != "")
{
string var = filebcs.substr(
string varString = filebcs.substr(
0, filebcs.find_last_of("."));
int len = var.length();
var = var.substr(len-1,len);
int len = varString.length();
varString = varString.substr(len-1, len);
cout << "Boundary condition from file:"
<< filebcs << endl;
std::vector<LibUtilities::FieldDefinitionsSharedPtr>
FieldDef;
std::vector<LibUtilities::
FieldDefinitionsSharedPtr> FieldDef;
std::vector<std::vector<NekDouble> > FieldData;
LibUtilities::FieldIO f(m_session->GetComm());
f.Import(filebcs,FieldDef, FieldData);
f.Import(filebcs, FieldDef, FieldData);
// copy FieldData into locExpList
locExpList->ExtractDataToCoeffs(
FieldDef[0], FieldData[0],
FieldDef[0]->m_fields[0],
locExpList->UpdateCoeffs());
locExpList->BwdTrans_IterPerExp(
locExpList->GetCoeffs(),
locExpList->UpdatePhys());
locExpList->FwdTrans_BndConstrained(
locExpList->GetPhys(),
locExpList->UpdateCoeffs());
......@@ -2181,10 +2190,11 @@ namespace Nektar
{
LibUtilities::Equation condition =
boost::static_pointer_cast<
SpatialDomains::DirichletBoundaryCondition
>(m_bndConditions[i])->m_dirichletCondition;
SpatialDomains::DirichletBoundaryCondition>
(m_bndConditions[i])->
m_dirichletCondition;
condition.Evaluate(x0,x1,x2,time,
condition.Evaluate(x0, x1, x2, time,
locExpList->UpdatePhys());
locExpList->FwdTrans_BndConstrained(
......@@ -2192,14 +2202,14 @@ namespace Nektar
locExpList->UpdateCoeffs());
}
}
else if(m_bndConditions[i]->GetBoundaryConditionType()
== SpatialDomains::eNeumann)
else if (m_bndConditions[i]->GetBoundaryConditionType()
== SpatialDomains::eNeumann)
{
string filebcs = boost::static_pointer_cast<
SpatialDomains::NeumannBoundaryCondition>(
m_bndConditions[i])->m_filename;
if(filebcs != "")
if (filebcs != "")
{
string var = filebcs.substr(
0, filebcs.find_last_of("."));
......@@ -2209,17 +2219,18 @@ namespace Nektar
cout << "Boundary condition from file: "
<< filebcs << endl;
std::vector<LibUtilities::FieldDefinitionsSharedPtr>
FieldDef;
std::vector<LibUtilities::
FieldDefinitionsSharedPtr> FieldDef;
std::vector<std::vector<NekDouble> > FieldData;
LibUtilities::FieldIO f(m_session->GetComm());
f.Import(filebcs,FieldDef, FieldData);
f.Import(filebcs, FieldDef, FieldData);
// copy FieldData into locExpList
locExpList->ExtractDataToCoeffs(
FieldDef[0], FieldData[0],
FieldDef[0]->m_fields[0],
locExpList->UpdateCoeffs());
locExpList->BwdTrans_IterPerExp(
locExpList->GetCoeffs(),
locExpList->UpdatePhys());
......@@ -2235,29 +2246,33 @@ namespace Nektar
}
*/
locExpList->IProductWRTBase(locExpList->GetPhys(),
locExpList->UpdateCoeffs());
locExpList->IProductWRTBase(
locExpList->GetPhys(),
locExpList->UpdateCoeffs());
}
else
{
LibUtilities::Equation condition =
boost::static_pointer_cast<
SpatialDomains::NeumannBoundaryCondition
>(m_bndConditions[i])->m_neumannCondition;
condition.Evaluate(x0,x1,x2,time,
SpatialDomains::NeumannBoundaryCondition>
(m_bndConditions[i])->
m_neumannCondition;
condition.Evaluate(x0, x1, x2, time,
locExpList->UpdatePhys());
locExpList->IProductWRTBase(locExpList->GetPhys(),
locExpList->UpdateCoeffs());
locExpList->IProductWRTBase(
locExpList->GetPhys(),
locExpList->UpdateCoeffs());
}
}
else if(m_bndConditions[i]->GetBoundaryConditionType()
== SpatialDomains::eRobin)
else if (m_bndConditions[i]->GetBoundaryConditionType()
== SpatialDomains::eRobin)
{
string filebcs = boost::static_pointer_cast<
SpatialDomains::RobinBoundaryCondition
>(m_bndConditions[i])->m_filename;
if(filebcs != "")
string filebcs = boost::static_pointer_cast<
SpatialDomains::RobinBoundaryCondition>
(m_bndConditions[i])->m_filename;
if (filebcs != "")
{
//Never tested!!!
string var = filebcs.substr(
......@@ -2265,11 +2280,11 @@ namespace Nektar
int len = var.length();
var = var.substr(len-1,len);
std::vector<LibUtilities::FieldDefinitionsSharedPtr>
FieldDef;
std::vector<std::vector<NekDouble> > FieldData;
std::vector<LibUtilities::
FieldDefinitionsSharedPtr> FieldDef;
std::vector<std::vector<NekDouble> > FieldData;
LibUtilities::FieldIO f(m_session->GetComm());
f.Import(filebcs,FieldDef, FieldData);
f.Import(filebcs, FieldDef, FieldData);
// copy FieldData into locExpList
locExpList->ExtractDataToCoeffs(
......@@ -2298,28 +2313,29 @@ namespace Nektar
{
LibUtilities::Equation condition =
boost::static_pointer_cast<
SpatialDomains::RobinBoundaryCondition
>(m_bndConditions[i])->m_robinFunction;
LibUtilities::Equation coeff =
SpatialDomains::RobinBoundaryCondition>
(m_bndConditions[i])->
m_robinFunction;
LibUtilities::Equation coeff =
boost::static_pointer_cast<
SpatialDomains::RobinBoundaryCondition
>(m_bndConditions[i])->m_robinPrimitiveCoeff;
condition.Evaluate(x0,x1,x2,time,
SpatialDomains::RobinBoundaryCondition>
(m_bndConditions[i])->
m_robinPrimitiveCoeff;
condition.Evaluate(x0, x1, x2, time,
locExpList->UpdatePhys());
locExpList->IProductWRTBase(
locExpList->GetPhys(),
locExpList->UpdateCoeffs());