Go to the documentation of this file.
46 displacementComponentLaplacianFvMotionSolver,
54 Foam::displacementComponentLaplacianFvMotionSolver::
55 displacementComponentLaplacianFvMotionSolver
67 "cellDisplacement" + cmptName_,
68 mesh.time().timeName(),
70 IOobject::READ_IF_PRESENT,
75 cellMotionBoundaryTypes<scalar>(pointDisplacement_.boundaryField())
77 pointLocation_(
nullptr),
80 coeffDict().
found(
"interpolation")
90 coeffDict().
found(
"frozenPointsZone")
91 ? fvMesh_.pointZones().findZoneID
93 coeffDict().get<word>(
"frozenPointsZone")
98 if (coeffDict().getOrDefault(
"applyPointLocation",
true))
107 fvMesh_.time().timeName(),
118 Info<<
"displacementComponentLaplacianFvMotionSolver :"
119 <<
" Read pointVectorField "
120 << pointLocation_().name()
121 <<
" to be used for boundary conditions on points."
123 <<
"Boundary conditions:"
124 << pointLocation_().boundaryField().types() <<
endl;
142 interpolationPtr_->interpolate
152 Info<<
"displacementComponentLaplacianFvMotionSolver : applying "
153 <<
" boundary conditions on " << pointLocation_().name()
154 <<
" to new point location."
160 pointLocation_().primitiveFieldRef() = fvMesh_.points();
162 pointLocation_().primitiveFieldRef().replace
165 points0_ + pointDisplacement_.primitiveField()
168 pointLocation_().correctBoundaryConditions();
171 if (frozenPointsZone_ != -1)
173 const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
177 label pointi = pz[i];
179 pointLocation_()[pointi][cmpt_] = points0_[pointi];
195 points0_ + pointDisplacement_.primitiveField()
199 if (frozenPointsZone_ != -1)
201 const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
205 label pointi = pz[i];
207 curPoints[pointi][cmpt_] = points0_[pointi];
211 twoDCorrectPoints(curPoints);
222 movePoints(fvMesh_.points());
224 diffusivityPtr_->correct();
225 pointDisplacement_.boundaryFieldRef().updateCoeffs();
236 *diffusivityPtr_->operator()(),
238 "laplacian(diffusivity,cellDisplacement)"
259 diffusivityPtr_.reset(
nullptr);
263 coeffDict().
lookup(
"diffusivity")
int debug
Static debugging option.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
vectorField pointField
pointField is a vectorField.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to database if not present.
SolverPerformance< Type > solveSegregatedOrCoupled(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Virtual base class for displacement motion solver.
virtual void updateMesh(const mapPolyMesh &)
Update topology.
Ostream & endl(Ostream &os)
Add newline and flush stream.
conserve primitiveFieldRef()+
Mesh consisting of general polyhedral cells.
~displacementComponentLaplacianFvMotionSolver()
Destructor.
#define forAll(list, i)
Loop across all elements in list.
fvScalarMatrix TEqn(fvm::ddt(T)+fvm::div(phi, T) - fvm::laplacian(alphaEff, T)==radiation->ST(rhoCpRef, T)+fvOptions(T))
Base class for fvMesh based motionSolvers.
virtual void updateMesh(const mapPolyMesh &)
Update local data for topology changes.
messageStream Info
Information stream (uses stdout - output is on the master only)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Calculate the matrix for the laplacian of the field.
static autoPtr< motionDiffusivity > New(const fvMesh &mesh, Istream &mdData)
Select null constructed.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
void replace(const direction, const UList< cmptType > &)
Replace a component field of the field.
Lookup type of boundary radiation properties.
const dimensionSet dimViscosity
Macros for easy insertion into run-time selection tables.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
const dictionary & solverDict() const
Return the solver dictionary taking into account finalIteration.
virtual void solve()
Solve for motion.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
defineTypeNameAndDebug(combustionModel, 0)