50Foam::pointToPointPlanarInterpolation::calcCoordinateSystem
59 <<
"Need at least three non-colinear points"
60 <<
" to be able to interpolate."
69 scalar maxDist = ROOTVSMALL;
87 <<
"Cannot find any point that is different from first point"
88 <<
p0 <<
". Are all your points coincident?"
104 e2.removeCollinear(e1);
106 scalar magE2 =
mag(e2);
118 <<
"Cannot find points that define a plane with a valid normal."
119 <<
nl <<
"Have so far points " <<
p0 <<
" and " << p1
120 <<
". Are all your points on a single line instead of a plane?"
127 <<
" Used points " <<
p0 <<
' ' <<
points[index1]
129 <<
" to define coordinate system with normal " <<
n <<
endl;
131 return coordSystem::cartesian
140void Foam::pointToPointPlanarInterpolation::calcWeights
161 <<
"Did not find a corresponding sourcePoint for every face"
165 nearestVertex_.setSize(destPoints.size());
166 nearestVertexWeight_.setSize(destPoints.size());
169 nearestVertex_[i][0] = destToSource[i];
170 nearestVertex_[i][1] = -1;
171 nearestVertex_[i][2] = -1;
173 nearestVertexWeight_[i][0] = 1.0;
174 nearestVertexWeight_[i][1] = 0.0;
175 nearestVertexWeight_[i][2] = 0.0;
182 label v0 = nearestVertex_[i][0];
184 Pout<<
"For location " << destPoints[i]
185 <<
" sampling vertex " << v0
186 <<
" at:" << sourcePoints[v0]
187 <<
" distance:" <<
mag(sourcePoints[v0]-destPoints[i])
191 OBJstream str(
"destToSource.obj");
192 Pout<<
"pointToPointPlanarInterpolation::calcWeights :"
193 <<
" Dumping lines from face centres to original points to "
198 label v0 = nearestVertex_[i][0];
199 str.write(
linePointRef(destPoints[i], sourcePoints[v0]));
205 auto tlocalVertices = referenceCS_.localPosition(sourcePoints);
206 auto& localVertices = tlocalVertices.ref();
208 const boundBox bb(localVertices,
true);
209 const point bbMid(bb.centre());
212 <<
" Perturbing points with " << perturb_
213 <<
" fraction of a random position inside " << bb
214 <<
" to break any ties on regular meshes." <<
nl
226 List<vector2D> localVertices2D(localVertices.size());
229 localVertices2D[i][0] = localVertices[i][0];
230 localVertices2D[i][1] = localVertices[i][1];
235 auto tlocalFaceCentres = referenceCS_.localPosition(destPoints);
236 const pointField& localFaceCentres = tlocalFaceCentres();
240 Pout<<
"pointToPointPlanarInterpolation::calcWeights :"
241 <<
" Dumping triangulated surface to triangulation.stl" <<
endl;
242 s.write(
"triangulation.stl");
244 OBJstream str(
"localFaceCentres.obj");
245 Pout<<
"pointToPointPlanarInterpolation::calcWeights :"
246 <<
" Dumping face centres to " << str.
name() <<
endl;
248 forAll(localFaceCentres, i)
250 str.write(localFaceCentres[i]);
267 Pout<<
"source:" << i <<
" at:" << sourcePoints[i]
268 <<
" 2d:" << localVertices[i]
272 OBJstream str(
"stencil.obj");
273 Pout<<
"pointToPointPlanarInterpolation::calcWeights :"
274 <<
" Dumping stencil to " << str.
name() <<
endl;
278 label v0 = nearestVertex_[i][0];
279 label v1 = nearestVertex_[i][1];
280 label v2 = nearestVertex_[i][2];
282 Pout<<
"For location " << destPoints[i]
283 <<
" 2d:" << localFaceCentres[i]
284 <<
" sampling vertices" <<
nl
286 <<
" at:" << sourcePoints[v0]
287 <<
" weight:" << nearestVertexWeight_[i][0] <<
nl;
289 str.write(
linePointRef(destPoints[i], sourcePoints[v0]));
294 <<
" at:" << sourcePoints[v1]
295 <<
" weight:" << nearestVertexWeight_[i][1] <<
nl;
296 str.write(
linePointRef(destPoints[i], sourcePoints[v1]));
301 <<
" at:" << sourcePoints[v2]
302 <<
" weight:" << nearestVertexWeight_[i][2] <<
nl;
303 str.write(
linePointRef(destPoints[i], sourcePoints[v2]));
319 const scalar perturb,
320 const bool nearestOnly
324 nearestOnly_(nearestOnly),
326 nPoints_(sourcePoints.size())
330 referenceCS_ = calcCoordinateSystem(sourcePoints);
332 calcWeights(sourcePoints, destPoints);
346 referenceCS_(referenceCS),
347 nPoints_(sourcePoints.size())
349 calcWeights(sourcePoints, destPoints);
355 const scalar perturb,
356 const bool nearestOnly,
358 const label sourceSize,
364 nearestOnly_(nearestOnly),
365 referenceCS_(referenceCS),
366 nPoints_(sourceSize),
367 nearestVertex_(nearestVertex),
368 nearestVertexWeight_(nearestVertexWeight)
398 names[i] = times[i].name();
407 const label startSampleTime,
408 const scalar timeVal,
413 lo = startSampleTime;
416 for (label i = startSampleTime+1; i < times.
size(); i++)
418 if (times[i].value() > timeVal)
439 if (lo < times.
size()-1)
449 Pout<<
"findTime : Found time " << timeVal <<
" after"
450 <<
" index:" << lo <<
" time:" << times[lo].value()
455 Pout<<
"findTime : Found time " << timeVal <<
" inbetween"
456 <<
" index:" << lo <<
" time:" << times[lo].value()
457 <<
" and index:" << hi <<
" time:" << times[hi].value()
A 1D vector of objects of type <T> with a fixed length <N>.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
virtual const fileName & name() const
Get the name of the stream.
Type position(const Type &start, const Type &end)
Return a sample on the interval [start,end].
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
void size(const label n)
Older name for setAddressableSize.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
A Cartesian coordinate system.
Base class for coordinate system specification, the default coordinate system type is cartesian .
Interpolates between two sets of unstructured points using 2D Delaunay triangulation....
static bool findTime(const instantList ×, const label startSampleTime, const scalar timeVal, label &lo, label &hi)
Helper: find time. Return true if successful.
autoPtr< pointToPointPlanarInterpolation > clone() const
Construct and return a clone.
static wordList timeNames(const instantList &)
Helper: extract words of times.
bool nearestOnly() const
Whether to use nearest point only (avoids triangulation, projection)
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.
const volScalarField & p0
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Determine correspondence between points. See below.
#define DebugInFunction
Report an information message using Foam::Info.
List< label > labelList
A List of labels.
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
line< point, const point & > linePointRef
A line using referred points.
vectorField pointField
pointField is a vectorField.
vector point
Point is a vector.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
bool matchPoints(const UList< point > &pts0, const UList< point > &pts1, const UList< scalar > &matchDistance, const bool verbose, labelList &from0To1, const point &origin=point::zero)
Determine correspondence between pointFields. Gets passed.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
errorManipArg< error, int > exit(error &err, const int errNo=1)
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.