Go to the documentation of this file.
49 solidBodyDisplacementLaplacianFvMotionSolver,
55 displacementMotionSolver,
56 solidBodyDisplacementLaplacianFvMotionSolver,
64 Foam::solidBodyDisplacementLaplacianFvMotionSolver::
65 solidBodyDisplacementLaplacianFvMotionSolver
79 mesh.time().timeName(),
81 IOobject::READ_IF_PRESENT,
86 cellMotionBoundaryTypes<vector>(pointDisplacement_.boundaryField())
88 pointLocation_(
nullptr),
91 coeffDict().
found(
"interpolation")
101 coeffDict().
found(
"frozenPointsZone")
102 ? fvMesh_.pointZones().findZoneID
104 coeffDict().get<word>(
"frozenPointsZone")
112 fvMesh_.time().timeName(),
120 Info<<
"solidBodyDisplacementLaplacianFvMotionSolver:" <<
nl
121 <<
" diffusivity : " << diffusivityPtr_().type() <<
nl
122 <<
" frozenPoints zone : " << frozenPointsZone_ <<
endl;
139 Info<<
"solidBodyDisplacementLaplacianFvMotionSolver :"
140 <<
" Read pointVectorField "
142 <<
" to be used for boundary conditions on points."
144 <<
"Boundary conditions:"
145 << pointLocation_().boundaryField().types() <<
endl;
151 Foam::solidBodyDisplacementLaplacianFvMotionSolver::
152 solidBodyDisplacementLaplacianFvMotionSolver
168 mesh.time().timeName(),
170 IOobject::READ_IF_PRESENT,
175 cellMotionBoundaryTypes<vector>(pointDisplacement_.boundaryField())
177 pointLocation_(
nullptr),
180 coeffDict().
found(
"interpolation")
190 coeffDict().
found(
"frozenPointsZone")
191 ? fvMesh_.pointZones().findZoneID
193 coeffDict().get<word>(
"frozenPointsZone")
201 fvMesh_.time().timeName(),
209 Info<<
"solidBodyDisplacementLaplacianFvMotionSolver:" <<
nl
210 <<
" diffusivity : " << diffusivityPtr_().type() <<
nl
211 <<
" frozenPoints zone : " << frozenPointsZone_ <<
endl;
228 Info<<
"solidBodyDisplacementLaplacianFvMotionSolver :"
229 <<
" Read pointVectorField "
231 <<
" to be used for boundary conditions on points."
233 <<
"Boundary conditions:"
234 << pointLocation_().boundaryField().types() <<
endl;
252 if (!diffusivityPtr_)
257 coeffDict().
lookup(
"diffusivity")
261 return *diffusivityPtr_;
268 interpolationPtr_->interpolate
284 Info<<
"solidBodyDisplacementLaplacianFvMotionSolver : applying "
285 <<
" boundary conditions on " << pointLocation_().name()
286 <<
" to new point location."
290 pointLocation_().primitiveFieldRef() =
292 + pointDisplacement_.internalField();
294 pointLocation_().correctBoundaryConditions();
297 if (frozenPointsZone_ != -1)
299 const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
303 pointLocation_()[pz[i]] = newPoints[pz[i]];
315 newPoints + pointDisplacement_.primitiveField()
320 if (frozenPointsZone_ != -1)
322 const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
326 curPoints[pz[i]] = newPoints[pz[i]];
330 twoDCorrectPoints(curPoints);
341 movePoints(fvMesh_.points());
343 diffusivity().correct();
344 pointDisplacement_.boundaryFieldRef().updateCoeffs();
353 *diffusivity().
operator()(),
355 "laplacian(diffusivity,cellDisplacement)"
376 diffusivityPtr_.clear();
int debug
Static debugging option.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
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 primitive field of type <T> with automated input and output.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
Abstract base class for cell-centre mesh motion diffusivity.
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 void updateMesh(const mapPolyMesh &)
Update topology.
Ostream & endl(Ostream &os)
Add newline and flush stream.
conserve primitiveFieldRef()+
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
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.
motionDiffusivity & diffusivity()
Return reference to the diffusivity field.
messageStream Info
Information stream (stdout output on master, null elsewhere)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
~solidBodyDisplacementLaplacianFvMotionSolver()
Destructor.
Calculate the matrix for the laplacian of the field.
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
static autoPtr< motionDiffusivity > New(const fvMesh &mesh, Istream &mdData)
Select null constructed.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Lookup type of boundary radiation properties.
const dimensionSet dimViscosity
Virtual base class for displacement motion solver.
void transformPoints(vectorField &, const septernion &, const vectorField &)
Transform given vectorField of coordinates with the given septernion.
Macros for easy insertion into run-time selection tables.
virtual void solve()
Solve for motion.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
pointField points0(pointIOField(IOobject("points", mesh.time().constant(), polyMesh::meshSubDir, mesh, IOobject::MUST_READ, IOobject::NO_WRITE, false)))
virtual void updateMesh(const mapPolyMesh &)
Update local data for topology changes.
const dictionary & solverDict() const
Return the solver dictionary taking into account finalIteration.
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)