Chombo + EB + MF  3.2
Public Types | Public Member Functions | Static Public Member Functions | List of all members
NWOViscousTensorOp Class Reference

#include <NWOViscousTensorOp.H>

Inheritance diagram for NWOViscousTensorOp:
Inheritance graph
[legend]

Public Types

enum  prolongationType { piecewiseConstant = 0, linearInterp, NUM_INTERP_TYPE }
 

Public Member Functions

virtual void diagonalScale (LevelData< FArrayBox > &a_rhs, bool a_kappaWeighted)
 
virtual void setAlphaAndBeta (const Real &a_alpha, const Real &a_beta)
 
virtual void diagonalScale (LevelData< FArrayBox > &a_rhs)
 
virtual void divideByIdentityCoef (LevelData< FArrayBox > &a_rhs)
 
virtual ~NWOViscousTensorOp ()
 
 NWOViscousTensorOp (const DisjointBoxLayout &a_grids, const DisjointBoxLayout &a_gridsFine, const DisjointBoxLayout &a_gridsCoar, const RefCountedPtr< LevelData< FluxBox > > &a_eta, const RefCountedPtr< LevelData< FluxBox > > &a_lambda, const RefCountedPtr< LevelData< FArrayBox > > &a_acoef, Real a_alpha, Real a_beta, int a_refToFine, int a_refToCoar, const ProblemDomain &a_domain, const Real &a_dxLevel, const Real &a_dxCoar, BCFunc a_bc, Real a_safety=VTOP_DEFAULT_SAFETY, Real a_relaxTolerance=0.0, int a_relaxMinIter=1000)
 
 NWOViscousTensorOp (const DisjointBoxLayout &a_grids, const DisjointBoxLayout &a_gridsFine, const DisjointBoxLayout &a_gridsCoar, const RefCountedPtr< LevelData< FluxBox > > &a_eta, const RefCountedPtr< LevelData< FluxBox > > &a_lambda, const RefCountedPtr< LevelData< FArrayBox > > &a_acoef, Real a_alpha, Real a_beta, int a_refToFine, int a_refToCoar, const ProblemDomain &a_domain, const Real &a_dxLevel, const Real &a_dxCoar, BCHolder &a_bc, Real a_safety=VTOP_DEFAULT_SAFETY, Real a_relaxTolerance=0.0, int a_relaxMinIter=1000)
 
void fillGrad (const LevelData< FArrayBox > &a_phiFine)
 
void cellGrad (FArrayBox &a_gradPhi, const FArrayBox &a_phi, const Box &a_grid)
 
void define (const DisjointBoxLayout &a_grids, const DisjointBoxLayout &a_gridsFine, const DisjointBoxLayout &a_gridsCoar, const RefCountedPtr< LevelData< FluxBox > > &a_eta, const RefCountedPtr< LevelData< FluxBox > > &a_lambda, const RefCountedPtr< LevelData< FArrayBox > > &a_acoef, Real a_alpha, Real a_beta, int a_refToFine, int a_refToCoar, const ProblemDomain &a_domain, const Real &a_dxLevel, const Real &a_dxCoar, BCHolder &a_bc, Real a_safety=VTOP_DEFAULT_SAFETY, Real a_relaxTolerance=0.0, int a_relaxMinIter=1000)
 
virtual void residual (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs, bool a_homogeneous=false)
 
virtual void preCond (LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_residual)
 
virtual void applyOp (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi, bool a_homogeneous=false)
 
virtual void create (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs)
 
virtual void createCoarsened (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs, const int &a_refRat)
 
