displacementComponentLaplacianFvMotionSolver.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-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 "motionDiffusivity.H"
28 #include "fvmLaplacian.H"
30 #include "mapPolyMesh.H"
31 #include "volPointInterpolation.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  defineTypeNameAndDebug(displacementComponentLaplacianFvMotionSolver, 0);
38 
40  (
41  motionSolver,
42  displacementComponentLaplacianFvMotionSolver,
43  dictionary
44  );
45 }
46 
47 
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49 
50 Foam::displacementComponentLaplacianFvMotionSolver::
51 displacementComponentLaplacianFvMotionSolver
52 (
53  const polyMesh& mesh,
54  const IOdictionary& dict
55 )
56 :
58  fvMotionSolverCore(mesh),
59  cellDisplacement_
60  (
61  IOobject
62  (
63  "cellDisplacement" + cmptName_,
64  mesh.time().timeName(),
65  mesh,
68  ),
69  fvMesh_,
71  (
72  "cellDisplacement",
73  pointDisplacement_.dimensions(),
74  0
75  ),
76  cellMotionBoundaryTypes<scalar>(pointDisplacement_.boundaryField())
77  ),
78  pointLocation_(NULL),
79  diffusivityPtr_
80  (
81  motionDiffusivity::New(fvMesh_, coeffDict().lookup("diffusivity"))
82  ),
83  frozenPointsZone_
84  (
85  coeffDict().found("frozenPointsZone")
86  ? fvMesh_.pointZones().findZoneID(coeffDict().lookup("frozenPointsZone"))
87  : -1
88  )
89 {
90  Switch applyPointLocation
91  (
92  coeffDict().lookupOrDefault
93  (
94  "applyPointLocation",
95  true
96  )
97  );
98 
99  if (applyPointLocation)
100  {
101  pointLocation_.reset
102  (
103  new pointVectorField
104  (
105  IOobject
106  (
107  "pointLocation",
108  fvMesh_.time().timeName(),
109  fvMesh_,
112  ),
113  pointMesh::New(fvMesh_)
114  )
115  );
116 
117  //if (debug)
118  {
119  Info<< "displacementComponentLaplacianFvMotionSolver :"
120  << " Read pointVectorField "
121  << pointLocation_().name()
122  << " to be used for boundary conditions on points."
123  << nl
124  << "Boundary conditions:"
125  << pointLocation_().boundaryField().types() << endl;
126  }
127  }
128 }
129 
130 
131 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
132 
135 {}
136 
137 
138 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
139 
142 {
144  (
145  cellDisplacement_,
146  pointDisplacement_
147  );
148 
149  if (pointLocation_.valid())
150  {
151  if (debug)
152  {
153  Info<< "displacementComponentLaplacianFvMotionSolver : applying "
154  << " boundary conditions on " << pointLocation_().name()
155  << " to new point location."
156  << endl;
157  }
158 
159  // Apply pointLocation_ b.c. to mesh points.
160 
161  pointLocation_().internalField() = fvMesh_.points();
162 
163  pointLocation_().internalField().replace
164  (
165  cmpt_,
166  points0_ + pointDisplacement_.internalField()
167  );
168 
169  pointLocation_().correctBoundaryConditions();
170 
171  // Implement frozen points
172  if (frozenPointsZone_ != -1)
173  {
174  const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
175 
176  forAll(pz, i)
177  {
178  label pointI = pz[i];
179 
180  pointLocation_()[pointI][cmpt_] = points0_[pointI];
181  }
182  }
183 
184  twoDCorrectPoints(pointLocation_().internalField());
185 
186  return tmp<pointField>(pointLocation_().internalField());
187  }
188  else
189  {
190  tmp<pointField> tcurPoints(new pointField(fvMesh_.points()));
191 
192  tcurPoints().replace
193  (
194  cmpt_,
195  points0_ + pointDisplacement_.internalField()
196  );
197 
198  // Implement frozen points
199  if (frozenPointsZone_ != -1)
200  {
201  const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
202 
203  forAll(pz, i)
204  {
205  label pointI = pz[i];
206 
207  tcurPoints()[pointI][cmpt_] = points0_[pointI];
208  }
209  }
210 
211  twoDCorrectPoints(tcurPoints());
212 
213  return tcurPoints;
214  }
215 }
216 
217 
219 {
220  // The points have moved so before interpolation update
221  // the motionSolver accordingly
222  movePoints(fvMesh_.points());
223 
224  diffusivityPtr_->correct();
225  pointDisplacement_.boundaryField().updateCoeffs();
226 
228  (
230  (
231  diffusivityPtr_->operator()(),
232  cellDisplacement_,
233  "laplacian(diffusivity,cellDisplacement)"
234  )
235  );
236 }
237 
238 
240 (
241  const mapPolyMesh& mpm
242 )
243 {
245 
246  // Update diffusivity. Note two stage to make sure old one is de-registered
247  // before creating/registering new one.
248  diffusivityPtr_.reset(NULL);
249  diffusivityPtr_ = motionDiffusivity::New
250  (
251  fvMesh_,
252  coeffDict().lookup("diffusivity")
253  );
254 }
255 
256 
257 // ************************************************************************* //
tmp< GeometricField< Type, pointPatchField, pointMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &) const
Interpolate volField using inverse distance weighting.
static autoPtr< motionDiffusivity > New(const fvMesh &mesh, Istream &mdData)
Select null constructed.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
static const pointMesh & New(const polyMesh &mesh)
A subset of mesh points. The labels of points in the zone can be obtained from the addressing() list...
Definition: pointZone.H:62
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
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:741
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
messageStream Info
virtual void updateMesh(const mapPolyMesh &)
Update local data for topology changes.
dynamicFvMesh & mesh
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
Namespace for OpenFOAM.
Base class for fvMesh based motionSolvers.
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
stressControl lookup("compactNormalStress") >> compactNormalStress
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
#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
Macros for easy insertion into run-time selection tables.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
bool found
Calculate the matrix for the laplacian of the field.
solverPerformance solve(fvMatrix< Type > &, const dictionary &)
Solve returning the solution statistics given convergence tolerance.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
const Time & time() const
Return time.
Virtual base class for displacement motion solver.
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
A class for managing temporary objects.
Definition: PtrList.H:118
conserve internalField()+
defineTypeNameAndDebug(combustionModel, 0)