pointFile.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-2016 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 "pointFile.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 defineTypeNameAndDebug(pointFile, 0);
37 addToRunTimeSelectionTable(initialPointsMethod, pointFile, 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  pointFileName_(detailsDict().lookup("pointFile")),
62  insideOutsideCheck_(detailsDict().lookup("insideOutsideCheck")),
63  randomiseInitialGrid_(detailsDict().lookup("randomiseInitialGrid")),
64  randomPerturbationCoeff_
65  (
66  readScalar(detailsDict().lookup("randomPerturbationCoeff"))
67  )
68 {
69  Info<< " Inside/Outside check is " << insideOutsideCheck_.asText()
70  << endl;
71 }
72 
73 
74 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
75 
76 List<Vb::Point> pointFile::initialPoints() const
77 {
79  (
80  IOobject
81  (
82  pointFileName_.name(),
83  time().timeName(),
84  time(),
87  )
88  );
89 
90  Info<< " Inserting points from file " << pointFileName_ << endl;
91 
92  if (points.empty())
93  {
95  << "Point file contain no points"
96  << exit(FatalError) << endl;
97  }
98 
99  if (Pstream::parRun())
100  {
101  // Testing filePath to see if the file originated in a processor
102  // directory, if so, assume that the points in each processor file
103  // are unique. They are unlikely to belong on the current
104  // processor as the background mesh is unlikely to be the same.
105 
106  const bool isParentFile = (points.objectPath() != points.filePath());
107 
108  if (!isParentFile)
109  {
111  }
112  else
113  {
114  // Otherwise, this is assumed to be points covering the whole
115  // domain, so filter the points to be only those on this processor
116  boolList procPt(decomposition().positionOnThisProcessor(points));
117 
118  List<boolList> allProcPt(Pstream::nProcs());
119 
120  allProcPt[Pstream::myProcNo()] = procPt;
121 
122  Pstream::gatherList(allProcPt);
123 
124  Pstream::scatterList(allProcPt);
125 
126  forAll(procPt, ptI)
127  {
128  bool foundAlready = false;
129 
130  forAll(allProcPt, proci)
131  {
132  // If a processor with a lower index has found this point
133  // to insert already, defer to it and don't insert.
134  if (foundAlready)
135  {
136  allProcPt[proci][ptI] = false;
137  }
138  else if (allProcPt[proci][ptI])
139  {
140  foundAlready = true;
141  }
142  }
143  }
144 
145  procPt = allProcPt[Pstream::myProcNo()];
146 
147  inplaceSubset(procPt, points);
148  }
149  }
150 
151  Field<bool> insidePoints(points.size(), true);
152 
153  if (insideOutsideCheck_)
154  {
156  (
157  points,
159  *sqr(cellShapeControls().cellSize(points))
160  );
161  }
162 
163  DynamicList<Vb::Point> initialPoints(insidePoints.size()/10);
164 
165  forAll(insidePoints, i)
166  {
167  if (insidePoints[i])
168  {
169  point& p = points[i];
170 
171  if (randomiseInitialGrid_)
172  {
173  p.x() += randomPerturbationCoeff_*(rndGen().scalar01() - 0.5);
174  p.y() += randomPerturbationCoeff_*(rndGen().scalar01() - 0.5);
175  p.z() += randomPerturbationCoeff_*(rndGen().scalar01() - 0.5);
176  }
177 
178  initialPoints.append(Vb::Point(p.x(), p.y(), p.z()));
179  }
180  }
181 
182  initialPoints.shrink();
183 
184  label nPointsRejected = points.size() - initialPoints.size();
185 
186  if (Pstream::parRun())
187  {
188  reduce(nPointsRejected, sumOp<label>());
189  }
190 
191  if (nPointsRejected)
192  {
193  Info<< " " << nPointsRejected << " points rejected from "
194  << pointFileName_.name() << endl;
195  }
196 
197  return initialPoints;
198 }
199 
200 
201 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
202 
203 } // End namespace Foam
204 
205 // ************************************************************************* //
void inplaceSubset(const UList< T > &select, const T &value, ListType &)
Inplace extract elements of List when select is a certain value.
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
#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.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
vectorIOField pointIOField
pointIOField is a vectorIOField.
Definition: pointIOField.H:42
dimensionedSymmTensor sqr(const dimensionedVector &dv)
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:417
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
Macros for easy insertion into run-time selection tables.
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
pointFile(const dictionary &initialPointsDict, const Time &runTime, Random &rndGen, const conformationSurfaces &geometryToConformTo, const cellShapeControl &cellShapeControls, const autoPtr< backgroundMeshDecomposition > &decomposition)
Construct from components.
stressControl lookup("compactNormalStress") >> compactNormalStress
const pointField & points
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)
insidePoints((1 2 3))
defineTypeNameAndDebug(combustionModel, 0)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:393
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:399
vector point
Point is a vector.
Definition: point.H:41
virtual List< Vb::Point > initialPoints() const =0
Return the initial points for the conformalVoronoiMesh.
const conformationSurfaces & geometryToConformTo() const
messageStream Info
const Time & time() const
autoPtr< mapDistribute > distributePoints(List< PointType > &points) const
Distribute supplied the points to the appropriate processor.
Field< bool > wellInside(const pointField &samplePts, const scalarField &testDistSqr) const
Check if point is inside surfaces to conform to by at least.
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
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.