Commit e2e0785e authored by Dave Moxey's avatar Dave Moxey

Fixing reading of initial conditions and L2 error for homogeneous fields.

More verbose output for initial conditions in EquationSystem.


git-svn-id: https://gforge.sci.utah.edu/svn/nektar/trunk@3918 305cdda6-5ce1-45b3-a98d-dfc68c8b3305
parent e3659ba5
......@@ -1758,7 +1758,7 @@ namespace Nektar
(*m_exp)[i]->SetPhys(m_phys+m_phys_offset[i]);
err = std::max(err,(*m_exp)[i]->Linf());
}
m_comm->AllReduce(err, LibUtilities::ReduceMax);
m_comm->GetRowComm()->AllReduce(err, LibUtilities::ReduceMax);
return err;
}
......@@ -1826,35 +1826,36 @@ namespace Nektar
errl2 = (*m_exp)[i]->L2();
err += errl2*errl2;
}
m_comm->AllReduce(err, LibUtilities::ReduceSum);
m_comm->GetRowComm()->AllReduce(err, LibUtilities::ReduceSum);
return sqrt(err);
}
Array<OneD, NekDouble> ExpList::v_HomogeneousEnergy (void)
Array<OneD, NekDouble> ExpList::v_HomogeneousEnergy (void)
{
ASSERTL0(false,
ASSERTL0(false,
"This method is not defined or valid for this class type");
Array<OneD, NekDouble> NoEnergy(1,0.0);
Array<OneD, NekDouble> NoEnergy(1,0.0);
return NoEnergy;
}
Array<OneD, unsigned int> ExpList::v_GetZIDs(void)
{
ASSERTL0(false,
Array<OneD, unsigned int> ExpList::v_GetZIDs(void)
{
ASSERTL0(false,
"This method is not defined or valid for this class type");
Array<OneD, unsigned int> NoModes(1);
Array<OneD, unsigned int> NoModes(1);
return NoModes;
}
return NoModes;
}
Array<OneD, unsigned int> ExpList::v_GetYIDs(void)
{
ASSERTL0(false,
Array<OneD, unsigned int> ExpList::v_GetYIDs(void)
{
ASSERTL0(false,
"This method is not defined or valid for this class type");
Array<OneD, unsigned int> NoModes(1);
Array<OneD, unsigned int> NoModes(1);
return NoModes;
}
return NoModes;
}
/**
......@@ -2050,11 +2051,11 @@ namespace Nektar
v_ExtractDataToCoeffs(fielddef,fielddata,field,m_coeffs);
}
//3D-Base Flow (implementation in the homogeneous classes)
void ExpList::v_ExtractDataToCoeffs(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, std::string &field, bool BaseFlow3D)
{
ASSERTL0(false, "This method is not defined or valid for this class type");
}
//3D-Base Flow (implementation in the homogeneous classes)
void ExpList::v_ExtractDataToCoeffs(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, std::string &field, bool BaseFlow3D)
{
ASSERTL0(false, "This method is not defined or valid for this class type");
}
void ExpList::v_ExtractDataToCoeffs(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, std::string &field, Array<OneD, NekDouble> &coeffs)
{
......
......@@ -395,25 +395,25 @@ namespace Nektar
for(int n = 0; n < m_planes.num_elements(); ++n)
{
errL2 = m_planes[n]->L2(soln + cnt);
cnt += m_planes[n]->GetTotPoints();
err += errL2*errL2*local_w[n]*m_lhom*0.5;
cnt += m_planes[n]->GetTotPoints();
err += errL2*errL2*local_w[n]*m_lhom*0.5;
}
m_comm->GetColumnComm()->AllReduce(err, LibUtilities::ReduceSum);
return sqrt(err);
}
NekDouble ExpList3DHomogeneous1D::v_L2(void)
{
NekDouble errL2,err = 0;
Array<OneD, const NekDouble> w = m_homogeneousBasis->GetW();
Array<OneD, NekDouble> local_w(m_planes.num_elements());
for(int n = 0; n < m_planes.num_elements(); n++)
{
local_w[n] = w[m_transposition->GetPlaneID(n)];
}
Array<OneD, NekDouble> local_w(m_planes.num_elements());
for(int n = 0; n < m_planes.num_elements(); n++)
{
local_w[n] = w[m_transposition->GetPlaneID(n)];
}
for(int n = 0; n < m_planes.num_elements(); ++n)
{
......
......@@ -634,7 +634,7 @@ namespace Nektar
{
int i,n;
int ncoeffs_per_plane = m_planes[0]->GetNcoeffs();
// Determine mapping from element ids to location in
// expansion list
map<int, int> ElmtID_to_ExpID;
......@@ -655,21 +655,31 @@ namespace Nektar
}
}
void ExpListHomogeneous1D::v_AppendFieldData(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata)
void ExpListHomogeneous1D::v_AppendFieldData(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata)
{
v_AppendFieldData(fielddef,fielddata,m_coeffs);
}
//Extract the data in fielddata into the m_coeff list
void ExpListHomogeneous1D::v_ExtractDataToCoeffs(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, std::string &field)
void ExpListHomogeneous1D::v_ExtractDataToCoeffs(
SpatialDomains::FieldDefinitionsSharedPtr &fielddef,
std::vector<NekDouble> &fielddata,
std::string &field,
Array<OneD, NekDouble> &coeffs)
{
int i,n;
int offset = 0;
int nzmodes;
int nzmodes;
int datalen = fielddata.size()/fielddef->m_fields.size();
int ncoeffs_per_plane = m_planes[0]->GetNcoeffs();
// Build map of plane IDs lying on this process.
std::map<int,int> homoZids;
for (i = 0; i < m_planes.num_elements(); ++i)
{
homoZids[m_transposition->GetPlaneID(i)] = i;
}
for(i = 0; i < fielddef->m_basis.size(); ++i)
{
if(fielddef->m_basis[i] == m_homogeneousBasis->GetBasisType())
......@@ -678,7 +688,7 @@ namespace Nektar
break;
}
}
ASSERTL1(i != fielddef->m_basis.size()," Failed to determine number of modes");
ASSERTL1(i != fielddef->m_basis.size(),"Failed to determine number of modes");
// Find data location according to field definition
for(i = 0; i < fielddef->m_fields.size(); ++i)
......@@ -689,11 +699,10 @@ namespace Nektar
}
offset += datalen;
}
ASSERTL0(i != fielddef->m_fields.size(),
"Field " + field + " not found in data file");
ASSERTL0(i!= fielddef->m_fields.size(),"Field not found in data file");
// Determine mapping from element ids to location in
// expansion list
// Determine mapping from element ids to location in expansion list.
map<int, int> ElmtID_to_ExpID;
for(i = 0; i < m_planes[0]->GetExpSize(); ++i)
{
......@@ -701,36 +710,44 @@ namespace Nektar
}
int modes_offset = 0;
int planes_offset = 0;
Array<OneD, NekDouble> coeff_tmp;
int planes_offset = 0;
Array<OneD, NekDouble> coeff_tmp;
for(i = 0; i < fielddef->m_elementIDs.size(); ++i)
{
int eid = ElmtID_to_ExpID[fielddef->m_elementIDs[i]];
int datalen = (*m_exp)[eid]->CalcNumberOfCoefficients(fielddef->m_numModes,modes_offset);
int datalen = (*m_exp)[eid]->CalcNumberOfCoefficients(
fielddef->m_numModes,modes_offset);
if(fielddef->m_uniOrder == true) // reset modes_offset to zero
{
modes_offset = 0;
}
for(n = 0; n < nzmodes; ++n)
for(n = 0; n < nzmodes; ++n, offset += datalen)
{
planes_offset = fielddef->m_homogeneousZIDs[n];
std::map<int,int>::iterator it = homoZids.find(
fielddef->m_homogeneousZIDs[n]);
// Check to make sure this mode number lies in this field.
if (it == homoZids.end())
{
continue;
}
planes_offset = it->second;
if(datalen == (*m_exp)[eid]->GetNcoeffs())
{
Vmath::Vcopy(datalen,&fielddata[offset],1,&m_coeffs[m_coeff_offset[eid]+ planes_offset*ncoeffs_per_plane],1);
Vmath::Vcopy(datalen,&fielddata[offset],1,&coeffs[m_coeff_offset[eid]+planes_offset*ncoeffs_per_plane],1);
}
else // unpack data to new order
{
(*m_exp)[eid]->ExtractDataToCoeffs(fielddata, offset, fielddef->m_numModes,modes_offset,coeff_tmp = m_coeffs + m_coeff_offset[eid] + planes_offset*ncoeffs_per_plane);
(*m_exp)[eid]->ExtractDataToCoeffs(fielddata, offset, fielddef->m_numModes,modes_offset,coeff_tmp = coeffs + m_coeff_offset[eid] + planes_offset*ncoeffs_per_plane);
}
offset += datalen;
}
}
}
//Extract the data in fielddata into the m_coeff list (for 2D files into 3D cases)
void ExpListHomogeneous1D::v_ExtractDataToCoeffs(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, std::string &field, bool BaseFlow3D)
{
......@@ -739,7 +756,6 @@ namespace Nektar
int nzmodes = m_homogeneousBasis->GetNumModes();
int datalen = fielddata.size()/fielddef->m_fields.size();
int ncoeffs_per_plane = m_planes[0]->GetNcoeffs();
// Find data location according to field definition
for(i = 0; i < fielddef->m_fields.size(); ++i)
......
......@@ -184,7 +184,7 @@ namespace Nektar
virtual void v_AppendFieldData(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, Array<OneD, NekDouble> &coeffs);
virtual void v_ExtractDataToCoeffs(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, std::string &field);
virtual void v_ExtractDataToCoeffs(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, std::string &field, Array<OneD, NekDouble> &coeffs);
virtual void v_ExtractDataToCoeffs(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, std::string &field,bool BaseFlow3D);
......
......@@ -778,13 +778,24 @@ namespace Nektar
LibUtilities::EquationSharedPtr ffunc
= m_session->GetFunction(pFunctionName, pFieldName);
if (m_comm->GetRank() == 0)
{
cout << "- Field " << pFieldName << ": "
<< ffunc->GetExpression() << endl;
}
ffunc->Evaluate(x0,x1,x2,pTime,pArray);
}
else if (vType == LibUtilities::eFunctionTypeFile)
{
std::string filename
= m_session->GetFunctionFilename(pFunctionName, pFieldName);
cout << pFunctionName << " from file: " << filename << endl;
if (m_comm->GetRank() == 0)
{
cout << "- Field " << pFieldName << ": "
<< "from file " << filename << endl;
}
std::vector<SpatialDomains::FieldDefinitionsSharedPtr> FieldDef;
std::vector<std::vector<NekDouble> > FieldData;
......
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