solidBodyMotionSolver.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) 2016-2019 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 
26 #include "solidBodyMotionSolver.H"
28 #include "transformField.H"
29 #include "cellZoneMesh.H"
30 #include "cellSet.H"
31 #include "boolList.H"
32 #include "syncTools.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(solidBodyMotionSolver, 0);
40  (
41  motionSolver,
42  solidBodyMotionSolver,
43  dictionary
44  );
45 }
46 
47 
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49 
51 (
52  const polyMesh& mesh,
53  const dictionary& dict
54 )
55 :
56  points0MotionSolver(mesh, dict, typeName),
57  SBMFPtr_(solidBodyMotionFunction::New(coeffDict(), mesh.time())),
58  pointIDs_(),
59  moveAllCells_(false)
60 {
61  word cellZoneName =
62  coeffDict().lookupOrDefault<word>("cellZone", "none");
63 
64  word cellSetName =
65  coeffDict().lookupOrDefault<word>("cellSet", "none");
66 
67  if ((cellZoneName != "none") && (cellSetName != "none"))
68  {
69  FatalIOErrorInFunction(coeffDict())
70  << "Either cellZone OR cellSet can be supplied, but not both. "
71  << "If neither is supplied, all cells will be included"
72  << exit(FatalIOError);
73  }
74 
75  labelList cellIDs;
76  if (cellZoneName != "none")
77  {
78  Info<< "Applying solid body motion to cellZone " << cellZoneName
79  << endl;
80 
81  label zoneID = mesh.cellZones().findZoneID(cellZoneName);
82 
83  if (zoneID == -1)
84  {
86  << "Unable to find cellZone " << cellZoneName
87  << ". Valid cellZones are:"
88  << mesh.cellZones().names()
89  << exit(FatalError);
90  }
91 
92  cellIDs = mesh.cellZones()[zoneID];
93  }
94 
95  if (cellSetName != "none")
96  {
97  Info<< "Applying solid body motion to cellSet " << cellSetName
98  << endl;
99 
100  cellSet set(mesh, cellSetName);
101 
102  cellIDs = set.toc();
103  }
104 
105  label nCells = returnReduce(cellIDs.size(), sumOp<label>());
106  moveAllCells_ = nCells == 0;
107 
108  if (moveAllCells_)
109  {
110  Info<< "Applying solid body motion to entire mesh" << endl;
111  }
112  else
113  {
114  // collect point IDs of points in cell zone
115 
116  boolList movePts(mesh.nPoints(), false);
117 
118  forAll(cellIDs, i)
119  {
120  label celli = cellIDs[i];
121  const cell& c = mesh.cells()[celli];
122  forAll(c, j)
123  {
124  const face& f = mesh.faces()[c[j]];
125  forAll(f, k)
126  {
127  label pointi = f[k];
128  movePts[pointi] = true;
129  }
130  }
131  }
132 
133  syncTools::syncPointList(mesh, movePts, orEqOp<bool>(), false);
134 
135  DynamicList<label> ptIDs(mesh.nPoints());
136  forAll(movePts, i)
137  {
138  if (movePts[i])
139  {
140  ptIDs.append(i);
141  }
142  }
143 
144  pointIDs_.transfer(ptIDs);
145  }
146 }
147 
148 
149 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
150 
152 {}
153 
154 
155 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
156 
158 {
159  if (moveAllCells_)
160  {
161  return transformPoints(SBMFPtr_().transformation(), points0_);
162  }
163  else
164  {
165  tmp<pointField> ttransformedPts(new pointField(mesh().points()));
166  pointField& transformedPts = ttransformedPts.ref();
167 
168  UIndirectList<point>(transformedPts, pointIDs_) = transformPoints
169  (
170  SBMFPtr_().transformation(),
171  pointField(points0_, pointIDs_)
172  );
173 
174  return ttransformedPts;
175  }
176 }
177 
178 
179 // ************************************************************************* //
solidBodyMotionSolver(const polyMesh &, const dictionary &)
Construct from mesh and dictionary.
#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
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
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:174
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
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.
Spatial transformation functions for primitive fields.
Foam::cellZoneMesh.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
dynamicFvMesh & mesh
const pointField & points
label findZoneID(const word &zoneName) const
Find zone index given a name.
Definition: ZoneMesh.C:341
A class for handling words, derived from string.
Definition: word.H:59
const cellZoneMesh & cellZones() const
Return cell zone mesh.
Definition: polyMesh.H:482
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
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
const Time & time() const
Return time.
defineTypeNameAndDebug(combustionModel, 0)
labelList f(nPoints)
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:56
A collection of cell labels.
Definition: cellSet.H:48
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)
Synchronize values on all mesh points.
wordList names() const
Return a list of zone names.
Definition: ZoneMesh.C:256
messageStream Info
label nPoints() const
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
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
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.