ResistivityOp Class Reference

#include <ResistivityOp.H>

Inheritance diagram for ResistivityOp:

Inheritance graph
[legend]
Collaboration diagram for ResistivityOp:

Collaboration graph
[legend]

List of all members.


Detailed Description

operator L is defined as L(B) = alpha B + beta*(divF) where F = div(eta(grad B - grad B^T) + eta I div B ) alpha and beta are constants. eta is a function of space. beta and eta are incorporated into the flux F when getFlux is called. This is always 3 variables. SpaceDim refers to only the derivatives that are meaningful in this bizzare application.

AMRLevelOp functions

RefCountedPtr< LevelData
< FluxBox > > 
m_eta
DisjointBoxLayout m_grids
Real m_alpha
Real m_beta
int m_refToCoar
int m_refToFine
BCFunc m_bc
Real m_dx
Real m_dxCrse
ProblemDomain m_domain
LevelData< FArrayBoxm_lambda
LevelData< FArrayBoxm_grad
LevelDataOps< FArrayBoxm_levelOps
Copier m_exchangeCopier
TensorCFInterp m_interpWithCoarser
LayoutData< CFIVSm_loCFIVS [SpaceDim]
LayoutData< CFIVSm_hiCFIVS [SpaceDim]
LayoutData< TensorFineStencilSetm_hiTanStencilSets [SpaceDim]
LayoutData< TensorFineStencilSetm_loTanStencilSets [SpaceDim]
Vector< IntVectm_colors
static const int s_nComp
static const int s_nGradComp
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)
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)
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 cellGrad (FArrayBox &a_gradPhi, const FArrayBox &a_phi, const Box &a_grid)
void cfinterp (const LevelData< FArrayBox > &a_phiFine, const LevelData< FArrayBox > &a_phiCoarse)
virtual void fillGrad (const LevelData< FArrayBox > &a_phiFine)
 These functions are part of the LevelTGA interface......
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)
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)
void divergenceCC (LevelData< FArrayBox > &a_div, const LevelData< FArrayBox > &a_phi, const LevelData< FArrayBox > *a_phiC)
 take cell centered divergence of the inputs.
void setLambda ()
 ResistivityOp ()
 weak construction is bad
 ResistivityOp (const ResistivityOp &a_opin)
void operator= (const ResistivityOp &a_opin)

Public Member Functions

virtual void setAlphaAndBeta (const Real &a_alpha, const Real &a_beta)
 set the constants in the equation
virtual void diagonalScale (LevelData< FArrayBox > &a_rhs)
virtual ~ResistivityOp ()
virtual void getFlux (FluxBox &a_flux, const LevelData< FArrayBox > &a_data, const Box &a_grid, const DataIndex &a_dit, Real a_scale)
 ResistivityOp (const DisjointBoxLayout &a_grids, const DisjointBoxLayout &a_gridsFine, const DisjointBoxLayout &a_gridsCoar, const RefCountedPtr< LevelData< FluxBox > > &a_eta, 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)
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 applyOpNoBoundary (LevelData< FArrayBox > &a_lhs, const LevelData< FArrayBox > &a_phi)
 apply operator without any boundary or coarse-fine boundary conditions and no finer level
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)
void getFlux (FArrayBox &a_flux, const FArrayBox &a_data, const FArrayBox &a_gradData, const FArrayBox &a_etaFace, 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
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 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)

Constructor & Destructor Documentation

virtual ResistivityOp::~ResistivityOp (  )  [inline, virtual]

ResistivityOp::ResistivityOp ( const DisjointBoxLayout a_grids,
const DisjointBoxLayout a_gridsFine,
const DisjointBoxLayout a_gridsCoar,
const RefCountedPtr< LevelData< FluxBox > > &  a_eta,
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 
)

ResistivityOp::ResistivityOp (  )  [inline, private]

weak construction is bad

References MayDay::Error().

ResistivityOp::ResistivityOp ( const ResistivityOp a_opin  )  [inline, private]

References MayDay::Error().


Member Function Documentation

virtual void ResistivityOp::setAlphaAndBeta ( const Real a_alpha,
const Real a_beta 
) [inline, virtual]

set the constants in the equation

Implements TGAHelmOp< LevelData< FArrayBox > >.

References m_alpha, m_beta, and setLambda().

