Main Page | Namespace List | Class Hierarchy | Alphabetical List | Compound List | File List | Compound Members | File Members

ParticleAMRNS Class Reference

#include <ParticleAMRNS.H>

Inheritance diagram for ParticleAMRNS:

Inheritance graph
[legend]
Collaboration diagram for ParticleAMRNS:

Collaboration graph
[legend]
List of all members.

Public Member Functions

 ParticleAMRNS ()
 ParticleAMRNS (AMRLevel *a_coarser_level_ptr, const Box &a_prob_domain, int a_level, int a_ref_ratio)
 ParticleAMRNS (AMRLevel *a_coarser_level_ptr, const ProblemDomain &a_prob_domain, int a_level, int a_ref_ratio)
virtual ~ParticleAMRNS ()
virtual AMRLevelmakeAMRLevel (AMRLevel *a_coarser_level_ptr, const Box &a_problem_domain, int a_level, int a_ref_ratio) const
virtual AMRLevelmakeAMRLevel (AMRLevel *a_coarser_level_ptr, const ProblemDomain &a_problem_domain, int a_level, int a_ref_ratio) const
virtual void define (AMRLevel *a_coarse_level_ptr, const Box &a_problem_domain, int a_level, int a_ref_ratio)
virtual void define (AMRLevel *a_coarse_level_ptr, const ProblemDomain &a_problem_domain, int a_level, int a_ref_ratio)
virtual Real advance ()
 advance by one timestep (returns what?)

virtual void postTimeStep ()
 things to do after a timestep

virtual void tagCells (IntVectSet &a_tags)
 create tags

virtual void tagCellsInit (IntVectSet &a_tags)
 create tags at initialization

virtual void regrid (const Vector< Box > &a_new_grids)
 regrid

virtual void postRegrid (int a_lBase)
 perform any post-regridding ops -- lBase is the finest unchanged level

virtual void initialGrid (const Vector< Box > &a_new_grids)
 initialize grids

virtual void initialData ()
 initialize data

virtual void postInitialize ()
 things do do after initialization

virtual Real computeDt ()
 compute dt

virtual Real computeInitialDt ()
 compute dt with initial data

virtual void CFL (Real a_cfl)
 set CFL number

virtual void particleCFL (Real a_particle_cfl)
 set CFL number for computing timestep based on particles

virtual void refinementThreshold (Real a_refine_threshold)
 set refinement threshold

void limitSolverCoarsening (bool a_limitSolverCoarsening)
 set solver parameter

bool finestLevel () const
 is this this finest level?

bool isEmpty () const
 is this an empty level?

void computeVorticity (LevelData< FArrayBox > &a_vorticity) const
void computeKineticEnergy (LevelData< FArrayBox > &a_energy) const
Real particleKineticEnergy () const
void dumpParticlesASCII (std::ostream &os) const
Real Dx () const
 access functions

Real Nu () const
 viscous coefficient

ParticleAMRNSfinerNSPtr () const
 encapsulates the finer-level-pointer casting

ParticleAMRNScrseNSPtr () const
 encapsulates the coarser-level pointer casting

LevelData< FArrayBox > & newVel ()
const LevelData< FArrayBox > & newVel () const
LevelData< FArrayBox > & oldVel ()
const LevelData< FArrayBox > & oldVel () const
void velocity (LevelData< FArrayBox > &a_vel, Real a_time) const
 returns velocity at time t

LevelData< BinFab< DragParticle > > & particles ()
 returns particles

const LevelData< BinFab< DragParticle > > & particles () const
 returns particles (const version)

int numParticles () const
 number of particles on this AMR level

LevelData< FArrayBox > & newLambda ()
const LevelData< FArrayBox > & newLambda () const
LevelData< FArrayBox > & oldLambda ()
const LevelData< FArrayBox > & oldLambda () const
void lambda (LevelData< FArrayBox > &a_lambda, Real a_time) const
 returns lambda at time t

LevelData< FArrayBox > & newScal (const int a_comp)
const LevelData< FArrayBox > & newScal (const int a_comp) const
LevelData< FArrayBox > & oldScal (const int a_comp)
const LevelData< FArrayBox > & oldScal (const int a_comp) const
void scalars (LevelData< FArrayBox > &a_scalars, const Real a_time, const int a_comp) const
 returns scalars at time t

void setPhysBC (const PhysBCUtil &a_BC)
 sets physical BC object

PhysBCUtilgetPhysBCPtr () const
 gets physical BC object


Protected Member Functions

DisjointBoxLayout loadBalance (const Vector< Box > &a_grids)
void levelSetup (const DisjointBoxLayout &level_domain)
void defineParticles (List< DragParticle > &a_particles) const
 define particle list for entire domain

void readParameters ()
 read parameters from parmParse database

void finestLevel (bool a_finest_level)
 set whether this is the finest level or not

