bodyCentredCubic.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) 2012-2015 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 "bodyCentredCubic.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 defineTypeNameAndDebug(bodyCentredCubic, 0);
37 addToRunTimeSelectionTable(initialPointsMethod, bodyCentredCubic, dictionary);
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
42 (
43  const dictionary& initialPointsDict,
44  const Time& runTime,
45  Random& rndGen,
46  const conformationSurfaces& geometryToConformTo,
47  const cellShapeControl& cellShapeControls,
48  const autoPtr<backgroundMeshDecomposition>& decomposition
49 )
50 :
51  initialPointsMethod
52  (
53  typeName,
54  initialPointsDict,
55  runTime,
56  rndGen,
57  geometryToConformTo,
58  cellShapeControls,
59  decomposition
60  ),
61  initialCellSize_(readScalar(detailsDict().lookup("initialCellSize"))),
62  randomiseInitialGrid_(detailsDict().lookup("randomiseInitialGrid")),
63  randomPerturbationCoeff_
64  (
65  readScalar(detailsDict().lookup("randomPerturbationCoeff"))
66  )
67 {}
68 
69 
70 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
71 
72 List<Vb::Point> bodyCentredCubic::initialPoints() const
73 {
74  boundBox bb;
75 
76  // Pick up the bounds of this processor, or the whole geometry, depending
77  // on whether this is a parallel run.
78  if (Pstream::parRun())
79  {
80  bb = decomposition().procBounds();
81  }
82  else
83  {
85  }
86 
87  scalar x0 = bb.min().x();
88  scalar xR = bb.max().x() - x0;
89  label ni = label(xR/initialCellSize_);
90 
91  scalar y0 = bb.min().y();
92  scalar yR = bb.max().y() - y0;
93  label nj = label(yR/initialCellSize_);
94 
95  scalar z0 = bb.min().z();
96  scalar zR = bb.max().z() - z0;
97  label nk = label(zR/initialCellSize_);
98 
99  vector delta(xR/ni, yR/nj, zR/nk);
100 
101  delta *= pow((1.0/2.0),-(1.0/3.0));
102 
103  scalar pert = randomPerturbationCoeff_*cmptMin(delta);
104 
105  DynamicList<Vb::Point> initialPoints(ni*nj*nk/10);
106 
107  for (label i = 0; i < ni; i++)
108  {
109  for (label j = 0; j < nj; j++)
110  {
111  // Generating, testing and adding points one line at a time to
112  // reduce the memory requirement for cases with bounding boxes that
113  // are very large in comparison to the volume to be filled
114 
115  label pI = 0;
116 
117  pointField points(2*nk);
118 
119  for (label k = 0; k < nk; k++)
120  {
121  point pA
122  (
123  x0 + i*delta.x(),
124  y0 + j*delta.y(),
125  z0 + k*delta.z()
126  );
127 
128  point pB = pA + 0.5*delta;
129 
130  if (randomiseInitialGrid_)
131  {
132  pA.x() += pert*(rndGen().scalar01() - 0.5);
133  pA.y() += pert*(rndGen().scalar01() - 0.5);
134  pA.z() += pert*(rndGen().scalar01() - 0.5);
135  }
136 
137  if (Pstream::parRun())
138  {
139  if (decomposition().positionOnThisProcessor(pA))
140  {
141  // Add this point in parallel only if this position is
142  // on this processor.
143  points[pI++] = pA;
144  }
145  }
146  else
147  {
148  points[pI++] = pA;
149  }
150 
151  if (randomiseInitialGrid_)
152  {
153  pB.x() += pert*(rndGen().scalar01() - 0.5);
154  pB.y() += pert*(rndGen().scalar01() - 0.5);
155  pB.z() += pert*(rndGen().scalar01() - 0.5);
156  }
157 
158  if (Pstream::parRun())
159  {
160  if (decomposition().positionOnThisProcessor(pB))
161  {
162  // Add this point in parallel only if this position is
163  // on this processor.
164  points[pI++] = pB;
165  }
166  }
167  else
168  {
169  points[pI++] = pB;
170  }
171  }
172 
173  points.setSize(pI);
174 
175  Field<bool> insidePoints =
177  (
178  points,
180  *sqr(cellShapeControls().cellSize(points))
181  );
182 
183  forAll(insidePoints, i)
184  {
185  if (insidePoints[i])
186  {
187  const point& p(points[i]);
188 
189  initialPoints.append(Vb::Point(p.x(), p.y(), p.z()));
190  }
191  }
192  }
193  }
194 
195  return initialPoints.shrink();
196 }
197 
198 
199 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
200 
201 } // End namespace Foam
202 
203 // ************************************************************************* //
cachedRandom rndGen(label(0),-1)
scalar delta
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const cellShapeControl & cellShapeControls() const
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
virtual List< Vb::Point > initialPoints() const
Return the initial points for the conformalVoronoiMesh.
const Cmpt & x() const
Definition: VectorI.H:75
const point & min() const
Minimum describing the bounding box.
Definition: boundBoxI.H:54
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar y0(const dimensionedScalar &ds)
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
label k
Boltzmann constant.
Macros for easy insertion into run-time selection tables.
bodyCentredCubic(const dictionary &initialPointsDict, const Time &runTime, Random &rndGen, const conformationSurfaces &geometryToConformTo, const cellShapeControl &cellShapeControls, const autoPtr< backgroundMeshDecomposition > &decomposition)
Construct from components.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
stressControl lookup("compactNormalStress") >> compactNormalStress
const pointField & points
virtual List< Vb::Point > initialPoints() const
Return the initial points for the conformalVoronoiMesh.
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:393
vector point
Point is a vector.
Definition: point.H:41
const conformationSurfaces & geometryToConformTo() const
volScalarField & p
Field< bool > wellInside(const pointField &samplePts, const scalarField &testDistSqr) const
Check if point is inside surfaces to conform to by at least.
void cmptMin(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
const treeBoundBox & globalBounds() const
Return the global bounds.
const treeBoundBox & procBounds() const
Return the boundBox of this processor.
const backgroundMeshDecomposition & decomposition() const
scalar scalar01()
Scalar [0..1] (so including 0,1)
Definition: Random.C:67
scalar minimumSurfaceDistanceCoeffSqr_
Only allow the placement of initial points that are within the.
Namespace for OpenFOAM.