Commit 0ba7ee30 authored by Dave Moxey's avatar Dave Moxey

Merge branch 'fix/GetBoundaryToElemtMap3DH1D-consistency' into 'master'

Fix/GetBoundaryToElemtMap3 dh1 d consistency

This MR changes the ordering from DisContField3DHomogeneous1D::GetBoundaryToElmtMap to make it consistent with m_bndCondExpansions (i.e. sorting by boundary, instead of sorting by plane). This allows the boundary expansion from 3DH1D to be treated in the same way as the non-homogeneous case, reducing the need for expansion-specific code.

All parts of the code assuming the old ordering were updated to take these changes into account. 

See merge request !526
parents f5bf354e 6c96c9d9
......@@ -350,14 +350,26 @@ namespace Nektar
// If this mesh (or partition) has no BCs, skip this step
if (MapSize > 0)
{
for(int i = 0; i < nplanes; i++)
int i ,j, n, cnt;
int cntPlane = 0;
for (cnt=n=0; n < m_bndCondExpansions.num_elements(); ++n)
{
for(int j = 0; j < MapSize; j++)
int planeExpSize = m_planes[0]
->GetBndCondExpansions()[n]
->GetExpSize();
for (i = 0; i < planeExpSize ; ++i, ++cntPlane)
{
ElmtID[j+i*MapSize] = ElmtID_tmp[j]+i*nel_per_plane;
EdgeID[j+i*MapSize] = EdgeID_tmp[j];
for(j = 0; j < nplanes; j++)
{
ElmtID[cnt+i+j*planeExpSize] =
ElmtID_tmp[cntPlane]+j*nel_per_plane;
EdgeID[cnt+i+j*planeExpSize] =
EdgeID_tmp[cntPlane];
}
}
cnt += m_bndCondExpansions[n]->GetExpSize();
}
m_BCtoElmMap = Array<OneD, int>(nplanes*MapSize);
m_BCtoEdgMap = Array<OneD, int>(nplanes*MapSize);
......
......@@ -97,6 +97,7 @@ IF( NEKTAR_SOLVER_INCNAVIERSTOKES )
ADD_NEKTAR_TEST(ChanFlow_3DH2D_FFT)
ADD_NEKTAR_TEST(CylFlow_MovBody)
ADD_NEKTAR_TEST(KovaFlow_3DH1D_MVM_FFTW_Consistency)
ADD_NEKTAR_TEST(KovaFlow_m8_short_HOBC_3D1H)
ENDIF (NEKTAR_USE_FFTW)
IF (NEKTAR_USE_MPI)
......
......@@ -473,7 +473,6 @@ namespace Nektar
m_session->LoadParameter("Delta_HighOrderBC",delta,1/20.0);
cnt = 0;
int count = 0;
for(int n = 0; n < m_PBndConds.num_elements(); ++n)
{
cnt_start = cnt;
......@@ -484,17 +483,9 @@ namespace Nektar
if (m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
{
int cnt_exp = 0; int cnt_plane = 0;
int veloffset = 0;
for(int i = 0; i < m_PBndExp[n]->GetExpSize(); ++i, cnt_exp++)
for(int i = 0; i < m_PBndExp[n]->GetExpSize(); ++i, cnt++)
{
// count the expansion list in each plane for e3DH1D case
if(cnt_exp == m_expsize_per_plane[n])
{
cnt_exp = 0; cnt_plane++;
}
int cnt = cnt_plane * m_totexps_per_plane + cnt_exp + count;
// find element and edge of this expansion.
Bc = boost::dynamic_pointer_cast<StdRegions::StdExpansion>
(m_PBndExp[n]->GetExp(i));
......@@ -532,16 +523,10 @@ namespace Nektar
m_PhyoutfVel[j][0]);
}
cnt_plane = 0; cnt_exp = 0;
cnt = cnt_start;
veloffset = 0;
for(int i = 0; i < m_PBndExp[n]->GetExpSize(); ++i, cnt_exp++)
for(int i = 0; i < m_PBndExp[n]->GetExpSize(); ++i, cnt++)
{
// count the expansion list for each plane for e3DH1D
if(cnt_exp == m_expsize_per_plane[n])
{
cnt_exp = 0; cnt_plane++;
}
cnt = cnt_plane * m_totexps_per_plane + cnt_exp + count;
int elmtid = m_pressureBCtoElmtID[cnt];
elmt = m_fields[0]->GetExp(elmtid);
......@@ -622,21 +607,10 @@ namespace Nektar
}
}
cnt = cnt_start;
veloffset = 0;
int cnt_exp = 0; int cnt_plane = 0; //only useful in e3DH1D case
for(int i = 0; i < m_PBndExp[n]->GetExpSize(); ++i,cnt++)
{
if (m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
{
// count the expansion list for e3DH1D
if(cnt_exp == m_expsize_per_plane[n])
{
cnt_exp = 0; cnt_plane++;
}
cnt = cnt_plane * m_totexps_per_plane + cnt_exp + count;
cnt_exp++;
}
// find element and edge of this expansion.
Bc = boost::dynamic_pointer_cast<StdRegions::StdExpansion>
(m_PBndExp[n]->GetExp(i));
......@@ -824,23 +798,9 @@ namespace Nektar
if(boost::iequals(UBndConds[j][n]->GetUserDefined(),"HOutflow"))
{
cnt = cnt_start;
int cnt_exp = 0; int cnt_plane = 0; //only useful in e3DH1D case
for(int i = 0; i < UBndExp[0][n]->GetExpSize();
++i, cnt++)
{
if(m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
{
// count the expansion order for e3DH1D
if(cnt_exp == m_expsize_per_plane[n])
{
cnt_exp = 0; cnt_plane++;
}
cnt = cnt_plane * m_totexps_per_plane + cnt_exp + count;
cnt_exp++;
}
Pbc = StdRegions::StdExpansionSharedPtr
(m_PBndExp[n]->GetExp(i));
Bc = StdRegions::StdExpansionSharedPtr
......@@ -879,10 +839,6 @@ namespace Nektar
else
{
cnt += m_PBndExp[n]->GetExpSize();
if(m_fields[0]->GetExpType() == MultiRegions::e3DH1D)
{
count += m_expsize_per_plane[n];
}
}
}
}
......@@ -1194,7 +1150,7 @@ namespace Nektar
m_negWavenumberSq = Array<OneD, NekDouble>(HBCnumber);
int exp_size, exp_size_per_plane;
int j=0;
int i, j, k, n;
int K;
NekDouble sign = -1.0;
int cnt = 0;
......@@ -1232,31 +1188,19 @@ namespace Nektar
m_npointsZ = m_session->GetParameter("HomModesZ");
}
Array<OneD, int> coeff_count(m_PBndConds.num_elements(),0);
Array<OneD, int> coeffPlaneOffset(m_PBndConds.num_elements(),0);
cnt = 0;
for(int n = 0 ; n < m_PBndConds.num_elements(); ++n)
{
coeffPlaneOffset[n] = cnt;
if(boost::iequals(m_PBndConds[n]->GetUserDefined(),"H"))
{
cnt += m_PBndExp[n]->GetNcoeffs();
}
}
cnt = 0;
for(int k = 0; k < num_planes; k++)
int coeff_count = 0;
for(n = 0, j= 0, cnt = 0; n < m_PBndConds.num_elements(); ++n)
{
K = planes[k]/2;
for(int n = 0 ; n < m_PBndConds.num_elements(); ++n)
exp_size = m_PBndExp[n]->GetExpSize();
exp_size_per_plane = exp_size/num_planes;
if(boost::iequals(m_PBndConds[n]->GetUserDefined(),"H"))
{
exp_size = m_PBndExp[n]->GetExpSize();
exp_size_per_plane = exp_size/num_planes;
if(boost::iequals(m_PBndConds[n]->GetUserDefined(),"H"))
for(k = 0; k < num_planes; k++)
{
for(int i = 0; i < exp_size_per_plane; ++i,cnt++)
K = planes[k]/2;
for(i = 0; i < exp_size_per_plane; ++i, ++j, ++cnt)
{
m_HBCdata[j].m_globalElmtID = m_pressureBCtoElmtID[cnt];
elmt = m_pressure->GetExp(m_HBCdata[j].m_globalElmtID);
......@@ -1265,8 +1209,8 @@ namespace Nektar
m_HBCdata[j].m_bndElmtOffset = i+k*exp_size_per_plane;
m_HBCdata[j].m_elmtTraceID = m_pressureBCtoTraceID[cnt];
m_HBCdata[j].m_bndryElmtID = n;
m_HBCdata[j].m_coeffOffset = coeffPlaneOffset[n] + coeff_count[n];
coeff_count[n] += elmt->GetEdgeNcoeffs(m_HBCdata[j].m_elmtTraceID);
m_HBCdata[j].m_coeffOffset = coeff_count;
coeff_count += elmt->GetEdgeNcoeffs(m_HBCdata[j].m_elmtTraceID);
if(m_SingleMode)
{
......@@ -1283,15 +1227,15 @@ namespace Nektar
m_wavenumber[j] = 2*M_PI*sign*(NekDouble(K))/m_LhomZ;
m_negWavenumberSq[j] = -1.0*m_wavenumber[j]*m_wavenumber[j];
}
int assElmtID;
if(k%2==0)
{
if(m_HalfMode)
{
assElmtID = m_HBCdata[j].m_globalElmtID;
}
else
{
......@@ -1302,18 +1246,16 @@ namespace Nektar
{
assElmtID = m_HBCdata[j].m_globalElmtID - num_elm_per_plane;
}
m_HBCdata[j].m_assPhysOffset = m_pressure->GetPhys_Offset(assElmtID);
j = j+1;
}
sign = -1.0*sign;
}
else // setting if just standard BC no High order
{
cnt += exp_size_per_plane;
}
}
sign = -1.0*sign;
else
{
cnt += exp_size;
}
}
}
break;
......
......@@ -175,14 +175,7 @@ void OutputFld::Process(po::variables_map &vm)
int cnt = 0;
for(int n = 0; n < Border; ++n)
{
if (m_f->m_fielddef[0]->m_numHomogeneousDir == 1)
{
cnt += BndExp[0][n]->GetPlane(0)->GetExpSize();
}
else
{
cnt += BndExp[0][n]->GetExpSize();
}
cnt += BndExp[0][n]->GetExpSize();
}
Array<OneD, NekDouble> tmp_array;
......@@ -194,27 +187,14 @@ void OutputFld::Process(po::variables_map &vm)
}
// setup phys arrays of normals
for(int j = 0; j < BndExp[0][Border]->GetExpSize();
++j)
for(int j=0; j < BndExp[0][Border]->GetExpSize(); ++j,++cnt)
{
// find element and face of this expansion.
int cnt2;
if (m_f->m_fielddef[0]->m_numHomogeneousDir == 1)
{
int exp_per_plane = BndExp[0][Border]->GetPlane(0)->
GetExpSize();
cnt2 = cnt + (j%exp_per_plane);
}
else
{
cnt2 = cnt+j;
}
int elmtid = BoundarytoElmtID[cnt2];
int elmtid = BoundarytoElmtID[cnt];
elmt = m_f->m_exp[0]->GetExp(elmtid);
//identify boundary of element looking at.
int boundary = BoundarytoTraceID[cnt2];
int boundary = BoundarytoTraceID[cnt];
// Dimension specific part
switch(normdim)
......
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