CG  Version 25
Public Types | Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | List of all members
MovingGrids Class Reference

#include <MovingGrids.h>

Collaboration diagram for MovingGrids:
Collaboration graph
[legend]

Public Types

enum  MovingGridOption {
  notMoving, rotate, shift, oscillate,
  scale, matrixMotion, rigidBody, deformingBody,
  userDefinedMovingGrid, numberOfMovingGridOptions
}
 

Public Member Functions

 MovingGrids (Parameters &parameters)
 This class is used to move grids.
 
 MovingGrids (const MovingGrids &mg)
 
virtual ~MovingGrids ()
 
virtual int correctGrids (const real t1, const real t2, GridFunction &cgf1, GridFunction &cgf2)
 Corrector step for moving grids.
 
int debug () const
 
int detectCollisions (GridFunction &cgf1)
 Detect collisions.
 
virtual int getBoundaryAcceleration (MappedGrid &c, realSerialArray &gtt, int grid, real t0, int option)
 
bool getCorrectionHasConverged ()
 Return true if the correction steps for moving grids have converged. This is usually only an issue for "light" bodies.
 
virtual int getGridVelocity (GridFunction &gf0, const real &tGV)
 
Integrate * getIntegrate () const
 Return a pointer to the Integrate object.
 
const RealArray & getMoveParameters () const
 
real getMaximumRelativeCorrection ()
 Return the maximum relative change in the moving grid correction scheme. This is usually only an issue for "light" bodies.
 
int getNumberOfDeformingBodies () const
 Return the number of deforming bodies.
 
DeformingBodyMotiongetDeformingBody (const int bodyNumber)
 Return the object that describes a deforming body.
 
int getNumberOfMatrixMotionBodies () const
 Return the number of matrix motion bodies.
 
MatrixMotiongetMatrixMotionBody (const int bodyNumber)
 
int getNumberOfRigidBodies () const
 Return the number of rigid bodies.
 
RigidBodyMotiongetRigidBody (const int bodyNumber)
 Return the object that describes a rigid body.
 
real getTimeStepForRigidBodies () const
 
int getUserDefinedBoundaryAcceleration (MappedGrid &mg, realSerialArray &gtt, int grid, real t0, int option, const int side, const int axis)
 
int getUserDefinedGridVelocity (GridFunction &gf0, const real &t0, const int grid)
 
int get (const GenericDataBase &dir, const aString &name)
 
virtual int gridAccelerationBC (const int &grid, const real &t0, MappedGrid &c, realMappedGridFunction &u, realMappedGridFunction &f, realMappedGridFunction &gridVelocity, realSerialArray &normal, const Index &I1, const Index &I2, const Index &I3, const Index &I1g, const Index &I2g, const Index &I3g)
 Add the grid acceleration in the normal direction, n.x_tt, to the function f f is normally the function that holds the rhs for the pressure eqn.
 
bool gridIsMoving (int grid) const
 
virtual int moveDeformingBodies (const real &t1, const real &t2, const real &t3, const real &dt0, GridFunction &cgf1, GridFunction &cgf2, GridFunction &cgf3)
 
virtual int moveGrids (const real &t1, const real &t2, const real &t3, const real &dt0, GridFunction &cgf1, GridFunction &cgf2, GridFunction &cgf3)
 Move grids to the new time.
 
MovingGridOption movingGridOption (int grid) const
 
aString movingGridOptionName (MovingGridOption option) const
 
bool isMovingGridProblem () const
 
int put (GenericDataBase &dir, const aString &name) const
 
virtual int rigidBodyMotion (const real &t1, const real &t2, const real &t3, const real &dt0, GridFunction &cgf1, GridFunction &cgf2, GridFunction &cgf3)
 
virtual int saveToShowFile () const
 Save moving grid info the the show file. We save a "sequence" for each rigid body.
 
int setIsMovingGridProblem (bool trueOrFalse=TRUE)
 
virtual int update (CompositeGrid &cg, GenericGraphicsInterface &gi)
 Define moving grid properties interactively.
 