virtual void ResistivityOp::diagonalScale ( LevelData< FArrayBox > &  a_rhs  )  [inline, virtual]

Set the diagonal scaling of the operator. If you are solving rho(x) dphi/dt = L(phi), this would mean multiply by rho. In EB applications, even for constant coefficients, it means to multiply by kappa.

Implements TGAHelmOp< LevelData< FArrayBox > >.

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

virtual void ResistivityOp::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 > >.

virtual void ResistivityOp::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 > >.

virtual void ResistivityOp::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 > >.

virtual void ResistivityOp::applyOpNoBoundary ( LevelData< FArrayBox > &  a_ans,
const LevelData< FArrayBox > &  a_phi 
) [virtual]

apply operator without any boundary or coarse-fine boundary conditions and no finer level

Implements TGAHelmOp< LevelData< FArrayBox > >.

virtual void ResistivityOp::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 > >.

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

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

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

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

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

virtual void ResistivityOp::assign ( LevelData< FArrayBox > &  a_lhs,
const LevelData< FArrayBox > &  a_rhs 
) [virtual]

Set a_lhs equal to a_rhs.

Implements LinearOp< LevelData< FArrayBox > >.

virtual Real ResistivityOp::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 > >.

virtual void ResistivityOp::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 > >.

virtual void ResistivityOp::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 > >.

virtual void ResistivityOp::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 > >.

virtual Real ResistivityOp::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 > >.

virtual void ResistivityOp::setToZero ( LevelData< FArrayBox > &  a_lhs  )  [virtual]

Set a_lhs to zero.

Implements LinearOp< LevelData< FArrayBox > >.

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

Use your relaxtion 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 > >.

virtual void ResistivityOp::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 > >.

virtual void ResistivityOp::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 > >.

virtual void ResistivityOp::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 > >.

virtual int ResistivityOp::refToCoarser (  )  [inline, virtual]

returns 1 when there are no coarser AMRLevelOp objects

Implements AMRLevelOp< LevelData< FArrayBox > >.

References m_refToCoar.

virtual void ResistivityOp::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 > >.

virtual void ResistivityOp::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 > >.

virtual void ResistivityOp::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 > >.

virtual void ResistivityOp::AMRRestrict ( LevelData< FArrayBox > &  a_resCoarse,
const LevelData< FArrayBox > &  a_residual,
const LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_coarseCorrection 
) [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 > >.

virtual void ResistivityOp::AMRProlong ( LevelData< FArrayBox > &  a_correction,
const LevelData< FArrayBox > &  a_coarseCorrection 
) [virtual]

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

Implements AMRLevelOp< LevelData< FArrayBox > >.

virtual void ResistivityOp::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 > >.

virtual Real ResistivityOp::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 > >.

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

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

does homogeneous coarse/fine interpolation for phi

void ResistivityOp::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)

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

void ResistivityOp::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 > >.

void ResistivityOp::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 > >.

virtual void ResistivityOp::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 > >.

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

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

virtual void ResistivityOp::fillGrad ( const LevelData< FArrayBox > &  a_phi  )  [virtual]

These functions are part of the LevelTGA interface......

Implements LevelTGAHelmOp< LevelData< FArrayBox >, FluxBox >.

void ResistivityOp::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 
)

void ResistivityOp::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 
)

void ResistivityOp::divergenceCC ( LevelData< FArrayBox > &  a_div,
const LevelData< FArrayBox > &  a_phi,
const LevelData< FArrayBox > *  a_phiC 
)

take cell centered divergence of the inputs.

not part of the operator evaluation. this is to test whether the divergence of b converges at h^2 as b does if b is analytically divergence-free.

void ResistivityOp::setLambda (  )  [protected]

Referenced by setAlphaAndBeta().

void ResistivityOp::operator= ( const ResistivityOp a_opin  )  [inline, private]

References MayDay::Error().


Member Data Documentation

Referenced by setAlphaAndBeta().

Referenced by setAlphaAndBeta().

int ResistivityOp::m_refToCoar [protected]

Referenced by refToCoarser().

int ResistivityOp::m_refToFine [protected]

const int ResistivityOp::s_nComp [static, protected]

const int ResistivityOp::s_nGradComp [static, protected]


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

Generated on Tue Apr 14 14:23:58 2009 for Chombo + EB by  doxygen 1.5.5