void swapOldAndNewStates ()
 move new-time state into old-time state -- also advances time by dt

void resetStates (const Real a_time)
 exchange old and new states, resets time to a_time

void computeAdvectionVelocities (LevelData< FluxBox > &a_adv_vel, LevelData< FArrayBox > &a_dragForce)
 compute advection velocities

void predictVelocities (LevelData< FArrayBox > &a_uDelU, LevelData< FluxBox > &a_advVel, LevelData< FArrayBox > &a_dragForce)
void computeUStar (LevelData< FArrayBox > &a_uStar, LevelData< FArrayBox > &a_dragForce)
void updateParticlePositions (LevelData< FArrayBox > &a_newDrag, LevelData< BinFab< DragParticle > > &a_particles, const LevelData< FArrayBox > &a_fluidVel, bool a_completeUpdate)
 move particles and compute new-time drag forcing

void updateParticleForces (LevelData< BinFab< DragParticle > > &a_particles)
 update particle forces

void interpVelToParticles (LevelData< BinFab< DragParticle > > &a_particles, const LevelData< FArrayBox > &a_velocity)
 update particle velocities by interpolating a_vel to particles

void advectScalar (LevelData< FArrayBox > &a_new_scal, LevelData< FArrayBox > &a_old_scal, LevelData< FluxBox > &a_adv_vel, LevelFluxRegister *a_crseFluxRegPtr, LevelFluxRegister &a_fineFluxReg, Real a_dt)
 do scalar advection -- assumes all BC's already done

void advectDiffuseScalar (LevelData< FArrayBox > &a_new_scal, LevelData< FArrayBox > &a_old_scal, LevelData< FluxBox > &a_adv_vel, const Real &a_diffusive_coeff, const LevelData< FArrayBox > *a_crse_scal_old, const LevelData< FArrayBox > *a_crse_scal_new, const Real a_old_crse_time, const Real a_new_crse_time, LevelFluxRegister *a_crseFluxRegPtr, LevelFluxRegister &a_fineFluxReg, DomainGhostBC &a_scalPhysBC, Real a_dt)
 advect/diffuse scalar -- assumes all BC's already done?

void initializeGlobalPressure ()
 manages pressure initialization after init or regrid

void initializeLevelPressure (Real a_currentTime, Real a_dtInit)
 does level-based pressure initalization

void smoothVelocityField (int lbase)
 smooths composite velocity field after regridding

void computeLapVel (LevelData< FArrayBox > &a_lapVel, LevelData< FArrayBox > &a_vel, const LevelData< FArrayBox > *a_crseVelPtr)
 compute Laplacian(velocity)

void computeLapScal (LevelData< FArrayBox > &a_lapScal, LevelData< FArrayBox > &a_scal, DomainGhostBC &a_physBC, const LevelData< FArrayBox > *a_crseScalPtr)
 compute Laplacian(scalars)

void fillVelocity (LevelData< FArrayBox > &a_vel, Real a_time)
 fill grown velocity field using piecewise-constant interp

void fillVelocity (LevelData< FArrayBox > &a_vel, Real a_time, int vel_comp, int dest_comp, int num_comp)
 fill grown velocity field using piecewise-constant interp

void fillLambda (LevelData< FArrayBox > &a_lambda, Real a_time) const
void fillScalars (LevelData< FArrayBox > &a_scal, Real a_time, const int a_comp) const
 fill interior and boundary cells (assume a_scal already sized)

virtual int numPlotComps () const
 number of components in plotfiles

virtual void getPlotData (LevelData< FArrayBox > &a_plot_data) const
 stuff plotfile data

void defineRegridAMRSolver (AMRSolver &a_solver, const Vector< DisjointBoxLayout > &a_grids, const Vector< ProblemDomain > &a_domains, const Vector< Real > &a_amrDx, const Vector< int > &a_refRatios, const int &a_lBase)
 defines (helmholtz) AMRSolver used in smoothing during regridding


Protected Attributes

LevelData< FArrayBox > * m_vel_old_ptr
 velocity vector at old time

LevelData< FArrayBox > * m_vel_new_ptr
 velocity vector at new time

LevelData< BinFab< DragParticle > > m_particles
 particles

LevelData< FArrayBox > * m_lambda_old_ptr
 freestream preservation qty at old time

LevelData< FArrayBox > * m_lambda_new_ptr
 freestream preservation qty at new time

Vector< LevelData< FArrayBox > * > m_scal_old
 scalars at old time

Vector< LevelData< FArrayBox > * > m_scal_new
 scalars at new time

Vector< Tuple< Real, SpaceDim > > m_dataPoints
 points where data is dumped out

Real m_dx
 grid spacing

Real m_cfl
 cfl number

Real m_particleCFL
 cfl number for particle updates

bool m_limitSolverCoarsening
 solver parameter -- passed on through when defining solvers