int updateToMatchGrid (CompositeGrid &cg)
 
int userDefinedMotion (const real &t1, const real &t2, const real &t3, const real &dt0, GridFunction &cgf1, GridFunction &cgf2, GridFunction &cgf3)
 
int userDefinedTransformMotion (const real &t1, const real &t2, const real &t3, const real &dt0, GridFunction &cgf1, GridFunction &cgf2, GridFunction &cgf3, const int grid)
 
int userDefinedGridAccelerationBC (const int &grid, const real &t0, MappedGrid &c, realMappedGridFunction &u, realMappedGridFunction &f, realMappedGridFunction &gridVelocity, realSerialArray &normal, const Index &I1, const Index &I2, const Index &I3, const Index &I1g, const Index &I2g, const Index &I3g)
 
int updateUserDefinedMotion (CompositeGrid &cg, GenericGraphicsInterface &gi)
 

Static Public Attributes

static int debug0 =0
 

Protected Member Functions

int initialize ()
 
int getRamp (real t, real rampInterval, real &ramp, real &rampSpeed, real &rampAcceleration)
 Return ramp value and derivatives.
 

Protected Attributes

Parametersparameters
 
IntegerArray gridsToMove
 
bool isInitialized
 
bool movingGridProblem
 
IntegerArray moveOption
 
IntegerArray movingGrid
 
RealArray moveParameters
 
int numberOfMatrixMotionBodies
 
MatrixMotion ** matrixMotionBody
 
int numberOfRigidBodies
 
RigidBodyMotion ** body
 
Integrate * integrate
 
bool useHybridGridsForSurfaceIntegrals
 
int rigidBodyInfoCount
 
RealArray rigidBodyInfo
 
RealArray rigidBodyInfoTime
 
int numberOfRigidBodyInfoNames
 
aString * rigidBodyInfoName
 
bool limitForces
 
real maximumAllowableForce
 
real maximumAllowableTorque
 
bool correctionHasConverged
 
real maximumRelativeCorrection
 
int numberOfDeformingBodies
 
DeformingBodyMotion ** deformingBodyList
 

Member Enumeration Documentation

Enumerator
notMoving 
rotate 
shift 
oscillate 
scale 
matrixMotion 
rigidBody 
deformingBody 
userDefinedMovingGrid 
numberOfMovingGridOptions 

Constructor & Destructor Documentation

MovingGrids::MovingGrids ( Parameters parameters_)
MovingGrids::MovingGrids ( const MovingGrids mg)
MovingGrids::~MovingGrids ( )
virtual

Member Function Documentation

int MovingGrids::correctGrids ( const real  t1,
const real  t2,
GridFunction cgf1,
GridFunction cgf2 
)
virtual

Corrector step for moving grids.

This function is called at the corrector step to update the moving grids. For example, in a predictor corrector type algorithm we may want to correct the forces and torques on bodies since the solution can depend on these (For INS the pressure BC depends on the acceleration on the boundary ).

Parameters
t1,cgf1(input) : solution at the old time
t2,cgf2(input) : solution at the new time (these are valid values)

References DeformingBodyMotion::correct(), correctionHasConverged, deformingBodyList, maximumRelativeCorrection, numberOfDeformingBodies, and rigidBodyMotion().

int MovingGrids::debug ( ) const
inline
int MovingGrids::detectCollisions ( GridFunction cgf1)
int MovingGrids::get ( const GenericDataBase &  dir,
const aString &  name 
)
int MovingGrids::getBoundaryAcceleration ( MappedGrid &  c,
realSerialArray &  gtt,
int  grid,
real  t0,
int  option 
)
virtual
bool MovingGrids::getCorrectionHasConverged ( )

Return true if the correction steps for moving grids have converged. This is usually only an issue for "light" bodies.

References correctionHasConverged.

DeformingBodyMotion & MovingGrids::getDeformingBody ( const int  bodyNumber)

Return the object that describes a deforming body.

Parameters
bodyNumber(input) : a body number starting from 0 and ending at getNumberOfDeformingBodies()-1.

