surfaceSlipDisplacementPointPatchVectorField.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-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 
28 #include "Time.H"
29 #include "transformField.H"
30 #include "fvMesh.H"
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 
40 template<>
41 const char*
42 NamedEnum<surfaceSlipDisplacementPointPatchVectorField::projectMode, 3>::
43 names[] =
44 {
45  "nearest",
46  "pointNormal",
47  "fixedNormal"
48 };
49 
51  surfaceSlipDisplacementPointPatchVectorField::projectModeNames_;
52 
53 
54 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
55 
56 void surfaceSlipDisplacementPointPatchVectorField::calcProjection
57 (
58  vectorField& displacement
59 ) const
60 {
61  const polyMesh& mesh = patch().boundaryMesh().mesh()();
62  const pointField& localPoints = patch().localPoints();
63  const labelList& meshPoints = patch().meshPoints();
64 
65  //const scalar deltaT = mesh.time().deltaTValue();
66 
67  // Construct large enough vector in direction of projectDir so
68  // we're guaranteed to hit something.
69 
70  //- Per point projection vector:
71  const scalar projectLen = mag(mesh.bounds().max()-mesh.bounds().min());
72 
73  // For case of fixed projection vector:
74  vector projectVec(0, 0, 0);
75  if (projectMode_ == FIXEDNORMAL)
76  {
77  vector n = projectDir_/mag(projectDir_);
78  projectVec = projectLen*n;
79  }
80 
81 
82  // Get fixed points (bit of a hack)
83  const pointZone* zonePtr = NULL;
84 
85  if (frozenPointsZone_.size() > 0)
86  {
87  const pointZoneMesh& pZones = mesh.pointZones();
88 
89  zonePtr = &pZones[frozenPointsZone_];
90 
91  Pout<< "surfaceSlipDisplacementPointPatchVectorField : Fixing all "
92  << zonePtr->size() << " points in pointZone " << zonePtr->name()
93  << endl;
94  }
95 
96  // Get the starting locations from the motionSolver
97  const pointField& points0 = mesh.lookupObject<displacementMotionSolver>
98  (
99  "dynamicMeshDict"
100  ).points0();
101 
102 
103  pointField start(meshPoints.size());
104  forAll(start, i)
105  {
106  start[i] = points0[meshPoints[i]] + displacement[i];
107  }
108 
109  label nNotProjected = 0;
110 
111  if (projectMode_ == NEAREST)
112  {
113  List<pointIndexHit> nearest;
114  labelList hitSurfaces;
116  (
117  start,
118  scalarField(start.size(), sqr(projectLen)),
119  hitSurfaces,
120  nearest
121  );
122 
123  forAll(nearest, i)
124  {
125  if (zonePtr && (zonePtr->whichPoint(meshPoints[i]) >= 0))
126  {
127  // Fixed point. Reset to point0 location.
128  displacement[i] = points0[meshPoints[i]] - localPoints[i];
129  }
130  else if (nearest[i].hit())
131  {
132  displacement[i] =
133  nearest[i].hitPoint()
134  - points0[meshPoints[i]];
135  }
136  else
137  {
138  nNotProjected++;
139 
140  if (debug)
141  {
142  Pout<< " point:" << meshPoints[i]
143  << " coord:" << localPoints[i]
144  << " did not find any surface within " << projectLen
145  << endl;
146  }
147  }
148  }
149  }
150  else
151  {
152  // Do tests on all points. Combine later on.
153 
154  // 1. Check if already on surface
155  List<pointIndexHit> nearest;
156  {
157  labelList nearestSurface;
159  (
160  start,
161  scalarField(start.size(), sqr(SMALL)),
162  nearestSurface,
163  nearest
164  );
165  }
166 
167  // 2. intersection. (combined later on with information from nearest
168  // above)
169  vectorField projectVecs(start.size(), projectVec);
170 
171  if (projectMode_ == POINTNORMAL)
172  {
173  projectVecs = projectLen*patch().pointNormals();
174  }
175 
176  // Knock out any wedge component
177  scalarField offset(start.size(), 0.0);
178  if (wedgePlane_ >= 0 && wedgePlane_ < vector::nComponents)
179  {
180  forAll(offset, i)
181  {
182  offset[i] = start[i][wedgePlane_];
183  start[i][wedgePlane_] = 0;
184  projectVecs[i][wedgePlane_] = 0;
185  }
186  }
187 
188  List<pointIndexHit> rightHit;
189  {
190  labelList rightSurf;
192  (
193  start,
194  start+projectVecs,
195  rightSurf,
196  rightHit
197  );
198  }
199 
200  List<pointIndexHit> leftHit;
201  {
202  labelList leftSurf;
204  (
205  start,
206  start-projectVecs,
207  leftSurf,
208  leftHit
209  );
210  }
211 
212  // 3. Choose either -fixed, nearest, right, left.
213  forAll(displacement, i)
214  {
215  if (zonePtr && (zonePtr->whichPoint(meshPoints[i]) >= 0))
216  {
217  // Fixed point. Reset to point0 location.
218  displacement[i] = points0[meshPoints[i]] - localPoints[i];
219  }
220  else if (nearest[i].hit())
221  {
222  // Found nearest.
223  displacement[i] =
224  nearest[i].hitPoint()
225  - points0[meshPoints[i]];
226  }
227  else
228  {
229  pointIndexHit interPt;
230 
231  if (rightHit[i].hit())
232  {
233  if (leftHit[i].hit())
234  {
235  if
236  (
237  magSqr(rightHit[i].hitPoint()-start[i])
238  < magSqr(leftHit[i].hitPoint()-start[i])
239  )
240  {
241  interPt = rightHit[i];
242  }
243  else
244  {
245  interPt = leftHit[i];
246  }
247  }
248  else
249  {
250  interPt = rightHit[i];
251  }
252  }
253  else
254  {
255  if (leftHit[i].hit())
256  {
257  interPt = leftHit[i];
258  }
259  }
260 
261 
262  if (interPt.hit())
263  {
264  if (wedgePlane_ >= 0 && wedgePlane_ < vector::nComponents)
265  {
266  interPt.rawPoint()[wedgePlane_] += offset[i];
267  }
268  displacement[i] = interPt.rawPoint()-points0[meshPoints[i]];
269  }
270  else
271  {
272  nNotProjected++;
273 
274  if (debug)
275  {
276  Pout<< " point:" << meshPoints[i]
277  << " coord:" << localPoints[i]
278  << " did not find any intersection between"
279  << " ray from " << start[i]-projectVecs[i]
280  << " to " << start[i]+projectVecs[i] << endl;
281  }
282  }
283  }
284  }
285  }
286 
287  reduce(nNotProjected, sumOp<label>());
288 
289  if (nNotProjected > 0)
290  {
291  Info<< "surfaceSlipDisplacement :"
292  << " on patch " << patch().name()
293  << " did not project " << nNotProjected
294  << " out of " << returnReduce(localPoints.size(), sumOp<label>())
295  << " points." << endl;
296  }
297 }
298 
299 
300 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
301 
304 (
305  const pointPatch& p,
307 )
308 :
309  pointPatchVectorField(p, iF),
310  projectMode_(NEAREST),
311  projectDir_(Zero),
312  wedgePlane_(-1)
313 {}
314 
315 
318 (
319  const pointPatch& p,
321  const dictionary& dict
322 )
323 :
324  pointPatchVectorField(p, iF, dict),
325  surfacesDict_(dict.subDict("geometry")),
326  projectMode_(projectModeNames_.read(dict.lookup("projectMode"))),
327  projectDir_(dict.lookup("projectDirection")),
328  wedgePlane_(dict.lookupOrDefault("wedgePlane", -1)),
329  frozenPointsZone_(dict.lookupOrDefault("frozenPointsZone", word::null))
330 {}
331 
332 
335 (
337  const pointPatch& p,
339  const pointPatchFieldMapper&
340 )
341 :
342  pointPatchVectorField(p, iF),
343  surfacesDict_(ppf.surfacesDict_),
344  projectMode_(ppf.projectMode_),
345  projectDir_(ppf.projectDir_),
346  wedgePlane_(ppf.wedgePlane_),
347  frozenPointsZone_(ppf.frozenPointsZone_)
348 {}
349 
350 
353 (
355 )
356 :
358  surfacesDict_(ppf.surfacesDict_),
359  projectMode_(ppf.projectMode_),
360  projectDir_(ppf.projectDir_),
361  wedgePlane_(ppf.wedgePlane_),
362  frozenPointsZone_(ppf.frozenPointsZone_)
363 {}
364 
365 
368 (
371 )
372 :
373  pointPatchVectorField(ppf, iF),
374  surfacesDict_(ppf.surfacesDict_),
375  projectMode_(ppf.projectMode_),
376  projectDir_(ppf.projectDir_),
377  wedgePlane_(ppf.wedgePlane_),
378  frozenPointsZone_(ppf.frozenPointsZone_)
379 {}
380 
381 
382 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
383 
384 const searchableSurfaces&
386 {
387  if (surfacesPtr_.empty())
388  {
389  surfacesPtr_.reset
390  (
392  (
393  IOobject
394  (
395  "abc", // dummy name
396  db().time().constant(), // directory
397  "triSurface", // instance
398  db().time(), // registry
401  ),
402  surfacesDict_,
403  true // use single region naming shortcut
404  )
405  );
406  }
407  return surfacesPtr_();
408 }
409 
410 
412 (
413  const Pstream::commsTypes commsType
414 )
415 {
416  vectorField displacement(this->patchInternalField());
417 
418  // Calculate displacement to project points onto surface
419  calcProjection(displacement);
420 
421  // Get internal field to insert values into
422  Field<vector>& iF = const_cast<Field<vector>&>(this->primitiveField());
423 
424  //setInInternalField(iF, motionU);
425  setInInternalField(iF, displacement);
426 
428 }
429 
430 
432 {
434  os.writeKeyword("geometry") << surfacesDict_
435  << token::END_STATEMENT << nl;
436  os.writeKeyword("projectMode") << projectModeNames_[projectMode_]
437  << token::END_STATEMENT << nl;
438  os.writeKeyword("projectDirection") << projectDir_
439  << token::END_STATEMENT << nl;
440  os.writeKeyword("wedgePlane") << wedgePlane_
441  << token::END_STATEMENT << nl;
442  if (frozenPointsZone_ != word::null)
443  {
444  os.writeKeyword("frozenPointsZone") << frozenPointsZone_
445  << token::END_STATEMENT << nl;
446  }
447 }
448 
449 
450 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
451 
453 (
456 );
457 
458 
459 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
460 
461 } // End namespace Foam
462 
463 // ************************************************************************* //
Virtual base class for displacement motion solver.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
Displacement follows a triSurface. Use in a displacementMotionSolver as a bc on the pointDisplacement...
const boundBox & bounds() const
Return mesh bounding box.
Definition: polyMesh.H:427
const point & min() const
Minimum describing the bounding box.
Definition: boundBoxI.H:54
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
commsTypes
Types of communications.
Definition: UPstream.H:64
bool hit() const
Is there a hit.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:633
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::blocking)
Update the patch field.
Foam::pointPatchFieldMapper.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
void setInInternalField(Field< Type1 > &iF, const Field< Type1 > &pF, const labelList &meshPoints) const
Given the internal field and a patch field,.
const point & max() const
Maximum describing the bounding box.
Definition: boundBoxI.H:60
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:53
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:52
Macros for easy insertion into run-time selection tables.
const searchableSurfaces & surfaces() const
Surface to follow. Demand loads surfaceNames.
const objectRegistry & db() const
Return local objectRegistry.
const pointMesh & mesh() const
Return the mesh reference.
Spatial transformation functions for primitive fields.
const word & name() const
Return name.
Definition: zone.H:150
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
static const direction nComponents
Number of components in this vector space.
Definition: VectorSpace.H:96
dynamicFvMesh & mesh
void findNearest(const pointField &, const scalarField &nearestDistSqr, labelList &surfaces, List< pointIndexHit > &) const
Find nearest. Return -1 (and a miss()) or surface and nearest.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static const word null
An empty word.
Definition: word.H:77
Container for searchableSurfaces.
const pointBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: pointPatch.H:139
static const zero Zero
Definition: zero.H:91
virtual const word & name() const =0
Return name.
const Point & rawPoint() const
Return point with no checking.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
void findAnyIntersection(const pointField &start, const pointField &end, labelList &surfaces, List< pointIndexHit > &) const
Find any intersection. Return hit point information and.
static const char nl
Definition: Ostream.H:262
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
label whichPoint(const label globalPointID) const
Helper function to re-direct to zone::localID(...)
Definition: pointZone.C:126
const pointZoneMesh & pointZones() const
Return point zone mesh.
Definition: polyMesh.H:457
pointPatchField< vector > pointPatchVectorField
virtual const vectorField & pointNormals() const =0
Return point normals.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
virtual const vectorField & localPoints() const =0
Return mesh points.
const pointPatch & patch() const
Return patch.
A subset of mesh points. The labels of points in the zone can be obtained from the addressing() list...
Definition: pointZone.H:62
Constant dispersed-phase particle diameter model.
Basic pointPatch represents a set of points from the mesh.
Definition: pointPatch.H:56
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
messageStream Info
const Field< Type > & primitiveField() const
Return internal field reference.
makePointPatchTypeField(pointPatchVectorField, solidBodyMotionDisplacementPointPatchVectorField)
dimensioned< scalar > mag(const dimensioned< Type > &)
label n
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
virtual void write(Ostream &) const
Write.
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::blocking)
Evaluate the patch field.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Enum read(Istream &) const
Read a word from Istream and return the corresponding.
Definition: NamedEnum.C:61
tmp< Field< Type > > patchInternalField() const
Return field created from appropriate internal field values.
surfaceSlipDisplacementPointPatchVectorField(const pointPatch &, const DimensionedField< vector, pointMesh > &)
Construct from patch and internal field.
Namespace for OpenFOAM.
virtual const labelList & meshPoints() const =0
Return mesh points.
IOporosityModelList pZones(mesh)
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451