Real m_dt_save
 dt from most recent timestep

int m_level_steps
 number of steps on this level

CoarseAverage m_coarse_average
 averaging from fine to coarse level for velocity

CoarseAverage m_coarse_average_lambda
 averaging from fine to coarse level for freestream preservation qty

CCProjector m_projection
 cell-centered projection object

ParticleProjector m_particle_projector
 particle drag force projector object

LevelFluxRegister m_flux_register
 momentum flux register

LevelFluxRegister m_lambda_flux_reg
 freestream preservation qty flux register

Vector< LevelFluxRegister * > m_scal_fluxreg_ptrs
 conserved scalar flux register

QuadCFInterp m_velCFInterp
 does quadratic coarse-fine interpolation for velocities

VelocityPoissonOp m_velocityPoissonOp
 PoissonOp for velocities.

PoissonOp m_scalarsPoissonOp
 PoissonOp for advected scalars (for diffusion).

LevelHelmholtzSolver m_viscous_solver
 levelsolver for viscous solves

Real m_refine_threshold
 refinement threshold for regridding

bool m_finest_level
 true if this is the finest extant level

bool m_is_empty
 true if this level is an empty level (no grids)

bool m_regrid_smoothing_done
 internal indicator whether post-regrid smoothing stuff has been done

PhysBCUtilm_physBCPtr
 physical BC object


Static Protected Attributes

const int s_num_vel_comps = SpaceDim
 number of components of m_state (vel comps + lambda)

int s_num_scal_comps
 number of scalars

const int s_max_scal_comps = SpaceDim
 maximum possible number of scalar components

const char * s_vel_names [s_num_vel_comps]
 names of components

const char * s_scal_names [s_max_scal_comps]
 names of scalars

Vector< Reals_scal_coeffs
 diffusion coefficients for scalar advection-diffusion

Vector< Reals_domLength
 size of physical domain

bool s_ppInit
 have we initialized static variables from parmParse?

Real s_init_shrink
 factor to shrink initial timestep

Real s_max_dt
 maximum allowable timestep

Real s_prescribedDt
bool s_project_initial_vel
 should we project initial velocity field (default is yes)

bool s_initialize_pressures
 should we initialize pressures after grid generation & regridding?

bool s_smooth_after_regrid
 do antidiffusive smoothing after regridding?

bool s_dump_smoothing_plots
 dump before-and-after plotfiles for post-regrid smoothing (for debug)

Real s_regrid_smoothing_coeff
 antidiffusive smoothing coefficient (multiplies nu*dtCrse)

int s_num_init_passes
 number of times to iterate when initializing pressures

bool s_reflux_momentum
 should we reflux momentum?

bool s_reflux_normal_momentum
 if we're refluxing momentum, should we reflux the normal component?

bool s_set_bogus_values
 debugging option -- initially set arrays to s_bogus_value

Real s_bogus_value
 debugging option -- value to which arrays are initially set

bool s_tag_vorticity
 tag on vorticity?

Real s_vort_factor
 factor for tagging -- tag if vorticity greater than factor*max_vort

int s_tags_grow
 amount by which to grow tagged regions before calling meshrefine

bool s_implicit_reflux
 do implicit refluxing for velocity?

bool s_applyFreestreamCorrection
 apply freestream preservation correction?

bool s_reflux_scal
 should we reflux scalars?

bool s_implicit_scal_reflux
 do implicit refluxing for advected/diffused scalars?

int s_viscous_solver_type
 viscous solver type: 0 = backwards euler, 1=Crank-Nicolson, 2=TGA

Real s_viscous_solver_tol
 viscous solver tolerance (defaults to 1e-7)

int s_viscous_num_smooth_up
 multigrid parameter for viscous solves

int s_viscous_num_smooth_down
 multigrid parameter for viscous solves

Real s_nu
 viscosity

bool s_do_particles
 if true, do particle stuff, if false, don't

bool s_read_particles_from_file
 read particle info from file?

string s_particle_file
 file containing particle info

Real s_particle_drag_coeff
 particle drag coefficient

RealVect s_particle_body_force
 particle body force (gravity)

Real s_delta_epsilon
 particle delta function parameter

bool s_use_image_particles
 if true, use image particles in particle projector near phys boundaries

int s_order_vel_interp
 order of velocity interpolation for particles

bool s_specifyInitialGrids
 do we specify initial grids in a gridFile?

string s_initialGridFile
 gridfile where initial grids are spec'd, if relevant

bool s_init_vel_from_vorticity
 initialize velocity field from vorticity field

Real s_backgroundVel
 background velocity (in the z-direction)

bool s_write_particles
 -- I/O options -- include particles in plotfile?

bool s_write_divergence
 include divergence in plotfile?

bool s_write_lambda
 include lambda in plotfile?