References deformingBodyList, numberOfDeformingBodies, OV_ABORT(), and printF().

Referenced by DomainSolver::getMovingGridOption().

int MovingGrids::getGridVelocity ( GridFunction gf0,
const real &  tGV 
)
virtual
Integrate * MovingGrids::getIntegrate ( ) const

Return a pointer to the Integrate object.

References integrate.

Referenced by rigidBodyBoundaryProjection().

MatrixMotion & MovingGrids::getMatrixMotionBody ( const int  bodyNumber)

/brief Return the object that describes a matrix motion body. /param bodyNumber (input) : a body number starting from 0 and ending at getNumberOfMatrixMotionBodies()-1.

References matrixMotionBody, numberOfMatrixMotionBodies, OV_ABORT(), and printF().

Referenced by DomainSolver::getMovingGridOption().

real MovingGrids::getMaximumRelativeCorrection ( )

Return the maximum relative change in the moving grid correction scheme. This is usually only an issue for "light" bodies.

References maximumRelativeCorrection.

Referenced by rigidBodyMotion().

const RealArray & MovingGrids::getMoveParameters ( ) const

References moveParameters.

int MovingGrids::getNumberOfDeformingBodies ( ) const

Return the number of deforming bodies.

References numberOfDeformingBodies.

Referenced by DomainSolver::buildMovingGridOptionsDialog(), and DomainSolver::getMovingGridOption().

int MovingGrids::getNumberOfMatrixMotionBodies ( ) const

Return the number of matrix motion bodies.

References numberOfMatrixMotionBodies.

Referenced by DomainSolver::buildMovingGridOptionsDialog(), and DomainSolver::getMovingGridOption().

int MovingGrids::getNumberOfRigidBodies ( ) const
int MovingGrids::getRamp ( real  t,
real  rampInterval,
real &  ramp,
real &  rampSpeed,
real &  rampAcceleration 
)
protected

Return ramp value and derivatives.

 cubic ramp [0,1] : ramp(t)=t*t*(3-2*t) = 3t^2 - 2t^3
                     r'=6t*(1.-t)
 ramp(1)=1

Referenced by getBoundaryAcceleration(), getGridVelocity(), gridAccelerationBC(), and moveGrids().

RigidBodyMotion & MovingGrids::getRigidBody ( const int  bodyNumber)

Return the object that describes a rigid body.

Parameters
bodyNumber(input) : a body number starting from 0 and ending at getNumberOfRigidBodies()-1.

References body, numberOfRigidBodies, OV_ABORT(), and printF().

Referenced by DomainSolver::getMovingGridOption(), and rigidBodyBoundaryProjection().

real MovingGrids::getTimeStepForRigidBodies ( ) const
int MovingGrids::getUserDefinedBoundaryAcceleration ( MappedGrid &  mg,
realSerialArray &  gtt,
int  grid,
real  t0,
int  option,
const int  side,
const int  axis 
)
int MovingGrids::getUserDefinedGridVelocity ( GridFunction gf0,
const real &  t0,
const int  grid 
)
int MovingGrids::gridAccelerationBC ( const int &  grid,
const real &  t0,
MappedGrid &  c,
realMappedGridFunction &  u,
realMappedGridFunction &  f,
realMappedGridFunction &  gridVelocity,
realSerialArray &  normal,
const Index &  J1,
const Index &  J2,
const Index &  J3,
const Index &  J1g,
const Index &  J2g,
const Index &  J3g 
)
virtual
bool MovingGrids::gridIsMoving ( int  grid) const

References grid.

int MovingGrids::initialize ( )
protected
bool MovingGrids::isMovingGridProblem ( ) const

References movingGridProblem.

int MovingGrids::moveDeformingBodies ( const real &  t1,
const real &  t2,
const real &  t3,
const real &  dt0,
GridFunction cgf1,
GridFunction cgf2,
GridFunction cgf3 
)
virtual
int MovingGrids::moveGrids ( const real &  t1,
const real &  t2,
const real &  t3,
const real &  dt0,
GridFunction cgf1,
GridFunction cgf2,
GridFunction cgf3 
)
virtual

