A global linear system. More...

#include <GlobalLinSysPETScFull.h>

Inheritance diagram for Nektar::MultiRegions::GlobalLinSysPETScFull:

Public Member Functions
	GlobalLinSysPETScFull (const GlobalLinSysKey &pLinSysKey, const std::weak_ptr< ExpList > &pExpList, const std::shared_ptr< AssemblyMap > &pLocToGloMap)
	Constructor for full direct matrix solve. More...

	~GlobalLinSysPETScFull () override=default

Public Member Functions inherited from Nektar::MultiRegions::GlobalLinSysPETSc
	GlobalLinSysPETSc (const GlobalLinSysKey &pKey, const std::weak_ptr< ExpList > &pExp, const std::shared_ptr< AssemblyMap > &pLocToGloMap)
	Constructor for full direct matrix solve. More...

	~GlobalLinSysPETSc () override
	Clean up PETSc objects. More...

Public Member Functions inherited from Nektar::MultiRegions::GlobalLinSys
	GlobalLinSys (const GlobalLinSysKey &pKey, const std::weak_ptr< ExpList > &pExpList, const std::shared_ptr< AssemblyMap > &pLocToGloMap)
	Constructor for full direct matrix solve. More...

virtual	~GlobalLinSys ()

const GlobalLinSysKey &	GetKey (void) const
	Returns the key associated with the system. More...

const std::weak_ptr< ExpList > &	GetLocMat (void) const

void	InitObject ()

void	Initialise (const std::shared_ptr< AssemblyMap > &pLocToGloMap)

void	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 using a specified local to global map. More...

std::shared_ptr< GlobalLinSys >	GetSharedThisPtr ()
	Returns a shared pointer to the current object. More...

int	GetNumBlocks ()

DNekScalMatSharedPtr	GetBlock (unsigned int n)

void	DropBlock (unsigned int n)

DNekScalBlkMatSharedPtr	GetStaticCondBlock (unsigned int n)

void	DropStaticCondBlock (unsigned int n)

void	SolveLinearSystem (const int pNumRows, const Array< OneD, const NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const AssemblyMapSharedPtr &locToGloMap, const int pNumDir=0)
	Solve the linear system for given input and output vectors. More...

Static Public Member Functions
static GlobalLinSysSharedPtr	create (const GlobalLinSysKey &pLinSysKey, const std::weak_ptr< ExpList > &pExpList, const std::shared_ptr< AssemblyMap > &pLocToGloMap)
	Creates an instance of this class. More...

Static Public Attributes
static std::string	className
	Name of class. More...

Protected Member Functions
void	v_Solve (const Array< OneD, const NekDouble > &in, Array< OneD, NekDouble > &out, const AssemblyMapSharedPtr &locToGloMap, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray) override
	Solve the linear system for given input and output vectors using a specified local to global map. More...

void	v_DoMatrixMultiply (const Array< OneD, const NekDouble > &input, Array< OneD, NekDouble > &output) override
	Apply matrix-vector multiplication using local approach and the assembly map. More...

Protected Member Functions inherited from Nektar::MultiRegions::GlobalLinSysPETSc
void	SetUpScatter ()
	Set up PETSc local (equivalent to Nektar++ global) and global (equivalent to universal) scatter maps. More...

void	SetUpMatVec (int nGlobal, int nDir)
	Construct PETSc matrix and vector handles. More...

void	SetUpSolver (NekDouble tolerance)
	Set up KSP solver object. More...

void	CalculateReordering (const Array< OneD, const int > &glo2uniMap, const Array< OneD, const int > &glo2unique, const AssemblyMapSharedPtr &pLocToGloMap)
	Calculate a reordering of universal IDs for PETSc. More...

void	v_SolveLinearSystem (const int pNumRows, const Array< OneD, const NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const AssemblyMapSharedPtr &locToGloMap, const int pNumDir) override
	Solve linear system using PETSc. More...

virtual void	v_DoMatrixMultiply (const Array< OneD, const NekDouble > &pInput, Array< OneD, NekDouble > &pOutput)=0

Protected Member Functions inherited from Nektar::MultiRegions::GlobalLinSys
virtual void	v_Solve (const Array< OneD, const NekDouble > &in, Array< OneD, NekDouble > &out, const AssemblyMapSharedPtr &locToGloMap, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)=0
	Solve a linear system based on mapping. More...

virtual void	v_SolveLinearSystem (const int pNumRows, const Array< OneD, const NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const AssemblyMapSharedPtr &locToGloMap, const int pNumDir)=0
	Solve a basic matrix system. More...

virtual void	v_InitObject ()