bool s_write_time_derivatives
 include du/dt and d(lambda)/dt in plotfile?

bool s_write_vorticity
 include vorticity in plotfile?

bool s_write_scalars
 include scalars in plotfile?

bool s_write_dScalar_dt
 include d(scalars)/dt in plotfile?

bool s_write_strains
 include strain-rates in plotfile? (2D only)

bool s_write_grad_eLambda
 include freestream preservation correcion in plotfile?

bool s_write_error
 include error (if there is an exact solution) in plotfile?

bool s_write_proc_ids
 include processor distribution visualization in plotfile?

bool s_compute_scal_err
 compute error and report norms after composite timesteps

bool s_write_initial_solve_rhs
 dump out initial vorticity->stream function RHS's into plotfiles

bool s_write_grids
 write out grids and processor assignments


Constructor & Destructor Documentation

ParticleAMRNS::ParticleAMRNS  ) 
 

ParticleAMRNS::ParticleAMRNS AMRLevel a_coarser_level_ptr,
const Box a_prob_domain,
int  a_level,
int  a_ref_ratio
 

ParticleAMRNS::ParticleAMRNS AMRLevel a_coarser_level_ptr,
const ProblemDomain a_prob_domain,
int  a_level,
int  a_ref_ratio
 

virtual ParticleAMRNS::~ParticleAMRNS  )  [virtual]
 


Member Function Documentation

virtual Real ParticleAMRNS::advance  )  [virtual]
 

advance by one timestep (returns what?)

Implements AMRLevel.

void ParticleAMRNS::advectDiffuseScalar LevelData< FArrayBox > &  a_new_scal,
LevelData< FArrayBox > &  a_old_scal,
LevelData< FluxBox > &  a_adv_vel,
const Real a_diffusive_coeff,
const LevelData< FArrayBox > *  a_crse_scal_old,
const LevelData< FArrayBox > *  a_crse_scal_new,
const Real  a_old_crse_time,
const Real  a_new_crse_time,
LevelFluxRegister a_crseFluxRegPtr,
LevelFluxRegister a_fineFluxReg,
DomainGhostBC a_scalPhysBC,
Real  a_dt
[protected]
 

advect/diffuse scalar -- assumes all BC's already done?

void ParticleAMRNS::advectScalar LevelData< FArrayBox > &  a_new_scal,
LevelData< FArrayBox > &  a_old_scal,
LevelData< FluxBox > &  a_adv_vel,
LevelFluxRegister a_crseFluxRegPtr,
LevelFluxRegister a_fineFluxReg,
Real  a_dt
[protected]
 

do scalar advection -- assumes all BC's already done

virtual void ParticleAMRNS::CFL Real  a_cfl  )  [virtual]
 

set CFL number

void ParticleAMRNS::computeAdvectionVelocities LevelData< FluxBox > &  a_adv_vel,
LevelData< FArrayBox > &  a_dragForce
[protected]
 

compute advection velocities

virtual Real ParticleAMRNS::computeDt  )  [virtual]
 

compute dt

Implements AMRLevel.

virtual Real ParticleAMRNS::computeInitialDt  )  [virtual]
 

compute dt with initial data

Implements AMRLevel.

void ParticleAMRNS::computeKineticEnergy LevelData< FArrayBox > &  a_energy  )  const
 

void ParticleAMRNS::computeLapScal LevelData< FArrayBox > &  a_lapScal,
LevelData< FArrayBox > &  a_scal,
DomainGhostBC a_physBC,
const LevelData< FArrayBox > *  a_crseScalPtr
[protected]
 

compute Laplacian(scalars)

void ParticleAMRNS::computeLapVel LevelData< FArrayBox > &  a_lapVel,
LevelData< FArrayBox > &  a_vel,
const LevelData< FArrayBox > *  a_crseVelPtr
[protected]
 

compute Laplacian(velocity)

void ParticleAMRNS::computeUStar LevelData< FArrayBox > &  a_uStar,
LevelData< FArrayBox > &  a_dragForce
[protected]
 

void ParticleAMRNS::computeVorticity LevelData< FArrayBox > &  a_vorticity  )  const
 

ParticleAMRNS* ParticleAMRNS::crseNSPtr  )  const
 

encapsulates the coarser-level pointer casting

virtual void ParticleAMRNS::define AMRLevel a_coarse_level_ptr,
const ProblemDomain a_problem_domain,
int  a_level,
int  a_ref_ratio
[virtual]
 

Defines this AMRLevel.

  • a_coarser_level_ptr (not modified): pointer to next coarser level object.
  • a_problem_domain (not modified): problem domain of this level.
  • a_level (not modified): index of this level. The base level is zero.
  • a_ref_ratio (not modified): the refinement ratio between this level and the next finer level.

Reimplemented from AMRLevel.

virtual void ParticleAMRNS::define AMRLevel a_coarse_level_ptr,
const Box a_problem_domain,
int  a_level,
int  a_ref_ratio
[virtual]
 

