componentDisplacementMotionSolver.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 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 
27 #include "mapPolyMesh.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33  defineTypeNameAndDebug(componentDisplacementMotionSolver, 0);
34 }
35 
36 
37 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
38 
39 Foam::direction Foam::componentDisplacementMotionSolver::cmpt
40 (
41  const word& cmptName
42 ) const
43 {
44  if (cmptName == "x")
45  {
46  return vector::X;
47  }
48  else if (cmptName == "y")
49  {
50  return vector::Y;
51  }
52  else if (cmptName == "z")
53  {
54  return vector::Z;
55  }
56  else
57  {
59  (
60  "componentDisplacementMotionSolver::"
61  "componentDisplacementMotionSolver"
62  "(const polyMesh& mesh, const IOdictionary&)"
63  ) << "Given component name " << cmptName << " should be x, y or z"
64  << exit(FatalError);
65 
66  return 0;
67  }
68 }
69 
70 
71 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
72 
73 Foam::componentDisplacementMotionSolver::componentDisplacementMotionSolver
74 (
75  const polyMesh& mesh,
76  const IOdictionary& dict,
77  const word& type
78 )
79 :
80  motionSolver(mesh, dict, type),
81  cmptName_(coeffDict().lookup("component")),
82  cmpt_(cmpt(cmptName_)),
83  points0_
84  (
86  (
87  IOobject
88  (
89  "points",
90  time().constant(),
92  mesh,
93  IOobject::MUST_READ,
95  false
96  )
97  ).component(cmpt_)
98  ),
99  pointDisplacement_
100  (
101  IOobject
102  (
103  "pointDisplacement" + cmptName_,
104  mesh.time().timeName(),
105  mesh,
108  ),
109  pointMesh::New(mesh)
110  )
111 {
112  if (points0_.size() != mesh.nPoints())
113  {
115  (
116  "componentDisplacementMotionSolver::"
117  "componentDisplacementMotionSolver\n"
118  "(\n"
119  " const polyMesh&,\n"
120  " const IOdictionary&\n"
121  ")"
122  ) << "Number of points in mesh " << mesh.nPoints()
123  << " differs from number of points " << points0_.size()
124  << " read from file "
125  <<
126  IOobject
127  (
128  "points",
129  mesh.time().constant(),
131  mesh,
134  false
135  ).filePath()
136  << exit(FatalError);
137  }
138 }
139 
140 
141 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
142 
144 {}
145 
146 
147 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
148 
150 {
151  // No local data to update
152 }
153 
154 
156 {
157  // pointMesh already updates pointFields.
158 
160 
161  // Map points0_. Bit special since we somehow have to come up with
162  // a sensible points0 position for introduced points.
163  // Find out scaling between points0 and current points
164 
165  // Get the new points either from the map or the mesh
166  const scalarField points
167  (
168  mpm.hasMotionPoints()
169  ? mpm.preMotionPoints().component(cmpt_)
170  : mesh().points().component(cmpt_)
171  );
172 
173  // Get extents of points0 and points and determine scale
174  const scalar scale =
175  (gMax(points0_)-gMin(points0_))
176  /(gMax(points)-gMin(points));
177 
178  scalarField newPoints0(mpm.pointMap().size());
179 
180  forAll(newPoints0, pointI)
181  {
182  label oldPointI = mpm.pointMap()[pointI];
183 
184  if (oldPointI >= 0)
185  {
186  label masterPointI = mpm.reversePointMap()[oldPointI];
187 
188  if (masterPointI == pointI)
189  {
190  newPoints0[pointI] = points0_[oldPointI];
191  }
192  else
193  {
194  // New point. Assume motion is scaling.
195  newPoints0[pointI] =
196  points0_[oldPointI]
197  + scale*(points[pointI]-points[masterPointI]);
198  }
199  }
200  else
201  {
203  (
204  "displacementLaplacianFvMotionSolver::updateMesh"
205  "(const mapPolyMesh& mpm)"
206  ) << "Cannot work out coordinates of introduced vertices."
207  << " New vertex " << pointI << " at coordinate "
208  << points[pointI] << exit(FatalError);
209  }
210  }
211  points0_.transfer(newPoints0);
212 }
213 
214 
215 // ************************************************************************* //
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:309
const pointField & points
unsigned char direction
Definition: direction.H:43
virtual void updateMesh(const mapPolyMesh &)=0
Update local data for topology changes.
Definition: motionSolver.C:178
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
Definition: Field.C:578
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
const labelList & reversePointMap() const
Reverse point map.
Definition: mapPolyMesh.H:465
static const pointMesh & New(const polyMesh &mesh)
bool hasMotionPoints() const
Has valid preMotionPoints?
Definition: mapPolyMesh.H:615
A class for handling words, derived from string.
Definition: word.H:59
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
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:741
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual void updateMesh(const mapPolyMesh &)
Update local data for topology changes.
dynamicFvMesh & mesh
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
Virtual base class for mesh motion solver.
Definition: motionSolver.H:53
Namespace for OpenFOAM.
const pointField & preMotionPoints() const
Pre-motion point positions.
Definition: mapPolyMesh.H:609
stressControl lookup("compactNormalStress") >> compactNormalStress
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1035
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
volScalarField & p
Definition: createFields.H:51
#define forAll(list, i)
Definition: UList.H:421
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
virtual void movePoints(const pointField &)
Update local data for geometry changes.
Type gMin(const FieldField< Field, Type > &f)
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
error FatalError
Constant dispersed-phase particle diameter model.
const labelList & pointMap() const
Old point map.
Definition: mapPolyMesh.H:392
const Time & time() const
Return time.
Type gMax(const FieldField< Field, Type > &f)
label nPoints() const
defineTypeNameAndDebug(combustionModel, 0)