virtual void	v_Initialise (const std::shared_ptr< AssemblyMap > &pLocToGloMap)

virtual int	v_GetNumBlocks ()
	Get the number of blocks in this system. More...

virtual DNekScalMatSharedPtr	v_GetBlock (unsigned int n)
	Retrieves the block matrix from n-th expansion using the matrix key provided by the m_linSysKey. More...

virtual void	v_DropBlock (unsigned int n)
	Releases the local block matrix from NekManager of n-th expansion using the matrix key provided by the m_linSysKey. More...

virtual DNekScalBlkMatSharedPtr	v_GetStaticCondBlock (unsigned int n)
	Retrieves a the static condensation block matrices from n-th expansion using the matrix key provided by the m_linSysKey. More...

virtual void	v_DropStaticCondBlock (unsigned int n)
	Releases the static condensation block matrices from NekManager of n-th expansion using the matrix key provided by the m_linSysKey. More...

PreconditionerSharedPtr	CreatePrecon (AssemblyMapSharedPtr asmMap)
	Create a preconditioner object from the parameters defined in the supplied assembly map. More...

Private Attributes
std::shared_ptr< AssemblyMap >	m_locToGloMap

Additional Inherited Members
Protected Attributes inherited from Nektar::MultiRegions::GlobalLinSysPETSc
Mat	m_matrix
	PETSc matrix object. More...

Vec	m_x
	PETSc vector objects used for local storage. More...

Vec	m_b

Vec	m_locVec

KSP	m_ksp
	KSP object that represents solver system. More...

PC	m_pc
	PCShell for preconditioner. More...

PETScMatMult	m_matMult
	Enumerator to select matrix multiplication type. More...

std::vector< int >	m_reorderedMap
	Reordering that takes universal IDs to a unique row in the PETSc matrix. More...

VecScatter	m_ctx
	PETSc scatter context that takes us between Nektar++ global ordering and PETSc vector ordering. More...

int	m_nLocal
	Number of unique degrees of freedom on this process. More...

PreconditionerSharedPtr	m_precon

Protected Attributes inherited from Nektar::MultiRegions::GlobalLinSys
const GlobalLinSysKey	m_linSysKey
	Key associated with this linear system. More...

const std::weak_ptr< ExpList >	m_expList
	Local Matrix System. More...

const std::map< int, RobinBCInfoSharedPtr >	m_robinBCInfo
	Robin boundary info. More...

bool	m_verbose

Detailed Description

A global linear system.

Definition at line 50 of file GlobalLinSysPETScFull.h.

Constructor & Destructor Documentation

◆ GlobalLinSysPETScFull()

Nektar::MultiRegions::GlobalLinSysPETScFull::GlobalLinSysPETScFull	(	const GlobalLinSysKey &	pLinSysKey,
		const std::weak_ptr< ExpList > &	pExpList,
		const std::shared_ptr< AssemblyMap > &	pLocToGloMap
	)

Constructor for full direct matrix solve.

Definition at line 57 of file GlobalLinSysPETScFull.cpp.

    : GlobalLinSys(pLinSysKey, pExp, pLocToGloMap),
      GlobalLinSysPETSc(pLinSysKey, pExp, pLocToGloMap)
{
    const int nDirDofs = pLocToGloMap->GetNumGlobalDirBndCoeffs();
 
    int i, j, n, cnt, gid1, gid2, loc_lda;
    NekDouble sign1, sign2, value;
    DNekScalMatSharedPtr loc_mat;
 
    // CALCULATE REORDERING MAPPING
    CalculateReordering(pLocToGloMap->GetGlobalToUniversalMap(),
                        pLocToGloMap->GetGlobalToUniversalMapUnique(),
                        pLocToGloMap);
 
    // SET UP VECTORS AND MATRIX
    SetUpMatVec(pLocToGloMap->GetNumGlobalCoeffs(), nDirDofs);
 
    // SET UP SCATTER OBJECTS
    SetUpScatter();
 
    // CONSTRUCT KSP OBJECT
    LibUtilities::SessionReaderSharedPtr pSession =
        m_expList.lock()->GetSession();
    string variable                    = pLocToGloMap->GetVariable();
    NekDouble IterativeSolverTolerance = NekConstants::kNekIterativeTol;
    if (pSession->DefinesGlobalSysSolnInfo(variable,
                                           "IterativeSolverTolerance"))
    {
        IterativeSolverTolerance = boost::lexical_cast<double>(
            pSession->GetGlobalSysSolnInfo(variable, "IterativeSolverTolerance")
                .c_str());
    }
    else if (pSession->DefinesParameter("IterativeSolverTolerance"))
    {
        IterativeSolverTolerance =
            pSession->GetParameter("IterativeSolverTolerance");
    }
    SetUpSolver(IterativeSolverTolerance);
 
    // POPULATE MATRIX
    for (n = cnt = 0; n < m_expList.lock()->GetNumElmts(); ++n)
    {
        loc_mat = GetBlock(n);
        loc_lda = loc_mat->GetRows();
 
        for (i = 0; i < loc_lda; ++i)
        {
            gid1  = pLocToGloMap->GetLocalToGlobalMap(cnt + i) - nDirDofs;
            sign1 = pLocToGloMap->GetLocalToGlobalSign(cnt + i);
            if (gid1 >= 0)
            {
                int gid1ro = m_reorderedMap[gid1];
                for (j = 0; j < loc_lda; ++j)
                {
                    gid2 =
                        pLocToGloMap->GetLocalToGlobalMap(cnt + j) - nDirDofs;
                    sign2 = pLocToGloMap->GetLocalToGlobalSign(cnt + j);
                    if (gid2 >= 0)
                    {
                        int gid2ro = m_reorderedMap[gid2];
                        value      = sign1 * sign2 * (*loc_mat)(i, j);
                        MatSetValue(m_matrix, gid1ro, gid2ro, value,
                                    ADD_VALUES);
                    }
                }
            }
        }
        cnt += loc_lda;
    }
 
    // ASSEMBLE MATRIX
    MatAssemblyBegin(m_matrix, MAT_FINAL_ASSEMBLY);
    MatAssemblyEnd(m_matrix, MAT_FINAL_ASSEMBLY);
}