Defines this AMRLevel.

  • a_coarser_level_ptr (not modified): pointer to next coarser level object.
  • a_problem_domain (not modified): problem domain of this level.
  • a_level (not modified): index of this level. The base level is zero.
  • a_ref_ratio (not modified): the refinement ratio between this level and the next finer level.

Reimplemented from AMRLevel.

void ParticleAMRNS::defineParticles List< DragParticle > &  a_particles  )  const [protected]
 

define particle list for entire domain

Defines all particles in domain; a_particles will then be distributed among the AMR levels

void ParticleAMRNS::defineRegridAMRSolver AMRSolver a_solver,
const Vector< DisjointBoxLayout > &  a_grids,
const Vector< ProblemDomain > &  a_domains,
const Vector< Real > &  a_amrDx,
const Vector< int > &  a_refRatios,
const int &  a_lBase
[protected]
 

defines (helmholtz) AMRSolver used in smoothing during regridding

void ParticleAMRNS::dumpParticlesASCII std::ostream &  os  )  const
 

Real ParticleAMRNS::Dx  )  const
 

access functions

void ParticleAMRNS::fillLambda LevelData< FArrayBox > &  a_lambda,
Real  a_time
const [protected]
 

resize temporary lambda for advection tracing computation and fill interior and boundary cells

void ParticleAMRNS::fillScalars LevelData< FArrayBox > &  a_scal,
Real  a_time,
const int  a_comp
const [protected]
 

fill interior and boundary cells (assume a_scal already sized)

void ParticleAMRNS::fillVelocity LevelData< FArrayBox > &  a_vel,
Real  a_time,
int  vel_comp,
int  dest_comp,
int  num_comp
[protected]
 

fill grown velocity field using piecewise-constant interp

void ParticleAMRNS::fillVelocity LevelData< FArrayBox > &  a_vel,
Real  a_time
[protected]
 

fill grown velocity field using piecewise-constant interp

ParticleAMRNS* ParticleAMRNS::finerNSPtr  )  const
 

encapsulates the finer-level-pointer casting

void ParticleAMRNS::finestLevel bool  a_finest_level  )  [protected]
 

set whether this is the finest level or not

bool ParticleAMRNS::finestLevel  )  const
 

is this this finest level?

PhysBCUtil* ParticleAMRNS::getPhysBCPtr  )  const
 

gets physical BC object

virtual void ParticleAMRNS::getPlotData LevelData< FArrayBox > &  a_plot_data  )  const [protected, virtual]
 

stuff plotfile data

virtual void ParticleAMRNS::initialData  )  [virtual]
 

initialize data

Implements AMRLevel.

virtual void ParticleAMRNS::initialGrid const Vector< Box > &  a_new_grids  )  [virtual]
 

initialize grids

Implements AMRLevel.

void ParticleAMRNS::initializeGlobalPressure  )  [protected]
 

manages pressure initialization after init or regrid

void ParticleAMRNS::initializeLevelPressure Real  a_currentTime,
Real  a_dtInit
[protected]
 

does level-based pressure initalization

void ParticleAMRNS::interpVelToParticles LevelData< BinFab< DragParticle > > &  a_particles,
const LevelData< FArrayBox > &  a_velocity
[protected]
 

update particle velocities by interpolating a_vel to particles

bool ParticleAMRNS::isEmpty  )  const
 

is this an empty level?

void ParticleAMRNS::lambda LevelData< FArrayBox > &  a_lambda,
Real  a_time
const
 

returns lambda at time t

void ParticleAMRNS::levelSetup const DisjointBoxLayout level_domain  )  [protected]
 

void ParticleAMRNS::limitSolverCoarsening bool  a_limitSolverCoarsening  ) 
 

set solver parameter

DisjointBoxLayout ParticleAMRNS::loadBalance const Vector< Box > &  a_grids  )  [protected]
 

virtual AMRLevel* ParticleAMRNS::makeAMRLevel AMRLevel a_coarser_level_ptr,
const ProblemDomain a_problem_domain,
int  a_level,
int  a_ref_ratio
const [virtual]
 

virtual AMRLevel* ParticleAMRNS::makeAMRLevel AMRLevel a_coarser_level_ptr,
const Box a_problem_domain,
int  a_level,
int  a_ref_ratio
const [virtual]
 

const LevelData<FArrayBox>& ParticleAMRNS::newLambda  )  const
 

LevelData<FArrayBox>& ParticleAMRNS::newLambda  ) 
 

const LevelData<FArrayBox>& ParticleAMRNS::newScal const int  a_comp  )  const
 

LevelData<FArrayBox>& ParticleAMRNS::newScal const int  a_comp  ) 
 

const LevelData<FArrayBox>& ParticleAMRNS::newVel  )  const
 

