Commit 3758077b authored by Dave Moxey's avatar Dave Moxey

Merge branch 'feature/varP-periodic' into 'master'

Variable P with periodic boundary conditions

This MR extends the variable P functionality to simulations with periodic boundary conditions.

This was a lot easier than I expected, but it seems to be working fine.

See merge request !658
parents 7798193c b4bd897a
......@@ -6,6 +6,8 @@ v4.4.0
**Library:**
- Add support for variable polynomial order for 3D simulations with continuous
Galerkin discretisation (!604)
- Add support for variable polynomial order with periodic boundary conditions
(!658)
**IncNavierStokesSolver:**
- Add ability to simulate additional scalar fields (!624)
......
......@@ -1305,7 +1305,7 @@ namespace Nektar
int p, q, numModes0, numModes1;
int cnt = 0;
int intDofCnt;
int meshVertId, meshEdgeId, meshFaceId;
int meshVertId, meshEdgeId, meshEdgeId2, meshFaceId, meshFaceId2;
int globalId;
int nEdgeInteriorCoeffs;
int firstNonDirGraphVertId;
......@@ -1395,6 +1395,32 @@ namespace Nektar
}
}
}
// Add non-local periodic dofs to the map
for (pIt = periodicEdges.begin(); pIt != periodicEdges.end(); ++pIt)
{
for (i = 0; i < pIt->second.size(); ++i)
{
meshEdgeId2 = pIt->second[i].id;
if (dofs[1].count(meshEdgeId2) == 0)
{
dofs[1][meshEdgeId2] = 1e6;
}
}
}
for (pIt = periodicFaces.begin(); pIt != periodicFaces.end(); ++pIt)
{
for (i = 0; i < pIt->second.size(); ++i)
{
meshFaceId2 = pIt->second[i].id;
if (faceModes[0].count(meshFaceId2) == 0)
{
faceModes[0][meshFaceId2] = 1e6;
faceModes[1][meshFaceId2] = 1e6;
}
}
}
// Now use information from all partitions to determine
// the correct size
map<int, int>::iterator dofIt, dofIt2;
......@@ -1413,6 +1439,19 @@ namespace Nektar
{
dofs[1][edgeId[i]] = (int) (edgeDof[i]+0.5);
}
// Periodic edges
for (pIt = periodicEdges.begin(); pIt != periodicEdges.end(); ++pIt)
{
meshEdgeId = pIt->first;
for (i = 0; i < pIt->second.size(); ++i)
{
meshEdgeId2 = pIt->second[i].id;
if (dofs[1][meshEdgeId2] < dofs[1][meshEdgeId])
{
dofs[1][meshEdgeId] = dofs[1][meshEdgeId2];
}
}
}
// faces
Array<OneD, long> faceId (faceModes[0].size());
Array<OneD, NekDouble> faceP (faceModes[0].size());
......@@ -1428,11 +1467,32 @@ namespace Nektar
Gs::Gather(faceP, Gs::gs_min, tmp2);
Gs::Gather(faceQ, Gs::gs_min, tmp2);
Gs::Finalise(tmp2);
int P, Q;
for (i=0; i < faceModes[0].size(); i++)
{
faceModes[0][faceId[i]] = (int) (faceP[i]+0.5);
faceModes[1][faceId[i]] = (int) (faceQ[i]+0.5);
}
// Periodic faces
for (pIt = periodicFaces.begin(); pIt != periodicFaces.end(); ++pIt)
{
meshFaceId = pIt->first;
for (i = 0; i < pIt->second.size(); ++i)
{
meshFaceId2 = pIt->second[i].id;
if (faceModes[0][meshFaceId2] < faceModes[0][meshFaceId])
{
faceModes[0][meshFaceId] = faceModes[0][meshFaceId2];
}
if (faceModes[1][meshFaceId2] < faceModes[1][meshFaceId])
{
faceModes[1][meshFaceId] = faceModes[1][meshFaceId2];
}
}
}
// Calculate number of dof in each face
int P, Q;
for (i=0; i < faceModes[0].size(); i++)
{
P = faceModes[0][faceId[i]];
Q = faceModes[1][faceId[i]];
if (faceType[faceId[i]] == LibUtilities::eQuadrilateral)
......@@ -1912,7 +1972,7 @@ namespace Nektar
if (m_staticCondLevel < (bottomUpGraph->GetNlevels()-1))
{
Array<OneD, int> vwgts_perm(
dofs[0].size() + dofs[1].size() + dofs[2].size()
graph[0].size() + graph[1].size() + graph[2].size()
- firstNonDirGraphVertId);
for (i = 0; i < locExpVector.size(); ++i)
......
......@@ -66,6 +66,7 @@ IF( NEKTAR_SOLVER_INCNAVIERSTOKES )
ADD_NEKTAR_TEST(KovaFlow_m8)
ADD_NEKTAR_TEST(KovaFlow_m8_short_HOBC)
ADD_NEKTAR_TEST(KovaFlow_varP)
ADD_NEKTAR_TEST(KovaFlow_varP_per)
#ADD_NEKTAR_TEST(KovaFlow_Oseen_m8)
ADD_NEKTAR_TEST_LENGTHY(KovaFlow_3DH1D_P5_20modes_MVM)
ADD_NEKTAR_TEST_LENGTHY(KovaFlow_3DH1D_P5_20modes_MVM_Deal)
......@@ -122,6 +123,7 @@ IF( NEKTAR_SOLVER_INCNAVIERSTOKES )
ADD_NEKTAR_TEST(ChanFlow_3DH1D_Parallel_mode2)
ADD_NEKTAR_TEST(ChanFlow_m3_par)
ADD_NEKTAR_TEST_LENGTHY(ChanFlow_m8_BodyForce_par)
ADD_NEKTAR_TEST(KovaFlow_varP_per_par)
ADD_NEKTAR_TEST_LENGTHY(Hex_channel_m8_par)
ADD_NEKTAR_TEST(Hex_channel_varP_par)
ADD_NEKTAR_TEST_LENGTHY(Pyr_channel_m6_par)
......
<?xml version="1.0" encoding="utf-8"?>
<test>
<description>Kovasznay Flow variable P, periodic BC</description>
<executable>IncNavierStokesSolver</executable>
<parameters>-I GlobalSysSoln=DirectMultiLevelStaticCond KovaFlow_varP_per.xml</parameters>
<files>
<file description="Session File">KovaFlow_varP_per.xml</file>
</files>
<metrics>
<metric type="L2" id="1">
<value variable="u" tolerance="1e-12">5.68743e-05</value>
<value variable="v" tolerance="1e-12">0.000123899</value>
<value variable="p" tolerance="1e-12">0.000569749</value>
</metric>
<metric type="Linf" id="2">
<value variable="u" tolerance="1e-12">0.00014278</value>
<value variable="v" tolerance="1e-12">0.000409628</value>
<value variable="p" tolerance="1e-12">0.00342738</value>
</metric>
</metrics>
</test>
<?xml version="1.0" encoding="utf-8" ?>
<NEKTAR xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
xsi:noNamespaceSchemaLocation="http://www.nektar.info/schema/nektar.xsd">
<EXPANSIONS>
<E COMPOSITE="C[7]" NUMMODES="7" FIELDS="u,v,p" TYPE="MODIFIED" />
<E COMPOSITE="C[8]" NUMMODES="6" FIELDS="u,v,p" TYPE="MODIFIED" />
<E COMPOSITE="C[9]" NUMMODES="8" FIELDS="u,v,p" TYPE="MODIFIED" />
<E COMPOSITE="C[10]" NUMMODES="6" FIELDS="u,v,p" TYPE="MODIFIED" />
<E COMPOSITE="C[11]" NUMMODES="9" FIELDS="u,v,p" TYPE="MODIFIED" />
</EXPANSIONS>
<CONDITIONS>
<SOLVERINFO>
<I PROPERTY="SolverType" VALUE="VelocityCorrectionScheme" />
<I PROPERTY="EQTYPE" VALUE="UnsteadyNavierStokes" />
<I PROPERTY="AdvectionForm" VALUE="Convective" />
<I PROPERTY="Projection" VALUE="Galerkin" />
<I PROPERTY="TimeIntegrationMethod" VALUE="IMEXOrder1" />
</SOLVERINFO>
<PARAMETERS>
<P> TimeStep = 0.001 </P>
<P> NumSteps = 10 </P>
<P> IO_CheckSteps = 1 </P>
<P> IO_InfoSteps = 1 </P>
<P> Kinvis = 0.025 </P>
</PARAMETERS>
<VARIABLES>
<V ID="0"> u </V>
<V ID="1"> v </V>
<V ID="2"> p </V>
</VARIABLES>
<BOUNDARYREGIONS>
<B ID="0"> C[1] </B>
<B ID="1"> C[2] </B>
<B ID="2"> C[3] </B>
<B ID="3"> C[4] </B>
</BOUNDARYREGIONS>
<BOUNDARYCONDITIONS>
<REGION REF="0">
<D VAR="u" VALUE="1-1.619099729265964*cos(2*PI*y)" />
<D VAR="v" VALUE="-0.248344108585656*sin(2*PI*y)" />
<N VAR="p" USERDEFINEDTYPE="H" VALUE="0" />
</REGION>
<REGION REF="1">
<D VAR="u" VALUE="1-0.381463333531742*cos(2*PI*y)" />
<D VAR="v" VALUE="-0.058510399212408*sin(2*PI*y)" />
<D VAR="p" VALUE="0.427242862585425" />
</REGION>
<REGION REF="2">
<P VAR="u" VALUE="[3]" />
<P VAR="v" VALUE="[3]" />
<P VAR="p" VALUE="[3]" />
</REGION>
<REGION REF="3">
<P VAR="u" VALUE="[2]" />
<P VAR="v" VALUE="[2]" />
<P VAR="p" VALUE="[2]" />
</REGION>
</BOUNDARYCONDITIONS>
<FUNCTION NAME="InitialConditions">
<E VAR="u" VALUE="(1-exp(-0.963740544195769*x)*cos(2*PI*y))" />
<E VAR="v"
VALUE="(-0.963740544195769/(2*PI))*exp(-0.963740544195769*x)*sin(2*PI*y)" />
<E VAR="p" VALUE="0.5*(1-exp(-2*0.963740544195769*x))" />
</FUNCTION>
<FUNCTION NAME="ExactSolution">
<E VAR="u" VALUE="(1-exp(-0.963740544195769*x)*cos(2*PI*y))" />
<E VAR="v"
VALUE="(-0.963740544195769/(2*PI))*exp(-0.963740544195769*x)*sin(2*PI*y)" />
<E VAR="p" VALUE="0.5*(1-exp(-2*0.963740544195769*x))" />
</FUNCTION>
</CONDITIONS>
<GEOMETRY DIM="2" SPACE="2">
<VERTEX COMPRESSED="B64Z-LittleEndian" BITSIZE="64">eJxjYEAGD/aj0gwMjFjlf9jDRJhQ5D/Yo8szY5VHmM+C1fwPcP2s2N0Hl2dDk18/MZnnf1aFLUyEHUpHO9a8Cq6u2YtuPweqfgz3caK5H2r+XpgIF5o8uvu4scoj/MeD1X5E+PGiut8WXZ4Pq36E+fxY5RHuE0CTdze//3OXWwLcf4JQ+vG+paKKtTm26OYLQemmsp+NV/7326CbL4zi/py9UPPh8QMANMBNrQAA</VERTEX>
<EDGE COMPRESSED="B64Z-LittleEndian" BITSIZE="64">eJx1kskSgjAQBUFckDUIigv6/5/pgelLVyWXV/3yMjWTpCj2dQjtQkvxoBzchlZickflx9CTuAk9i6l3UT6F1uJrRhupz7di5unE1Ou1Dw/iSX1PyiX5twyPYnTOKHUX+dzjXczcD/XPP1jFvOtTTN8v+bz3W8wcH/n8g03MHF/5rJ+Y+f7dGwRp</EDGE>
<ELEMENT>
<Q COMPRESSED="B64Z-LittleEndian" BITSIZE="64">eJx1zskOgyAABmFsRaq0bl3t4vs/ppfh4CTl8iVDAn8I+1PhAY9/eo0RG92XnvCErd4tvcOMZ/1zwR4H7YrqI056p1Gf8aq9Sf2Gd+0t+x/4xJf2tuoLvrW3U//gV3uz+g9X7d0A9rcC/gAA</Q>
</ELEMENT>
<COMPOSITE>
<C ID="7"> Q[0-2] </C> <!-- Domain -->
<C ID="8"> Q[3-5] </C> <!-- Domain -->
<C ID="9"> Q[6-7] </C> <!-- Domain -->
<C ID="10"> Q[8-9] </C> <!-- Domain -->
<C ID="11"> Q[10-11] </C> <!-- Domain -->
<C ID="1"> E[23,25,27,29] </C>
<C ID="2"> E[3,6,9,12] </C>
<C ID="3"> E[11,21,30] </C>
<C ID="4"> E[0,13,22] </C>
</COMPOSITE>
<DOMAIN> C[7-11] </DOMAIN>
</GEOMETRY>
</NEKTAR>
<?xml version="1.0" encoding="utf-8"?>
<test>
<description>Kovasznay Flow variable P, periodic BC (parallel)</description>
<executable>IncNavierStokesSolver</executable>
<parameters>--use-metis -I GlobalSysSoln=XxtMultiLevelStaticCond KovaFlow_varP_per.xml</parameters>
<processes> 3 </processes>
<files>
<file description="Session File">KovaFlow_varP_per.xml</file>
</files>
<metrics>
<metric type="L2" id="1">
<value variable="u" tolerance="1e-12">4.74673e-05</value>
<value variable="v" tolerance="1e-12">0.000201486</value>
<value variable="p" tolerance="1e-12">0.000565444</value>
</metric>
<metric type="Linf" id="2">
<value variable="u" tolerance="1e-12">0.00011373</value>
<value variable="v" tolerance="1e-12"> 0.000382783</value>
<value variable="p" tolerance="1e-12">0.00344565</value>
</metric>
</metrics>
</test>
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