regionModel1D.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-2013 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 "regionModel1D.H"
27 
28 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 namespace regionModels
33 {
34  defineTypeNameAndDebug(regionModel1D, 0);
35 }
36 }
37 
38 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
39 
40 void Foam::regionModels::regionModel1D::constructMeshObjects()
41 {
42 
43  nMagSfPtr_.reset
44  (
46  (
47  IOobject
48  (
49  "nMagSf",
50  time().timeName(),
51  regionMesh(),
54  ),
55  regionMesh(),
56  dimensionedScalar("zero", dimArea, 0.0)
57  )
58  );
59 }
60 
61 
62 void Foam::regionModels::regionModel1D::initialise()
63 {
64  if (debug)
65  {
66  Pout<< "regionModel1D::initialise()" << endl;
67  }
68 
69  // Calculate boundaryFaceFaces and boundaryFaceCells
70 
71  DynamicList<label> faceIDs;
72  DynamicList<label> cellIDs;
73 
74  label localPyrolysisFaceI = 0;
75 
76  const polyBoundaryMesh& rbm = regionMesh().boundaryMesh();
77 
78  forAll(intCoupledPatchIDs_, i)
79  {
80  const label patchI = intCoupledPatchIDs_[i];
81  const polyPatch& ppCoupled = rbm[patchI];
82  forAll(ppCoupled, localFaceI)
83  {
84  label faceI = ppCoupled.start() + localFaceI;
85  label cellI = -1;
86  label nFaces = 0;
87  label nCells = 0;
88  do
89  {
90  label ownCellI = regionMesh().faceOwner()[faceI];
91  if (ownCellI != cellI)
92  {
93  cellI = ownCellI;
94  }
95  else
96  {
97  cellI = regionMesh().faceNeighbour()[faceI];
98  }
99  nCells++;
100  cellIDs.append(cellI);
101  const cell& cFaces = regionMesh().cells()[cellI];
102  faceI = cFaces.opposingFaceLabel(faceI, regionMesh().faces());
103  faceIDs.append(faceI);
104  nFaces++;
105  } while (regionMesh().isInternalFace(faceI));
106 
107  boundaryFaceOppositeFace_[localPyrolysisFaceI] = faceI;
108  faceIDs.remove(); //remove boundary face.
109  nFaces--;
110 
111  boundaryFaceFaces_[localPyrolysisFaceI].transfer(faceIDs);
112  boundaryFaceCells_[localPyrolysisFaceI].transfer(cellIDs);
113 
114  localPyrolysisFaceI++;
115  nLayers_ = nCells;
116  }
117  }
118 
119  boundaryFaceOppositeFace_.setSize(localPyrolysisFaceI);
120  boundaryFaceFaces_.setSize(localPyrolysisFaceI);
121  boundaryFaceCells_.setSize(localPyrolysisFaceI);
122 
123  surfaceScalarField& nMagSf = nMagSfPtr_();
124 
125  localPyrolysisFaceI = 0;
126  forAll(intCoupledPatchIDs_, i)
127  {
128  const label patchI = intCoupledPatchIDs_[i];
129  const polyPatch& ppCoupled = rbm[patchI];
130  const vectorField& pNormals = ppCoupled.faceNormals();
131  nMagSf.boundaryField()[patchI] =
132  regionMesh().Sf().boundaryField()[patchI] & pNormals;
133  forAll(pNormals, localFaceI)
134  {
135  const vector& n = pNormals[localFaceI];
136  const labelList& faces = boundaryFaceFaces_[localPyrolysisFaceI++];
137  forAll (faces, faceI)
138  {
139  const label faceID = faces[faceI];
140  nMagSf[faceID] = regionMesh().Sf()[faceID] & n;
141  }
142  }
143  }
144 }
145 
146 
147 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
148 
150 {
151  if (regionModel::read())
152  {
153  moveMesh_.readIfPresent("moveMesh", coeffs_);
154 
155  return true;
156  }
157  else
158  {
159  return false;
160  }
161 }
162 
163 
165 {
166  if (regionModel::read(dict))
167  {
168  moveMesh_.readIfPresent("moveMesh", coeffs_);
169 
170  return true;
171  }
172  else
173  {
174  return false;
175  }
176 }
177 
178 
180 (
181  const scalarList& deltaV,
182  const scalar minDelta
183 )
184 {
185  tmp<labelField> tcellMoveMap(new labelField(regionMesh().nCells(), 0));
186  labelField& cellMoveMap = tcellMoveMap();
187 
188  if (!moveMesh_)
189  {
190  return cellMoveMap;
191  }
192 
193  pointField oldPoints = regionMesh().points();
194  pointField newPoints = oldPoints;
195 
196  const polyBoundaryMesh& bm = regionMesh().boundaryMesh();
197 
198  label totalFaceId = 0;
199  forAll(intCoupledPatchIDs_, localPatchI)
200  {
201  label patchI = intCoupledPatchIDs_[localPatchI];
202  const polyPatch pp = bm[patchI];
203  const vectorField& cf = regionMesh().Cf().boundaryField()[patchI];
204 
205  forAll(pp, patchFaceI)
206  {
207  const labelList& faces = boundaryFaceFaces_[totalFaceId];
208  const labelList& cells = boundaryFaceCells_[totalFaceId];
209 
210  const vector n = pp.faceNormals()[patchFaceI];
211  const vector sf = pp.faceAreas()[patchFaceI];
212 
213  List<point> oldCf(faces.size() + 1);
214  oldCf[0] = cf[patchFaceI];
215  forAll(faces, i)
216  {
217  oldCf[i + 1] = regionMesh().faceCentres()[faces[i]];
218  }
219 
220  vector newDelta = vector::zero;
221  point nbrCf = oldCf[0];
222 
223  forAll(faces, i)
224  {
225  const label faceI = faces[i];
226  const label cellI = cells[i];
227 
228  const face f = regionMesh().faces()[faceI];
229 
230  newDelta += (deltaV[cellI]/mag(sf))*n;
231 
232  vector localDelta = vector::zero;
233  forAll(f, pti)
234  {
235  const label pointI = f[pti];
236 
237  if
238  (
239  mag((nbrCf - (oldPoints[pointI] + newDelta)) & n)
240  > minDelta
241  )
242  {
243  newPoints[pointI] = oldPoints[pointI] + newDelta;
244  localDelta = newDelta;
245  cellMoveMap[cellI] = 1;
246  }
247  }
248  nbrCf = oldCf[i + 1] + localDelta;
249  }
250  // Modify boundary
251  const label bFaceI = boundaryFaceOppositeFace_[totalFaceId];
252  const face f = regionMesh().faces()[bFaceI];
253  const label cellI = cells[cells.size() - 1];
254  newDelta += (deltaV[cellI]/mag(sf))*n;
255  forAll(f, pti)
256  {
257  const label pointI = f[pti];
258  if
259  (
260  mag((nbrCf - (oldPoints[pointI] + newDelta)) & n)
261  > minDelta
262  )
263  {
264  newPoints[pointI] = oldPoints[pointI] + newDelta;
265  cellMoveMap[cellI] = 1;
266  }
267  }
268  totalFaceId ++;
269  }
270  }
271  // Move points
272  regionMesh().movePoints(newPoints);
273 
274  return tcellMoveMap;
275 }
276 
277 
278 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
279 
280 Foam::regionModels::regionModel1D::regionModel1D
281 (
282  const fvMesh& mesh,
283  const word& regionType
284 )
285 :
286  regionModel(mesh, regionType),
287  boundaryFaceFaces_(),
288  boundaryFaceCells_(),
289  boundaryFaceOppositeFace_(),
290  nLayers_(0),
291  nMagSfPtr_(NULL),
292  moveMesh_(false)
293 {}
294 
295 
296 Foam::regionModels::regionModel1D::regionModel1D
297 (
298  const fvMesh& mesh,
299  const word& regionType,
300  const word& modelName,
301  bool readFields
302 )
303 :
304  regionModel(mesh, regionType, modelName, false),
305  boundaryFaceFaces_(regionMesh().nCells()),
306  boundaryFaceCells_(regionMesh().nCells()),
307  boundaryFaceOppositeFace_(regionMesh().nCells()),
308  nLayers_(0),
309  nMagSfPtr_(NULL),
310  moveMesh_(true)
311 {
312  if (active_)
313  {
314  constructMeshObjects();
315  initialise();
316 
317  if (readFields)
318  {
319  read();
320  }
321  }
322 }
323 
324 
325 Foam::regionModels::regionModel1D::regionModel1D
326 (
327  const fvMesh& mesh,
328  const word& regionType,
329  const word& modelName,
330  const dictionary& dict,
331  bool readFields
332 )
333 :
334  regionModel(mesh, regionType, modelName, dict, readFields),
335  boundaryFaceFaces_(regionMesh().nCells()),
336  boundaryFaceCells_(regionMesh().nCells()),
337  boundaryFaceOppositeFace_(regionMesh().nCells()),
338  nLayers_(0),
339  nMagSfPtr_(NULL),
340  moveMesh_(false)
341 {
342  if (active_)
343  {
344  constructMeshObjects();
345  initialise();
346 
347  if (readFields)
348  {
349  read(dict);
350  }
351  }
352 }
353 
354 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
355 
357 {}
358 
359 
360 // ************************************************************************* //
tmp< labelField > moveMesh(const scalarList &deltaV, const scalar minDelta=0.0)
Move mesh points according to change in cell volumes.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const vectorField::subField faceAreas() const
Return face normals.
Definition: polyPatch.C:290
dimensioned< scalar > mag(const dimensioned< Type > &)
labelList f(nPoints)
virtual ~regionModel1D()
Destructor.
virtual bool read()
Read control parameters from dictionary.
Field< label > labelField
Specialisation of Field<T> for label.
Definition: labelField.H:49
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
A class for handling words, derived from string.
Definition: word.H:59
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
defineTypeNameAndDebug(regionModel, 0)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
volScalarField sf(fieldObject, mesh)
Namespace for OpenFOAM.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
label n
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Foam::polyBoundaryMesh.
Field< vector > vectorField
Specialisation of Field<T> for vector.
#define forAll(list, i)
Definition: UList.H:421
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
const Field< PointType > & faceNormals() const
Return face normals for patch.
const cellShapeList & cells
List< label > labelList
A List of labels.
Definition: labelList.H:56
static const Vector zero
Definition: Vector.H:80
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
bool read(const char *, int32_t &)
Definition: int32IO.C:87
virtual bool read()
Read control parameters from dictionary.
Definition: regionModel.C:171
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
Definition: PtrList.H:118
word timeName
Definition: getTimeIndex.H:3
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53