#include <PreconCfsBRJ.h>

Inheritance diagram for Nektar::PreconCfsBRJ:

Public Member Functions
	PreconCfsBRJ (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const LibUtilities::SessionReaderSharedPtr &pSession, const LibUtilities::CommSharedPtr &vComm)

	~PreconCfsBRJ ()

virtual bool	v_UpdatePreconMatCheck (const Array< OneD, const NekDouble > &res, const NekDouble dtLambda) override

Public Member Functions inherited from Nektar::PreconCfsOp
	PreconCfsOp (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const LibUtilities::SessionReaderSharedPtr &pSession, const LibUtilities::CommSharedPtr &vComm)

virtual	~PreconCfsOp ()

void	SetOperators (const NekPreconCfsOperators &in)

Public Member Functions inherited from Nektar::PreconCfs
	PreconCfs (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const LibUtilities::SessionReaderSharedPtr &pSession, const LibUtilities::CommSharedPtr &vComm)

virtual	~PreconCfs ()

void	DoPreconCfs (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const bool &flag)

void	BuildPreconCfs (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, const Array< OneD, NekDouble >> &intmp, const NekDouble time, const NekDouble lambda)

void	InitObject ()

bool	UpdatePreconMatCheck (const Array< OneD, const NekDouble > &res, const NekDouble dtLambda)

Static Public Member Functions
static PreconCfsOpSharedPtr	create (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const LibUtilities::SessionReaderSharedPtr &pSession, const LibUtilities::CommSharedPtr &vComm)
	Creates an instance of this class. More...

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

Protected Member Functions
virtual void	v_InitObject () override

Protected Member Functions inherited from Nektar::PreconCfs
void	DoNullPrecon (const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const bool &flag)

Protected Attributes
int	m_PreconItsStep

int	m_BRJRelaxParam

Array< OneD, Array< OneD, SNekBlkMatSharedPtr > >	m_PreconMatVarsSingle

unsigned int	m_max_nblocks

unsigned int	m_max_nElmtDof

std::vector< simd< NekSingle >, tinysimd::allocator< simd< NekSingle > > >	m_sBlkDiagMat

std::vector< int >	m_inputIdx

Array< OneD, SNekBlkMatSharedPtr >	m_TraceJacSingle

TensorOfArray4D< NekSingle >	m_TraceJacArraySingle

Array< OneD, SNekBlkMatSharedPtr >	m_TraceJacDerivSingle

TensorOfArray4D< NekSingle >	m_TraceJacDerivArraySingle

Array< OneD, Array< OneD, NekSingle > >	m_TraceJacDerivSignSingle

TensorOfArray5D< NekSingle >	m_TraceIPSymJacArraySingle

PrecType	m_PreconMatStorage

Protected Attributes inherited from Nektar::PreconCfsOp
NekPreconCfsOperators	m_operator

Protected Attributes inherited from Nektar::PreconCfs
LibUtilities::CommSharedPtr	m_Comm

bool	m_verbose

int	m_spacedim

int	m_PreconMatFreezNumb

int	m_PreconTimesCounter

NekDouble	m_DtLambdaPreconMat = -1.0

NekDouble	m_BndEvaluateTime

bool	m_CalcPreconMatFlag = false

Private Member Functions
virtual void	v_DoPreconCfs (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, NekDouble > &pInput, Array< OneD, NekDouble > &pOutput, const bool &flag) override

virtual void	v_BuildPreconCfs (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, const Array< OneD, NekDouble >> &intmp, const NekDouble time, const NekDouble lambda) override

void	PreconBlkDiag (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)

template<typename DataType >
void	MinusOffDiag2Rhs (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const size_t nvariables, const size_t nCoeffs, const Array< OneD, const Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray, bool flagUpdateDervFlux, Array< OneD, Array< OneD, NekDouble >> &FwdFluxDeriv, Array< OneD, Array< OneD, NekDouble >> &BwdFluxDeriv, TensorOfArray3D< NekDouble > &qfield, TensorOfArray3D< NekDouble > &wspTrace, Array< OneD, Array< OneD, DataType >> &wspTraceDataType, const TensorOfArray4D< DataType > &TraceJacArray, const TensorOfArray4D< DataType > &TraceJacDerivArray, const Array< OneD, const Array< OneD, DataType >> &TraceJacDerivSign, const TensorOfArray5D< DataType > &TraceIPSymJacArray)

template<typename TypeNekBlkMatSharedPtr >
void	AllocatePreconBlkDiagCoeff (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, Array< OneD, Array< OneD, TypeNekBlkMatSharedPtr >> &gmtxarray, const int &nscale=1)