void reflux (const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, LevelData< FArrayBox > &residual, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
 
virtual void getFlux (FluxBox &a_flux, const LevelData< FArrayBox > &a_data, const Box &a_grid, const DataIndex &a_dit, Real a_scale)
 
void applyOpNoBoundary (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi)
 
void getFlux (FArrayBox &a_flux, const FArrayBox &a_data, const FArrayBox &a_etaFace, const FArrayBox &a_lambdaFace, const Box &a_facebox, int a_dir, int ref=1)
 
void computeOperatorNoBCs (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi)
 utility function which computes operator after all bc's have been set More...
 
virtual void assign (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs)
 
virtual Real dotProduct (const LevelData< FArrayBox > &a_1, const LevelData< FArrayBox > &a_2)
 
virtual void incr (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_x, Real a_scale)
 
virtual void axby (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_x, const LevelData< FArrayBox > &a_y, Real a, Real b)
 
virtual void scale (LevelData< FArrayBox > &a_lhs, const Real &a_scale)
 
virtual Real norm (const LevelData< FArrayBox > &a_x, int a_ord)
 
virtual Real dx () const
 
virtual Real dxCrse () const
 
virtual void setToZero (LevelData< FArrayBox > &a_x)
 
MGLevelOp functions
virtual void relax (LevelData< FArrayBox > &a_e, const LevelData< FArrayBox > &a_residual, int iterations)
 
virtual void createCoarser (LevelData< FArrayBox > &a_coarse, const LevelData< FArrayBox > &a_fine, bool ghosted)
 
virtual void restrictResidual (LevelData< FArrayBox > &a_resCoarse, LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_rhsFine)
 
virtual void prolongIncrement (LevelData< FArrayBox > &a_phiThisLevel, const LevelData< FArrayBox > &a_correctCoarse)
 
- Public Member Functions inherited from LevelTGAHelmOp< LevelData< FArrayBox >, FluxBox >
 LevelTGAHelmOp ()
 
 LevelTGAHelmOp (bool a_isTimeIndependent)
 
virtual ~LevelTGAHelmOp ()
 Destructor. More...
 
- Public Member Functions inherited from TGAHelmOp< LevelData< FArrayBox > >
 TGAHelmOp ()
 
 TGAHelmOp (bool a_isTimeDependent)
 
virtual ~TGAHelmOp ()
 Destructor. More...
 
virtual void kappaScale (LevelData< FArrayBox > &a_rhs)
 for eb only. kappa weight the rhs but do not multiply by the identity coefficient More...
 
virtual void setTime (Real a_time)
 
virtual void setTime (Real a_oldTime, Real a_mu, Real a_dt)
 
bool isTimeDependent () const
 Returns true if the operator is time-dependent, false otherwise. More...
 
virtual LevelData< FArrayBox > & identityCoef ()
 Allows access to the identity coefficient data for the operator. More...
 
- Public Member Functions inherited from AMRLevelOp< LevelData< FArrayBox > >
virtual void dumpAMR (Vector< LevelData< FArrayBox > * > &a_data, string name)
 
virtual void dumpLevel (LevelData< FArrayBox > &a_data, string name)
 
 AMRLevelOp ()
 Constructor. More...
 
virtual void dumpStuff (Vector< LevelData< FArrayBox > * > data, string filename)
 
virtual ~AMRLevelOp ()
 Destructor. More...
 
virtual void buildCopier (Copier &a_copier, const LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs)
 
virtual void assignCopier (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs, const Copier &a_copier)
 
virtual void zeroCovered (LevelData< FArrayBox > &a_lhs, LevelData< FArrayBox > &a_rhs, const Copier &a_copier)
 
virtual Real localMaxNorm (const LevelData< FArrayBox > &a_phi)
 
virtual void AMRProlongS (LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection, LevelData< FArrayBox > &a_temp, const Copier &a_copier)
 
virtual void AMRProlongS_2 (LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection, LevelData< FArrayBox > &a_temp, const Copier &a_copier, const Copier &a_cornerCopier, const AMRLevelOp< LevelData< FArrayBox > > *a_crsOp)
 
virtual void AMRRestrictS (LevelData< FArrayBox > &a_resCoarse, const LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection, LevelData< FArrayBox > &scratch, bool a_skip_res=false)
 
virtual unsigned int orderOfAccuracy (void) const
 
virtual void enforceCFConsistency (LevelData< FArrayBox > &a_coarseCorrection, const LevelData< FArrayBox > &a_correction)
 This routine is for operators with orderOfAccuracy()>2. More...
 
- Public Member Functions inherited from MGLevelOp< LevelData< FArrayBox > >
 MGLevelOp ()
 Constructor. More...
 
virtual ~MGLevelOp ()
 Destructor. More...
 
virtual void relaxNF (LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > *a_phiCoarse, const LevelData< FArrayBox > &a_rhs, int a_iterations)
 specialized no-fine relax function, useful for full-multigrid schemes, defaults to regular relax More...
 
virtual void restrictResidual (LevelData< FArrayBox > &a_resCoarse, LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > *a_phiCoarse, const LevelData< FArrayBox > &a_rhsFine, bool homogeneous)
 full-multigrid version of restrictResidual, useful for FAS-type schemes. Defaults to standard restriction More...
 
virtual void restrictR (LevelData< FArrayBox > &a_phiCoarse, const LevelData< FArrayBox > &a_phiFine)
 simple restriction function More...
 
void addObserver (MGLevelOpObserver< LevelData< FArrayBox > > *a_observer)
 
virtual void applyOpMg (LevelData< FArrayBox > &a_lhs, LevelData< FArrayBox > &a_phi, LevelData< FArrayBox > *a_phiCoarse, bool a_homogeneous)
 Apply an operator. More...
 
virtual void residualNF (LevelData< FArrayBox > &a_lhs, LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > *a_phiCoarse, const LevelData< FArrayBox > &a_rhs, bool a_homogeneous=false)
 
void removeObserver (MGLevelOpObserver< LevelData< FArrayBox > > *a_observer)
 
void addCoarserObserver (MGLevelOp< LevelData< FArrayBox > > *a_operator, int a_coarseningFactor)
 
void notifyObserversOfChange ()
 This should be called whenever the operator's data is updated. More...
 
virtual void finerOperatorChanged (const MGLevelOp< LevelData< FArrayBox > > &a_operator, int a_coarseningFactor)
 
int numObservers () const
 Returns the number of objects observing this operator. More...
 
- Public Member Functions inherited from LinearOp< LevelData< FArrayBox > >
virtual ~LinearOp ()
 
virtual void clear (LevelData< FArrayBox > &a_lhs)
 
virtual void assignLocal (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_rhs)
 
virtual void mDotProduct (const LevelData< FArrayBox > &a_1, const int a_sz, const LevelData< FArrayBox > a_2[], Real a_mdots[])
 
virtual void write (const LevelData< FArrayBox > *a, const char *filename)
 
- Public Member Functions inherited from MGLevelOpObserver< LevelData< FArrayBox > >
 MGLevelOpObserver ()
 Base level Constructor. Called by all subclasses. More...
 
virtual ~MGLevelOpObserver ()
 Destructor. More...
 
virtual void operatorChanged (const MGLevelOp< LevelData< FArrayBox > > &a_operator)
 
void setObservee (MGLevelOp< LevelData< FArrayBox > > *a_observee)
 
void clearObservee ()
 

Static Public Member Functions

static void loHiCenterFace (Box &a_loBox, int &a_hasLo, Box &a_hiBox, int &a_hasHi, Box &a_centerBox, const ProblemDomain &a_eblg, const Box &a_inBox, const int &a_dir)
 
static void getFaceDivAndGrad (FArrayBox &a_faceDiv, FArrayBox &a_faceGrad, const FArrayBox &a_data, const FArrayBox &a_gradData, const ProblemDomain &a_domain, const Box &a_faceBox, const int &a_faceDir, const Real a_dx)
 
static void getFlux (FArrayBox &a_flux, const FArrayBox &a_data, const FArrayBox &a_etaFace, const FArrayBox &a_lamFace, const Box &a_faceBox, const ProblemDomain &a_domain, const Real &a_dx, const Real &a_beta, int a_dir, int a_ref)
 

AMRLevelOp functions

static int s_prolongType
 prolongation type More...
 
RefCountedPtr< LevelData< FluxBox > > m_eta
 
RefCountedPtr< LevelData< FluxBox > > m_lambda
 
RefCountedPtr< LevelData< FArrayBox > > m_acoef
 
Real m_alpha
 
Real m_beta
 
int m_refToCoar
 
int m_refToFine
 
BCHolder m_bc
 
DisjointBoxLayout m_grids
 
DisjointBoxLayout m_gridsCoar
 
DisjointBoxLayout m_gridsFine
 
int m_ncomp
 
Real m_dx
 
Real m_dxCrse
 
ProblemDomain m_domain
 
Real m_safety
 
Real m_relaxTolerance
 
int m_relaxMinIter
 
LevelDataOps< FArrayBoxm_levelOps
 
Copier m_exchangeCopier
 
NWOQuadCFInterp m_interpWithCoarser
 
Vector< IntVectm_colors
 
LevelData< FArrayBoxm_relaxCoef
 
virtual int refToCoarser ()
 
virtual void AMRResidual (LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_phiCoarse, const LevelData< FArrayBox > &a_rhs, bool a_homogeneousPhysBC, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
 
virtual void AMRResidualNC (LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_rhs, bool a_homogeneousPhysBC, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
 
virtual void AMRResidualNF (LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_phiCoarse, const LevelData< FArrayBox > &a_rhs, bool a_homogeneousPhysBC)
 
virtual void AMRRestrict (LevelData< FArrayBox > &a_resCoarse, const LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection, bool a_skip_res=false)
 
virtual void AMRProlong (LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection)
 
virtual void AMRUpdateResidual (LevelData< FArrayBox > &a_residual, const LevelData< FArrayBox > &a_correction, const LevelData< FArrayBox > &a_coarseCorrection)
 
virtual Real AMRNorm (const LevelData< FArrayBox > &a_coarseResid, const LevelData< FArrayBox > &a_fineResid, const int &a_refRat, const int &a_ord)
 
virtual void outputLevel (LevelData< FArrayBox > &a_rhs, string &a_name)
 
virtual void outputAMR (Vector< LevelData< FArrayBox > * > &a_rhs, string &a_name)
 
void homogeneousCFInterp (LevelData< FArrayBox > &a_phif)
 
void homogeneousCFInterpPhi (LevelData< FArrayBox > &a_phif, const DataIndex &a_datInd, int a_idir, Side::LoHiSide a_hiorlo)
 
void homogeneousCFInterpTanGrad (LevelData< FArrayBox > &a_tanGrad, const LevelData< FArrayBox > &a_phi, const DataIndex &a_datInd, int a_idir, Side::LoHiSide a_hiorlo)
 
void interpOnIVSHomo (LevelData< FArrayBox > &a_phif, const DataIndex &a_datInd, const int a_idir, const Side::LoHiSide a_hiorlo, const IntVectSet &a_interpIVS)
 
void AMROperator (LevelData< FArrayBox > &a_LofPhi, const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_phiCoarse, bool a_homogeneousDomBC, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
 
void AMROperatorNF (LevelData< FArrayBox > &a_LofPhi, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > &a_phiCoarse, bool a_homogeneousBC)
 
virtual void AMROperatorNC (LevelData< FArrayBox > &a_LofPhi, const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phi, bool a_homogeneousBC, AMRLevelOp< LevelData< FArrayBox > > *a_finerOp)
 
void cfinterp (const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phiCoarse)
 
RefCountedPtr< LevelData< FluxBox > > getEta () const
 access function More...
 
RefCountedPtr< LevelData< FluxBox > > getLambda () const
 
RefCountedPtr< LevelData< FArrayBox > > getACoef () const
 
Real getAlpha () const
 
Real getBeta () const
 
void defineRelCoef ()
 
 NWOViscousTensorOp ()
 weak construction is bad More...
 
 NWOViscousTensorOp (const NWOViscousTensorOp &a_opin)
 
void operator= (const NWOViscousTensorOp &a_opin)
 

Detailed Description

operator is alpha a I + (divF) = alpha*acoef I + beta*div(eta(grad B + grad BT) + lambda I div B ) beta, lambda, and eta are incorporated into the flux F

Member Enumeration Documentation

◆ prolongationType

Enumerator
piecewiseConstant 
linearInterp 
NUM_INTERP_TYPE 

Constructor & Destructor Documentation

◆ ~NWOViscousTensorOp()

virtual NWOViscousTensorOp::~NWOViscousTensorOp ( )
inlinevirtual

◆ NWOViscousTensorOp() [1/4]

NWOViscousTensorOp::NWOViscousTensorOp ( const DisjointBoxLayout a_grids,
const DisjointBoxLayout a_gridsFine,
const DisjointBoxLayout a_gridsCoar,
const RefCountedPtr< LevelData< FluxBox > > &  a_eta,
const RefCountedPtr< LevelData< FluxBox > > &  a_lambda,
const RefCountedPtr< LevelData< FArrayBox > > &  a_acoef,
Real  a_alpha,
Real  a_beta,
int  a_refToFine,
int  a_refToCoar,
const ProblemDomain a_domain,
const Real a_dxLevel,
const Real a_dxCoar,
BCFunc  a_bc,
Real  a_safety = VTOP_DEFAULT_SAFETY,
Real  a_relaxTolerance = 0.0,
int  a_relaxMinIter = 1000 
)

◆ NWOViscousTensorOp() [2/4]

NWOViscousTensorOp::NWOViscousTensorOp ( const DisjointBoxLayout a_grids,
const DisjointBoxLayout a_gridsFine,
const DisjointBoxLayout a_gridsCoar,
const RefCountedPtr< LevelData< FluxBox > > &  a_eta,
const RefCountedPtr< LevelData< FluxBox > > &  a_lambda,
const RefCountedPtr< LevelData< FArrayBox > > &  a_acoef,
Real  a_alpha,
Real  a_beta,
int  a_refToFine,
int  a_refToCoar,
const ProblemDomain a_domain,
const Real a_dxLevel,
const Real a_dxCoar,
BCHolder a_bc,
Real  a_safety = VTOP_DEFAULT_SAFETY,
Real  a_relaxTolerance = 0.0,
int  a_relaxMinIter = 1000 
)

◆ NWOViscousTensorOp() [3/4]

NWOViscousTensorOp::NWOViscousTensorOp ( )
inlineprivate

weak construction is bad

References MayDay::Error().

Referenced by ~NWOViscousTensorOp().

◆ NWOViscousTensorOp() [4/4]

NWOViscousTensorOp::NWOViscousTensorOp ( const NWOViscousTensorOp a_opin)
inlineprivate

References MayDay::Error().

Member Function Documentation

◆ diagonalScale() [1/2]

virtual void NWOViscousTensorOp::diagonalScale ( LevelData< FArrayBox > &  a_rhs,
bool  a_kappaWeighted 
)
inlinevirtual

Scales the right hand side of the Helmholtz equation by the identity term in the operator. If you are solving rho(x) dphi/dt = L(phi), this means multiply by rho (or kappa * rho in embedded boundary calculations.

Parameters
a_rhsThe right hand side of the equation to be scaled.
a_kappaWeightedIf set to true, a_rhs will be scaled by the volume fraction in addition to the identity term.

Reimplemented from TGAHelmOp< LevelData< FArrayBox > >.

◆ setAlphaAndBeta()

virtual void NWOViscousTensorOp::setAlphaAndBeta ( const Real a_alpha,
const Real a_beta 
)
inlinevirtual

Sets the scaling constants in the Helmholtz equation.

Parameters
a_alphaThe scaling constant that multiplies the identity term.
a_betaThe scaling constant that multiplies the derivative term.

Implements TGAHelmOp< LevelData< FArrayBox > >.

References defineRelCoef(), m_alpha, and m_beta.

◆ diagonalScale() [2/2]

virtual void NWOViscousTensorOp::diagonalScale ( LevelData< FArrayBox > &  a_rhs)
inlinevirtual

Scales the right hand side of the Helmholtz equation by the identity term in the operator. This version assumes volume fraction weighting.

Parameters
a_rhsThe right hand side of the equation to be scaled.

Reimplemented from TGAHelmOp< LevelData< FArrayBox > >.

References BoxLayout::dataIterator(), LevelData< T >::disjointBoxLayout(), m_acoef, LayoutIterator::ok(), and SpaceDim.

◆ divideByIdentityCoef()

virtual void NWOViscousTensorOp::divideByIdentityCoef ( LevelData< FArrayBox > &  a_rhs)
inlinevirtual

Divides the right hand side of the Helmholtz equation by the identity coefficient rho(x) in the equation rho(x) dphi/dt = L(phi).

Parameters
a_rhsThe right hand side of the equation to be scaled.

Implements TGAHelmOp< LevelData< FArrayBox > >.

References BoxLayout::dataIterator(), LevelData< T >::disjointBoxLayout(), m_acoef, LayoutIterator::ok(), and SpaceDim.

◆ fillGrad()

void NWOViscousTensorOp::fillGrad ( const LevelData< FArrayBox > &  a_phiFine)
inlinevirtual

Fillgrad here is a placeholder that does not do anything. The whole fillGrad paradigm is an artifact of the how old ViscousTensorOp worked. Sometimes performance tuning eventually forces you to do the right thing.

Implements LevelTGAHelmOp< LevelData< FArrayBox >, FluxBox >.

References applyOp(), cellGrad(), create(), createCoarsened(), define(), getFaceDivAndGrad(), loHiCenterFace(), preCond(), reflux(), residual(), and VTOP_DEFAULT_SAFETY.

◆ loHiCenterFace()

static void NWOViscousTensorOp::loHiCenterFace ( Box a_loBox,
int &  a_hasLo,
Box a_hiBox,
int &  a_hasHi,
Box a_centerBox,
const ProblemDomain a_eblg,
const Box a_inBox,
const int &  a_dir 
)
static

Referenced by fillGrad().

◆ cellGrad()

void NWOViscousTensorOp::cellGrad ( FArrayBox a_gradPhi,
const FArrayBox a_phi,
const Box a_grid 
)

Referenced by fillGrad().

◆ getFaceDivAndGrad()

static void NWOViscousTensorOp::getFaceDivAndGrad ( FArrayBox a_faceDiv,
FArrayBox a_faceGrad,
const FArrayBox a_data,
const FArrayBox a_gradData,
const ProblemDomain a_domain,
const Box a_faceBox,
const int &  a_faceDir,
const Real  a_dx 
)
static

Referenced by fillGrad().

◆ define()

void NWOViscousTensorOp::define ( const DisjointBoxLayout a_grids,
const DisjointBoxLayout a_gridsFine,
const DisjointBoxLayout a_gridsCoar,
const RefCountedPtr< LevelData< FluxBox > > &  a_eta,
const RefCountedPtr< LevelData< FluxBox > > &  a_lambda,
const RefCountedPtr< LevelData< FArrayBox > > &  a_acoef,
Real  a_alpha,
Real  a_beta,
int  a_refToFine,
int  a_refToCoar,
const ProblemDomain a_domain,
const Real a_dxLevel,
const Real a_dxCoar,
BCHolder a_bc,
Real  a_safety = VTOP_DEFAULT_SAFETY,
Real  a_relaxTolerance = 0.0,
int  a_relaxMinIter = 1000 
)

◆ residual()

virtual void NWOViscousTensorOp::residual ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_rhs,
bool  a_homogeneous = false 
)
virtual

Say you are solving L(phi) = rhs. Make a_lhs = L(a_phi) - a_rhs. If a_homogeneous is true, evaluate the operator using homogeneous boundary conditions.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by fillGrad().

◆ preCond()

virtual void NWOViscousTensorOp::preCond ( LevelData< FArrayBox > &  a_cor,
const LevelData< FArrayBox > &  a_residual 
)
virtual

Given the current state of the residual the correction, apply your preconditioner to a_cor.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by fillGrad().

◆ applyOp()

virtual void NWOViscousTensorOp::applyOp ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_phi,
bool  a_homogeneous = false 
)
virtual

In the context of solving L(phi) = rhs, set a_lhs = L(a_phi). If a_homogeneous is true, evaluate the operator using homogeneous boundary conditions.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by fillGrad().

◆ create()

virtual void NWOViscousTensorOp::create ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_rhs 
)
virtual

Creat data holder a_lhs that mirrors a_rhs. You do not need to copy the data of a_rhs, just make a holder the same size.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by fillGrad().

◆ createCoarsened()

virtual void NWOViscousTensorOp::createCoarsened ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_rhs,
const int &  a_refRat 
)
virtual

◆ reflux()

void NWOViscousTensorOp::reflux ( const LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_phi,
LevelData< FArrayBox > &  residual,
AMRLevelOp< LevelData< FArrayBox > > *  a_finerOp 
)

Referenced by fillGrad().

◆ getFlux() [1/3]

virtual void NWOViscousTensorOp::getFlux ( FluxBox a_flux,
const LevelData< FArrayBox > &  a_data,
const Box a_grid,
const DataIndex a_dit,
Real  a_scale 
)
inlinevirtual

◆ applyOpNoBoundary()

void NWOViscousTensorOp::applyOpNoBoundary ( LevelData< FArrayBox > &  a_ans,
const LevelData< FArrayBox > &  a_phi 
)
inlinevirtual

Apply the differential operator without any boundary or coarse-fine boundary conditions and no finer level

Parameters
a_ansThe result of the application of the operator to a_phi.
a_phiThe data to which the operator is applied.

Implements TGAHelmOp< LevelData< FArrayBox > >.

References assign(), axby(), computeOperatorNoBCs(), dotProduct(), getFlux(), incr(), norm(), and scale().

◆ getFlux() [2/3]

void NWOViscousTensorOp::getFlux ( FArrayBox a_flux,
const FArrayBox a_data,
const FArrayBox a_etaFace,
const FArrayBox a_lambdaFace,
const Box a_facebox,
int  a_dir,
int  ref = 1 
)

◆ getFlux() [3/3]

static void NWOViscousTensorOp::getFlux ( FArrayBox a_flux,
const FArrayBox a_data,
const FArrayBox a_etaFace,
const FArrayBox a_lamFace,
const Box a_faceBox,
const ProblemDomain a_domain,
const Real a_dx,
const Real a_beta,
int  a_dir,
int  a_ref 
)
static

◆ computeOperatorNoBCs()

void NWOViscousTensorOp::computeOperatorNoBCs ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_phi 
)

