multiSolidBodyMotionFvMesh.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2012 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 "volFields.H"
29 #include "transformField.H"
30 #include "cellZoneMesh.H"
31 #include "boolList.H"
32 #include "syncTools.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(multiSolidBodyMotionFvMesh, 0);
40  (
41  dynamicFvMesh,
42  multiSolidBodyMotionFvMesh,
43  IOobject
44  );
45 }
46 
47 
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49 
50 Foam::multiSolidBodyMotionFvMesh::multiSolidBodyMotionFvMesh(const IOobject& io)
51 :
52  dynamicFvMesh(io),
53  dynamicMeshCoeffs_
54  (
56  (
57  IOobject
58  (
59  "dynamicMeshDict",
60  io.time().constant(),
61  *this,
62  IOobject::MUST_READ_IF_MODIFIED,
63  IOobject::NO_WRITE,
64  false
65  )
66  ).subDict(typeName + "Coeffs")
67  ),
68  undisplacedPoints_
69  (
70  IOobject
71  (
72  "points",
73  io.time().constant(),
74  meshSubDir,
75  *this,
76  IOobject::MUST_READ,
77  IOobject::NO_WRITE,
78  false
79  )
80  )
81 {
82  if (undisplacedPoints_.size() != nPoints())
83  {
85  (
86  "multiSolidBodyMotionFvMesh::multiSolidBodyMotionFvMesh"
87  "(const IOobject&)",
88  dynamicMeshCoeffs_
89  ) << "Read " << undisplacedPoints_.size()
90  << " undisplaced points from " << undisplacedPoints_.objectPath()
91  << " but the current mesh has " << nPoints()
92  << exit(FatalIOError);
93  }
94 
95 
96  zoneIDs_.setSize(dynamicMeshCoeffs_.size());
97  SBMFs_.setSize(dynamicMeshCoeffs_.size());
98  pointIDs_.setSize(dynamicMeshCoeffs_.size());
99  label zoneI = 0;
100 
101  forAllConstIter(dictionary, dynamicMeshCoeffs_, iter)
102  {
103  if (iter().isDict())
104  {
105  zoneIDs_[zoneI] = cellZones().findZoneID(iter().keyword());
106 
107  if (zoneIDs_[zoneI] == -1)
108  {
110  (
111  "multiSolidBodyMotionFvMesh::"
112  "multiSolidBodyMotionFvMesh(const IOobject&)",
113  dynamicMeshCoeffs_
114  ) << "Cannot find cellZone named " << iter().keyword()
115  << ". Valid zones are " << cellZones().names()
116  << exit(FatalIOError);
117  }
118 
119  const dictionary& subDict = iter().dict();
120 
121  SBMFs_.set
122  (
123  zoneI,
124  solidBodyMotionFunction::New(subDict, io.time())
125  );
126 
127  // Collect points of cell zone.
128  const cellZone& cz = cellZones()[zoneIDs_[zoneI]];
129 
130  boolList movePts(nPoints(), false);
131 
132  forAll(cz, i)
133  {
134  label cellI = cz[i];
135  const cell& c = cells()[cellI];
136  forAll(c, j)
137  {
138  const face& f = faces()[c[j]];
139  forAll(f, k)
140  {
141  label pointI = f[k];
142  movePts[pointI] = true;
143  }
144  }
145  }
146 
147  syncTools::syncPointList(*this, movePts, orEqOp<bool>(), false);
148 
149  DynamicList<label> ptIDs(nPoints());
150  forAll(movePts, i)
151  {
152  if (movePts[i])
153  {
154  ptIDs.append(i);
155  }
156  }
157 
158  pointIDs_[zoneI].transfer(ptIDs);
159 
160  Info<< "Applying solid body motion " << SBMFs_[zoneI].type()
161  << " to " << pointIDs_[zoneI].size() << " points of cellZone "
162  << iter().keyword() << endl;
163 
164  zoneI++;
165  }
166  }
167  zoneIDs_.setSize(zoneI);
168  SBMFs_.setSize(zoneI);
169  pointIDs_.setSize(zoneI);
170 }
171 
172 
173 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
174 
176 {}
177 
178 
179 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
180 
182 {
183  static bool hasWarned = false;
184 
185  pointField transformedPts(undisplacedPoints_);
186 
187  forAll(zoneIDs_, i)
188  {
189  const labelList& zonePoints = pointIDs_[i];
190 
191  UIndirectList<point>(transformedPts, zonePoints) =
192  transform
193  (
194  SBMFs_[i].transformation(),
195  pointField(transformedPts, zonePoints)
196  );
197  }
198 
199  fvMesh::movePoints(transformedPts);
200 
201  if (foundObject<volVectorField>("U"))
202  {
203  const_cast<volVectorField&>(lookupObject<volVectorField>("U"))
205  }
206  else if (!hasWarned)
207  {
208  hasWarned = true;
209 
210  WarningIn("multiSolidBodyMotionFvMesh::update()")
211  << "Did not find volVectorField U."
212  << " Not updating U boundary conditions." << endl;
213  }
214 
215  return true;
216 }
217 
218 
219 // ************************************************************************* //
const cellZoneMesh & cellZones() const
Return cell zone mesh.
Definition: polyMesh.H:469
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:390
Foam::cellZoneMesh.
A subset of mesh cells.
Definition: cellZone.H:61
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:56
labelList f(nPoints)
static autoPtr< solidBodyMotionFunction > New(const dictionary &SBMFCoeffs, const Time &runTime)
Select constructed from the SBMFCoeffs dictionary and Time.
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
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:310
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const cellList & cells() const
messageStream Info
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:638
label size() const
Return number of elements in list.
Definition: DLListBaseI.H:77
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
Namespace for OpenFOAM.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
void setSize(const label)
Reset size of List.
Definition: List.C:318
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOerror FatalIOError
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
#define WarningIn(functionName)
Report a warning using Foam::Warning.
A List with indirect addressing.
Definition: fvMatrix.H:106
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:465
#define forAll(list, i)
Definition: UList.H:421
label k
Boltzmann constant.
Spatial transformation functions for primitive fields.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
fileName objectPath() const
Return complete path + object name.
Definition: IOobject.H:363
Macros for easy insertion into run-time selection tables.
static void syncPointList(const polyMesh &, List< T > &, const CombineOp &cop, const T &nullValue, const TransformOp &top)
Synchronize values on all mesh points.
bool isDict(const word &) const
Check if entry is a sub-dictionary.
Definition: dictionary.C:607
void correctBoundaryConditions()
Correct boundary field.
Constant dispersed-phase particle diameter model.
virtual bool update()
Update the mesh for both mesh motion and topology change.
const Time & time() const
Return time.
Definition: IOobject.C:251
label findZoneID(const word &zoneName) const
Find zone index given a name.
Definition: ZoneMesh.C:354
const dimensionedScalar c
Speed of light in a vacuum.
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1060
wordList names() const
Return a list of zone names.
Definition: ZoneMesh.C:269
Abstract base class for geometry and/or topology changing fvMesh.
Definition: dynamicFvMesh.H:51
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
label nPoints() const
bool movePoints()
Do what is neccessary if the mesh has moved.
defineTypeNameAndDebug(combustionModel, 0)