void	AllocateSIMDPreconBlkMatDiag (const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
	This function creates the matrix structure for the block diagonal operator. It organizes the way that the matrix. More...

void	AllocateNekBlkMatDig (SNekBlkMatSharedPtr &mat, const Array< OneD, unsigned int > nrow, const Array< OneD, unsigned int > ncol)

Friends
class	MemoryManager< PreconCfsBRJ >

Detailed Description

Block Relaxed(weighted) Jacobi iterative (BRJ) Preconditioner for CFS

Solves a linear system using iterative methods.

Definition at line 49 of file PreconCfsBRJ.h.

Constructor & Destructor Documentation

◆ PreconCfsBRJ()

Nektar::PreconCfsBRJ::PreconCfsBRJ	(	const Array< OneD, MultiRegions::ExpListSharedPtr > &	pFields,
		const LibUtilities::SessionReaderSharedPtr &	pSession,
		const LibUtilities::CommSharedPtr &	vComm
	)

Definition at line 53 of file PreconCfsBRJ.cpp.

     : PreconCfsOp(pFields, pSession, vComm)
 {
     pSession->LoadParameter("PreconItsStep", m_PreconItsStep, 7);
     pSession->LoadParameter("BRJRelaxParam", m_BRJRelaxParam, 1.0);
  
     size_t nvariables  = pFields.size();
     m_PreconMatStorage = eDiagonal;
  
     m_PreconMatVarsSingle =
         Array<OneD, Array<OneD, SNekBlkMatSharedPtr>>(nvariables);
     for (size_t i = 0; i < nvariables; i++)
     {
         m_PreconMatVarsSingle[i] = Array<OneD, SNekBlkMatSharedPtr>(nvariables);
     }
     AllocatePreconBlkDiagCoeff(pFields, m_PreconMatVarsSingle);
  
     AllocateSIMDPreconBlkMatDiag(pFields);
 }

References AllocatePreconBlkDiagCoeff(), AllocateSIMDPreconBlkMatDiag(), Nektar::eDiagonal, m_BRJRelaxParam, m_PreconItsStep, m_PreconMatStorage, and m_PreconMatVarsSingle.

◆ ~PreconCfsBRJ()

Nektar::PreconCfsBRJ::~PreconCfsBRJ ( )

inline

Definition at line 71 of file PreconCfsBRJ.h.

71 {};

Member Function Documentation

◆ AllocateNekBlkMatDig()

void Nektar::PreconCfsBRJ::AllocateNekBlkMatDig	(	SNekBlkMatSharedPtr &	mat,
		const Array< OneD, unsigned int >	nrow,
		const Array< OneD, unsigned int >	ncol
	)

inlineprivate

Definition at line 241 of file PreconCfsBRJ.h.

     {
         mat =
             MemoryManager<SNekBlkMat>::AllocateSharedPtr(nrow, ncol, eDIAGONAL);
         SNekMatSharedPtr loc_matNvar;
         for (size_t nelm = 0; nelm < nrow.size(); ++nelm)
         {
             int nrowsVars = nrow[nelm];
             int ncolsVars = ncol[nelm];
  
             loc_matNvar = MemoryManager<SNekMat>::AllocateSharedPtr(
                 nrowsVars, ncolsVars, 0.0);
             mat->SetBlock(nelm, nelm, loc_matNvar);
         }
     }

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and Nektar::eDIAGONAL.

Referenced by AllocatePreconBlkDiagCoeff(), and v_BuildPreconCfs().

◆ AllocatePreconBlkDiagCoeff()

template<typename TypeNekBlkMatSharedPtr >

void Nektar::PreconCfsBRJ::AllocatePreconBlkDiagCoeff	(	const Array< OneD, MultiRegions::ExpListSharedPtr > &	pFields,
		Array< OneD, Array< OneD, TypeNekBlkMatSharedPtr >> &	gmtxarray,
		const int &	nscale = `1`
	)

private

Definition at line 502 of file PreconCfsBRJ.cpp.

 {
  
     size_t nvars  = pFields.size();
     size_t nelmts = pFields[0]->GetNumElmts();
     size_t nelmtcoef;
     Array<OneD, unsigned int> nelmtmatdim(nelmts);
     for (size_t i = 0; i < nelmts; i++)
     {
         nelmtcoef      = pFields[0]->GetExp(i)->GetNcoeffs();
         nelmtmatdim[i] = nelmtcoef * nscale;
     }
  
     for (size_t i = 0; i < nvars; i++)
     {
         for (size_t j = 0; j < nvars; j++)
         {
             AllocateNekBlkMatDig(gmtxarray[i][j], nelmtmatdim, nelmtmatdim);
         }
     }
 }

References AllocateNekBlkMatDig().

Referenced by PreconCfsBRJ().

◆ AllocateSIMDPreconBlkMatDiag()

void Nektar::PreconCfsBRJ::AllocateSIMDPreconBlkMatDiag ( const Array< OneD, MultiRegions::ExpListSharedPtr > & pFields )

inlineprivate

This function creates the matrix structure for the block diagonal operator. It organizes the way that the matrix.

by SIMD instructions. The degrees of freedom of each element are reorganized, so they are placed in the correct location to perfom SIMD instructions.

Definition at line 143 of file PreconCfsBRJ.h.

     {
         using vec_t = simd<NekSingle>;
  
         int TotMatLen         = 0;
         int TotLen            = 0;
         const auto nTotElmt   = pFields[0]->GetNumElmts();
         const auto nvariables = pFields.size();
         const auto vecwidth   = vec_t::width;
  
         m_max_nblocks  = 0;
         m_max_nElmtDof = 0;
  
         for (int ne = 0; ne < nTotElmt; ne++)
         {
             const auto nElmtDof = pFields[0]->GetNcoeffs(ne) * nvariables;
             const auto nblocks  = nElmtDof / vecwidth;
             unsigned int totblocks =
                 (nElmtDof % vecwidth) ? nblocks + 1 : nblocks;
  
             m_max_nblocks =
                 (m_max_nblocks > totblocks) ? m_max_nblocks : totblocks;
             m_max_nElmtDof =
                 (m_max_nElmtDof > nElmtDof) ? m_max_nElmtDof : nElmtDof;
  
             TotLen += totblocks * vecwidth;
             TotMatLen += nElmtDof * totblocks;
         }
  
         m_sBlkDiagMat.resize(TotMatLen);
         m_inputIdx.resize(TotLen);
  
         // generate a index list of vector width aware mapping from
         // local coeff storage over all variables to elemental storage
         // over variables
         unsigned int ncoeffs = pFields[0]->GetNcoeffs();
         for (int ne = 0, cnt1 = 0; ne < nTotElmt; ne++)
         {
             const auto nElmtCoeff  = pFields[0]->GetNcoeffs(ne);
             const auto nElmtDof    = nElmtCoeff * nvariables;
             const auto nblocks     = nElmtDof / vecwidth;
             const auto nCoefOffset = pFields[0]->GetCoeff_Offset(ne);
             int i                  = 0;
             int i0                 = 0;
             int inOffset, j;
  
             for (int m = 0; m < nvariables; m++)
             {
                 inOffset = m * ncoeffs + nCoefOffset;
  
                 if (m && (vecwidth - i0 < nElmtCoeff))
                 {
                     // May need to add entries from later variables to
                     // remainder of last variable if the vector width
                     // was not exact multiple of number of elemental
                     // coeffs
                     for (i = 0; i0 < vecwidth; ++i, ++i0)
                     {
                         m_inputIdx[cnt1++] = inOffset + i;
                     }
                 }
                 else
                 {
                     i = 0;
                 }
  
                 // load up other vectors in variable that fit into vector
                 // width
                 for (j = 0; (j + 1) * vecwidth < nElmtCoeff - i; ++j)
                 {
                     for (i0 = 0; i0 < vecwidth; ++i0)
                     {
                         m_inputIdx[cnt1++] = inOffset + i + j * vecwidth + i0;
                     }
                 }
  
                 // load up any residual data for this variable
                 for (i0 = 0, j = j * vecwidth; j < nElmtCoeff - i; ++j, ++i0)
                 {
                     m_inputIdx[cnt1++] = inOffset + i + j;
                 }
             }
  
             const auto endwidth = nElmtDof - nblocks * vecwidth;
  
             // fill out rest of index to match vector width with last entry
             if (endwidth)
             {
                 for (i0 = endwidth; i0 < vecwidth; ++i0)
                 {
                     m_inputIdx[cnt1++] = inOffset + i + j - 1;
                 }
             }
             ASSERTL1(cnt1 <= TotLen, "m_inputIdx over extended");
         }
     }

References ASSERTL1.

Referenced by PreconCfsBRJ().

◆ create()

static PreconCfsOpSharedPtr Nektar::PreconCfsBRJ::create	(	const Array< OneD, MultiRegions::ExpListSharedPtr > &	pFields,
		const LibUtilities::SessionReaderSharedPtr &	pSession,
		const LibUtilities::CommSharedPtr &	vComm
	)

inlinestatic

Creates an instance of this class.

Definition at line 55 of file PreconCfsBRJ.h.

     {
         PreconCfsOpSharedPtr p = MemoryManager<PreconCfsBRJ>::AllocateSharedPtr(
             pFields, pSession, vComm);
         return p;
     }

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and CellMLToNektar.cellml_metadata::p.

◆ MinusOffDiag2Rhs()

template<typename DataType >

void Nektar::PreconCfsBRJ::MinusOffDiag2Rhs	(	const Array< OneD, MultiRegions::ExpListSharedPtr > &	pFields,
		const size_t	nvariables,
		const size_t	nCoeffs,
		const Array< OneD, const Array< OneD, NekDouble >> &	inarray,
		Array< OneD, Array< OneD, NekDouble >> &	outarray,
		bool	flagUpdateDervFlux,
		Array< OneD, Array< OneD, NekDouble >> &	FwdFluxDeriv,
		Array< OneD, Array< OneD, NekDouble >> &	BwdFluxDeriv,
		TensorOfArray3D< NekDouble > &	qfield,
		TensorOfArray3D< NekDouble > &	wspTrace,
		Array< OneD, Array< OneD, DataType >> &	wspTraceDataType,
		const TensorOfArray4D< DataType > &	TraceJacArray,
		const TensorOfArray4D< DataType > &	TraceJacDerivArray,
		const Array< OneD, const Array< OneD, DataType >> &	TraceJacDerivSign,
		const TensorOfArray5D< DataType > &	TraceIPSymJacArray
	)

private

Definition at line 374 of file PreconCfsBRJ.cpp.

 {
     boost::ignore_unused(flagUpdateDervFlux, qfield, TraceJacDerivArray,
                          TraceJacDerivSign, FwdFluxDeriv, BwdFluxDeriv,
                          TraceIPSymJacArray);
  
     size_t nTracePts = pFields[0]->GetTrace()->GetNpoints();
     size_t npoints   = pFields[0]->GetNpoints();
     size_t nDim      = m_spacedim;
  
     Array<OneD, Array<OneD, NekDouble>> outpnts(nvariables);
     for (size_t i = 0; i < nvariables; i++)
     {
         outpnts[i] = Array<OneD, NekDouble>(npoints, 0.0);
         pFields[i]->BwdTrans(outarray[i], outpnts[i]);
     }
  
     // Store forwards/backwards space along trace space
     Array<OneD, Array<OneD, NekDouble>> Fwd;
     Array<OneD, Array<OneD, NekDouble>> Bwd;
     Array<OneD, Array<OneD, NekDouble>> FwdFlux;
     Array<OneD, Array<OneD, NekDouble>> BwdFlux;
     TensorOfArray3D<NekDouble> numDerivBwd(nDim);
     TensorOfArray3D<NekDouble> numDerivFwd(nDim);
     size_t indexwspTrace = 0;
     Fwd                  = wspTrace[indexwspTrace], indexwspTrace++;
     Bwd                  = wspTrace[indexwspTrace], indexwspTrace++;
     FwdFlux              = wspTrace[indexwspTrace], indexwspTrace++;
     BwdFlux              = wspTrace[indexwspTrace], indexwspTrace++;
  
     LibUtilities::Timer timer;
     for (size_t i = 0; i < nvariables; ++i)
     {
         timer.Start();
         pFields[i]->GetFwdBwdTracePhys(outpnts[i], Fwd[i], Bwd[i]);
         timer.Stop();
         timer.AccumulateRegion("ExpList::GetFwdBwdTracePhys", 10);
     }
  
     size_t indexwspTraceDataType = 0;
     Array<OneD, Array<OneD, DataType>> Fwdarray(nvariables);
     for (size_t m = 0; m < nvariables; ++m)
     {
         Fwdarray[m] = wspTraceDataType[indexwspTraceDataType],
         indexwspTraceDataType++;
     }
     Array<OneD, DataType> Fwdreslt;
     Fwdreslt = wspTraceDataType[indexwspTraceDataType], indexwspTraceDataType++;
  
     for (size_t m = 0; m < nvariables; ++m)
     {
         for (size_t i = 0; i < nTracePts; ++i)
         {
             Fwdarray[m][i] = DataType(Fwd[m][i]);
         }
     }
     for (size_t m = 0; m < nvariables; ++m)
     {
         Vmath::Zero(nTracePts, Fwdreslt, 1);
         for (size_t n = 0; n < nvariables; ++n)
         {
             for (size_t p = 0; p < nTracePts; ++p)
             {
                 Fwdreslt[p] += TraceJacArray[0][m][n][p] * Fwdarray[n][p];
             }
         }
         for (size_t i = 0; i < nTracePts; ++i)
         {
             FwdFlux[m][i] = NekDouble(Fwdreslt[i]);
         }
     }
  
     for (size_t m = 0; m < nvariables; ++m)
     {
         for (size_t i = 0; i < nTracePts; ++i)
         {
             Fwdarray[m][i] = DataType(Bwd[m][i]);
         }
     }
     for (size_t m = 0; m < nvariables; ++m)
     {
         Vmath::Zero(nTracePts, Fwdreslt, 1);
         for (size_t n = 0; n < nvariables; ++n)
         {
             for (size_t p = 0; p < nTracePts; ++p)
             {
                 Fwdreslt[p] += TraceJacArray[1][m][n][p] * Fwdarray[n][p];
             }
         }
         for (size_t i = 0; i < nTracePts; ++i)
         {
             BwdFlux[m][i] = NekDouble(Fwdreslt[i]);
         }
     }
  
     for (size_t i = 0; i < nvariables; ++i)
     {
         Vmath::Fill(nCoeffs, 0.0, outarray[i], 1);
         timer.Start();
         pFields[i]->AddTraceIntegralToOffDiag(FwdFlux[i], BwdFlux[i],
                                               outarray[i]);
         timer.Stop();
         timer.AccumulateRegion("ExpList::AddTraceIntegralToOffDiag", 10);
     }
  
     for (size_t i = 0; i < nvariables; ++i)
     {
         for (size_t p = 0; p < nCoeffs; ++p)
         {
             outarray[i][p] =
                 -m_DtLambdaPreconMat * outarray[i][p] + inarray[i][p];
         }
     }
 }

References Nektar::LibUtilities::Timer::AccumulateRegion(), Vmath::Fill(), Nektar::PreconCfs::m_DtLambdaPreconMat, Nektar::PreconCfs::m_spacedim, CellMLToNektar.cellml_metadata::p, Nektar::LibUtilities::Timer::Start(), Nektar::LibUtilities::Timer::Stop(), and Vmath::Zero().

Referenced by v_DoPreconCfs().

◆ PreconBlkDiag()

void Nektar::PreconCfsBRJ::PreconBlkDiag	(	const Array< OneD, MultiRegions::ExpListSharedPtr > &	pFields,
		const Array< OneD, NekDouble > &	inarray,
		Array< OneD, NekDouble > &	outarray
	)

private

Definition at line 297 of file PreconCfsBRJ.cpp.

 {
     unsigned int nvariables = pFields.size();
  
     int nTotElmt = pFields[0]->GetNumElmts();
  
     using vec_t         = simd<NekSingle>;
     const auto vecwidth = vec_t::width;
  
     // vectorized matrix multiply
     std::vector<vec_t, tinysimd::allocator<vec_t>> Sinarray(m_max_nblocks);
     std::vector<vec_t, tinysimd::allocator<vec_t>> Soutarray(m_max_nElmtDof);
     // std::vector<vec_t, tinysimd::allocator<vec_t>> tmp;
  
     alignas(vec_t::alignment) std::array<NekSingle, vec_t::width> tmp;
  
     for (int ne = 0, cnt = 0, icnt = 0, icnt1 = 0; ne < nTotElmt; ne++)
     {
         const auto nElmtCoef = pFields[0]->GetNcoeffs(ne);
         const auto nElmtDof  = nElmtCoef * nvariables;
         const auto nblocks   = (nElmtDof % vecwidth) ? nElmtDof / vecwidth + 1
                                                      : nElmtDof / vecwidth;
  
         // gather data into blocks - could probably be done with a
         // gather call? can be replaced with a gather op when working
         for (int j = 0; j < nblocks; ++j, icnt += vecwidth)
         {
             for (int i = 0; i < vecwidth; ++i)
             {
                 tmp[i] = inarray[m_inputIdx[icnt + i]];
             }
  
             Sinarray[j].load(tmp.data());
         }
  
         // Do matrix multiply
         // first block just needs multiplying
         vec_t in = Sinarray[0];
         for (int i = 0; i < nElmtDof; ++i)
         {
             Soutarray[i] = m_sBlkDiagMat[cnt++] * in;
         }
  
         // rest of blocks are multiply add operations;
         for (int n = 1; n < nblocks; ++n)
         {
             in = Sinarray[n];
             for (int i = 0; i < nElmtDof; ++i)
             {
                 Soutarray[i].fma(m_sBlkDiagMat[cnt++], in);
             }
         }
  
         // get block aligned index for this expansion
         NekSingle val;
         for (int i = 0; i < nElmtDof; ++i)
         {
             // Get hold of data
             Soutarray[i].store(tmp.data());
  
             // Sum vector width
             val = tmp[0];
             for (int j = 1; j < vecwidth; ++j)
             {
                 val += tmp[j];
             }
             // put data into outarray
             outarray[m_inputIdx[icnt1 + i]] = NekDouble(val);
         }
  
         icnt1 += nblocks * vecwidth;
     }
 }

References m_inputIdx, m_max_nblocks, m_max_nElmtDof, and m_sBlkDiagMat.

Referenced by v_DoPreconCfs().

◆ v_BuildPreconCfs()

void Nektar::PreconCfsBRJ::v_BuildPreconCfs	(	const Array< OneD, MultiRegions::ExpListSharedPtr > &	pFields,
		const Array< OneD, const Array< OneD, NekDouble >> &	intmp,
		const NekDouble	time,
		const NekDouble	lambda
	)

overrideprivatevirtual

Copy array into column major blocks of vector width

Reimplemented from Nektar::PreconCfsOp.

Definition at line 188 of file PreconCfsBRJ.cpp.

 {
     boost::ignore_unused(pFields);
  
     if (0 < m_PreconItsStep)
     {
         SNekBlkMatSharedPtr PreconMatSingle;
         using vec_t    = simd<NekSingle>;
         int nvariables = pFields.size();
         int nelmts     = pFields[0]->GetNumElmts();
         Array<OneD, unsigned int> matdim(nelmts);
         for (int i = 0; i < nelmts; i++)
         {
             matdim[i] = pFields[0]->GetExp(i)->GetNcoeffs() * nvariables;
         }
         AllocateNekBlkMatDig(PreconMatSingle, matdim, matdim);
  
         m_operator.DoCalcPreconMatBRJCoeff(
             intmp, m_PreconMatVarsSingle, PreconMatSingle, m_TraceJacSingle,
             m_TraceJacDerivSingle, m_TraceJacDerivSignSingle,
             m_TraceJacArraySingle, m_TraceJacDerivArraySingle,
             m_TraceIPSymJacArraySingle);
  
         if (m_verbose && m_Comm->GetRank() == 0)
         {
             cout << "     ## CalcuPreconMat " << endl;
         }
  
         // copy matrix to simd layout
         // load matrix
         int cnt             = 0;
         const auto vecwidth = vec_t::width;
  
         alignas(vec_t::alignment) std::array<NekSingle, vec_t::width> tmp;
  
         for (int ne = 0; ne < nelmts; ne++)
         {
             const auto nElmtDof = matdim[ne];
             const auto nblocks  = nElmtDof / vecwidth;
  
             const NekSingle *mmat =
                 PreconMatSingle->GetBlockPtr(ne, ne)->GetRawPtr();
             /// Copy array into column major blocks of vector width
             for (int i1 = 0; i1 < nblocks; ++i1)
             {
                 for (int j = 0; j < nElmtDof; ++j)
                 {
                     for (int i = 0; i < vecwidth; ++i)
                     {
                         tmp[i] = mmat[j + (i1 * vecwidth + i) * nElmtDof];
                     }
                     // store contiguous vec_t array.
                     m_sBlkDiagMat[cnt++].load(tmp.data());
                 }
             }
  
             const auto endwidth = nElmtDof - nblocks * vecwidth;
  
             // end rows that do not fit into vector widths
             if (endwidth)
             {
                 for (int j = 0; j < nElmtDof; ++j)
                 {
                     for (int i = 0; i < endwidth; ++i)
                     {
                         tmp[i] = mmat[j + (nblocks * vecwidth + i) * nElmtDof];
                     }
  
                     for (int i = endwidth; i < vecwidth; ++i)
                     {
                         tmp[i] = 0.0;
                     }
                     m_sBlkDiagMat[cnt++].load(tmp.data());
                 }
             }
         }
     }
  
     m_BndEvaluateTime   = time;
     m_DtLambdaPreconMat = lambda;
  
     m_CalcPreconMatFlag  = false;
     m_PreconTimesCounter = 1;
 }

References AllocateNekBlkMatDig(), Nektar::NekPreconCfsOperators::DoCalcPreconMatBRJCoeff(), Nektar::PreconCfs::m_BndEvaluateTime, Nektar::PreconCfs::m_CalcPreconMatFlag, Nektar::PreconCfs::m_Comm, Nektar::PreconCfs::m_DtLambdaPreconMat, Nektar::PreconCfsOp::m_operator, m_PreconItsStep, m_PreconMatVarsSingle, Nektar::PreconCfs::m_PreconTimesCounter, m_sBlkDiagMat, m_TraceIPSymJacArraySingle, m_TraceJacArraySingle, m_TraceJacDerivArraySingle, m_TraceJacDerivSignSingle, m_TraceJacDerivSingle, m_TraceJacSingle, and Nektar::PreconCfs::m_verbose.

◆ v_DoPreconCfs()

void Nektar::PreconCfsBRJ::v_DoPreconCfs	(	const Array< OneD, MultiRegions::ExpListSharedPtr > &	pFields,
		const Array< OneD, NekDouble > &	pInput,
		Array< OneD, NekDouble > &	pOutput,
		const bool &	flag
	)

overrideprivatevirtual

Reimplemented from Nektar::PreconCfsOp.

Definition at line 81 of file PreconCfsBRJ.cpp.

 {
     boost::ignore_unused(flag);
  
     size_t nBRJIterTot = m_PreconItsStep;
     if (0 == nBRJIterTot)
     {
         DoNullPrecon(inarray, outarray, flag);
     }
     else
     {
         const NekDouble BRJParam   = m_BRJRelaxParam;
         const NekDouble OmBRJParam = 1.0 - BRJParam;
  
         size_t nvariables = pFields.size();
         size_t npoints    = pFields[0]->GetNcoeffs();
         size_t ntotpnt    = inarray.size();
  
         ASSERTL0(nvariables * npoints == ntotpnt,
                  "nvariables*npoints!=ntotpnt in PreconCoeff");
  
         Array<OneD, NekDouble> rhs(ntotpnt);
  
         Array<OneD, NekDouble> outN(ntotpnt);
         Array<OneD, NekDouble> outTmp(ntotpnt);
         Array<OneD, Array<OneD, NekDouble>> rhs2d(nvariables);
         Array<OneD, Array<OneD, NekDouble>> out_2d(nvariables);
         Array<OneD, Array<OneD, NekDouble>> outTmp_2d(nvariables);
         for (size_t m = 0; m < nvariables; m++)
         {
             size_t moffset = m * npoints;
             rhs2d[m]       = rhs + moffset;
             out_2d[m]      = outarray + moffset;
             outTmp_2d[m]   = outTmp + moffset;
             pFields[m]->MultiplyByMassMatrix(inarray + moffset, rhs2d[m]);
         }
  
         size_t nphysic   = pFields[0]->GetNpoints();
         size_t nTracePts = pFields[0]->GetTrace()->GetNpoints();
         TensorOfArray3D<NekDouble> qfield(m_spacedim);
         for (int i = 0; i < m_spacedim; i++)
         {
             qfield[i] = Array<OneD, Array<OneD, NekDouble>>(nvariables);
             for (size_t j = 0; j < nvariables; j++)
             {
                 qfield[i][j] = Array<OneD, NekDouble>(nphysic);
             }
         }
         size_t ntmpTrace = 4 + 2 * m_spacedim;
         TensorOfArray3D<NekDouble> tmpTrace(ntmpTrace);
         for (size_t i = 0; i < ntmpTrace; i++)
         {
             tmpTrace[i] = Array<OneD, Array<OneD, NekDouble>>(nvariables);
             for (size_t j = 0; j < nvariables; j++)
             {
                 tmpTrace[i][j] = Array<OneD, NekDouble>(nTracePts);
             }
         }
         Array<OneD, Array<OneD, NekDouble>> FwdFluxDeriv(nvariables);
         Array<OneD, Array<OneD, NekDouble>> BwdFluxDeriv(nvariables);
         for (size_t j = 0; j < nvariables; j++)
         {
             FwdFluxDeriv[j] = Array<OneD, NekDouble>(nTracePts);
             BwdFluxDeriv[j] = Array<OneD, NekDouble>(nTracePts);
         }
  
         bool flagUpdateDervFlux = false;
  
         const size_t nwspTraceDataType = nvariables + 1;
         Array<OneD, Array<OneD, NekSingle>> wspTraceDataType(nwspTraceDataType);
         for (size_t m = 0; m < nwspTraceDataType; m++)
         {
             wspTraceDataType[m] = Array<OneD, NekSingle>(nTracePts);
         }
  
         LibUtilities::Timer timer;
         timer.Start();
         PreconBlkDiag(pFields, rhs, outarray);
         timer.Stop();
         timer.AccumulateRegion("PreconCfsBRJ::PreconBlkDiag", 2);
  
         for (size_t nrelax = 0; nrelax < nBRJIterTot - 1; nrelax++)
         {
             Vmath::Smul(ntotpnt, OmBRJParam, outarray, 1, outN, 1);
  
             timer.Start();
             MinusOffDiag2Rhs(
                 pFields, nvariables, npoints, rhs2d, out_2d, flagUpdateDervFlux,
                 FwdFluxDeriv, BwdFluxDeriv, qfield, tmpTrace, wspTraceDataType,
                 m_TraceJacArraySingle, m_TraceJacDerivArraySingle,
                 m_TraceJacDerivSignSingle, m_TraceIPSymJacArraySingle);
             timer.Stop();
             timer.AccumulateRegion("PreconCfsBRJ::MinusOffDiag2Rhs", 2);
  
             timer.Start();
             PreconBlkDiag(pFields, outarray, outTmp);
             timer.Stop();
             timer.AccumulateRegion("PreconCfsBRJ::PreconBlkDiag", 2);
  
             Vmath::Svtvp(ntotpnt, BRJParam, outTmp, 1, outN, 1, outarray, 1);
         }
     }
 }

References Nektar::LibUtilities::Timer::AccumulateRegion(), ASSERTL0, Nektar::PreconCfs::DoNullPrecon(), m_BRJRelaxParam, m_PreconItsStep, Nektar::PreconCfs::m_spacedim, m_TraceIPSymJacArraySingle, m_TraceJacArraySingle, m_TraceJacDerivArraySingle, m_TraceJacDerivSignSingle, MinusOffDiag2Rhs(), PreconBlkDiag(), Vmath::Smul(), Nektar::LibUtilities::Timer::Start(), Nektar::LibUtilities::Timer::Stop(), and Vmath::Svtvp().

◆ v_InitObject()

void Nektar::PreconCfsBRJ::v_InitObject ( )

overrideprotectedvirtual

Reimplemented from Nektar::PreconCfsOp.

Definition at line 76 of file PreconCfsBRJ.cpp.

 {
     PreconCfsOp::v_InitObject();
 }

References Nektar::PreconCfsOp::v_InitObject().

◆ v_UpdatePreconMatCheck()

bool Nektar::PreconCfsBRJ::v_UpdatePreconMatCheck	(	const Array< OneD, const NekDouble > &	res,
		const NekDouble	dtLambda
	)

overridevirtual

Reimplemented from Nektar::PreconCfs.

Definition at line 276 of file PreconCfsBRJ.cpp.

 {
     boost::ignore_unused(res);
  
     bool flag = false;
  
     if (m_CalcPreconMatFlag || (m_DtLambdaPreconMat != dtLambda))
     {
         flag = true;
     }
  
     if (m_PreconMatFreezNumb < m_PreconTimesCounter)
     {
         flag = true;
     }
  
     m_CalcPreconMatFlag = flag;
     return flag;
 }

References Nektar::PreconCfs::m_CalcPreconMatFlag, Nektar::PreconCfs::m_DtLambdaPreconMat, Nektar::PreconCfs::m_PreconMatFreezNumb, and Nektar::PreconCfs::m_PreconTimesCounter.

Friends And Related Function Documentation

◆ MemoryManager< PreconCfsBRJ >

friend class MemoryManager< PreconCfsBRJ >

friend

Definition at line 1 of file PreconCfsBRJ.h.

Member Data Documentation

◆ className

std::string Nektar::PreconCfsBRJ::className

static

Initial value:

=
    GetPreconCfsOpFactory().RegisterCreatorFunction(
        "PreconCfsBRJ", PreconCfsBRJ::create,
        "Block Relaxed Jacobi Preconditioner for CFS.")

Name of the class.

Definition at line 66 of file PreconCfsBRJ.h.

◆ m_BRJRelaxParam

int Nektar::PreconCfsBRJ::m_BRJRelaxParam

protected

Definition at line 78 of file PreconCfsBRJ.h.

Referenced by PreconCfsBRJ(), and v_DoPreconCfs().

◆ m_inputIdx

std::vector<int> Nektar::PreconCfsBRJ::m_inputIdx

protected

Definition at line 86 of file PreconCfsBRJ.h.

Referenced by PreconBlkDiag().

◆ m_max_nblocks

unsigned int Nektar::PreconCfsBRJ::m_max_nblocks

protected

Definition at line 82 of file PreconCfsBRJ.h.

Referenced by PreconBlkDiag().

◆ m_max_nElmtDof

unsigned int Nektar::PreconCfsBRJ::m_max_nElmtDof

protected

Definition at line 83 of file PreconCfsBRJ.h.

Referenced by PreconBlkDiag().

◆ m_PreconItsStep

int Nektar::PreconCfsBRJ::m_PreconItsStep

protected

Definition at line 77 of file PreconCfsBRJ.h.

Referenced by PreconCfsBRJ(), v_BuildPreconCfs(), and v_DoPreconCfs().

◆ m_PreconMatStorage

PrecType Nektar::PreconCfsBRJ::m_PreconMatStorage

protected

Definition at line 95 of file PreconCfsBRJ.h.

Referenced by PreconCfsBRJ().

◆ m_PreconMatVarsSingle

Array<OneD, Array<OneD, SNekBlkMatSharedPtr> > Nektar::PreconCfsBRJ::m_PreconMatVarsSingle

protected

Definition at line 80 of file PreconCfsBRJ.h.

Referenced by PreconCfsBRJ(), and v_BuildPreconCfs().

◆ m_sBlkDiagMat

std::vector<simd<NekSingle>, tinysimd::allocator<simd<NekSingle> > > Nektar::PreconCfsBRJ::m_sBlkDiagMat

protected

Definition at line 85 of file PreconCfsBRJ.h.

Referenced by PreconBlkDiag(), and v_BuildPreconCfs().

◆ m_TraceIPSymJacArraySingle

TensorOfArray5D<NekSingle> Nektar::PreconCfsBRJ::m_TraceIPSymJacArraySingle

protected

Definition at line 93 of file PreconCfsBRJ.h.

Referenced by v_BuildPreconCfs(), and v_DoPreconCfs().

◆ m_TraceJacArraySingle

TensorOfArray4D<NekSingle> Nektar::PreconCfsBRJ::m_TraceJacArraySingle

protected

Definition at line 89 of file PreconCfsBRJ.h.

Referenced by v_BuildPreconCfs(), and v_DoPreconCfs().

◆ m_TraceJacDerivArraySingle

TensorOfArray4D<NekSingle> Nektar::PreconCfsBRJ::m_TraceJacDerivArraySingle

protected

Definition at line 91 of file PreconCfsBRJ.h.

Referenced by v_BuildPreconCfs(), and v_DoPreconCfs().

◆ m_TraceJacDerivSignSingle

Array<OneD, Array<OneD, NekSingle> > Nektar::PreconCfsBRJ::m_TraceJacDerivSignSingle

protected

Definition at line 92 of file PreconCfsBRJ.h.

Referenced by v_BuildPreconCfs(), and v_DoPreconCfs().

◆ m_TraceJacDerivSingle

Array<OneD, SNekBlkMatSharedPtr> Nektar::PreconCfsBRJ::m_TraceJacDerivSingle

protected

Definition at line 90 of file PreconCfsBRJ.h.

Referenced by v_BuildPreconCfs().

◆ m_TraceJacSingle

Array<OneD, SNekBlkMatSharedPtr> Nektar::PreconCfsBRJ::m_TraceJacSingle

protected

Definition at line 88 of file PreconCfsBRJ.h.

Referenced by v_BuildPreconCfs().

Public Member Functions

Static Public Member Functions

Static Public Attributes

Protected Member Functions

Protected Attributes

Private Member Functions

Friends

Detailed Description

Constructor & Destructor Documentation

◆ PreconCfsBRJ()

◆ ~PreconCfsBRJ()

Member Function Documentation

◆ AllocateNekBlkMatDig()

◆ AllocatePreconBlkDiagCoeff()

◆ AllocateSIMDPreconBlkMatDiag()

◆ create()

◆ MinusOffDiag2Rhs()

◆ PreconBlkDiag()

◆ v_BuildPreconCfs()

◆ v_DoPreconCfs()

◆ v_InitObject()

◆ v_UpdatePreconMatCheck()

Friends And Related Function Documentation

◆ MemoryManager< PreconCfsBRJ >

Member Data Documentation

◆ className

◆ m_BRJRelaxParam

◆ m_inputIdx

◆ m_max_nblocks

◆ m_max_nElmtDof

◆ m_PreconItsStep

◆ m_PreconMatStorage

◆ m_PreconMatVarsSingle

◆ m_sBlkDiagMat

◆ m_TraceIPSymJacArraySingle

◆ m_TraceJacArraySingle

◆ m_TraceJacDerivArraySingle

◆ m_TraceJacDerivSignSingle

◆ m_TraceJacDerivSingle

◆ m_TraceJacSingle