Merge branch 'feature/ADRSolverProjection' into 'master'

Add projection equation system to ADRSolver

This MR adds a projection operator to the ADRSolver, which can be used to impose C^0 connectivity on an input forcing field and is useful in e.g. the creation of boundary conditions on complex geometries.

See merge request !675

CHANGELOG.md

@@ -19,6 +19,9 @@ v4.4.0
 - Adjust CMake dependencies to reduce compile time (!671)
 - Homogeneous1D dealiasing improvements (!622)
 
+**ADRSolver:**
+- Add a projection equation system for C^0 projections (!675)
+
 **APESolver:**
 - Use a continuous basefield projection and revert to constant c formulation (!664)
 - Added ability to compute CFL number (!664)

docs/user-guide/solvers/adr.tex

@@ -24,6 +24,8 @@ see the table below.
 \textbf{Equation to solve} & \textbf{EquationType} & \textbf{Dimensions}   &
 \textbf{Projections} \\
 \midrule
+$u = f$ &
+    \inltt{Projection} & All & Continuous/Discontinuous \\
 $\nabla^2 u = 0$ & 
     \inltt{Laplace} & All &  Continuous/Discontinuous	\\
 $\nabla^2 u  =  f$ & 

solvers/ADRSolver/CMakeLists.txt

