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-2022 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 "pointPatchDist.H"
29 #include "pointConstraints.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  defineTypeNameAndDebug(interpolatingSolidBodyMotionSolver, 0);
36 
38  (
39  motionSolver,
40  interpolatingSolidBodyMotionSolver,
41  dictionary
42  );
43 }
44 
45 
46 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 
49 (
50  const word& name,
51  const polyMesh& mesh,
52  const dictionary& dict
53 )
54 :
55  points0MotionSolver(name, mesh, dict, typeName),
56  SBMFPtr_(solidBodyMotionFunction::New(coeffDict(), mesh.time())),
57  patches_(wordReList(coeffDict().lookup("patches"))),
58  patchSet_(mesh.boundaryMesh().patchSet(patches_)),
59  CofG_(coeffDict().lookup("CofG")),
60  di_(coeffDict().lookup<scalar>("innerDistance")),
61  do_(coeffDict().lookup<scalar>("outerDistance")),
62  scale_
63  (
64  IOobject
65  (
66  "motionScale",
67  mesh.time().timeName(),
68  mesh,
71  false
72  ),
73  pointMesh::New(mesh),
75  )
76 {
77  // Calculate scaling factor everywhere
78 
79  {
80  const pointMesh& pMesh = pointMesh::New(mesh);
81 
82  pointPatchDist pDist(pMesh, patchSet_, points0());
83 
84  // Scaling: 1 up to di then linear down to 0 at do away from patches
85  scale_.primitiveFieldRef() =
86  min
87  (
88  max
89  (
90  (do_ - pDist.primitiveField())/(do_ - di_),
91  scalar(0)
92  ),
93  scalar(1)
94  );
95 
96  // Convert the scale function to a cosine
97  scale_.primitiveFieldRef() =
98  min
99  (
100  max
101  (
102  0.5
103  - 0.5
104  *cos(scale_.primitiveField()
106  scalar(0)
107  ),
108  scalar(1)
109  );
110 
111  pointConstraints::New(pMesh).constrain(scale_);
112  scale_.write();
113  }
114 }
115 
116 
117 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
118 
120 {}
121 
122 
123 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
124 
127 {
128  const pointField& points0 = this->points0();
129  const septernion s = SBMFPtr_().transformation();
130 
131  tmp<pointField> tpoints(new pointField(points0));
132  pointField& points = tpoints.ref();
133 
134  forAll(points, pointi)
135  {
136  // Move non-stationary points
137  if (scale_[pointi] > small)
138  {
139  // Use solid-body motion where scale = 1
140  if (scale_[pointi] > 1 - small)
141  {
142  points[pointi] =
143  CofG_ + s.transformPoint(points0[pointi] - CofG_);
144  }
145  // Slerp septernion interpolation
146  else
147  {
148  const septernion ss(slerp(septernion::I, s, scale_[pointi]));
149 
150  points[pointi] =
151  CofG_ + ss.transformPoint(points0[pointi] - CofG_);
152  }
153  }
154  }
155 
156  return tpoints;
157 }
158 
159 
160 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:453
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
const dimensionSet dimless
Septernion class used to perform translations and rotations in 3D space.
Definition: septernion.H:65
interpolatingSolidBodyMotionSolver(const word &name, const polyMesh &, const dictionary &dict)
Construct from polyMesh and dictionary.
fvMesh & mesh
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool usePatchGroups=true) const
Return the set of patch IDs corresponding to the given names.
Macros for easy insertion into run-time selection tables.
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
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.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
stressControl lookup("compactNormalStress") >> compactNormalStress
dimensionedScalar cos(const dimensionedScalar &ds)
const pointField & points
static word timeName(const scalar, const int precision=curPrecision_)
Return time name of given scalar time.
Definition: Time.C:666
A class for handling words, derived from string.
Definition: word.H:59
static autoPtr< solidBodyMotionFunction > New(const dictionary &SBMFCoeffs, const Time &runTime)
Select constructed from the SBMFCoeffs dictionary and Time.
Calculation of distance to nearest patch for all points.
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
static const septernion I
Definition: septernion.H:83
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
const Time & time() const
Return time.
defineTypeNameAndDebug(combustionModel, 0)
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Virtual base class for displacement motion solvers.
quaternion slerp(const quaternion &qa, const quaternion &qb, const scalar t)
Spherical linear interpolation of quaternions.
Definition: quaternion.C:55
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
List< wordRe > wordReList
A List of wordRe (word or regular expression)
Definition: wordReList.H:50
A class for managing temporary objects.
Definition: PtrList.H:53
vector transformPoint(const vector &v) const
Transform the given coordinate point.
Definition: septernionI.H:82
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
void constrain(GeometricField< Type, pointPatchField, pointMesh > &pf, const bool overrideValue=false) const
Apply boundary conditions (single-patch constraints) and.
Namespace for OpenFOAM.