Go to the documentation of this file.
34 #include "surfaceInterpolate.H"
52 surfaceAlignedSBRStressFvMotionSolver,
60 Foam::surfaceAlignedSBRStressFvMotionSolver::
61 surfaceAlignedSBRStressFvMotionSolver
68 surfaceNames_(coeffDict().
lookup(
"surfaces")),
69 surfaceMesh_(surfaceNames_.size()),
75 fvMesh_.time().timeName(),
83 maxAng_(coeffDict().getOrDefault<scalar>(
"maxAng", 80)),
84 minAng_(coeffDict().getOrDefault<scalar>(
"minAng", 20)),
85 accFactor_(coeffDict().getOrDefault<scalar>(
"accFactor", 1)),
86 smoothFactor_(coeffDict().getOrDefault<scalar>(
"smoothFactor", 0.9)),
87 nNonOrthogonalCorr_(coeffDict().get<label>(
"nNonOrthogonalCorr")),
88 pointDisplacement_(pointDisplacement()),
94 fvMesh_.time().timeName(),
96 IOobject::READ_IF_PRESENT,
102 minSigmaDiff_(coeffDict().getOrDefault<scalar>(
"minSigmaDiff", 1
e-4))
114 mesh.time().constant(),
136 void Foam::surfaceAlignedSBRStressFvMotionSolver::calculateCellRot()
138 cellRot_.primitiveFieldRef() =
Zero;
139 pointDisplacement_.primitiveFieldRef() =
Zero;
148 for (label faceI = 0; faceI < fvMesh_.nInternalFaces(); faceI++)
150 start[faceI] = Cc[fvMesh_.faceOwner()[faceI]];
151 end[faceI] = Cc[fvMesh_.faceNeighbour()[faceI]];
154 DynamicList<label> hitCells;
156 forAll(surfaceMesh_, surfi)
158 List<pointIndexHit> hit(start.size());
159 surfaceMesh_[surfi].findLineAny(start,
end, hit);
163 const vectorField& nf = surfaceMesh_[surfi].faceNormals();
169 DynamicList<label> cellsHit;
173 if (hit[facei].hit())
176 const vector hitPoint = hit[facei].hitPoint();
178 if (fvMesh_.isInternalFace(facei))
180 const vector cCOne = Cc[fvMesh_.faceOwner()[facei]];
181 const vector cCTwo = Cc[fvMesh_.faceNeighbour()[facei]];
183 if (
mag(cCOne - hitPoint) <
mag(cCTwo - hitPoint))
185 rotCellId = fvMesh_.faceOwner()[facei];
189 rotCellId = fvMesh_.faceNeighbour()[facei];
195 if (isA<processorPolyPatch>(pbm[patchi]))
197 const point& ownCc = Cc[fvMesh_.faceOwner()[facei]];
199 const vector cCentreOne = ownCc - hitPoint;
202 refCast<const processorPolyPatch>(pbm[patchi])
203 .neighbFaceCellCentres()[facei];
205 const vector cCentreTwo = nbrCc - hitPoint;
207 if (cCentreOne < cCentreTwo)
209 rotCellId = fvMesh_.faceOwner()[facei];
222 const labelList& cFaces = fvMesh_.cells()[rotCellId];
224 scalar cosMax(-GREAT);
229 mag(nf[hit[facei].index()] & nSfMesh[cFaces[
k]]);
248 cellRot_[rotCellId] =
249 nSfMesh[
faceId]^nf[hit[facei].index()];
251 const scalar magRot =
mag(cellRot_[rotCellId]);
255 const scalar theta =
::asin(magRot);
256 quaternion q(cellRot_[rotCellId]/magRot, theta);
259 fvMesh_.cellPoints(rotCellId);
263 const label pointId = cPoints[j];
265 pointsCount[pointId]++;
268 fvMesh_.points()[pointId];
270 pointDisplacement_[pointId] +=
271 (
R & (pointPos - hitPoint))
272 - (pointPos - hitPoint);
280 vectorField& pd = pointDisplacement_.primitiveFieldRef();
284 point /= pointsCount[pointi];
294 this->movePoints(fvMesh_.points());
298 diffusivity().correct();
310 fvMesh_.time().timeName(),
317 cellMotionBoundaryTypes<vector>
319 pointDisplacement().boundaryField()
330 Ud[i] = pointInter.
interpolate(pointDisplacement_);
355 fvMesh_.time().timeName(),
367 mu.primitiveFieldRef() = (1.0/V);
378 const scalar diffSigmaD =
379 gSum(
mag(sigmaD_.oldTime().primitiveField()))
382 if (
mag(diffSigmaD) > minSigmaDiff_)
384 sigmaD_ = magNewSigmaD;
391 pointDisplacement_.boundaryFieldRef().updateCoeffs();
397 for (
int nonOrth=0; nonOrth<=nNonOrthogonalCorr_; nonOrth++)
405 "laplacian(diffusivity,cellDisplacement)"
425 DEqn.solveSegregatedOrCoupled(DEqn.solverDict());
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
List< label > labelList
A List of labels.
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,...
tmp< GeometricField< cmptType, PatchField, GeoMesh > > component(const direction) const
Return a component of the field.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
void replace(const direction d, const GeometricField< cmptType, PatchField, GeoMesh > &gcf)
Replace specified field component with content from another field.
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
const dimensionedScalar mu
Atomic mass unit.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensioned< symmTensor > dimensionedSymmTensor
Dimensioned tensor obtained from generic dimensioned type.
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.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Calculate the divergence of the given field.
Unit conversion functions.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
IOoject and searching on triSurface.
Type gSum(const FieldField< Field, Type > &f)
Type interpolate(const GeometricField< Type, pointPatchField, pointMesh > &psip) const
Interpolate field.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
Field< vector > vectorField
Specialisation of Field<T> for vector.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
virtual void solve()
Solve for motion.
#define R(A, B, C, D, E, F, K, M)
Field< label > labelField
Specialisation of Field<T> for label.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Calculate the matrix for the laplacian of the field.
~surfaceAlignedSBRStressFvMotionSolver()
Destructor.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
Container to calculate weights for interpolating directly from vertices of cell using Mean Value Coor...
Mesh motion solver for an fvMesh. Based on solving the cell-centre solid-body rotation stress equatio...
Lookup type of boundary radiation properties.
const dimensionSet dimViscosity
dimensionedScalar lambda("lambda", dimTime/sqr(dimLength), laminarTransport)
Macros for easy insertion into run-time selection tables.
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Vector< scalar > vector
A scalar version of the templated Vector.
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
tmp< GeometricField< Type, PatchField, GeoMesh > > T() const
Return transpose (only if it is a tensor field)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label k
Boltzmann constant.
const dimensionedScalar e
Elementary charge.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
vector point
Point is a vector.
Graphite solid properties.
defineTypeNameAndDebug(combustionModel, 0)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
Tensor< scalar > tensor
Tensor of scalars, i.e. Tensor<scalar>.
dimensionedScalar asin(const dimensionedScalar &ds)
Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistr...
static const Identity< scalar > I
const dimensionSet dimless
Dimensionless.
void smooth(volScalarField &field, const scalar coeff)
dimensionedScalar cos(const dimensionedScalar &ds)