References Nektar::MultiRegions::GlobalLinSysPETSc::CalculateReordering(), Nektar::MultiRegions::GlobalLinSys::GetBlock(), Nektar::NekConstants::kNekIterativeTol, Nektar::MultiRegions::GlobalLinSys::m_expList, Nektar::MultiRegions::GlobalLinSysPETSc::m_matrix, Nektar::MultiRegions::GlobalLinSysPETSc::m_reorderedMap, Nektar::MultiRegions::GlobalLinSysPETSc::SetUpMatVec(), Nektar::MultiRegions::GlobalLinSysPETSc::SetUpScatter(), and Nektar::MultiRegions::GlobalLinSysPETSc::SetUpSolver().

◆ ~GlobalLinSysPETScFull()

Nektar::MultiRegions::GlobalLinSysPETScFull::~GlobalLinSysPETScFull ( )

overridedefault

Member Function Documentation

◆ create()

static GlobalLinSysSharedPtr Nektar::MultiRegions::GlobalLinSysPETScFull::create	(	const GlobalLinSysKey &	pLinSysKey,
		const std::weak_ptr< ExpList > &	pExpList,
		const std::shared_ptr< AssemblyMap > &	pLocToGloMap
	)

inlinestatic

Creates an instance of this class.

Definition at line 54 of file GlobalLinSysPETScFull.h.

    {
        return MemoryManager<GlobalLinSysPETScFull>::AllocateSharedPtr(
            pLinSysKey, pExpList, pLocToGloMap);
    }

References Nektar::MemoryManager< DataType >::AllocateSharedPtr().

◆ v_DoMatrixMultiply()

void Nektar::MultiRegions::GlobalLinSysPETScFull::v_DoMatrixMultiply	(	const Array< OneD, const NekDouble > &	input,
		Array< OneD, NekDouble > &	output
	)

overrideprotectedvirtual

Apply matrix-vector multiplication using local approach and the assembly map.

Parameters

input	Vector input.
output	Result of multiplication.

Implements Nektar::MultiRegions::GlobalLinSysPETSc.

Definition at line 219 of file GlobalLinSysPETScFull.cpp.

{
    std::shared_ptr<MultiRegions::ExpList> expList = m_expList.lock();
 
    int nLocDofs = m_locToGloMap->GetNumLocalCoeffs();
 
    Array<OneD, NekDouble> tmp(nLocDofs);
    Array<OneD, NekDouble> tmp1(nLocDofs);
 
    m_locToGloMap->GlobalToLocal(input, tmp);
 
    // Perform matrix-vector operation A*d_i
    expList->GeneralMatrixOp(m_linSysKey, tmp, tmp1);
 
    m_locToGloMap->Assemble(tmp1, output);
}

References Nektar::MultiRegions::GlobalLinSys::m_expList, Nektar::MultiRegions::GlobalLinSys::m_linSysKey, and m_locToGloMap.

◆ v_Solve()

