Helmholtz operator using LocalRegions implementation. More...

Inheritance diagram for Nektar::Collections::Helmholtz_IterPerExp:

Public Member Functions
	~Helmholtz_IterPerExp () final=default

void	operator() (const Array< OneD, const NekDouble > &input, Array< OneD, NekDouble > &output, Array< OneD, NekDouble > &output1, Array< OneD, NekDouble > &output2, Array< OneD, NekDouble > &wsp) final
	Perform operation. More...

void	operator() (int dir, const Array< OneD, const NekDouble > &input, Array< OneD, NekDouble > &output, Array< OneD, NekDouble > &wsp) final

void	CheckFactors (StdRegions::FactorMap factors, int coll_phys_offset) override
	Check the validity of supplied constant factors. More...

Public Member Functions inherited from Nektar::Collections::Operator
	Operator (std::vector< StdRegions::StdExpansionSharedPtr > pCollExp, std::shared_ptr< CoalescedGeomData > GeomData, StdRegions::FactorMap factors)
	Constructor. More...

virtual	~Operator ()=default

virtual COLLECTIONS_EXPORT void	operator() (const Array< OneD, const NekDouble > &input, Array< OneD, NekDouble > &output0, Array< OneD, NekDouble > &output1, Array< OneD, NekDouble > &output2, Array< OneD, NekDouble > &wsp=NullNekDouble1DArray)=0
	Perform operation. More...

virtual COLLECTIONS_EXPORT void	operator() (int dir, const Array< OneD, const NekDouble > &input, Array< OneD, NekDouble > &output, Array< OneD, NekDouble > &wsp=NullNekDouble1DArray)=0

virtual COLLECTIONS_EXPORT void	CheckFactors (StdRegions::FactorMap factors, int coll_phys_offset)=0
	Check the validity of the supplied factor map. More...

unsigned int	GetWspSize ()
	Get the size of the required workspace. More...

unsigned int	GetNumElmt ()
	Get number of elements. More...

StdRegions::StdExpansionSharedPtr	GetExpSharedPtr ()
	Get expansion pointer. More...

unsigned int	GetInputSize ()

unsigned int	GetOutputSize ()

Protected Attributes
Array< TwoD, const NekDouble >	m_derivFac

Array< OneD, const NekDouble >	m_jac

int	m_dim

int	m_coordim

StdRegions::FactorMap	m_factors

NekDouble	m_lambda

bool	m_HasVarCoeffDiff

Array< OneD, Array< OneD, NekDouble > >	m_diff

const StdRegions::ConstFactorType	m_factorCoeffDef [3][3]

Protected Attributes inherited from Nektar::Collections::Operator
bool	m_isDeformed

StdRegions::StdExpansionSharedPtr	m_stdExp

unsigned int	m_numElmt
	number of elements that the operator is applied on More...

unsigned int	m_nqe

unsigned int	m_wspSize

unsigned int	m_inputSize
	number of modes or quadrature points that are passed as input to an operator More...

unsigned int	m_outputSize
	number of modes or quadrature points that are taken as output from an operator More...

Private Member Functions
	Helmholtz_IterPerExp (vector< StdRegions::StdExpansionSharedPtr > pCollExp, CoalescedGeomDataSharedPtr pGeomData, StdRegions::FactorMap factors)

Private Member Functions inherited from Nektar::Collections::Helmholtz_Helper
	Helmholtz_Helper ()

Additional Inherited Members
Protected Member Functions inherited from Nektar::Collections::Operator
virtual int	v_GetInputSize ()
	This purely virtual function needs to be set-up for every operator inside Collections. It is responsible for returning the size of input collection, that the operator is applied on either in physical or coefficient space. More...

virtual int	v_GetOutputSize ()
	This purely virtual function needs to be set-up for every operator inside Collections. It is responsible for returning the output size either in physical or coefficient space of an operator inside Collections. More...

Detailed Description

Helmholtz operator using LocalRegions implementation.

Definition at line 173 of file library/Collections/Helmholtz.cpp.

Constructor & Destructor Documentation