LevelData<FArrayBox>& ParticleAMRNS::newVel  ) 
 

Real ParticleAMRNS::Nu  )  const
 

viscous coefficient

int ParticleAMRNS::numParticles  )  const
 

number of particles on this AMR level

virtual int ParticleAMRNS::numPlotComps  )  const [protected, virtual]
 

number of components in plotfiles

const LevelData<FArrayBox>& ParticleAMRNS::oldLambda  )  const
 

LevelData<FArrayBox>& ParticleAMRNS::oldLambda  ) 
 

const LevelData<FArrayBox>& ParticleAMRNS::oldScal const int  a_comp  )  const
 

LevelData<FArrayBox>& ParticleAMRNS::oldScal const int  a_comp  ) 
 

const LevelData<FArrayBox>& ParticleAMRNS::oldVel  )  const
 

LevelData<FArrayBox>& ParticleAMRNS::oldVel  ) 
 

virtual void ParticleAMRNS::particleCFL Real  a_particle_cfl  )  [virtual]
 

set CFL number for computing timestep based on particles

Real ParticleAMRNS::particleKineticEnergy  )  const
 

const LevelData<BinFab<DragParticle> >& ParticleAMRNS::particles  )  const
 

returns particles (const version)

LevelData<BinFab<DragParticle> >& ParticleAMRNS::particles  ) 
 

returns particles

virtual void ParticleAMRNS::postInitialize  )  [virtual]
 

things do do after initialization

Implements AMRLevel.

virtual void ParticleAMRNS::postRegrid int  a_lBase  )  [virtual]
 

perform any post-regridding ops -- lBase is the finest unchanged level

Reimplemented from AMRLevel.

virtual void ParticleAMRNS::postTimeStep  )  [virtual]
 

things to do after a timestep

Implements AMRLevel.

void ParticleAMRNS::predictVelocities LevelData< FArrayBox > &  a_uDelU,
LevelData< FluxBox > &  a_advVel,
LevelData< FArrayBox > &  a_dragForce
[protected]
 

void ParticleAMRNS::readParameters  )  [protected]
 

read parameters from parmParse database

virtual void ParticleAMRNS::refinementThreshold Real  a_refine_threshold  )  [virtual]
 

set refinement threshold

virtual void ParticleAMRNS::regrid const Vector< Box > &  a_new_grids  )  [virtual]
 

regrid

Implements AMRLevel.

void ParticleAMRNS::resetStates const Real  a_time  )  [protected]
 

exchange old and new states, resets time to a_time

void ParticleAMRNS::scalars LevelData< FArrayBox > &  a_scalars,
const Real  a_time,
const int  a_comp
const
 

returns scalars at time t

void ParticleAMRNS::setPhysBC const PhysBCUtil a_BC  ) 
 

sets physical BC object

void ParticleAMRNS::smoothVelocityField int  lbase  )  [protected]
 

smooths composite velocity field after regridding

void ParticleAMRNS::swapOldAndNewStates  )  [protected]
 

move new-time state into old-time state -- also advances time by dt

virtual void ParticleAMRNS::tagCells IntVectSet a_tags  )  [virtual]
 

create tags

Implements AMRLevel.

virtual void ParticleAMRNS::tagCellsInit IntVectSet a_tags  )  [virtual]
 

create tags at initialization

Implements AMRLevel.

void ParticleAMRNS::updateParticleForces LevelData< BinFab< DragParticle > > &  a_particles  )  [protected]
 

update particle forces

void ParticleAMRNS::updateParticlePositions LevelData< FArrayBox > &  a_newDrag,
LevelData< BinFab< DragParticle > > &  a_particles,
const LevelData< FArrayBox > &  a_fluidVel,
bool  a_completeUpdate
[protected]
 

move particles and compute new-time drag forcing

This is a little complicated, in order to avoid creating a second particle list. We first move particles to a predicted position P^* (where P refers to the particle state). Then, we compute an estimate of the particle forcing term at the new time and place it in a_newDrag. Then, if completeUpdate is true, we go back and correct the new particle state to be formally second-order accurate.

void ParticleAMRNS::velocity LevelData< FArrayBox > &  a_vel,
Real  a_time
const
 

returns velocity at time t


Member Data Documentation

Real ParticleAMRNS::m_cfl [protected]
 

cfl number

CoarseAverage ParticleAMRNS::m_coarse_average [protected]
 

averaging from fine to coarse level for velocity

CoarseAverage ParticleAMRNS::m_coarse_average_lambda [protected]
 

averaging from fine to coarse level for freestream preservation qty

Vector<Tuple<Real, SpaceDim> > ParticleAMRNS::m_dataPoints [protected]
 

points where data is dumped out

Real ParticleAMRNS::m_dt_save [protected]
 

dt from most recent timestep

Real ParticleAMRNS::m_dx [protected]
 