utility function which computes operator after all bc's have been set

Referenced by applyOpNoBoundary().

◆ assign()

virtual void NWOViscousTensorOp::assign ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_rhs 
)
virtual

Set a_lhs equal to a_rhs.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by applyOpNoBoundary().

◆ dotProduct()

virtual Real NWOViscousTensorOp::dotProduct ( const LevelData< FArrayBox > &  a_1,
const LevelData< FArrayBox > &  a_2 
)
virtual

Compute and return the dot product of a_1 and a_2. In most contexts, this means return the sum over all data points of a_1*a_2.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by applyOpNoBoundary().

◆ incr()

virtual void NWOViscousTensorOp::incr ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_x,
Real  a_scale 
)
virtual

Increment by scaled amount (a_lhs += a_scale*a_x).

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by applyOpNoBoundary().

◆ axby()

virtual void NWOViscousTensorOp::axby ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_x,
const LevelData< FArrayBox > &  a_y,
Real  a_a,
Real  a_b 
)
virtual

Set input to a scaled sum (a_lhs = a_a*a_x + a_b*a_y).

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by applyOpNoBoundary().

◆ scale()

virtual void NWOViscousTensorOp::scale ( LevelData< FArrayBox > &  a_lhs,
const Real a_scale 
)
virtual

Multiply the input by a given scale (a_lhs *= a_scale).

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by applyOpNoBoundary().