Move grids to the new time.

The basic idea is to advance the solution from grid g1 at time t1 to grid g3 at time t3 using the grid velocity from grid g2 at time t2:

  • g3(t3) <- g1(t1) + (t3-t1)*d(g2(t2))/dt Steps:
  • detect collisions
  • compute forces and moments on rigid and deforming bodies
  • compute grid velocity of g2
  • move grids
  • compute grid velocity of g3
Parameters
t1,cgf1(input) : grid and solution at time t1
t2,cgf2(input) : grid and solution at time t2. The grid velocity, cgf2.gridVelocity, is computed and used from this time.
t3(input) : new time
cgf3(output): holds new grid at time t3. The grid velocity, cgf3.gridVelocity, is computed at the. end of the step. NOTE: cgf3 must be different from cgf1.

References amplitude, assert(), axis, body, GridFunction::cg, cg, d, Parameters::dbase, debug(), deformingBody, deformingBodyList, detectCollisions(), display(), getGridVelocity(), RigidBodyMotion::getInitialConditions(), MatrixMotion::getMotion(), RigidBodyMotion::getPosition(), getRamp(), RigidBodyMotion::getRotationMatrix(), grid, isInitialized, mapPointer, matrixMotion, matrixMotionBody, moveDeformingBodies(), moveOption, moveParameters, movingGrid, movingGridProblem, notMoving, numberOfDeformingBodies, numberOfMatrixMotionBodies, numberOfRigidBodies, GridFunction::numberOfTransformMappings, omega, oscillate, OV_ABORT(), parameters, printF(), r, Parameters::readInitialConditionFromRestartFile, Parameters::readInitialConditionFromShowFile, DeformingBodyMotion::regenerateComponentGrids(), rigidBody, rigidBodyMotion(), rotate, scale, shift, sx, sy, sz, GridFunction::t, GridFunction::transform, userDefinedMotion(), userDefinedMovingGrid, userDefinedTransformMotion(), and x3.

MovingGrids::MovingGridOption MovingGrids::movingGridOption ( int  grid) const
aString MovingGrids::movingGridOptionName ( MovingGridOption  option) const
int MovingGrids::put ( GenericDataBase &  dir,
const aString &  name 
) const

References dir.

int MovingGrids::rigidBodyMotion ( const real &  t1,
const real &  t2,
const real &  t3,
const real &  dt0,
GridFunction cgf1,
GridFunction cgf2,
GridFunction cgf3 
)
virtual
int MovingGrids::saveToShowFile ( ) const
virtual

Save moving grid info the the show file. We save a "sequence" for each rigid body.

References all, Parameters::dbase, debug(), N(), numberOfRigidBodies, parameters, printF(), rigidBodyInfo, rigidBodyInfoCount, rigidBodyInfoName, and rigidBodyInfoTime.

Referenced by DomainSolver::advance(), Cgsm::saveShow(), DomainSolver::saveShow(), and Cgmp::solve().

int MovingGrids::setIsMovingGridProblem ( bool  trueOrFalse = TRUE)

References movingGridProblem.

int MovingGrids::update ( CompositeGrid &  cg,
GenericGraphicsInterface &  gi 
)
virtual

Define moving grid properties interactively.

This routine is used to define how bodies move. A moving body may be constructed with multiple grids (e.g. a sphere may be covered with 2 or three patches). There are two main steps to define a motion of a body

  • define the motion of a body
  • define which grids belong to the body The types of motions are
  • "matrix motions" define predetermined rigid motions and are given by x(t) = R(t) x(0) + g(t) where R(t) is a "rotation" (or scaling) matrix and g(t) is a translation.
  • "rigid body motions" define motions of rigid bodies whose position depend on the forces and torques imposed by the fluid.
  • "deforming body motions" define bodies that can move and deform (either specified deformations or deformations that depend on the fluid).
  • "user defined motion" are motions defined by the user in the file UserDefinedMotion.C