grid spacing

bool ParticleAMRNS::m_finest_level [protected]
 

true if this is the finest extant level

LevelFluxRegister ParticleAMRNS::m_flux_register [protected]
 

momentum flux register

bool ParticleAMRNS::m_is_empty [protected]
 

true if this level is an empty level (no grids)

LevelFluxRegister ParticleAMRNS::m_lambda_flux_reg [protected]
 

freestream preservation qty flux register

LevelData<FArrayBox>* ParticleAMRNS::m_lambda_new_ptr [protected]
 

freestream preservation qty at new time

LevelData<FArrayBox>* ParticleAMRNS::m_lambda_old_ptr [protected]
 

freestream preservation qty at old time

int ParticleAMRNS::m_level_steps [protected]
 

number of steps on this level

bool ParticleAMRNS::m_limitSolverCoarsening [protected]
 

solver parameter -- passed on through when defining solvers

ParticleProjector ParticleAMRNS::m_particle_projector [protected]
 

particle drag force projector object

Real ParticleAMRNS::m_particleCFL [protected]
 

cfl number for particle updates

LevelData<BinFab<DragParticle> > ParticleAMRNS::m_particles [protected]
 

particles

PhysBCUtil* ParticleAMRNS::m_physBCPtr [protected]
 

physical BC object

CCProjector ParticleAMRNS::m_projection [protected]
 

cell-centered projection object

Real ParticleAMRNS::m_refine_threshold [protected]
 

refinement threshold for regridding

bool ParticleAMRNS::m_regrid_smoothing_done [protected]
 

internal indicator whether post-regrid smoothing stuff has been done

Vector<LevelFluxRegister*> ParticleAMRNS::m_scal_fluxreg_ptrs [protected]
 

conserved scalar flux register

Vector<LevelData<FArrayBox>*> ParticleAMRNS::m_scal_new [protected]
 

scalars at new time

Vector<LevelData<FArrayBox>*> ParticleAMRNS::m_scal_old [protected]
 

scalars at old time

PoissonOp ParticleAMRNS::m_scalarsPoissonOp [protected]
 

PoissonOp for advected scalars (for diffusion).

LevelData<FArrayBox>* ParticleAMRNS::m_vel_new_ptr [protected]
 

velocity vector at new time

LevelData<FArrayBox>* ParticleAMRNS::m_vel_old_ptr [protected]
 

velocity vector at old time

QuadCFInterp ParticleAMRNS::m_velCFInterp [protected]
 

does quadratic coarse-fine interpolation for velocities

VelocityPoissonOp ParticleAMRNS::m_velocityPoissonOp [protected]
 

PoissonOp for velocities.

LevelHelmholtzSolver ParticleAMRNS::m_viscous_solver [protected]
 

levelsolver for viscous solves

bool ParticleAMRNS::s_applyFreestreamCorrection [static, protected]
 

apply freestream preservation correction?

Real ParticleAMRNS::s_backgroundVel [static, protected]
 

background velocity (in the z-direction)

Real ParticleAMRNS::s_bogus_value [static, protected]
 

debugging option -- value to which arrays are initially set

bool ParticleAMRNS::s_compute_scal_err [static, protected]
 

compute error and report norms after composite timesteps

Real ParticleAMRNS::s_delta_epsilon [static, protected]
 

particle delta function parameter

bool ParticleAMRNS::s_do_particles [static, protected]
 

if true, do particle stuff, if false, don't

Vector<Real> ParticleAMRNS::s_domLength [static, protected]
 

size of physical domain

bool ParticleAMRNS::s_dump_smoothing_plots [static, protected]
 

dump before-and-after plotfiles for post-regrid smoothing (for debug)

bool ParticleAMRNS::s_implicit_reflux [static, protected]
 

do implicit refluxing for velocity?

bool ParticleAMRNS::s_implicit_scal_reflux [static, protected]
 

do implicit refluxing for advected/diffused scalars?

Real ParticleAMRNS::s_init_shrink [static, protected]
 

factor to shrink initial timestep

bool ParticleAMRNS::s_init_vel_from_vorticity [static, protected]
 

initialize velocity field from vorticity field

string ParticleAMRNS::s_initialGridFile [static, protected]
 

gridfile where initial grids are spec'd, if relevant

bool ParticleAMRNS::s_initialize_pressures [static, protected]
 

should we initialize pressures after grid generation & regridding?

Real ParticleAMRNS::s_max_dt [static, protected]
 

maximum allowable timestep

const int ParticleAMRNS::s_max_scal_comps = SpaceDim [static, protected]
 

maximum possible number of scalar components

Real ParticleAMRNS::s_nu [static, protected]
 

viscosity

int ParticleAMRNS::s_num_init_passes [static, protected]
 

number of times to iterate when initializing pressures

int ParticleAMRNS::s_num_scal_comps [static, protected]
 