◆ ~Helmholtz_IterPerExp()

Nektar::Collections::Helmholtz_IterPerExp::~Helmholtz_IterPerExp ( )

finaldefault

◆ Helmholtz_IterPerExp()

Nektar::Collections::Helmholtz_IterPerExp::Helmholtz_IterPerExp	(	vector< StdRegions::StdExpansionSharedPtr >	pCollExp,
		CoalescedGeomDataSharedPtr	pGeomData,
		StdRegions::FactorMap	factors
	)

inlineprivate

Definition at line 497 of file library/Collections/Helmholtz.cpp.

        : Operator(pCollExp, pGeomData, factors), Helmholtz_Helper()
    {
        m_dim     = pCollExp[0]->GetShapeDimension();
        m_coordim = pCollExp[0]->GetCoordim();
        int nqtot = m_stdExp->GetTotPoints();
 
        m_derivFac = pGeomData->GetDerivFactors(pCollExp);
        m_jac      = pGeomData->GetJac(pCollExp);
        m_wspSize  = (2 * m_coordim + 1) * nqtot;
 
        m_lambda          = 1.0;
        m_HasVarCoeffDiff = false;
        m_factors         = StdRegions::NullFactorMap;
        this->CheckFactors(factors, 0);
    }

References CheckFactors(), m_coordim, m_derivFac, m_dim, m_factors, m_HasVarCoeffDiff, m_jac, m_lambda, Nektar::Collections::Operator::m_stdExp, Nektar::Collections::Operator::m_wspSize, and Nektar::StdRegions::NullFactorMap.

Member Function Documentation

◆ CheckFactors()

void Nektar::Collections::Helmholtz_IterPerExp::CheckFactors	(	StdRegions::FactorMap	factors,
		int	coll_phys_offset
	)

inlineoverridevirtual

Check the validity of supplied constant factors.

Parameters

factors	Map of factors
coll_phys_offset	Unused

Implements Nektar::Collections::Operator.

Definition at line 421 of file library/Collections/Helmholtz.cpp.

    {
        // If match previous factors, nothing to do.
        if (m_factors == factors)
        {
            return;
        }
 
        m_factors = factors;
 
        // Check Lambda constant of Helmholtz operator
        auto x = factors.find(StdRegions::eFactorLambda);
        ASSERTL1(
            x != factors.end(),
            "Constant factor not defined: " +
                std::string(
                    StdRegions::ConstFactorTypeMap[StdRegions::eFactorLambda]));
        m_lambda = x->second;
 
        // If varcoeffs not supplied, nothing else to do.
        m_HasVarCoeffDiff = false;
        auto d            = factors.find(StdRegions::eFactorCoeffD00);
        if (d == factors.end())
        {
            return;
        }
 
        m_diff = Array<OneD, Array<OneD, NekDouble>>(m_coordim);
        for (int i = 0; i < m_coordim; ++i)
        {
            m_diff[i] = Array<OneD, NekDouble>(m_coordim, 0.0);
        }
 
        for (int i = 0; i < m_coordim; ++i)
        {
            d = factors.find(m_factorCoeffDef[i][i]);
            if (d != factors.end())
            {
                m_diff[i][i] = d->second;
            }
            else
            {
                m_diff[i][i] = 1.0;
            }
 
            for (int j = i + 1; j < m_coordim; ++j)
            {
                d = factors.find(m_factorCoeffDef[i][j]);
                if (d != factors.end())
                {
                    m_diff[i][j] = m_diff[j][i] = d->second;
                }
            }
        }
        m_HasVarCoeffDiff = true;
    }

References ASSERTL1, Nektar::StdRegions::ConstFactorTypeMap, Nektar::UnitTests::d(), Nektar::StdRegions::eFactorCoeffD00, Nektar::StdRegions::eFactorLambda, Nektar::VarcoeffHashingTest::factors, m_coordim, m_diff, m_factorCoeffDef, m_factors, m_HasVarCoeffDiff, and m_lambda.

Referenced by Helmholtz_IterPerExp().

◆ operator()() [1/2]