◆ norm()

virtual Real NWOViscousTensorOp::norm ( const LevelData< FArrayBox > &  a_rhs,
int  a_ord 
)
virtual

Return the norm of a_rhs. a_ord == 0 max norm, a_ord == 1 sum(abs(a_rhs)), else, L(a_ord) norm.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by applyOpNoBoundary().

◆ dx()

virtual Real NWOViscousTensorOp::dx ( ) const
inlinevirtual

Return dx at this level of refinement

Reimplemented from LinearOp< LevelData< FArrayBox > >.

References m_dx.

◆ dxCrse()

virtual Real NWOViscousTensorOp::dxCrse ( ) const
inlinevirtual

◆ setToZero()

virtual void NWOViscousTensorOp::setToZero ( LevelData< FArrayBox > &  a_lhs)
virtual

Set a_lhs to zero.

Implements LinearOp< LevelData< FArrayBox > >.

Referenced by dxCrse().

◆ relax()

virtual void NWOViscousTensorOp::relax ( LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_residual,
int  a_iterations 
)
virtual

Use your relaxation operator to remove the high frequency wave numbers from the correction so that it may be averaged to a coarser refinement. A point relaxtion scheme, for example takes the form a_correction -= lambda*(L(a_correction) - a_residual).

Implements MGLevelOp< LevelData< FArrayBox > >.

