Go to the documentation of this file.
48 displacementLaplacianFvMotionSolver,
54 displacementMotionSolver,
55 displacementLaplacianFvMotionSolver,
63 Foam::displacementLaplacianFvMotionSolver::displacementLaplacianFvMotionSolver
76 mesh.time().timeName(),
78 IOobject::READ_IF_PRESENT,
83 cellMotionBoundaryTypes<vector>(pointDisplacement_.boundaryField())
85 pointLocation_(
nullptr),
88 coeffDict().
found(
"interpolation")
98 coeffDict().
found(
"frozenPointsZone")
99 ? fvMesh_.pointZones().findZoneID
101 coeffDict().get<word>(
"frozenPointsZone")
109 fvMesh_.time().timeName(),
117 Info<<
"displacementLaplacianFvMotionSolver:" <<
nl
118 <<
" diffusivity : " << diffusivityPtr_().type() <<
nl
119 <<
" frozenPoints zone : " << frozenPointsZone_ <<
endl;
136 Info<<
"displacementLaplacianFvMotionSolver :"
137 <<
" Read pointVectorField "
139 <<
" to be used for boundary conditions on points."
141 <<
"Boundary conditions:"
142 << pointLocation_().boundaryField().types() <<
endl;
148 Foam::displacementLaplacianFvMotionSolver::
149 displacementLaplacianFvMotionSolver
164 mesh.time().timeName(),
166 IOobject::READ_IF_PRESENT,
171 cellMotionBoundaryTypes<vector>(pointDisplacement_.boundaryField())
173 pointLocation_(
nullptr),
176 coeffDict().
found(
"interpolation")
186 coeffDict().
found(
"frozenPointsZone")
187 ? fvMesh_.pointZones().findZoneID
189 coeffDict().get<word>(
"frozenPointsZone")
197 fvMesh_.time().timeName(),
205 Info<<
"displacementLaplacianFvMotionSolver:" <<
nl
206 <<
" diffusivity : " << diffusivityPtr_().type() <<
nl
207 <<
" frozenPoints zone : " << frozenPointsZone_ <<
endl;
224 Info<<
"displacementLaplacianFvMotionSolver :"
225 <<
" Read pointVectorField "
227 <<
" to be used for boundary conditions on points."
229 <<
"Boundary conditions:"
230 << pointLocation_().boundaryField().types() <<
endl;
248 if (!diffusivityPtr_)
253 coeffDict().
lookup(
"diffusivity")
257 return *diffusivityPtr_;
264 interpolationPtr_->interpolate
274 Info<<
"displacementLaplacianFvMotionSolver : applying "
275 <<
" boundary conditions on " << pointLocation_().name()
276 <<
" to new point location."
280 pointLocation_().primitiveFieldRef() =
282 + pointDisplacement_.primitiveField();
284 pointLocation_().correctBoundaryConditions();
287 if (frozenPointsZone_ != -1)
289 const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
293 pointLocation_()[pz[i]] =
points0()[pz[i]];
305 points0() + pointDisplacement_.primitiveField()
310 if (frozenPointsZone_ != -1)
312 const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
316 curPoints[pz[i]] =
points0()[pz[i]];
320 twoDCorrectPoints(curPoints);
331 movePoints(fvMesh_.points());
333 diffusivity().correct();
334 pointDisplacement_.boundaryFieldRef().updateCoeffs();
345 *diffusivity().
operator()(),
347 "laplacian(diffusivity,cellDisplacement)"
368 diffusivityPtr_.clear();
int debug
Static debugging option.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
~displacementLaplacianFvMotionSolver()
Destructor.
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.
virtual void solve()
Solve for motion.
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.
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual void updateMesh(const mapPolyMesh &)
Update topology.
conserve primitiveFieldRef()+
motionDiffusivity & diffusivity()
Return reference to the diffusivity field.
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.
messageStream Info
Information stream (stdout output on master, null elsewhere)
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.
Lookup type of boundary radiation properties.
const dimensionSet dimViscosity
Virtual base class for displacement motion solver.
Macros for easy insertion into run-time selection tables.
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
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)