Commit 8f64a550 authored by Dave Moxey's avatar Dave Moxey
Browse files

Fix for reading of restart files in parallel.

parent 17b697cb
......@@ -2105,7 +2105,7 @@ namespace Nektar
void ExpList::v_ExtractDataToCoeffs(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, std::string &field, Array<OneD, NekDouble> &coeffs)
{
int i;
int i, cnt;
int offset = 0;
int datalen = fielddata.size()/fielddef->m_fields.size();
......@@ -2121,44 +2121,108 @@ namespace Nektar
if(i == fielddef->m_fields.size())
{
cerr << "Field (" << field << ") not found in data file. Setting it to zero. " << endl;
cerr << "Field (" << field << ") not found in data file; "
<< "Setting it to zero. " << endl;
Vmath::Zero(coeffs.num_elements(),coeffs,1);
}
else
{
// 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;
// loop in reverse order so that in case where using
// a Homogeneous expansion it sets geometry ids to
// first part of m_exp list. Otherwise will set to
// second (complex) expansion
// loop in reverse order so that in case where using a
// Homogeneous expansion it sets geometry ids to first part of
// m_exp list. Otherwise will set to second (complex) expansion
for(i = (*m_exp).size()-1; i >=0; --i)
{
ElmtID_to_ExpID[(*m_exp)[i]->GetGeom()->GetGlobalID()] = i;
}
// Communicate number of modes between processors.
LibUtilities::CommSharedPtr vComm =
m_session->GetComm()->GetRowComm();
int n = vComm->GetSize();
int p = vComm->GetRank();
Array<OneD, int> numEls(n, 0);
numEls[p] = m_exp->size();
vComm->AllReduce(numEls, LibUtilities::ReduceSum);
int totEls = Vmath::Vsum(n, numEls, 1);
Array<OneD, int> elOffsets(n, 0);
elOffsets[0] = 0;
for (i = 1; i < n; ++i)
{
elOffsets[i] = elOffsets[i-1] + numEls[i-1];
}
Array<OneD, int> coeffsPerEl (totEls, 0.0);
Array<OneD, int> elmtGlobalIds(totEls, 0.0);
// Determine number of coefficients in each local (to this
// partition) element.
int modes_offset = 0;
for(cnt = i = 0; i < fielddef->m_elementIDs.size(); ++i)
{
if (ElmtID_to_ExpID.count(fielddef->m_elementIDs[i]) == 0)
{
continue;
}
int eid = ElmtID_to_ExpID[fielddef->m_elementIDs[i]];
int datalen = (*m_exp)[eid]->CalcNumberOfCoefficients(
fielddef->m_numModes,modes_offset);
if(fielddef->m_uniOrder == true)
{
modes_offset = 0;
}
elmtGlobalIds[cnt+elOffsets[p]] = fielddef->m_elementIDs[i];
coeffsPerEl [cnt+elOffsets[p]] = datalen;
cnt++;
}
vComm->AllReduce(coeffsPerEl, LibUtilities::ReduceSum);
vComm->AllReduce(elmtGlobalIds, LibUtilities::ReduceSum);
map<int,int> coeffsElmtMap;
for (i = 0; i < totEls; ++i)
{
coeffsElmtMap[elmtGlobalIds[i]] = coeffsPerEl[i];
}
for(i = 0; i < fielddef->m_elementIDs.size(); ++i)
{
if (ElmtID_to_ExpID.count(fielddef->m_elementIDs[i]) == 0)
{
ASSERTL0(coeffsElmtMap.count(fielddef->m_elementIDs[i])
> 0, "Couldn't find element!");
offset += coeffsElmtMap[fielddef->m_elementIDs[i]];
continue;
}
int eid = ElmtID_to_ExpID[fielddef->m_elementIDs[i]];
int datalen = (*m_exp)[eid]->CalcNumberOfCoefficients(fielddef->m_numModes,modes_offset);
if(fielddef->m_uniOrder == true) // reset modes_offset to zero
int datalen = coeffsElmtMap[fielddef->m_elementIDs[i]];
if(fielddef->m_uniOrder == true)
{
modes_offset = 0;
}
// copy data if of same length as expansion
// Copy data if it is the same length as expansion.
if(datalen == (*m_exp)[eid]->GetNcoeffs())
{
Vmath::Vcopy(datalen,&fielddata[offset],1,&coeffs[m_coeff_offset[eid]],1);
Vmath::Vcopy(datalen, &fielddata[offset], 1,
&coeffs[m_coeff_offset[eid]], 1);
}
else // unpack data to new order
{
(*m_exp)[eid]->ExtractDataToCoeffs(&fielddata[offset], fielddef->m_numModes, modes_offset, &coeffs[m_coeff_offset[eid]]);
(*m_exp)[eid]->ExtractDataToCoeffs(
&fielddata[offset], fielddef->m_numModes,
modes_offset, &coeffs[m_coeff_offset[eid]]);
}
offset += datalen;
}
}
......
Supports Markdown
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