void Nektar::Collections::Helmholtz_IterPerExp::operator()	(	const Array< OneD, const NekDouble > &	input,
		Array< OneD, NekDouble > &	output0,
		Array< OneD, NekDouble > &	output1,
		Array< OneD, NekDouble > &	output2,
		Array< OneD, NekDouble > &	wsp
	)

inlinefinalvirtual

Perform operation.

Implements Nektar::Collections::Operator.

Definition at line 180 of file library/Collections/Helmholtz.cpp.

    {
        const int nCoeffs = m_stdExp->GetNcoeffs();
        const int nPhys   = m_stdExp->GetTotPoints();
 
        ASSERTL1(input.size() >= m_numElmt * nCoeffs,
                 "input array size is insufficient");
        ASSERTL1(output.size() >= m_numElmt * nCoeffs,
                 "output array size is insufficient");
 
        Array<OneD, NekDouble> tmpphys, t1;
        Array<OneD, Array<OneD, NekDouble>> dtmp(3);
        Array<OneD, Array<OneD, NekDouble>> tmp(3);
 
        tmpphys = wsp;
        for (int i = 1; i < m_coordim + 1; ++i)
        {
            dtmp[i - 1] = wsp + i * nPhys;
            tmp[i - 1]  = wsp + (i + m_coordim) * nPhys;
        }
 
        for (int i = 0; i < m_numElmt; ++i)
        {
            m_stdExp->BwdTrans(input + i * nCoeffs, tmpphys);
 
            // local derivative
            m_stdExp->PhysDeriv(tmpphys, dtmp[0], dtmp[1], dtmp[2]);
 
            // determine mass matrix term
            if (m_isDeformed)
            {
                Vmath::Vmul(nPhys, m_jac + i * nPhys, 1, tmpphys, 1, tmpphys,
                            1);
            }
            else
            {
                Vmath::Smul(nPhys, m_jac[i], tmpphys, 1, tmpphys, 1);
            }
 
            m_stdExp->IProductWRTBase(tmpphys, t1 = output + i * nCoeffs);
            Vmath::Smul(nCoeffs, m_lambda, output + i * nCoeffs, 1,
                        t1 = output + i * nCoeffs, 1);
 
            if (m_isDeformed)
            {
                // calculate full derivative
                for (int j = 0; j < m_coordim; ++j)
                {
                    Vmath::Vmul(nPhys,
                                m_derivFac[j * m_dim].origin() + i * nPhys, 1,
                                &dtmp[0][0], 1, &tmp[j][0], 1);
 
                    for (int k = 1; k < m_dim; ++k)
                    {
                        Vmath::Vvtvp(
                            nPhys,
                            m_derivFac[j * m_dim + k].origin() + i * nPhys, 1,
                            &dtmp[k][0], 1, &tmp[j][0], 1, &tmp[j][0], 1);
                    }
                }
 
                if (m_HasVarCoeffDiff)
                {
                    // calculate dtmp[i] = dx/dxi sum_j diff[0][j] tmp[j]
                    //                   + dy/dxi sum_j diff[1][j] tmp[j]
                    //                   + dz/dxi sum_j diff[2][j] tmp[j]
 
                    // First term
                    Vmath::Smul(nPhys, m_diff[0][0], &tmp[0][0], 1, &tmpphys[0],
                                1);
                    for (int l = 1; l < m_coordim; ++l)
                    {
                        Vmath::Svtvp(nPhys, m_diff[0][l], &tmp[l][0], 1,
                                     &tmpphys[0], 1, &tmpphys[0], 1);
                    }
 
                    for (int j = 0; j < m_dim; ++j)
                    {
                        Vmath::Vmul(nPhys, m_derivFac[j].origin() + i * nPhys,
                                    1, &tmpphys[0], 1, &dtmp[j][0], 1);
                    }
 
                    // Second and third terms
                    for (int k = 1; k < m_coordim; ++k)
                    {
                        Vmath::Smul(nPhys, m_diff[k][0], &tmp[0][0], 1,
                                    &tmpphys[0], 1);
                        for (int l = 1; l < m_coordim; ++l)
                        {
                            Vmath::Svtvp(nPhys, m_diff[k][l], &tmp[l][0], 1,
                                         &tmpphys[0], 1, &tmpphys[0], 1);
                        }
 
                        for (int j = 0; j < m_dim; ++j)
                        {
                            Vmath::Vvtvp(nPhys,
                                         m_derivFac[j + k * m_dim].origin() +
                                             i * nPhys,
                                         1, &tmpphys[0], 1, &dtmp[j][0], 1,
                                         &dtmp[j][0], 1);
                        }
                    }
                }
                else
                {
                    // calculate dx/dxi tmp[0] + dy/dxi tmp[1]
                    //                         + dz/dxi tmp[2]
                    for (int j = 0; j < m_dim; ++j)
                    {
                        Vmath::Vmul(nPhys, m_derivFac[j].origin() + i * nPhys,
                                    1, &tmp[0][0], 1, &dtmp[j][0], 1);
 
                        for (int k = 1; k < m_coordim; ++k)
                        {
                            Vmath::Vvtvp(nPhys,
                                         m_derivFac[j + k * m_dim].origin() +
                                             i * nPhys,
                                         1, &tmp[k][0], 1, &dtmp[j][0], 1,
                                         &dtmp[j][0], 1);
                        }
                    }
                }
 
                // calculate Iproduct WRT Std Deriv
                for (int j = 0; j < m_dim; ++j)
                {
                    // multiply by Jacobian
                    Vmath::Vmul(nPhys, m_jac + i * nPhys, 1, dtmp[j], 1,
                                dtmp[j], 1);
 
                    m_stdExp->IProductWRTDerivBase(j, dtmp[j], tmp[0]);
                    Vmath::Vadd(nCoeffs, tmp[0], 1, output + i * nCoeffs, 1,
                                t1 = output + i * nCoeffs, 1);
                }
            }
            else
            {
                // calculate full derivative
                for (int j = 0; j < m_coordim; ++j)
                {
                    Vmath::Smul(nPhys, m_derivFac[j * m_dim][i], &dtmp[0][0], 1,
                                &tmp[j][0], 1);
 
                    for (int k = 1; k < m_dim; ++k)
                    {
                        Vmath::Svtvp(nPhys, m_derivFac[j * m_dim + k][i],
                                     &dtmp[k][0], 1, &tmp[j][0], 1, &tmp[j][0],
                                     1);
                    }
                }
 
                if (m_HasVarCoeffDiff)
                {
                    // calculate dtmp[i] = dx/dxi sum_j diff[0][j] tmp[j]
                    //                   + dy/dxi sum_j diff[1][j] tmp[j]
                    //                   + dz/dxi sum_j diff[2][j] tmp[j]
 
                    // First term
                    Vmath::Smul(nPhys, m_diff[0][0], &tmp[0][0], 1, &tmpphys[0],
                                1);
                    for (int l = 1; l < m_coordim; ++l)
                    {
                        Vmath::Svtvp(nPhys, m_diff[0][l], &tmp[l][0], 1,
                                     &tmpphys[0], 1, &tmpphys[0], 1);
                    }
 
                    for (int j = 0; j < m_dim; ++j)
                    {
                        Vmath::Smul(nPhys, m_derivFac[j][i], &tmpphys[0], 1,
                                    &dtmp[j][0], 1);
                    }
 
                    // Second and third terms
                    for (int k = 1; k < m_coordim; ++k)
                    {
                        Vmath::Smul(nPhys, m_diff[k][0], &tmp[0][0], 1,
                                    &tmpphys[0], 1);
                        for (int l = 1; l < m_coordim; ++l)
                        {
                            Vmath::Svtvp(nPhys, m_diff[k][l], &tmp[l][0], 1,
                                         &tmpphys[0], 1, &tmpphys[0], 1);
                        }
 
                        for (int j = 0; j < m_dim; ++j)
                        {
                            Vmath::Svtvp(nPhys, m_derivFac[j + k * m_dim][i],
                                         &tmpphys[0], 1, &dtmp[j][0], 1,
                                         &dtmp[j][0], 1);
                        }
                    }
                }
                else
                {
                    // calculate dx/dxi tmp[0] + dy/dxi tmp[2]
                    //                         + dz/dxi tmp[3]
                    for (int j = 0; j < m_dim; ++j)
                    {
                        Vmath::Smul(nPhys, m_derivFac[j][i], &tmp[0][0], 1,
                                    &dtmp[j][0], 1);
 
                        for (int k = 1; k < m_coordim; ++k)
                        {
                            Vmath::Svtvp(nPhys, m_derivFac[j + k * m_dim][i],
                                         &tmp[k][0], 1, &dtmp[j][0], 1,
                                         &dtmp[j][0], 1);
                        }
                    }
                }
 
                // calculate Iproduct WRT Std Deriv
                for (int j = 0; j < m_dim; ++j)
                {
                    // multiply by Jacobian
                    Vmath::Smul(nPhys, m_jac[i], dtmp[j], 1, dtmp[j], 1);
 
                    m_stdExp->IProductWRTDerivBase(j, dtmp[j], tmp[0]);
                    Vmath::Vadd(nCoeffs, tmp[0], 1, output + i * nCoeffs, 1,
                                t1 = output + i * nCoeffs, 1);
                }
            }
        }
    }