References all, amplitude, assert(), axis, body, boundaryCondition(), c, Parameters::dbase, debug(), debug0, DeformingBodyMotion::defineBody(), deformingBody, deformingBodyList, dir, e, f, RigidBodyMotion::getDensity(), grid, includeGhost, DeformingBodyMotion::initialize(), initialize(), DeformingBodyMotion::initializePast(), integrate, isInitialized, limitForces, matrixMotion, matrixMotionBody, maximumAllowableForce, maximumAllowableTorque, moveOption, moveParameters, movingGrid, movingGridOption(), notMoving, numberOfDeformingBodies, numberOfMatrixMotionBodies, numberOfRigidBodies, numberOfRigidBodyInfoNames, ok, omega, oscillate, parameters, printF(), rigidBody, rigidBodyInfo, rigidBodyInfoTime, rotate, scale, RigidBodyMotion::setInitialCentreOfMass(), RigidBodyMotion::setInitialConditions(), RigidBodyMotion::setMass(), shift, side, sx, sy, sz, tol, u, MatrixMotion::update(), DeformingBodyMotion::update(), RigidBodyMotion::update(), updateUserDefinedMotion(), useHybridGridsForSurfaceIntegrals, userDefinedMovingGrid, and x.

int MovingGrids::updateToMatchGrid ( CompositeGrid &  cg)

References all, grid, moveOption, moveParameters, movingGrid, and R.

int MovingGrids::updateUserDefinedMotion ( CompositeGrid &  cg,
GenericGraphicsInterface &  gi 
)

References a, linearMotion, omega, rampMotion, sinusoidalMotion, and ta.

Referenced by update().

int MovingGrids::userDefinedGridAccelerationBC ( const int &  grid,
const real &  t0,
MappedGrid &  c,
realMappedGridFunction &  u,
realMappedGridFunction &  f,
realMappedGridFunction &  gridVelocity,
realSerialArray &  normal,
const Index &  I1,
const Index &  I2,
const Index &  I3,
const Index &  I1g,
const Index &  I2g,
const Index &  I3g 
)
int MovingGrids::userDefinedMotion ( const real &  t1,
const real &  t2,
const real &  t3,
const real &  dt0,
GridFunction cgf1,
GridFunction cgf2,
GridFunction cgf3 
)

Referenced by moveGrids().

int MovingGrids::userDefinedTransformMotion ( const real &  t1,
const real &  t2,
const real &  t3,
const real &  dt0,
GridFunction cgf1,
GridFunction cgf2,
GridFunction cgf3,
const int  grid 
)

Member Data Documentation

RigidBodyMotion** MovingGrids::body
protected
bool MovingGrids::correctionHasConverged
protected
int MovingGrids::debug0 =0
static

Referenced by debug(), and update().

DeformingBodyMotion** MovingGrids::deformingBodyList
protected
IntegerArray MovingGrids::gridsToMove
protected

Referenced by get().

Integrate* MovingGrids::integrate
protected
bool MovingGrids::isInitialized
protected
bool MovingGrids::limitForces
protected
MatrixMotion** MovingGrids::matrixMotionBody
protected
real MovingGrids::maximumAllowableForce
protected
real MovingGrids::maximumAllowableTorque
protected
real MovingGrids::maximumRelativeCorrection
protected
IntegerArray MovingGrids::moveOption
protected
RealArray MovingGrids::moveParameters
protected
IntegerArray MovingGrids::movingGrid
protected
bool MovingGrids::movingGridProblem
protected
int MovingGrids::numberOfDeformingBodies
protected
int MovingGrids::numberOfMatrixMotionBodies
protected
int MovingGrids::numberOfRigidBodies
protected
int MovingGrids::numberOfRigidBodyInfoNames
protected
Parameters& MovingGrids::parameters
protected
RealArray MovingGrids::rigidBodyInfo
protected
int MovingGrids::rigidBodyInfoCount
protected
aString* MovingGrids::rigidBodyInfoName
protected
RealArray MovingGrids::rigidBodyInfoTime
protected
bool MovingGrids::useHybridGridsForSurfaceIntegrals
protected

Referenced by MovingGrids(), and update().


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