Commit c0e7ecfb authored by David Moxey's avatar David Moxey
Browse files

Add new inheritance structure to remove duplicate code inside staticcond classes

parent c60ce231
......@@ -34,6 +34,7 @@ GlobalLinSysDirectStaticCond.cpp
GlobalLinSysIterative.cpp
GlobalLinSysIterativeFull.cpp
GlobalLinSysIterativeStaticCond.cpp
GlobalLinSysStaticCond.cpp
GlobalMatrix.cpp
GlobalMatrixKey.cpp
GlobalOptimizationParameters.cpp
......@@ -74,6 +75,7 @@ GlobalLinSysDirectStaticCond.h
GlobalLinSysIterative.h
GlobalLinSysIterativeFull.h
GlobalLinSysIterativeStaticCond.h
GlobalLinSysStaticCond.h
GlobalMatrix.h
GlobalMatrixKey.h
GlobalOptimizationParameters.h
......
......@@ -1617,7 +1617,7 @@ namespace Nektar
int nq = (*m_exp)[expansion]->GetTotPoints();
// printing the fields of that zone
outfile << " <DataArray type=\"Float32\" Name=\""
outfile << " <DataArray type=\"Float64\" Name=\""
<< var << "\">" << endl;
outfile << " ";
const Array<OneD, NekDouble> phys = m_phys + m_phys_offset[expansion];
......
......@@ -72,22 +72,12 @@ namespace Nektar
const AssemblyMapSharedPtr &pLocToGloMap,
const int pNumDir)
{
DNekVec Vin(pInput.num_elements(),pInput);
// DNekVec Vout(pOutput.num_elements(),pOutput,eWrapper);
ASSERTL1(pInput.num_elements() <= pOutput.num_elements(),
"output array must be at least as long as input array");
DNekVec Vout(pInput.num_elements(),pOutput,eWrapper);
m_linSys->Solve(Vin,Vout);
}
const int nHomDofs = pNumRows - pNumDir;
/// Solve the linear system for given input and output vectors
/// using a specified local to global map.
void GlobalLinSysDirect::v_Solve( const Array<OneD, const NekDouble> &in,
Array<OneD, NekDouble> &out,
const AssemblyMapSharedPtr &pLocToGloMap,
const Array<OneD, const NekDouble> &pDirForcing)
{
ASSERTL0(false, "Not implemented for this GlobalLinSys type.");
DNekVec Vin (nHomDofs, pInput + pNumDir);
DNekVec Vout(nHomDofs, pOutput + pNumDir, eWrapper);
m_linSys->Solve(Vin, Vout);
}
}
}
......@@ -46,44 +46,26 @@ namespace Nektar
class ExpList;
/// A global linear system.
class GlobalLinSysDirect : public GlobalLinSys
class GlobalLinSysDirect : virtual public GlobalLinSys
{
public:
/// Default constructor
// MULTI_REGIONS_EXPORT GlobalLinSysDirect(void);
/// Constructor for full direct matrix solve.
MULTI_REGIONS_EXPORT GlobalLinSysDirect(
const GlobalLinSysKey &pKey,
const boost::weak_ptr<ExpList> &pExp,
const boost::shared_ptr<AssemblyMap>
&pLocToGloMap);
public:
MULTI_REGIONS_EXPORT GlobalLinSysDirect(
const GlobalLinSysKey &pKey,
const boost::weak_ptr<ExpList> &pExp,
const boost::shared_ptr<AssemblyMap> &pLocToGloMap);
MULTI_REGIONS_EXPORT virtual ~GlobalLinSysDirect();
protected:
/// Basic linear system object.
DNekLinSysSharedPtr m_linSys;
/// Solve the linear system for given input and output vectors
/// using a specified local to global map.
virtual void v_Solve(
const Array<OneD, const NekDouble> &in,
Array<OneD, NekDouble> &out,
const AssemblyMapSharedPtr &locToGloMap,
const Array<OneD, const NekDouble> &dirForcing = NullNekDouble1DArray);
/// Solve the linear system for given input and output vectors.
virtual void v_SolveLinearSystem(
const int pNumRows,
const Array<OneD,const NekDouble> &pInput,
Array<OneD, NekDouble> &pOutput,
const AssemblyMapSharedPtr &locToGloMap,
const int pNumDir = 0);
private:
protected:
/// Basic linear system object.
DNekLinSysSharedPtr m_linSys;
/// Solve the linear system for given input and output vectors.
virtual void v_SolveLinearSystem(
const int pNumRows,
const Array<OneD,const NekDouble> &pInput,
Array<OneD, NekDouble> &pOutput,
const AssemblyMapSharedPtr &locToGloMap,
const int pNumDir = 0);
};
}
}
......
......@@ -81,7 +81,8 @@ namespace Nektar
const boost::weak_ptr<ExpList> &pExp,
const boost::shared_ptr<AssemblyMap>
&pLocToGloMap)
: GlobalLinSysDirect(pLinSysKey, pExp, pLocToGloMap)
: GlobalLinSys(pLinSysKey, pExp, pLocToGloMap),
GlobalLinSysDirect(pLinSysKey, pExp, pLocToGloMap)
{
ASSERTL1(m_linSysKey.GetGlobalSysSolnType()==eDirectFullMatrix,
......@@ -112,7 +113,7 @@ namespace Nektar
bool dirForcCalculated = (bool) pDirForcing.num_elements();
int nDirDofs = pLocToGloMap->GetNumGlobalDirBndCoeffs();
int nGlobDofs = pLocToGloMap->GetNumGlobalCoeffs();
Array<OneD, NekDouble> tmp(nGlobDofs), tmp2;
Array<OneD, NekDouble> tmp(nGlobDofs);
if(nDirDofs)
{
......@@ -126,8 +127,7 @@ namespace Nektar
}
else
{
// Calculate the dirichlet forcing and substract it
// from the rhs
// Calculate Dirichlet forcing and subtract it from the rhs
m_expList.lock()->GeneralMatrixOp(
m_linSysKey, pOutput, tmp, eGlobal);
......@@ -137,8 +137,7 @@ namespace Nektar
tmp.get(), 1);
}
SolveLinearSystem(nGlobDofs, tmp + nDirDofs,
tmp2 = pOutput + nDirDofs,
SolveLinearSystem(nGlobDofs, tmp, pOutput,
pLocToGloMap, nDirDofs);
}
else
......@@ -172,10 +171,10 @@ namespace Nektar
switch(m_linSysKey.GetMatrixType())
{
// case for all symmetric matices
case StdRegions::eMass:
case StdRegions::eHelmholtz:
case StdRegions::eLaplacian:
case StdRegions::eHybridDGHelmBndLam:
case StdRegions::eMass:
case StdRegions::eHelmholtz:
case StdRegions::eLaplacian:
case StdRegions::eHybridDGHelmBndLam:
{
if( (2*(bwidth+1)) < rows)
{
......@@ -192,18 +191,18 @@ namespace Nektar
::AllocateSharedPtr(rows, cols, zero,
matStorage);
}
break;
}
break;
case StdRegions::eLinearAdvectionReaction:
case StdRegions::eLinearAdvectionDiffusionReaction:
case StdRegions::eLinearAdvectionReaction:
case StdRegions::eLinearAdvectionDiffusionReaction:
{
matStorage = eFULL;
Gmat = MemoryManager<DNekMat>
::AllocateSharedPtr(rows, cols, zero,
matStorage);
break;
}
break;
default:
default:
{
NEKERROR(ErrorUtil::efatal, "Add MatrixType to switch "
"statement");
......
......@@ -36,6 +36,7 @@
#define NEKTAR_LIB_MULTIREGIONS_GLOBALLINSYSDIRECTSTATICCOND_H
#include <MultiRegions/GlobalLinSysDirect.h>
#include <MultiRegions/GlobalLinSysStaticCond.h>
#include <MultiRegions/MultiRegionsDeclspec.h>
namespace Nektar
......@@ -50,18 +51,21 @@ namespace Nektar
GlobalLinSysDirectStaticCondSharedPtr;
/// A global linear system.
class GlobalLinSysDirectStaticCond : public GlobalLinSysDirect
class GlobalLinSysDirectStaticCond : virtual public GlobalLinSysDirect,
virtual public GlobalLinSysStaticCond
{
public:
/// Creates an instance of this class
static GlobalLinSysSharedPtr create(
const GlobalLinSysKey &pLinSysKey,
const boost::weak_ptr<ExpList> &pExpList,
const boost::shared_ptr<AssemblyMap>
&pLocToGloMap)
const GlobalLinSysKey &pLinSysKey,
const boost::weak_ptr<ExpList> &pExpList,
const boost::shared_ptr<AssemblyMap> &pLocToGloMap)
{
return MemoryManager<GlobalLinSysDirectStaticCond>
GlobalLinSysDirectStaticCondSharedPtr ret =
MemoryManager<GlobalLinSysDirectStaticCond>
::AllocateSharedPtr(pLinSysKey, pExpList, pLocToGloMap);
ret->InitObject();
return ret;
}
/// Name of class
......@@ -70,74 +74,39 @@ namespace Nektar
/// Constructor for full direct matrix solve.
MULTI_REGIONS_EXPORT GlobalLinSysDirectStaticCond(
const GlobalLinSysKey &mkey,
const boost::weak_ptr<ExpList> &pExpList,
const boost::shared_ptr<AssemblyMap>
&locToGloMap);
const GlobalLinSysKey &mkey,
const boost::weak_ptr<ExpList> &pExpList,
const boost::shared_ptr<AssemblyMap> &locToGloMap);
/// Constructor for full direct matrix solve.
MULTI_REGIONS_EXPORT GlobalLinSysDirectStaticCond(
const GlobalLinSysKey &mkey,
const boost::weak_ptr<ExpList> &pExpList,
const DNekScalBlkMatSharedPtr pSchurCompl,
const DNekScalBlkMatSharedPtr pBinvD,
const DNekScalBlkMatSharedPtr pC,
const DNekScalBlkMatSharedPtr pInvD,
const boost::shared_ptr<AssemblyMap>
&locToGloMap);
// MULTI_REGIONS_EXPORT GlobalLinSysDirectStaticCond(
// const DNekScalBlkMatSharedPtr pSchurCompl,
// const DNekScalBlkMatSharedPtr pBinvD,
// const DNekScalBlkMatSharedPtr pC,
// const DNekScalBlkMatSharedPtr pInvD,
// const boost::shared_ptr<AssemblyMap>
// &pLocToGloMap);
const GlobalLinSysKey &mkey,
const boost::weak_ptr<ExpList> &pExpList,
const DNekScalBlkMatSharedPtr pSchurCompl,
const DNekScalBlkMatSharedPtr pBinvD,
const DNekScalBlkMatSharedPtr pC,
const DNekScalBlkMatSharedPtr pInvD,
const boost::shared_ptr<AssemblyMap> &locToGloMap);
MULTI_REGIONS_EXPORT virtual ~GlobalLinSysDirectStaticCond();
private:
/// Schur complement for Direct Static Condensation.
GlobalLinSysDirectStaticCondSharedPtr m_recursiveSchurCompl;
/// Block matrices at this level
DNekScalBlkMatSharedPtr m_schurCompl;
DNekScalBlkMatSharedPtr m_BinvD;
DNekScalBlkMatSharedPtr m_C;
DNekScalBlkMatSharedPtr m_invD;
/// Solve the linear system for given input and output vectors
/// using a specified local to global map.
virtual void v_Solve(
const Array<OneD, const NekDouble> &in,
Array<OneD, NekDouble> &out,
const AssemblyMapSharedPtr &locToGloMap,
const Array<OneD, const NekDouble> &dirForcing
= NullNekDouble1DArray);
/// Initialise this object
void Initialise(
const boost::shared_ptr<AssemblyMap>& locToGloMap,
MatrixStorage matStorage);
protected:
virtual void v_AssembleSchurComplement(
boost::shared_ptr<AssemblyMap> pLocToGloMap);
virtual GlobalLinSysStaticCondSharedPtr v_Recurse(
const GlobalLinSysKey &mkey,
const boost::weak_ptr<ExpList> &pExpList,
const DNekScalBlkMatSharedPtr pSchurCompl,
const DNekScalBlkMatSharedPtr pBinvD,
const DNekScalBlkMatSharedPtr pC,
const DNekScalBlkMatSharedPtr pInvD,
const boost::shared_ptr<AssemblyMap> &l2gMap);
private:
/// Matrix Storage type for known matrices
MatrixStorage DetermineMatrixStorage(
const boost::shared_ptr<AssemblyMap>& locToGloMap);
/// Set up the storage for the Schur complement or the top level
/// of the multi-level Schur complement.
void SetupTopLevel(
const boost::shared_ptr<AssemblyMap>& locToGloMap);
/// Assemble the Schur complement matrix.
void AssembleSchurComplement(
const boost::shared_ptr<AssemblyMap>& locToGloMap,
const MatrixStorage matStorage);
///
void ConstructNextLevelCondensedSystem(
const boost::shared_ptr<AssemblyMap>& locToGloMap);
};
};
}
}
......
......@@ -166,7 +166,7 @@ namespace Nektar
// from the rhs
m_expList.lock()->GeneralMatrixOp(
m_linSysKey, pOutput, tmp, eGlobal);
Vmath::Vsub(nGlobDofs,
pInput.get(), 1,
tmp.get(), 1,
......
This diff is collapsed.
///////////////////////////////////////////////////////////////////////////////
//
// File: GlobalLinSysStaticCond.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: A collection of routines common to statically condensed systems.
//
///////////////////////////////////////////////////////////////////////////////
#ifndef NEKTAR_LIB_MULTIREGIONS_GLOBALLINSYSSTATICCOND_H
#define NEKTAR_LIB_MULTIREGIONS_GLOBALLINSYSSTATICCOND_H
#include <MultiRegions/GlobalMatrix.h>
#include <MultiRegions/GlobalLinSysIterative.h>
#include <LibUtilities/LinearAlgebra/SparseMatrixFwd.hpp>
namespace Nektar
{
namespace MultiRegions
{
// Forward declarations
class ExpList;
class GlobalLinSysStaticCond;
typedef boost::shared_ptr<GlobalLinSysStaticCond>
GlobalLinSysStaticCondSharedPtr;
/// A global linear system.
class GlobalLinSysStaticCond : virtual public GlobalLinSys
{
public:
/// Constructor for full direct matrix solve.
GlobalLinSysStaticCond(
const GlobalLinSysKey &mkey,
const boost::weak_ptr<ExpList> &pExpList,
const boost::shared_ptr<AssemblyMap> &locToGloMap);
virtual ~GlobalLinSysStaticCond();
protected:
virtual DNekScalBlkMatSharedPtr v_PreSolve(
int scLevel,
NekVector<NekDouble> &F_GlobBnd)
{
return m_schurCompl;
}
virtual void v_BasisTransform(
Array<OneD, NekDouble>& pInOut,
int offset)
{
}
virtual void v_BasisInvTransform(
Array<OneD, NekDouble>& pInOut,
int offset)
{
}
virtual void v_AssembleSchurComplement(
boost::shared_ptr<AssemblyMap> pLoctoGloMap)
{
}
virtual int v_GetNumBlocks();
virtual DNekScalBlkMatSharedPtr v_GetStaticCondBlock(unsigned int n)
{
return GlobalLinSys::GetStaticCondBlock(n);
}
virtual GlobalLinSysStaticCondSharedPtr v_Recurse(
const GlobalLinSysKey &mkey,
const boost::weak_ptr<ExpList> &pExpList,
const DNekScalBlkMatSharedPtr pSchurCompl,
const DNekScalBlkMatSharedPtr pBinvD,
const DNekScalBlkMatSharedPtr pC,
const DNekScalBlkMatSharedPtr pInvD,
const boost::shared_ptr<AssemblyMap> &locToGloMap) = 0;
/// Schur complement for Direct Static Condensation.
GlobalLinSysStaticCondSharedPtr m_recursiveSchurCompl;
/// Block Schur complement matrix.
DNekScalBlkMatSharedPtr m_schurCompl;
/// Block \f$ BD^{-1} \f$ matrix.
DNekScalBlkMatSharedPtr m_BinvD;
/// Block \f$ C \f$ matrix.
DNekScalBlkMatSharedPtr m_C;
/// Block \f$ D^{-1} \f$ matrix.
DNekScalBlkMatSharedPtr m_invD;
/// Local to global map.
boost::shared_ptr<AssemblyMap> m_locToGloMap;
/// Workspace array for matrix multiplication
Array<OneD, NekDouble> m_wsp;
/// Solve the linear system for given input and output vectors
/// using a specified local to global map.
virtual void v_Solve(
const Array<OneD, const NekDouble> &in,
Array<OneD, NekDouble> &out,
const AssemblyMapSharedPtr &locToGloMap,
const Array<OneD, const NekDouble> &dirForcing
= NullNekDouble1DArray);
virtual void v_InitObject();
/// Initialise this object
void Initialise(
const boost::shared_ptr<AssemblyMap>& locToGloMap);
/// Set up the storage for the Schur complement or the top level
/// of the multi-level Schur complement.
void SetupTopLevel(
const boost::shared_ptr<AssemblyMap>& locToGloMap);
///
void ConstructNextLevelCondensedSystem(
const boost::shared_ptr<AssemblyMap>& locToGloMap);
};
}
}
#endif
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