References ASSERTL1, m_coordim, m_derivFac, m_diff, m_dim, m_HasVarCoeffDiff, Nektar::Collections::Operator::m_isDeformed, m_jac, m_lambda, Nektar::Collections::Operator::m_numElmt, Nektar::Collections::Operator::m_stdExp, Nektar::MovementTests::origin, Vmath::Smul(), Vmath::Svtvp(), Vmath::Vadd(), Vmath::Vmul(), and Vmath::Vvtvp().

◆ operator()() [2/2]

void Nektar::Collections::Helmholtz_IterPerExp::operator()	(	int	dir,
		const Array< OneD, const NekDouble > &	input,
		Array< OneD, NekDouble > &	output,
		Array< OneD, NekDouble > &	wsp
	)

inlinefinalvirtual

Implements Nektar::Collections::Operator.

Definition at line 407 of file library/Collections/Helmholtz.cpp.

    {
        NEKERROR(ErrorUtil::efatal, "Not valid for this operator.");
    }

References Nektar::ErrorUtil::efatal, and NEKERROR.

Member Data Documentation

◆ m_coordim

int Nektar::Collections::Helmholtz_IterPerExp::m_coordim

protected

Definition at line 483 of file library/Collections/Helmholtz.cpp.

Referenced by CheckFactors(), Helmholtz_IterPerExp(), and operator()().

