FreeFOAM The Cross-Platform CFD Toolkit
mixerFvMesh.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) 1991-2010 OpenCFD Ltd.
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 "mixerFvMesh.H"
27 #include <OpenFOAM/Time.H>
28 #include <meshTools/regionSplit.H>
31 #include <OpenFOAM/mapPolyMesh.H>
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  defineTypeNameAndDebug(mixerFvMesh, 0);
38 
39  addToRunTimeSelectionTable(topoChangerFvMesh, mixerFvMesh, IOobject);
40 }
41 
42 
43 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
44 
45 void Foam::mixerFvMesh::addZonesAndModifiers()
46 {
47  // Add zones and modifiers for motion action
48 
49  if
50  (
51  pointZones().size()
52  || faceZones().size()
53  || cellZones().size()
54  || topoChanger_.size()
55  )
56  {
57  Info<< "void mixerFvMesh::addZonesAndModifiers() : "
58  << "Zones and modifiers already present. Skipping."
59  << endl;
60 
61  return;
62  }
63 
64  Info<< "Time = " << time().timeName() << endl
65  << "Adding zones and modifiers to the mesh" << endl;
66 
67  // Add zones
68  List<pointZone*> pz(1);
69 
70  // Add an empty zone for cut points
71 
72  pz[0] = new pointZone
73  (
74  "cutPointZone",
75  labelList(0),
76  0,
77  pointZones()
78  );
79 
80 
81  // Do face zones for slider
82 
83  List<faceZone*> fz(3);
84 
85  // Inner slider
86  const word innerSliderName(motionDict_.subDict("slider").lookup("inside"));
87  const polyPatch& innerSlider =
88  boundaryMesh()[boundaryMesh().findPatchID(innerSliderName)];
89 
90  labelList isf(innerSlider.size());
91 
92  forAll (isf, i)
93  {
94  isf[i] = innerSlider.start() + i;
95  }
96 
97  fz[0] = new faceZone
98  (
99  "insideSliderZone",
100  isf,
101  boolList(innerSlider.size(), false),
102  0,
103  faceZones()
104  );
105 
106  // Outer slider
107  const word outerSliderName(motionDict_.subDict("slider").lookup("outside"));
108  const polyPatch& outerSlider =
109  boundaryMesh()[boundaryMesh().findPatchID(outerSliderName)];
110 
111  labelList osf(outerSlider.size());
112 
113  forAll (osf, i)
114  {
115  osf[i] = outerSlider.start() + i;
116  }
117 
118  fz[1] = new faceZone
119  (
120  "outsideSliderZone",
121  osf,
122  boolList(outerSlider.size(), false),
123  1,
124  faceZones()
125  );
126 
127  // Add empty zone for cut faces
128  fz[2] = new faceZone
129  (
130  "cutFaceZone",
131  labelList(0),
132  boolList(0, false),
133  2,
134  faceZones()
135  );
136 
137  List<cellZone*> cz(1);
138 
139  // Mark every cell with its topological region
140  regionSplit rs(*this);
141 
142  // Get the region of the cell containing the origin.
143  label originRegion = rs[findNearestCell(cs().origin())];
144 
145  labelList movingCells(nCells());
146  label nMovingCells = 0;
147 
148  forAll(rs, cellI)
149  {
150  if (rs[cellI] == originRegion)
151  {
152  movingCells[nMovingCells] = cellI;
153  nMovingCells++;
154  }
155  }
156 
157  movingCells.setSize(nMovingCells);
158  Info << "Number of cells in the moving region: " << nMovingCells << endl;
159 
160  cz[0] = new cellZone
161  (
162  "movingCells",
163  movingCells,
164  0,
165  cellZones()
166  );
167 
168  Info << "Adding point, face and cell zones" << endl;
169  addZones(pz, fz, cz);
170 
171  // Add a topology modifier
172  Info << "Adding topology modifiers" << endl;
175  (
176  0,
177  new slidingInterface
178  (
179  "mixerSlider",
180  0,
181  topoChanger_,
182  outerSliderName + "Zone",
183  innerSliderName + "Zone",
184  "cutPointZone",
185  "cutFaceZone",
186  outerSliderName,
187  innerSliderName,
189  )
190  );
192 
193  write();
194 }
195 
196 
197 void Foam::mixerFvMesh::calcMovingMasks() const
198 {
199  if (debug)
200  {
201  Info<< "void mixerFvMesh::calcMovingMasks() const : "
202  << "Calculating point and cell masks"
203  << endl;
204  }
205 
206  if (movingPointsMaskPtr_)
207  {
208  FatalErrorIn("void mixerFvMesh::calcMovingMasks() const")
209  << "point mask already calculated"
210  << abort(FatalError);
211  }
212 
213  // Set the point mask
214  movingPointsMaskPtr_ = new scalarField(points().size(), 0);
215  scalarField& movingPointsMask = *movingPointsMaskPtr_;
216 
217  const cellList& c = cells();
218  const faceList& f = faces();
219 
220  const labelList& cellAddr =
221  cellZones()[cellZones().findZoneID("movingCells")];
222 
223  forAll (cellAddr, cellI)
224  {
225  const cell& curCell = c[cellAddr[cellI]];
226 
227  forAll (curCell, faceI)
228  {
229  // Mark all the points as moving
230  const face& curFace = f[curCell[faceI]];
231 
232  forAll (curFace, pointI)
233  {
234  movingPointsMask[curFace[pointI]] = 1;
235  }
236  }
237  }
238 
239  const word innerSliderZoneName
240  (
241  word(motionDict_.subDict("slider").lookup("inside"))
242  + "Zone"
243  );
244 
245  const labelList& innerSliderAddr =
246  faceZones()[faceZones().findZoneID(innerSliderZoneName)];
247 
248  forAll (innerSliderAddr, faceI)
249  {
250  const face& curFace = f[innerSliderAddr[faceI]];
251 
252  forAll (curFace, pointI)
253  {
254  movingPointsMask[curFace[pointI]] = 1;
255  }
256  }
257 
258  const word outerSliderZoneName
259  (
260  word(motionDict_.subDict("slider").lookup("outside"))
261  + "Zone"
262  );
263 
264  const labelList& outerSliderAddr =
265  faceZones()[faceZones().findZoneID(outerSliderZoneName)];
266 
267  forAll (outerSliderAddr, faceI)
268  {
269  const face& curFace = f[outerSliderAddr[faceI]];
270 
271  forAll (curFace, pointI)
272  {
273  movingPointsMask[curFace[pointI]] = 0;
274  }
275  }
276 }
277 
278 
279 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
280 
281 // Construct from components
282 Foam::mixerFvMesh::mixerFvMesh
283 (
284  const IOobject& io
285 )
286 :
287  topoChangerFvMesh(io),
288  motionDict_
289  (
291  (
292  IOobject
293  (
294  "dynamicMeshDict",
295  time().constant(),
296  *this,
299  )
300  ).subDict(typeName + "Coeffs")
301  ),
302  csPtr_
303  (
305  (
306  "coordinateSystem",
307  motionDict_.subDict("coordinateSystem")
308  )
309  ),
310  rpm_(readScalar(motionDict_.lookup("rpm"))),
311  movingPointsMaskPtr_(NULL)
312 {
313  addZonesAndModifiers();
314 
315  Info<< "Mixer mesh:" << nl
316  << " origin: " << cs().origin() << nl
317  << " axis: " << cs().axis() << nl
318  << " rpm: " << rpm_ << endl;
319 }
320 
321 
322 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
323 
325 {
326  deleteDemandDrivenData(movingPointsMaskPtr_);
327 }
328 
329 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
330 
331 // Return moving points mask. Moving points marked with 1
332 const Foam::scalarField& Foam::mixerFvMesh::movingPointsMask() const
333 {
334  if (!movingPointsMaskPtr_)
335  {
336  calcMovingMasks();
337  }
338 
339  return *movingPointsMaskPtr_;
340 }
341 
342 
344 {
345  // Rotational speed needs to be converted from rpm
346  movePoints
347  (
348  csPtr_->globalPosition
349  (
350  csPtr_->localPosition(points())
351  + vector(0, rpm_*360.0*time().deltaT().value()/60.0, 0)
352  *movingPointsMask()
353  )
354  );
355 
356  // Make changes. Use inflation (so put new points in topoChangeMap)
357  autoPtr<mapPolyMesh> topoChangeMap = topoChanger_.changeMesh(true);
358 
359  if (topoChangeMap.valid())
360  {
361  if (debug)
362  {
363  Info << "Mesh topology is changing" << endl;
364  }
365 
366  deleteDemandDrivenData(movingPointsMaskPtr_);
367  }
368 
369  movePoints
370  (
371  csPtr_->globalPosition
372  (
373  csPtr_->localPosition(oldPoints())
374  + vector(0, rpm_*360.0*time().deltaT().value()/60.0, 0)
375  *movingPointsMask()
376  )
377  );
378 
379  return true;
380 }
381 
382 
383 // ************************ vim: set sw=4 sts=4 et: ************************ //