interpolatingSolidBodyMotionSolver.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2018-2024 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
28 #include "pointDist.H"
29 #include "pointConstraints.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
36 
38  (
42  );
43 }
44 
45 
46 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
47 
48 void Foam::interpolatingSolidBodyMotionSolver::calcScale()
49 {
50  const pointMesh& pMesh = pointMesh::New(mesh());
51 
52  const pointDist pDist(pMesh, patchSet_, points0(), do_);
53 
54  // Scaling: 1 up to di then linear down to 0 at do away from patches
55  scale_.primitiveFieldRef() =
56  min
57  (
58  max
59  (
60  (do_ - pDist.primitiveField())/(do_ - di_),
61  scalar(0)
62  ),
63  scalar(1)
64  );
65 
66  // Convert the scale function to a cosine
67  scale_.primitiveFieldRef() =
68  min
69  (
70  max
71  (
72  0.5
73  - 0.5
74  *cos(scale_.primitiveField()
76  scalar(0)
77  ),
78  scalar(1)
79  );
80 
81  pointConstraints::New(pMesh).constrain(scale_);
82  scale_.write();
83 }
84 
85 
86 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
87 
89 (
90  const word& name,
91  const polyMesh& mesh,
92  const dictionary& dict
93 )
94 :
95  points0MotionSolver(name, mesh, dict, typeName),
96  SBMFPtr_(solidBodyMotionFunction::New(coeffDict(), mesh.time())),
97  patches_(wordReList(coeffDict().lookup("patches"))),
98  patchSet_(mesh.boundaryMesh().patchSet(patches_)),
99  CofG_(coeffDict().lookup("CofG")),
100  di_(coeffDict().lookup<scalar>("innerDistance")),
101  do_(coeffDict().lookup<scalar>("outerDistance")),
102  scale_
103  (
104  IOobject
105  (
106  "motionScale",
107  mesh.time().name(),
108  mesh,
109  IOobject::NO_READ,
110  IOobject::NO_WRITE
111  ),
112  pointMesh::New(mesh),
114  )
115 {
116  // Calculate scaling factor everywhere
117  calcScale();
118 }
119 
120 
121 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
122 
124 {}
125 
126 
127 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
128 
131 {
132  const pointField& points0 = this->points0();
133  const septernion s = SBMFPtr_().transformation();
134 
135  tmp<pointField> tpoints(new pointField(points0));
136  pointField& points = tpoints.ref();
137 
138  forAll(points, pointi)
139  {
140  // Move non-stationary points
141  if (scale_[pointi] > small)
142  {
143  // Use solid-body motion where scale = 1
144  if (scale_[pointi] > 1 - small)
145  {
146  points[pointi] =
147  CofG_ + s.transformPoint(points0[pointi] - CofG_);
148  }
149  // Slerp septernion interpolation
150  else
151  {
152  const septernion ss(slerp(septernion::I, s, scale_[pointi]));
153 
154  points[pointi] =
155  CofG_ + ss.transformPoint(points0[pointi] - CofG_);
156  }
157  }
158  }
159 
160  return tpoints;
161 }
162 
163 
165 (
166  const polyTopoChangeMap& map
167 )
168 {
169  // Pending implementation of the inverse transformation of points0
171 }
172 
173 
175 {
177 
178  // scale is resized by the meshToMesh mapper
179  scale_ = Zero;
180  calcScale();
181 }
182 
183 
184 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
static pointMesh & New(const word &name, const polyMesh &mesh)
Construct and return the named DemandDrivenMeshObject.
Internal::FieldType & primitiveFieldRef()
Return a reference to the primitive field.
const Internal::FieldType & primitiveField() const
Return a const-reference to the primitive field.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Solid-body motion of the mesh specified by a run-time selectable motion function. Applies SLERP inter...
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
virtual void topoChange(const polyTopoChangeMap &)
Update local data for topology changes.
interpolatingSolidBodyMotionSolver(const word &name, const polyMesh &, const dictionary &dict)
Construct from polyMesh and dictionary.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
Virtual base class for mesh motion solver.
Definition: motionSolver.H:57
const polyMesh & mesh() const
Return reference to mesh.
Definition: motionSolver.H:137
void constrain(PointField< Type > &pf, const bool overrideValue=false) const
Apply boundary conditions (single-patch constraints) and.
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:53
Virtual base class for displacement motion solvers.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
pointField & points0()
Return reference to the reference field.
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:51
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
virtual bool write(const bool write=true) const
Write using setting from DB.
Septernion class used to perform translations and rotations in 3D space.
Definition: septernion.H:66
vector transformPoint(const vector &v) const
Transform the given coordinate point.
Definition: septernionI.H:82
static const septernion I
Definition: septernion.H:83
Base class for defining solid-body motions.
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
A class for handling words, derived from string.
Definition: word.H:62
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:381
const pointField & points
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
addToRunTimeSelectionTable(polyPatch, mergedCyclicPolyPatch, word)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
const dimensionSet dimless
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
defineTypeNameAndDebug(combustionModel, 0)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
quaternion slerp(const quaternion &qa, const quaternion &qb, const scalar t)
Spherical linear interpolation of quaternions.
Definition: quaternion.C:55
dimensionedScalar cos(const dimensionedScalar &ds)
dictionary dict