◆ m_derivFac

Array<TwoD, const NekDouble> Nektar::Collections::Helmholtz_IterPerExp::m_derivFac

protected

Definition at line 480 of file library/Collections/Helmholtz.cpp.

Referenced by Helmholtz_IterPerExp(), and operator()().

◆ m_diff

Array<OneD, Array<OneD, NekDouble> > Nektar::Collections::Helmholtz_IterPerExp::m_diff

protected

Definition at line 487 of file library/Collections/Helmholtz.cpp.

Referenced by CheckFactors(), and operator()().

◆ m_dim

int Nektar::Collections::Helmholtz_IterPerExp::m_dim

protected

Definition at line 482 of file library/Collections/Helmholtz.cpp.

Referenced by Helmholtz_IterPerExp(), and operator()().

◆ m_factorCoeffDef

const StdRegions::ConstFactorType Nektar::Collections::Helmholtz_IterPerExp::m_factorCoeffDef[3][3]

protected

Initial value:

= {
        {StdRegions::eFactorCoeffD00, StdRegions::eFactorCoeffD01,
         StdRegions::eFactorCoeffD02},
        {StdRegions::eFactorCoeffD01, StdRegions::eFactorCoeffD11,
         StdRegions::eFactorCoeffD12},
        {StdRegions::eFactorCoeffD02, StdRegions::eFactorCoeffD12,
         StdRegions::eFactorCoeffD22}}