Go to the documentation of this file.
50 elasticityMotionSolver,
65 if (isA<fixedValuePointPatchVectorField>(pointBCs))
68 refCast<fixedValuePointPatchVectorField>(pointBCs);
69 fixedValueBCs == fixedValueBCs/scalar(
nSteps_);
81 if (isA<fixedValueFvPatchVectorField>(bField))
96 Foam::elasticityMotionSolver::elasticityMotionSolver
107 refCast<const fvMesh>(
mesh)
122 fixedValuePointPatchVectorField::typeName
136 pointMotionU_.boundaryField().types()
140 coeffDict().found(
"interpolation")
156 zeroGradientFvPatchScalarField::typeName
158 exponent_(this->coeffDict().
get<scalar>(
"exponent")),
159 nSteps_(this->coeffDict().
get<label>(
"steps")),
160 nIters_(this->coeffDict().
get<label>(
"iters")),
161 tolerance_(this->coeffDict().
get<scalar>(
"tolerance"))
182 setBoundaryConditions();
185 for (label istep = 0; istep < nSteps_; ++istep)
191 E_.primitiveFieldRef() = 1./
pow(vols, exponent_);
192 E_.correctBoundaryConditions();
194 for (label iter = 0; iter < nIters_; ++iter)
196 Info<<
"Iteration " << iter <<
endl;
197 cellMotionU_.storePrevIter();
205 scalar residual =
mag(dEqn.
solve().initialResidual());
206 cellMotionU_.relax();
209 fvMesh_.time().printExecutionTime(
Info);
212 if (residual < tolerance_)
214 Info<<
"\n***Reached mesh movement convergence limit for step "
216 <<
" iteration " << iter <<
"***\n\n";
222 interpolationPtr_->interpolate(cellMotionU_, pointMotionU_);
225 mesh().
points() + pointMotionU_.primitiveFieldRef()
229 fvMesh_.movePoints(newPoints);
230 fvMesh_.checkMesh(
true);
234 Info<<
" Writing new mesh points " <<
endl;
240 mesh().pointsInstance(),
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...
virtual const pointField & points() const
Return raw points.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
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)
virtual void solve()
Solve for motion.
static const pointMesh & New(const polyMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
static word timeName(const scalar t, const int precision=precision_)
Calculate the divergence of the given field.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
virtual tmp< pointField > curPoints() const
Return point location. Mesh is actually moved in solve()
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Mesh consisting of general polyhedral cells.
Calculate the matrix for the divergence of the given field and flux.
#define forAll(list, i)
Loop across all elements in list.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
A special matrix type and solver, designed for finite volume solutions of scalar equations.
volVectorField cellMotionU_
messageStream Info
Information stream (stdout output on master, null elsewhere)
A patch is a list of labels that address the faces in the global face list.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Calculate the matrix for the laplacian of the field.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Lookup type of boundary radiation properties.
Macros for easy insertion into run-time selection tables.
Mesh data needed to do the Finite Volume discretisation.
const scalarField & cellVolumes() const
virtual void movePoints(const pointField &)
Update local data for geometry changes.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Base class for interpolation of cell displacement fields, generated by fvMotionSolvers,...
Virtual base class for mesh motion solver.
const std::string patch
OpenFOAM patch number as a std::string.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
A special matrix type and solver, designed for finite volume solutions of scalar equations....
const polyMesh & mesh() const
Return reference to mesh.
SolverPerformance< Type > solve(const dictionary &)
Solve returning the solution statistics.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
const Time & time() const
Return the top-level database.
static const Vector< scalar > zero
virtual void updateMesh(const mapPolyMesh &)
Update the mesh corresponding to given map.
void setBoundaryConditions()
Set boundary conditions of cellMotionU based on pointMotionU.
static autoPtr< motionInterpolation > New(const fvMesh &mesh)
Select default.
defineTypeNameAndDebug(combustionModel, 0)
pointVectorField pointMotionU_
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const dimensionSet dimless
Dimensionless.
label nSteps_
Intermediate steps to solve the PDEs.