35#include "surfaceInterpolate.H"
61Foam::surfaceAlignedSBRStressFvMotionSolver::
62surfaceAlignedSBRStressFvMotionSolver
69 surfaceNames_(coeffDict().
lookup(
"surfaces")),
70 surfaceMesh_(surfaceNames_.size()),
84 maxAng_(coeffDict().getOrDefault<scalar>(
"maxAng", 80)),
85 minAng_(coeffDict().getOrDefault<scalar>(
"minAng", 20)),
86 accFactor_(coeffDict().getOrDefault<scalar>(
"accFactor", 1)),
87 smoothFactor_(coeffDict().getOrDefault<scalar>(
"smoothFactor", 0.9)),
88 nNonOrthogonalCorr_(coeffDict().get<label>(
"nNonOrthogonalCorr")),
89 pointDisplacement_(pointDisplacement()),
103 minSigmaDiff_(coeffDict().getOrDefault<scalar>(
"minSigmaDiff", 1
e-4))
137void Foam::surfaceAlignedSBRStressFvMotionSolver::calculateCellRot()
139 cellRot_.primitiveFieldRef() =
Zero;
140 pointDisplacement_.primitiveFieldRef() =
Zero;
149 for (label faceI = 0; faceI < fvMesh_.nInternalFaces(); faceI++)
151 start[faceI] = Cc[fvMesh_.faceOwner()[faceI]];
152 end[faceI] = Cc[fvMesh_.faceNeighbour()[faceI]];
155 DynamicList<label> hitCells;
157 forAll(surfaceMesh_, surfi)
159 List<pointIndexHit> hit(start.size());
160 surfaceMesh_[surfi].findLineAny(start, end, hit);
164 const vectorField& nf = surfaceMesh_[surfi].faceNormals();
170 DynamicList<label> cellsHit;
174 if (hit[facei].hit())
177 const vector hitPoint = hit[facei].hitPoint();
179 if (fvMesh_.isInternalFace(facei))
181 const vector cCOne = Cc[fvMesh_.faceOwner()[facei]];
182 const vector cCTwo = Cc[fvMesh_.faceNeighbour()[facei]];
184 if (
mag(cCOne - hitPoint) <
mag(cCTwo - hitPoint))
186 rotCellId = fvMesh_.faceOwner()[facei];
190 rotCellId = fvMesh_.faceNeighbour()[facei];
196 if (isA<processorPolyPatch>(pbm[patchi]))
198 const point& ownCc = Cc[fvMesh_.faceOwner()[facei]];
200 const vector cCentreOne = ownCc - hitPoint;
203 refCast<const processorPolyPatch>(pbm[patchi])
204 .neighbFaceCellCentres()[facei];
206 const vector cCentreTwo = nbrCc - hitPoint;
208 if (cCentreOne < cCentreTwo)
210 rotCellId = fvMesh_.faceOwner()[facei];
223 const labelList& cFaces = fvMesh_.cells()[rotCellId];
225 scalar cosMax(-GREAT);
230 mag(nf[hit[facei].index()] & nSfMesh[cFaces[
k]]);
249 cellRot_[rotCellId] =
250 nSfMesh[
faceId]^nf[hit[facei].index()];
252 const scalar magRot =
mag(cellRot_[rotCellId]);
256 const scalar theta =
::asin(magRot);
257 quaternion q(cellRot_[rotCellId]/magRot, theta);
260 fvMesh_.cellPoints(rotCellId);
264 const label pointId = cPoints[j];
266 pointsCount[pointId]++;
269 fvMesh_.points()[pointId];
271 pointDisplacement_[pointId] +=
272 (
R & (pointPos - hitPoint))
273 - (pointPos - hitPoint);
281 vectorField& pd = pointDisplacement_.primitiveFieldRef();
285 point /= pointsCount[pointi];
295 this->movePoints(fvMesh_.points());
299 diffusivity().correct();
311 fvMesh_.time().timeName(),
318 cellMotionBoundaryTypes<vector>
320 pointDisplacement().boundaryField()
331 Ud[i] = pointInter.
interpolate(pointDisplacement_);
356 fvMesh_.time().timeName(),
379 const scalar diffSigmaD =
380 gSum(
mag(sigmaD_.oldTime().primitiveField()))
383 if (
mag(diffSigmaD) > minSigmaDiff_)
385 sigmaD_ = magNewSigmaD;
392 pointDisplacement_.boundaryFieldRef().updateCoeffs();
398 for (
int nonOrth=0; nonOrth<=nNonOrthogonalCorr_; nonOrth++)
406 "laplacian(diffusivity,cellDisplacement)"
#define R(A, B, C, D, E, F, K, M)
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Graphite solid properties.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
tmp< GeometricField< Type, PatchField, GeoMesh > > T() const
Return transpose (only if it is a tensor field)
void replace(const direction d, const GeometricField< cmptType, PatchField, GeoMesh > &gcf)
Replace specified field component with content from another field.
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
tmp< GeometricField< cmptType, PatchField, GeoMesh > > component(const direction) const
Return a component of the field.
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,...
const word & constant() const
Return constant name.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Mesh motion solver for an fvMesh. Based on solving the cell-centre solid-body rotation stress equatio...
A special matrix type and solver, designed for finite volume solutions of scalar equations....
SolverPerformance< Type > solveSegregatedOrCoupled(const dictionary &)
Solve segregated or coupled returning the solution statistics.
const dictionary & solverDict() const
Return the solver dictionary taking into account finalIteration.
const Time & time() const
Return the top-level database.
const fvMesh & mesh() const
Return reference to the fvMesh to be moved.
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
Virtual base class for mesh motion solver.
Container to calculate weights for interpolating directly from vertices of cell using Mean Value Coor...
Type interpolate(const GeometricField< Type, pointPatchField, pointMesh > &psip) const
Interpolate field.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
Mesh consisting of general polyhedral cells.
Lookup type of boundary radiation properties.
Mesh motion solver for an fvMesh. Based on solving the cell-centre solid-body rotation stress equatio...
~surfaceAlignedSBRStressFvMotionSolver()
Destructor.
virtual void solve()
Solve for motion.
Tensor of scalars, i.e. Tensor<scalar>.
A class for managing temporary objects.
IOoject and searching on triSurface.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
const volScalarField & mu
Calculate the divergence of the given field.
Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistr...
Calculate the matrix for the laplacian of the field.
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.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
void smooth(volScalarField &field, const scalar coeff)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
const dimensionSet dimViscosity
dimensionedScalar asin(const dimensionedScalar &ds)
Type gSum(const FieldField< Field, Type > &f)
const dimensionSet dimless
Dimensionless.
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
List< label > labelList
A List of labels.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
static const Identity< scalar > I
vector point
Point is a vector.
constexpr scalar degToRad() noexcept
Multiplication factor for degrees to radians conversion.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Field< vector > vectorField
Specialisation of Field<T> for vector.
Field< label > labelField
Specialisation of Field<T> for label.
static constexpr const zero Zero
Global zero (0)
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar lambda("lambda", dimTime/sqr(dimLength), laminarTransport)
#define forAll(list, i)
Loop across all elements in list.
Unit conversion functions.