Referenced by dxCrse().

◆ createCoarser()

virtual void NWOViscousTensorOp::createCoarser ( LevelData< FArrayBox > &  a_coarse,
const LevelData< FArrayBox > &  a_fine,
bool  ghosted 
)
virtual

Create a coarsened (by two) version of the input data. This does not include averaging the data. So if a_fine is over a Box of (0, 0, 0) (63, 63, 63), a_fine should be over a Box (0, 0, 0) (31, 31, 31).

Implements MGLevelOp< LevelData< FArrayBox > >.

Referenced by dxCrse().

◆ restrictResidual()

virtual void NWOViscousTensorOp::restrictResidual ( LevelData< FArrayBox > &  a_resCoarse,
LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_rhsFine 
)
virtual

calculate restricted residual a_resCoarse[2h] = I[h->2h] (rhsFine[h] - L[h](phiFine[h])

Implements MGLevelOp< LevelData< FArrayBox > >.

Referenced by dxCrse().

◆ prolongIncrement()

virtual void NWOViscousTensorOp::prolongIncrement ( LevelData< FArrayBox > &  a_phiThisLevel,
const LevelData< FArrayBox > &  a_correctCoarse 
)
virtual

correct the fine solution based on coarse correction a_phiThisLevel += I[2h->h](a_correctCoarse)

Implements MGLevelOp< LevelData< FArrayBox > >.

Referenced by dxCrse().

◆ refToCoarser()

virtual int NWOViscousTensorOp::refToCoarser ( )
inlinevirtual

◆ AMRResidual()

virtual void NWOViscousTensorOp::AMRResidual ( LevelData< FArrayBox > &  a_residual,
const LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_phiCoarse,
const LevelData< FArrayBox > &  a_rhs,
bool  a_homogeneousPhysBC,
AMRLevelOp< LevelData< FArrayBox > > *  a_finerOp 
)
virtual

a_residual = a_rhs - L(a_phi, a_phiFine, a_phiCoarse)

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRResidualNC()

virtual void NWOViscousTensorOp::AMRResidualNC ( LevelData< FArrayBox > &  a_residual,
const LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_rhs,
bool  a_homogeneousPhysBC,
AMRLevelOp< LevelData< FArrayBox > > *  a_finerOp 
)
virtual

residual assuming no more coarser AMR levels

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRResidualNF()

virtual void NWOViscousTensorOp::AMRResidualNF ( LevelData< FArrayBox > &  a_residual,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_phiCoarse,
const LevelData< FArrayBox > &  a_rhs,
bool  a_homogeneousPhysBC 
)
virtual

a_residual = a_rhs - L(a_phi, a_phiCoarse)

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRRestrict()

virtual void NWOViscousTensorOp::AMRRestrict ( LevelData< FArrayBox > &  a_resCoarse,
const LevelData< FArrayBox > &  a_residual,
const LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_coarseCorrection,
bool  a_skip_res = false 
)
virtual

a_resCoarse = I[h-2h]( a_residual - L(a_correction, a_coarseCorrection)) it is assumed that a_resCoarse has already been filled in with the coarse version of AMRResidualNF and that this operation is free to overwrite in the overlap regions.

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRProlong()

virtual void NWOViscousTensorOp::AMRProlong ( LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_coarseCorrection 
)
virtual

a_correction += I[h->h](a_coarseCorrection)

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRUpdateResidual()

virtual void NWOViscousTensorOp::AMRUpdateResidual ( LevelData< FArrayBox > &  a_residual,
const LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_coarseCorrection 
)
virtual

a_residual = a_residual - L(a_correction, a_coarseCorrection)

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ AMRNorm()

virtual Real NWOViscousTensorOp::AMRNorm ( const LevelData< FArrayBox > &  a_coarseResid,
const LevelData< FArrayBox > &  a_fineResid,
const int &  a_refRat,
const int &  a_ord 
)
virtual

compute norm over all cells on coarse not covered by finer

Reimplemented from AMRLevelOp< LevelData< FArrayBox > >.

Referenced by refToCoarser().

◆ outputLevel()

virtual void NWOViscousTensorOp::outputLevel ( LevelData< FArrayBox > &  a_rhs,
string &  a_name 
)
inlinevirtual

Reimplemented from AMRLevelOp< LevelData< FArrayBox > >.

References writeLevelname().

Referenced by outputAMR().

◆ outputAMR()

virtual void NWOViscousTensorOp::outputAMR ( Vector< LevelData< FArrayBox > * > &  a_rhs,
string &  a_name 
)
inlinevirtual

◆ homogeneousCFInterp()

void NWOViscousTensorOp::homogeneousCFInterp ( LevelData< FArrayBox > &  a_phif)

Referenced by outputAMR().

◆ homogeneousCFInterpPhi()

void NWOViscousTensorOp::homogeneousCFInterpPhi ( LevelData< FArrayBox > &  a_phif,
const DataIndex a_datInd,
int  a_idir,
Side::LoHiSide  a_hiorlo 
)

does homogeneous coarse/fine interpolation for phi

Referenced by outputAMR().

◆ homogeneousCFInterpTanGrad()

void NWOViscousTensorOp::homogeneousCFInterpTanGrad ( LevelData< FArrayBox > &  a_tanGrad,
const LevelData< FArrayBox > &  a_phi,
const DataIndex a_datInd,
int  a_idir,
Side::LoHiSide  a_hiorlo 
)

does homogeneous coarse/fine interpolation for tangential gradient (needs phi ghost cells to be filled in first, so should call homogeneousCFInterpPhi first)

Referenced by outputAMR().

◆ interpOnIVSHomo()

void NWOViscousTensorOp::interpOnIVSHomo ( LevelData< FArrayBox > &  a_phif,
const DataIndex a_datInd,
const int  a_idir,
const Side::LoHiSide  a_hiorlo,
const IntVectSet a_interpIVS 
)

Referenced by outputAMR().

◆ AMROperator()

void NWOViscousTensorOp::AMROperator ( LevelData< FArrayBox > &  a_LofPhi,
const LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_phiCoarse,
bool  a_homogeneousDomBC,
AMRLevelOp< LevelData< FArrayBox > > *  a_finerOp 
)
virtual

Apply the AMR operator, including coarse-fine matching

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by outputAMR().

◆ AMROperatorNF()

void NWOViscousTensorOp::AMROperatorNF ( LevelData< FArrayBox > &  a_LofPhi,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > &  a_phiCoarse,
bool  a_homogeneousBC 
)
virtual

Apply the AMR operator, including coarse-fine matching. assume no finer AMR level

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by outputAMR().

◆ AMROperatorNC()

virtual void NWOViscousTensorOp::AMROperatorNC ( LevelData< FArrayBox > &  a_LofPhi,
const LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_phi,
bool  a_homogeneousBC,
AMRLevelOp< LevelData< FArrayBox > > *  a_finerOp 
)
virtual

Apply the AMR operator, including coarse-fine matching assume no coarser AMR level

Implements AMRLevelOp< LevelData< FArrayBox > >.

Referenced by outputAMR().

◆ cfinterp()

void NWOViscousTensorOp::cfinterp ( const LevelData< FArrayBox > &  a_phiFine,
const LevelData< FArrayBox > &  a_phiCoarse 
)

Referenced by outputAMR().

◆ getEta()

RefCountedPtr<LevelData<FluxBox> > NWOViscousTensorOp::getEta ( ) const
inline

access function

References m_eta.

◆ getLambda()

RefCountedPtr<LevelData<FluxBox> > NWOViscousTensorOp::getLambda ( ) const
inline

References m_lambda.

◆ getACoef()

RefCountedPtr<LevelData<FArrayBox> > NWOViscousTensorOp::getACoef ( ) const
inline

References m_acoef.

◆ getAlpha()

Real NWOViscousTensorOp::getAlpha ( ) const
inline

References m_alpha.

◆ getBeta()

Real NWOViscousTensorOp::getBeta ( ) const
inline

References defineRelCoef(), and m_beta.

◆ defineRelCoef()

void NWOViscousTensorOp::defineRelCoef ( )
protected

Referenced by getBeta(), and setAlphaAndBeta().

◆ operator=()

void NWOViscousTensorOp::operator= ( const NWOViscousTensorOp a_opin)
inlineprivate

References MayDay::Error().

Member Data Documentation

◆ s_prolongType

int NWOViscousTensorOp::s_prolongType
static

prolongation type

make this public/static to make it easy to change, while also ensuring that it remains constent across all instances

◆ m_eta

RefCountedPtr<LevelData<FluxBox> > NWOViscousTensorOp::m_eta
protected

◆ m_lambda

RefCountedPtr<LevelData<FluxBox> > NWOViscousTensorOp::m_lambda
protected

◆ m_acoef

RefCountedPtr<LevelData<FArrayBox> > NWOViscousTensorOp::m_acoef
protected

◆ m_alpha

Real NWOViscousTensorOp::m_alpha
protected

Referenced by getAlpha(), and setAlphaAndBeta().

◆ m_beta

Real NWOViscousTensorOp::m_beta
protected

Referenced by getBeta(), and setAlphaAndBeta().

◆ m_refToCoar

int NWOViscousTensorOp::m_refToCoar
protected

Referenced by refToCoarser().

◆ m_refToFine

int NWOViscousTensorOp::m_refToFine
protected

◆ m_bc

BCHolder NWOViscousTensorOp::m_bc
protected

◆ m_grids

DisjointBoxLayout NWOViscousTensorOp::m_grids
protected

◆ m_gridsCoar

DisjointBoxLayout NWOViscousTensorOp::m_gridsCoar
protected

◆ m_gridsFine

DisjointBoxLayout NWOViscousTensorOp::m_gridsFine
protected

◆ m_ncomp

int NWOViscousTensorOp::m_ncomp
protected

◆ m_dx

Real NWOViscousTensorOp::m_dx
protected

Referenced by dx().

◆ m_dxCrse

Real NWOViscousTensorOp::m_dxCrse
protected

Referenced by dxCrse().

◆ m_domain

ProblemDomain NWOViscousTensorOp::m_domain
protected

Referenced by getFlux().

◆ m_relaxCoef

LevelData<FArrayBox> NWOViscousTensorOp::m_relaxCoef

◆ m_safety

Real NWOViscousTensorOp::m_safety
protected

◆ m_relaxTolerance

Real NWOViscousTensorOp::m_relaxTolerance
protected

◆ m_relaxMinIter

int NWOViscousTensorOp::m_relaxMinIter
protected

◆ m_levelOps

LevelDataOps<FArrayBox> NWOViscousTensorOp::m_levelOps
protected

◆ m_exchangeCopier

Copier NWOViscousTensorOp::m_exchangeCopier
protected

◆ m_interpWithCoarser

NWOQuadCFInterp NWOViscousTensorOp::m_interpWithCoarser
protected

◆ m_colors

Vector<IntVect> NWOViscousTensorOp::m_colors
protected

The documentation for this class was generated from the following file: