uniformGrid.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) 2012-2015 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "uniformGrid.H"
30
31// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32
33namespace Foam
34{
35 defineTypeNameAndDebug(uniformGrid, 0);
36 addToRunTimeSelectionTable(initialPointsMethod, uniformGrid, dictionary);
37}
38
39
40// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41
43(
44 const dictionary& initialPointsDict,
45 const Time& runTime,
46 Random& rndGen,
47 const conformationSurfaces& geometryToConformTo,
48 const cellShapeControl& cellShapeControls,
49 const autoPtr<backgroundMeshDecomposition>& decomposition
50)
51:
52 initialPointsMethod
53 (
54 typeName,
55 initialPointsDict,
56 runTime,
57 rndGen,
58 geometryToConformTo,
59 cellShapeControls,
60 decomposition
61 ),
62 initialCellSize_(detailsDict().get<scalar>("initialCellSize")),
63 randomiseInitialGrid_(detailsDict().get<Switch>("randomiseInitialGrid")),
64 randomPerturbationCoeff_
65 (
66 detailsDict().get<scalar>("randomPerturbationCoeff")
67 )
68{}
69
70
71// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
72
74{
75 boundBox bb;
76
77 // Pick up the bounds of this processor, or the whole geometry, depending
78 // on whether this is a parallel run.
79 if (Pstream::parRun())
80 {
82 }
83 else
84 {
86 }
87
88 scalar x0 = bb.min().x();
89 scalar xR = bb.max().x() - x0;
90 label ni = label(xR/initialCellSize_);
91
92 scalar y0 = bb.min().y();
93 scalar yR = bb.max().y() - y0;
94 label nj = label(yR/initialCellSize_);
95
96 scalar z0 = bb.min().z();
97 scalar zR = bb.max().z() - z0;
98 label nk = label(zR/initialCellSize_);
99
100 vector delta(xR/ni, yR/nj, zR/nk);
101
102 delta *= pow((1.0),-(1.0/3.0));
103
104 scalar pert = randomPerturbationCoeff_*cmptMin(delta);
105
106 // Initialise points list
107 DynamicList<Vb::Point> initialPoints(scalar(ni)*nj*nk/10);
108
109 for (label i = 0; i < ni; i++)
110 {
111 for (label j = 0; j < nj; j++)
112 {
113 // Generating, testing and adding points one line at a time to
114 // reduce the memory requirement for cases with bounding boxes that
115 // are very large in comparison to the volume to be filled
116
117 label pI = 0;
118
119 pointField points(nk);
120
121 for (label k = 0; k < nk; k++)
122 {
123 point p
124 (
125 x0 + (i + 0.5)*delta.x(),
126 y0 + (j + 0.5)*delta.y(),
127 z0 + (k + 0.5)*delta.z()
128 );
129
130 if (randomiseInitialGrid_)
131 {
132 p.x() += pert*(rndGen().sample01<scalar>() - 0.5);
133 p.y() += pert*(rndGen().sample01<scalar>() - 0.5);
134 p.z() += pert*(rndGen().sample01<scalar>() - 0.5);
135 }
136
137 if
138 (
140 && !decomposition().positionOnThisProcessor(p)
141 )
142 {
143 // Skip this point if, in parallel, this position is not on
144 // this processor.
145 continue;
146 }
147
148 points[pI++] = p;
149 }
150
151 points.setSize(pI);
152
153 Field<bool> insidePoints =
155 (
156 points,
158 *sqr
159 (
160 cellShapeControls().cellSize(points)
161 )
162 );
163
165 {
166 if (insidePoints[i])
167 {
168 const point& p(points[i]);
169
170 initialPoints.append(Vb::Point(p.x(), p.y(), p.z()));
171 }
172 }
173 }
174 }
175
176 return initialPoints.shrink();
177}
178
179
180// ************************************************************************* //
label k
scalar delta
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:77
void setSize(const label n)
Alias for resize()
Definition: List.H:218
Type sample01()
Return a sample whose components lie in the range [0,1].
static bool & parRun() noexcept
Test if this a parallel run.
Definition: UPstream.H:433
const treeBoundBox & procBounds() const
Return the boundBox of this processor.
const point & min() const
Minimum describing the bounding box.
Definition: boundBoxI.H:91
Field< bool > wellInside(const pointField &samplePts, const scalarField &testDistSqr) const
Check if point is inside surfaces to conform to by at least.
const treeBoundBox & globalBounds() const
Return the global bounds.
scalar minimumSurfaceDistanceCoeffSqr_
Only allow the placement of initial points that are within the.
const conformationSurfaces & geometryToConformTo() const
const cellShapeControl & cellShapeControls() const
const backgroundMeshDecomposition & decomposition() const
Generate a uniform grid of points inside the surfaces to be conformed to of the conformalVoronoiMesh.
Definition: uniformGrid.H:55
virtual List< Vb::Point > initialPoints() const
Return the initial points for the conformalVoronoiMesh.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
volScalarField & p
engineTime & runTime
const pointField & points
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
Namespace for OpenFOAM.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar y0(const dimensionedScalar &ds)
void cmptMin(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:44
vector point
Point is a vector.
Definition: point.H:43
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
insidePoints((1 2 3))
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Random rndGen
Definition: createFields.H:23