multiSolidBodyMotionSolver.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) 2011-2021 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 "transformField.H"
29 #include "meshCellZones.H"
30 #include "boolList.H"
31 #include "syncTools.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  defineTypeNameAndDebug(multiSolidBodyMotionSolver, 0);
39  (
40  motionSolver,
41  multiSolidBodyMotionSolver,
42  dictionary
43  );
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
50 (
51  const polyMesh& mesh,
52  const dictionary& dict
53 )
54 :
55  points0MotionSolver(mesh, dict, typeName)
56 {
57  zoneIDs_.setSize(coeffDict().size());
58  SBMFs_.setSize(coeffDict().size());
59  pointIDs_.setSize(coeffDict().size());
60  label zonei = 0;
61 
62  forAllConstIter(dictionary, coeffDict(), iter)
63  {
64  if (iter().isDict())
65  {
66  zoneIDs_[zonei] = mesh.cellZones().findZoneID(iter().keyword());
67 
68  if (zoneIDs_[zonei] == -1)
69  {
71  (
72  coeffDict()
73  ) << "Cannot find cellZone named " << iter().keyword()
74  << ". Valid zones are " << mesh.cellZones().names()
75  << exit(FatalIOError);
76  }
77 
78  const dictionary& subDict = iter().dict();
79 
80  SBMFs_.set
81  (
82  zonei,
83  solidBodyMotionFunction::New(subDict, mesh.time())
84  );
85 
86  // Collect points of cell zone.
87  const cellZone& cz = mesh.cellZones()[zoneIDs_[zonei]];
88 
89  boolList movePts(mesh.nPoints(), false);
90 
91  forAll(cz, i)
92  {
93  label celli = cz[i];
94  const cell& c = mesh.cells()[celli];
95  forAll(c, j)
96  {
97  const face& f = mesh.faces()[c[j]];
98  forAll(f, k)
99  {
100  label pointi = f[k];
101  movePts[pointi] = true;
102  }
103  }
104  }
105 
106  syncTools::syncPointList(mesh, movePts, orEqOp<bool>(), false);
107 
108  DynamicList<label> ptIDs(mesh.nPoints());
109  forAll(movePts, i)
110  {
111  if (movePts[i])
112  {
113  ptIDs.append(i);
114  }
115  }
116 
117  pointIDs_[zonei].transfer(ptIDs);
118 
119  Info<< "Applying solid body motion " << SBMFs_[zonei].type()
120  << " to " << pointIDs_[zonei].size() << " points of cellZone "
121  << iter().keyword() << endl;
122 
123  zonei++;
124  }
125  }
126  zoneIDs_.setSize(zonei);
127  SBMFs_.setSize(zonei);
128  pointIDs_.setSize(zonei);
129 }
130 
131 
132 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
133 
135 {}
136 
137 
138 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
139 
141 {
142  tmp<pointField> ttransformedPts(new pointField(mesh().points()));
143  pointField& transformedPts = ttransformedPts.ref();
144 
145  forAll(zoneIDs_, i)
146  {
147  const labelList& zonePoints = pointIDs_[i];
148 
149  UIndirectList<point>(transformedPts, zonePoints) = transformPoints
150  (
151  SBMFs_[i].transformation(),
152  pointField(points0_, zonePoints)
153  );
154  }
155 
156  return ttransformedPts;
157 }
158 
159 
160 // ************************************************************************* //
label findZoneID(const word &zoneName) const
Find zone index given a name.
Definition: MeshZones.C:341
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
const meshCellZones & cellZones() const
Return cell zones.
Definition: polyMesh.H:482
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
void transformPoints(vectorField &, const septernion &, const vectorField &)
Transform given vectorField of coordinates with the given septernion.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const cellList & cells() const
label k
Boltzmann constant.
const dimensionedScalar c
Speed of light in a vacuum.
Macros for easy insertion into run-time selection tables.
wordList names() const
Return a list of zone names.
Definition: MeshZones.C:256
Spatial transformation functions for primitive fields.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
dynamicFvMesh & mesh
const pointField & points
static autoPtr< solidBodyMotionFunction > New(const dictionary &SBMFCoeffs, const Time &runTime)
Select constructed from the SBMFCoeffs dictionary and Time.
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1156
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
const Time & time() const
Return time.
defineTypeNameAndDebug(combustionModel, 0)
labelList f(nPoints)
A subset of mesh cells.
Definition: cellZone.H:61
Foam::meshCellZones.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:335
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:56
Virtual base class for displacement motion solvers.
A List with indirect addressing.
Definition: fvMatrix.H:106
static void syncPointList(const polyMesh &, List< T > &, const CombineOp &cop, const T &nullValue, const TransformOp &top)
Synchronise values on all mesh points.
void set(entry *)
Assign a new entry, overwrite any existing entry.
Definition: dictionary.C:1271
messageStream Info
label nPoints() const
multiSolidBodyMotionSolver(const polyMesh &, const dictionary &)
Construct from mesh and dictionary.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A class for managing temporary objects.
Definition: PtrList.H:53
Namespace for OpenFOAM.
IOerror FatalIOError