linearValveLayersFvMesh.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2017 OpenFOAM Foundation
9  Copyright (C) 2018-2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
30 #include "Time.H"
31 #include "slidingInterface.H"
32 #include "layerAdditionRemoval.H"
33 #include "pointField.H"
34 #include "mapPolyMesh.H"
35 #include "polyTopoChange.H"
37 
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 
40 namespace Foam
41 {
42  defineTypeNameAndDebug(linearValveLayersFvMesh, 0);
44  (
45  topoChangerFvMesh,
46  linearValveLayersFvMesh,
47  IOobject
48  );
49 }
50 
51 
52 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
53 
54 void Foam::linearValveLayersFvMesh::addZonesAndModifiers()
55 {
56  // Add zones and modifiers for motion action
57 
58  if
59  (
60  pointZones().size()
61  || faceZones().size()
62  || cellZones().size()
63  || topoChanger_.size()
64  )
65  {
67  << "Zones and modifiers already present. Skipping."
68  << endl;
69 
70  return;
71  }
72 
73  Info<< "Time = " << time().timeName() << endl
74  << "Adding zones and modifiers to the mesh" << endl;
75 
76  // Add zones
77  List<pointZone*> pz(1);
78  List<faceZone*> fz(4);
79  List<cellZone*> cz(0);
80 
81 
82  // An empty zone for cut points
83  pz[0] = new pointZone("cutPointZone", 0, pointZones());
84 
85 
86  // Do face zones for slider
87 
88  // Inner slider
89  const word innerSliderName
90  (
91  motionDict_.subDict("slider").get<word>("inside")
92  );
93  const polyPatch& innerSlider = boundaryMesh()[innerSliderName];
94 
95  fz[0] = new faceZone
96  (
97  "insideSliderZone",
98  identity(innerSlider.range()),
99  false, // none are flipped
100  0,
101  faceZones()
102  );
103 
104  // Outer slider
105  const word outerSliderName
106  (
107  motionDict_.subDict("slider").get<word>("outside")
108  );
109  const polyPatch& outerSlider = boundaryMesh()[outerSliderName];
110 
111  fz[1] = new faceZone
112  (
113  "outsideSliderZone",
114  identity(outsideSlider.range()),
115  false, // none are flipped
116  1,
117  faceZones()
118  );
119 
120  // An empty zone for cut faces
121  fz[2] = new faceZone("cutFaceZone", 2, faceZones());
122 
123  // Add face zone for layer addition
124  const word layerPatchName
125  (
126  motionDict_.subDict("layer").get<word>("patch")
127  );
128 
129  const polyPatch& layerPatch = boundaryMesh()[layerPatchName];
130 
131  fz[3] = new faceZone
132  (
133  "valveLayerZone",
134  identity(layerPatch.range()),
135  lpf,
136  true, // all are flipped
137  0,
138  faceZones()
139  );
140 
141 
142  Info<< "Adding point and face zones" << endl;
143  addZones(pz, fz, cz);
144 
145  // Add a topology modifier
146 
147  List<polyMeshModifier*> tm(2);
148 
149  tm[0] = new slidingInterface
150  (
151  "valveSlider",
152  0,
153  topoChanger_,
154  outerSliderName + "Zone",
155  innerSliderName + "Zone",
156  "cutPointZone",
157  "cutFaceZone",
158  outerSliderName,
159  innerSliderName,
161  true // Attach-detach action
162  );
163 
164  tm[1] =
165  new layerAdditionRemoval
166  (
167  "valveLayer",
168  1,
169  topoChanger_,
170  "valveLayerZone",
171  motionDict_.subDict("layer").get<scalar>("minThickness"),
172  motionDict_.subDict("layer").get<scalar>("maxThickness")
173  );
174 
175 
176  Info<< "Adding topology modifiers" << endl;
177  addTopologyModifiers(tm);
178 
179  // Write mesh
180  write();
181 }
182 
183 
184 void Foam::linearValveLayersFvMesh::makeLayersLive()
185 {
186  const polyTopoChanger& topoChanges = topoChanger_;
187 
188  // Enable layering
189  forAll(topoChanges, modI)
190  {
191  if (isA<layerAdditionRemoval>(topoChanges[modI]))
192  {
193  topoChanges[modI].enable();
194  }
195  else if (isA<slidingInterface>(topoChanges[modI]))
196  {
197  topoChanges[modI].disable();
198  }
199  else
200  {
202  << "Don't know what to do with mesh modifier "
203  << modI << " of type " << topoChanges[modI].type()
204  << abort(FatalError);
205  }
206  }
207 }
208 
209 
210 void Foam::linearValveLayersFvMesh::makeSlidersLive()
211 {
212  const polyTopoChanger& topoChanges = topoChanger_;
213 
214  // Enable sliding interface
215  forAll(topoChanges, modI)
216  {
217  if (isA<layerAdditionRemoval>(topoChanges[modI]))
218  {
219  topoChanges[modI].disable();
220  }
221  else if (isA<slidingInterface>(topoChanges[modI]))
222  {
223  topoChanges[modI].enable();
224  }
225  else
226  {
228  << "Don't know what to do with mesh modifier "
229  << modI << " of type " << topoChanges[modI].type()
230  << abort(FatalError);
231  }
232  }
233 }
234 
235 
236 bool Foam::linearValveLayersFvMesh::attached() const
237 {
238  const polyTopoChanger& topoChanges = topoChanger_;
239 
240  bool result = false;
241 
242  forAll(topoChanges, modI)
243  {
244  if (isA<slidingInterface>(topoChanges[modI]))
245  {
246  result =
247  result
248  || refCast<const slidingInterface>(topoChanges[modI]).attached();
249  }
250  }
251 
252  // Check thal all sliders are in sync (debug only)
253  forAll(topoChanges, modI)
254  {
255  if (isA<slidingInterface>(topoChanges[modI]))
256  {
257  if
258  (
259  result
260  != refCast<const slidingInterface>(topoChanges[modI]).attached()
261  )
262  {
264  << "Slider " << modI << " named "
265  << topoChanges[modI].name()
266  << " out of sync: Should be" << result
267  << abort(FatalError);
268  }
269  }
270  }
271 
272  return result;
273 }
274 
275 
276 Foam::tmp<Foam::pointField> Foam::linearValveLayersFvMesh::newPoints() const
277 {
278  tmp<pointField> tnewPoints
279  (
280  new pointField(points())
281  );
282 
283  pointField& np = tnewPoints();
284 
285  const word layerPatchName
286  (
287  motionDict_.subDict("layer").get<word>("patch")
288  );
289 
290  const polyPatch& layerPatch = boundaryMesh()[layerPatchName];
291 
292  const labelList& patchPoints = layerPatch.meshPoints();
293 
294  const vector vel
295  (
296  motionDict_.get<vector>("pistonVelocity")
297  );
298 
299  forAll(patchPoints, ppI)
300  {
301  np[patchPoints[ppI]] += vel*time().deltaTValue();
302  }
303 
304  return tnewPoints;
305 }
306 
307 
308 
309 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
310 
311 // Construct from components
312 Foam::linearValveLayersFvMesh::linearValveLayersFvMesh(const IOobject& io)
313 :
314  topoChangerFvMesh(io),
315  motionDict_
316  (
318  (
319  IOobject
320  (
321  "dynamicMeshDict",
322  time().constant(),
323  *this,
324  IOobject::MUST_READ_IF_MODIFIED,
325  IOobject::NO_WRITE,
326  false
327  )
328  ).optionalSubDict(typeName + "Coeffs")
329  )
330 {
331  addZonesAndModifiers();
332 }
333 
334 
335 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
336 
338 {}
339 
340 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
341 
343 {
344  // Detaching the interface
345  if (attached())
346  {
347  Info<< "Decoupling sliding interfaces" << endl;
348  makeSlidersLive();
349 
350  // Changing topology
351  resetMorph();
352  setMorphTimeIndex(3*time().timeIndex());
353  updateMesh();
354  }
355  else
356  {
357  Info<< "Sliding interfaces decoupled" << endl;
358  }
359 
360  // Perform layer action and mesh motion
361  makeLayersLive();
362 
363  // Changing topology
364  resetMorph();
365  setMorphTimeIndex(3*time().timeIndex() + 1);
366  updateMesh();
367 
368  if (topoChangeMap)
369  {
370  if (topoChangeMap().hasMotionPoints())
371  {
372  Info<< "Topology change; executing pre-motion" << endl;
373  movePoints(topoChangeMap().preMotionPoints());
374  }
375  }
376 
377  // Move points
378  movePoints(newPoints());
379 
380  // Attach the interface
381  Info<< "Coupling sliding interfaces" << endl;
382  makeSlidersLive();
383 
384  // Changing topology
385  resetMorph();
386  setMorphTimeIndex(3*time().timeIndex() + 2);
387  updateMesh();
388 
389  //Info<< "Moving points post slider attach" << endl;
390  //const pointField p = allPoints();
391  //movePoints(p);
392 
393  Info<< "Sliding interfaces coupled: " << attached() << endl;
394 }
395 
396 
397 // ************************************************************************* //
Foam::addToRunTimeSelectionTable
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
Foam::labelList
List< label > labelList
A List of labels.
Definition: List.H:67
Foam::pointField
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:44
Foam::IOdictionary
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:54
Foam::HashTable< regIOobject * >::size
label size() const noexcept
The number of elements in table.
Definition: HashTableI.H:52
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
InfoInFunction
#define InfoInFunction
Report an information message using Foam::Info.
Definition: messageStream.H:350
Foam::fvMesh::write
virtual bool write(const bool valid=true) const
Write mesh using IO settings from time.
Definition: fvMesh.C:1041
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::linearValveLayersFvMesh::~linearValveLayersFvMesh
virtual ~linearValveLayersFvMesh()
Destructor.
Definition: linearValveLayersFvMesh.C:337
mapPolyMesh.H
polyTopoChange.H
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
Foam::polyMesh::boundaryMesh
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:444
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::dictionary::get
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:107
layerAdditionRemoval.H
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::polyMesh::cellZones
const cellZoneMesh & cellZones() const noexcept
Return cell zone mesh.
Definition: polyMesh.H:492
Foam::dictionary::subDict
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
Definition: dictionary.C:460
Foam::polyMesh::faceZones
const faceZoneMesh & faceZones() const noexcept
Return face zone mesh.
Definition: polyMesh.H:486
linearValveLayersFvMesh.H
Foam::FatalError
error FatalError
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::vector
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:51
pointField.H
Time.H
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
slidingInterface.H
Foam::linearValveLayersFvMesh::update
virtual bool update()
Update the mesh for both mesh motion and topology change.
Definition: linearValveLayersFvMesh.C:342
Foam::topoChangerFvMesh
Abstract base class for a topology changing fvMesh.
Definition: topoChangerFvMesh.H:53
points
const pointField & points
Definition: gmvOutputHeader.H:1
Foam::identity
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
Definition: labelList.C:38
Foam::polyMesh::pointZones
const pointZoneMesh & pointZones() const noexcept
Return point zone mesh.
Definition: polyMesh.H:480
Foam::fvMesh::time
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:280
timeIndex
label timeIndex
Definition: getTimeIndex.H:30
constant
constant condensation/saturation model.
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
Foam::slidingInterface::INTEGRAL
Definition: slidingInterface.H:85
Foam::polyMesh::addZones
void addZones(const List< pointZone * > &pz, const List< faceZone * > &fz, const List< cellZone * > &cz)
Add mesh zones.
Definition: polyMesh.C:999
Foam::topoChangerFvMesh::topoChanger_
polyTopoChanger topoChanger_
Definition: topoChangerFvMesh.H:70