@@ -5,9 +5,10 @@ CMAKE_DEPENDENT_OPTION(NEKTAR_SOLVER_ADR
 IF( NEKTAR_SOLVER_ADR )
     SET(ADRSolverSource 
        ./ADRSolver.cpp 
-       ./EquationSystems/Laplace.cpp
-       ./EquationSystems/Poisson.cpp       
        ./EquationSystems/Helmholtz.cpp
+       ./EquationSystems/Laplace.cpp
+       ./EquationSystems/Poisson.cpp
+       ./EquationSystems/Projection.cpp
        ./EquationSystems/SteadyAdvectionDiffusion.cpp
        ./EquationSystems/SteadyAdvectionDiffusionReaction.cpp
        ./EquationSystems/UnsteadyAdvection.cpp
@@ -131,6 +132,9 @@ IF( NEKTAR_SOLVER_ADR )
     ADD_NEKTAR_TEST_LENGTHY(ExDiffusion_2D_LFRSD_deformed)
     ADD_NEKTAR_TEST_LENGTHY(ExDiffusion_2D_LFRSD_deformed_GAUSS)
 
+    # C^0 projection
+    ADD_NEKTAR_TEST(Projection2D)
+
     ADD_NEKTAR_TEST(Helmholtz1D_8modes_DG)
     ADD_NEKTAR_TEST(Helmholtz1D_8modes)
     ADD_NEKTAR_TEST(Helmholtz1D_8nodes)

solvers/ADRSolver/EquationSystems/Projection.cpp

+///////////////////////////////////////////////////////////////////////////////
+//
+// File Projection.cpp
+//
+// For more information, please see: http://www.nektar.info
+//
+// The MIT License
+//
+// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
+// Department of Aeronautics, Imperial College London (UK), and Scientific
+// Computing and Imaging Institute, University of Utah (USA).
+//
+// License for the specific language governing rights and limitations under
+// Permission is hereby granted, free of charge, to any person obtaining a
+// copy of this software and associated documentation files (the "Software"),
+// to deal in the Software without restriction, including without limitation
+// the rights to use, copy, modify, merge, publish, distribute, sublicense,
+// and/or sell copies of the Software, and to permit persons to whom the
+// Software is furnished to do so, subject to the following conditions:
+//
+// The above copyright notice and this permission notice shall be included
+// in all copies or substantial portions of the Software.
+//
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
+// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
+// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+// DEALINGS IN THE SOFTWARE.
+//
+// Description: Projection solve routines
+//
+///////////////////////////////////////////////////////////////////////////////
+
+#include <ADRSolver/EquationSystems/Projection.h>
+
+using namespace std;
+
+namespace Nektar
+{
+string Projection::className =
+    GetEquationSystemFactory().RegisterCreatorFunction("Projection",
+                                                       Projection::create);
+Projection::Projection(const LibUtilities::SessionReaderSharedPtr &pSession)
+    : EquationSystem(pSession)
+{
+}
+
+void Projection::v_InitObject()
+{
+    EquationSystem::v_InitObject();
+
+    EvaluateFunction(m_session->GetVariables(), m_fields, "Forcing");
+}
+
+Projection::~Projection()
+{
+}
+
+void Projection::v_DoSolve()
+{
+    for (int i = 0; i < m_fields.num_elements(); ++i)
+    {
+        // Zero field so initial conditions are zero
+        Vmath::Zero(m_fields[i]->GetNcoeffs(), m_fields[i]->UpdateCoeffs(), 1);
+        m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
+                              m_fields[i]->UpdateCoeffs());
+        m_fields[i]->SetPhysState(false);
+    }
+}
+
+void Projection::v_GenerateSummary(SolverUtils::SummaryList &s)
+{
+    EquationSystem::SessionSummary(s);
+    for (int i = 0; i < m_fields.num_elements(); ++i)
+    {
+        stringstream name;
+        name << "Forcing func [" << i << "]";
+        SolverUtils::AddSummaryItem(
+            s, name.str(),
+            m_session->GetFunction("Forcing", i)->GetExpression());
+    }
+}
+}

solvers/ADRSolver/EquationSystems/Projection.h

+///////////////////////////////////////////////////////////////////////////////
+//
+// File Projection.h
+//
+// For more information, please see: http://www.nektar.info
+//
+// The MIT License
+//
+// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
+// Department of Aeronautics, Imperial College London (UK), and Scientific
+// Computing and Imaging Institute, University of Utah (USA).
+//
+// License for the specific language governing rights and limitations under
+// Permission is hereby granted, free of charge, to any person obtaining a
+// copy of this software and associated documentation files (the "Software"),
+// to deal in the Software without restriction, including without limitation
+// the rights to use, copy, modify, merge, publish, distribute, sublicense,
+// and/or sell copies of the Software, and to permit persons to whom the
+// Software is furnished to do so, subject to the following conditions:
+//
+// The above copyright notice and this permission notice shall be included
+// in all copies or substantial portions of the Software.
+//
+// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
+// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
+// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+// DEALINGS IN THE SOFTWARE.
+//
+// Description: Projection solve routines
+//
+///////////////////////////////////////////////////////////////////////////////
+
+#ifndef NEKTAR_SOLVERS_ADRSOLVER_EQUATIONSYSTEMS_PROJECTION_H
+#define NEKTAR_SOLVERS_ADRSOLVER_EQUATIONSYSTEMS_PROJECTION_H
+
+#include <SolverUtils/EquationSystem.h>
+
+using namespace Nektar::SolverUtils;
+
+namespace Nektar
+{
+class Projection : public EquationSystem
+{
+public:
+    /// Class may only be instantiated through the MemoryManager.
+    friend class MemoryManager<Projection>;
+
+    /// Creates an instance of this class
+    static EquationSystemSharedPtr create(
+        const LibUtilities::SessionReaderSharedPtr &pSession)
+    {
+        EquationSystemSharedPtr p =
+            MemoryManager<Projection>::AllocateSharedPtr(pSession);
+        p->InitObject();
+        return p;
+    }
+
+    /// Name of class
+    static std::string className;
+
+    virtual ~Projection();
+
+protected:
+    Projection(const LibUtilities::SessionReaderSharedPtr &pSession);
+
+    virtual void v_InitObject();
+    virtual void v_GenerateSummary(SolverUtils::SummaryList &s);
+    virtual void v_DoSolve();
+
+private:
+};
+}
+
+#endif

solvers/ADRSolver/Tests/Projection2D.tst

+<?xml version="1.0" encoding="utf-8" ?>
+<test>
+    <description>ADRSolver C0 projection</description>
+    <executable>ADRSolver</executable>
+    <parameters>Projection2D.xml</parameters>
+    <files>
+        <file description="Session File">Projection2D.xml</file>
+    </files>
+    <metrics>
+        <metric type="L2" id="0">
+            <value variable="u" tolerance="1e-12">1.19001e-16</value>
+        </metric>
+        <metric type="Linf" id="1">
+            <value variable="u" tolerance="1e-12">4.16334e-16</value>
+        </metric>
+    </metrics>
+</test>

solvers/ADRSolver/Tests/Projection2D.xml

+<?xml version="1.0" encoding="utf-8" ?>
+<NEKTAR>
+    <!--
+        2 quad setup at linear order, 2 in y-direction, 1 in x-direction.
+        Set interior vertex to be 1, all other vertices to be 0 so that we get
+        a linear blend from 1 -> 0.
+    -->
+    <GEOMETRY DIM="2" SPACE="3">
+        <VERTEX COMPRESSED="B64Z-LittleEndian" BITSIZE="64">eJxjYMAOZs0EgZ32jCiiH+zR5Zmwyj+wh8kzYzUdIc+CVf8HuDwrVv0IeQCHGiOT</VERTEX>
+        <EDGE COMPRESSED="B64Z-LittleEndian" BITSIZE="64">eJxjYEAFjDhoJhw0Mw4aBljQ1MP4rDj4bGh8mHkADjAAOQAA</EDGE>
+        <ELEMENT>
+            <Q COMPRESSED="B64Z-LittleEndian" BITSIZE="64">eJxjYEAFjFCaCUoz4xBngdKsUJoNSgMAAtgAGQAA</Q>
+        </ELEMENT>
+        <COMPOSITE>
+            <C ID="0"> Q[0-1] </C>
+            <C ID="1"> E[0] </C>
+            <C ID="2"> E[1,4] </C>
+            <C ID="3"> E[5] </C>
+            <C ID="4"> E[3,6] </C>
+        </COMPOSITE>
+        <DOMAIN> C[0] </DOMAIN>
+    </GEOMETRY>
+    <EXPANSIONS>
+        <E COMPOSITE="C[0]"
+           NUMMODES="2,2" BASISTYPE="Modified_A,Modified_A"
+           NUMPOINTS="5,5" POINTSTYPE="GaussLobattoLegendre,GaussLobattoLegendre"
+           FIELDS="u" />
+    </EXPANSIONS>
+    <CONDITIONS>
+        <SOLVERINFO>
+            <I PROPERTY="EQTYPE" VALUE="Projection" />
+            <I PROPERTY="Projection" VALUE="Continuous" />
+        </SOLVERINFO>
+
+        <VARIABLES>
+            <V ID="0"> u </V>
+        </VARIABLES>
+
+        <BOUNDARYREGIONS>
+            <!-- LHS edge -->
+            <B ID="0"> C[4] </B>
+            <!-- other edges -->
+            <B ID="1"> C[1-3] </B>
+        </BOUNDARYREGIONS>
+
+        <BOUNDARYCONDITIONS>
+            <!-- This point is first to set LHS vertex mode to 1 -->
+            <REGION REF="0">
+                <D VAR="u" VALUE="1" />
+            </REGION>
+            <!-- All other points 0 -->
+            <REGION REF="1">
+                <D VAR="u" VALUE="0" />
+            </REGION>
+        </BOUNDARYCONDITIONS>
+
+        <!-- Not important since no interior modes -->
+        <FUNCTION NAME="Forcing">
+            <E VAR="u" VALUE="0" />
+        </FUNCTION>
+
+        <!-- Exact solution is a hat mode with linear blend in x-direction -->
+        <FUNCTION NAME="ExactSolution">
+            <E VAR="u" VALUE="(1.0-abs(2*(y-0.5)))*(1-x)" />
+        </FUNCTION>
+    </CONDITIONS>
+</NEKTAR>