number of scalars

const int ParticleAMRNS::s_num_vel_comps = SpaceDim [static, protected]
 

number of components of m_state (vel comps + lambda)

int ParticleAMRNS::s_order_vel_interp [static, protected]
 

order of velocity interpolation for particles

RealVect ParticleAMRNS::s_particle_body_force [static, protected]
 

particle body force (gravity)

Real ParticleAMRNS::s_particle_drag_coeff [static, protected]
 

particle drag coefficient

string ParticleAMRNS::s_particle_file [static, protected]
 

file containing particle info

bool ParticleAMRNS::s_ppInit [static, protected]
 

have we initialized static variables from parmParse?

Real ParticleAMRNS::s_prescribedDt [static, protected]
 

bool ParticleAMRNS::s_project_initial_vel [static, protected]
 

should we project initial velocity field (default is yes)

bool ParticleAMRNS::s_read_particles_from_file [static, protected]
 

read particle info from file?

bool ParticleAMRNS::s_reflux_momentum [static, protected]
 

should we reflux momentum?

bool ParticleAMRNS::s_reflux_normal_momentum [static, protected]
 

if we're refluxing momentum, should we reflux the normal component?

bool ParticleAMRNS::s_reflux_scal [static, protected]
 

should we reflux scalars?

Real ParticleAMRNS::s_regrid_smoothing_coeff [static, protected]
 

antidiffusive smoothing coefficient (multiplies nu*dtCrse)

Vector<Real> ParticleAMRNS::s_scal_coeffs [static, protected]
 

diffusion coefficients for scalar advection-diffusion

const char* ParticleAMRNS::s_scal_names[s_max_scal_comps] [static, protected]
 

names of scalars

bool ParticleAMRNS::s_set_bogus_values [static, protected]
 

debugging option -- initially set arrays to s_bogus_value

bool ParticleAMRNS::s_smooth_after_regrid [static, protected]
 

do antidiffusive smoothing after regridding?

bool ParticleAMRNS::s_specifyInitialGrids [static, protected]
 

do we specify initial grids in a gridFile?

bool ParticleAMRNS::s_tag_vorticity [static, protected]
 

tag on vorticity?

int ParticleAMRNS::s_tags_grow [static, protected]
 

amount by which to grow tagged regions before calling meshrefine

bool ParticleAMRNS::s_use_image_particles [static, protected]
 

if true, use image particles in particle projector near phys boundaries

const char* ParticleAMRNS::s_vel_names[s_num_vel_comps] [static, protected]
 

names of components

int ParticleAMRNS::s_viscous_num_smooth_down [static, protected]
 

multigrid parameter for viscous solves

int ParticleAMRNS::s_viscous_num_smooth_up [static, protected]
 

multigrid parameter for viscous solves

Real ParticleAMRNS::s_viscous_solver_tol [static, protected]
 

viscous solver tolerance (defaults to 1e-7)

int ParticleAMRNS::s_viscous_solver_type [static, protected]
 

viscous solver type: 0 = backwards euler, 1=Crank-Nicolson, 2=TGA

Real ParticleAMRNS::s_vort_factor [static, protected]
 

factor for tagging -- tag if vorticity greater than factor*max_vort

bool ParticleAMRNS::s_write_divergence [static, protected]
 

include divergence in plotfile?

bool ParticleAMRNS::s_write_dScalar_dt [static, protected]
 

include d(scalars)/dt in plotfile?

bool ParticleAMRNS::s_write_error [static, protected]
 

include error (if there is an exact solution) in plotfile?

bool ParticleAMRNS::s_write_grad_eLambda [static, protected]
 

include freestream preservation correcion in plotfile?

bool ParticleAMRNS::s_write_grids [static, protected]
 

write out grids and processor assignments

bool ParticleAMRNS::s_write_initial_solve_rhs [static, protected]
 

dump out initial vorticity->stream function RHS's into plotfiles

bool ParticleAMRNS::s_write_lambda [static, protected]
 

include lambda in plotfile?

bool ParticleAMRNS::s_write_particles [static, protected]
 

-- I/O options -- include particles in plotfile?

bool ParticleAMRNS::s_write_proc_ids [static, protected]
 

include processor distribution visualization in plotfile?

bool ParticleAMRNS::s_write_scalars [static, protected]
 

include scalars in plotfile?

bool ParticleAMRNS::s_write_strains [static, protected]
 

include strain-rates in plotfile? (2D only)

bool ParticleAMRNS::s_write_time_derivatives [static, protected]
 

include du/dt and d(lambda)/dt in plotfile?

bool ParticleAMRNS::s_write_vorticity [static, protected]
 

include vorticity in plotfile?


The documentation for this class was generated from the following file:
Generated on Wed Jun 2 13:58:53 2004 for Chombo&INSwithParticles by doxygen 1.3.2