void Nektar::MultiRegions::GlobalLinSysPETScFull::v_Solve	(	const Array< OneD, const NekDouble > &	pLocInput,
		Array< OneD, NekDouble > &	pLocOutput,
		const AssemblyMapSharedPtr &	pLocToGloMap,
		const Array< OneD, const NekDouble > &	pDirForcing = `NullNekDouble1DArray`
	)

overrideprotectedvirtual

Solve the linear system for given input and output vectors using a specified local to global map.

Solve the linear system using a full global matrix system.

Implements Nektar::MultiRegions::GlobalLinSys.

Definition at line 138 of file GlobalLinSysPETScFull.cpp.

{
    m_locToGloMap = pLocToGloMap;
 
    bool dirForcCalculated = (bool)pDirForcing.size();
    int nDirDofs           = pLocToGloMap->GetNumGlobalDirBndCoeffs();
    int nGlobDofs          = pLocToGloMap->GetNumGlobalCoeffs();
    int nLocDofs           = pLocToGloMap->GetNumLocalCoeffs();
 
    int nDirTotal                                  = nDirDofs;
    std::shared_ptr<MultiRegions::ExpList> expList = m_expList.lock();
    expList->GetComm()->GetRowComm()->AllReduce(nDirTotal,
                                                LibUtilities::ReduceSum);
 
    if (nDirTotal)
    {
        Array<OneD, NekDouble> rhs(nLocDofs);
 
        // Calculate the Dirichlet forcing
        if (dirForcCalculated)
        {
            // Assume pDirForcing is in local space
            ASSERTL0(
                pDirForcing.size() >= nLocDofs,
                "DirForcing is not of sufficient size. Is it in local space?");
            Vmath::Vsub(nLocDofs, pLocInput, 1, pDirForcing, 1, rhs, 1);
        }
        else
        {
            // Calculate initial condition and Dirichlet forcing and subtract it
            // from the rhs
            expList->GeneralMatrixOp(m_linSysKey, pLocOutput, rhs);
 
            // Iterate over all the elements computing Robin BCs where
            // necessary
            for (auto &r : m_robinBCInfo) // add robin mass matrix
            {
                RobinBCInfoSharedPtr rBC;
                Array<OneD, NekDouble> rhsloc;
 
                int n      = r.first;
                int offset = expList->GetCoeff_Offset(n);
 
                LocalRegions::ExpansionSharedPtr vExp = expList->GetExp(n);
                // Add local matrix contribution
                for (rBC = r.second; rBC; rBC = rBC->next)
                {
                    vExp->AddRobinTraceContribution(
                        rBC->m_robinID, rBC->m_robinPrimitiveCoeffs,
                        pLocOutput + offset, rhsloc = rhs + offset);
                }
            }
            Vmath::Vsub(nLocDofs, pLocInput, 1, rhs, 1, rhs, 1);
        }
 
        Array<OneD, NekDouble> diff(nLocDofs);
 
        // Solve for perturbation from initial guess in pOutput
        SolveLinearSystem(nGlobDofs, rhs, diff, pLocToGloMap, nDirDofs);
 
        // Add back initial and boundary condition
        Vmath::Vadd(nLocDofs, diff, 1, pLocOutput, 1, pLocOutput, 1);
    }
    else
    {
        SolveLinearSystem(nGlobDofs, pLocInput, pLocOutput, pLocToGloMap,
                          nDirDofs);
    }
}

References ASSERTL0, Nektar::MultiRegions::GlobalLinSys::m_expList, Nektar::MultiRegions::GlobalLinSys::m_linSysKey, m_locToGloMap, Nektar::MultiRegions::GlobalLinSys::m_robinBCInfo, Nektar::LibUtilities::ReduceSum, Nektar::MultiRegions::GlobalLinSys::SolveLinearSystem(), Vmath::Vadd(), and Vmath::Vsub().

Member Data Documentation

◆ className

string Nektar::MultiRegions::GlobalLinSysPETScFull::className

static

Initial value:

=
    GetGlobalLinSysFactory().RegisterCreatorFunction(
        "PETScFull", GlobalLinSysPETScFull::create, "PETSc Full Matrix.")

Name of class.

Registers the class with the Factory.

Definition at line 64 of file GlobalLinSysPETScFull.h.

◆ m_locToGloMap

std::shared_ptr<AssemblyMap> Nektar::MultiRegions::GlobalLinSysPETScFull::m_locToGloMap

private

Definition at line 87 of file GlobalLinSysPETScFull.h.

Referenced by v_DoMatrixMultiply(), and v_Solve().

Public Member Functions

Static Public Member Functions

Static Public Attributes

Protected Member Functions

Private Attributes

Additional Inherited Members