displacementComponentLaplacianFvMotionSolver.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 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 
52 (
53  const polyMesh& mesh,
54  const dictionary& dict
55 )
56 :
58  fvMotionSolver(mesh),
59  cellDisplacement_
60  (
61  IOobject
62  (
63  "cellDisplacement" + cmptName_,
64  mesh.time().timeName(),
65  mesh,
68  ),
69  fvMesh_,
70  dimensionedScalar(pointDisplacement_.dimensions(), 0),
71  cellMotionBoundaryTypes<scalar>(pointDisplacement_.boundaryField())
72  ),
73  pointLocation_(nullptr),
74  diffusivityPtr_
75  (
76  motionDiffusivity::New(fvMesh_, coeffDict().lookup("diffusivity"))
77  ),
78  frozenPointsZone_
79  (
80  coeffDict().found("frozenPointsZone")
81  ? fvMesh_.pointZones().findZoneID(coeffDict().lookup("frozenPointsZone"))
82  : -1
83  )
84 {
85  Switch applyPointLocation
86  (
87  coeffDict().lookupOrDefault
88  (
89  "applyPointLocation",
90  true
91  )
92  );
93 
94  if (applyPointLocation)
95  {
96  pointLocation_.reset
97  (
99  (
100  IOobject
101  (
102  "pointLocation",
103  fvMesh_.time().timeName(),
104  fvMesh_,
107  ),
108  pointMesh::New(fvMesh_)
109  )
110  );
111 
112  // if (debug)
113  {
114  Info<< "displacementComponentLaplacianFvMotionSolver :"
115  << " Read pointVectorField "
116  << pointLocation_().name()
117  << " to be used for boundary conditions on points."
118  << nl
119  << "Boundary conditions:"
120  << pointLocation_().boundaryField().types() << endl;
121  }
122  }
123 }
124 
125 
126 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
127 
130 {}
131 
132 
133 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
134 
137 {
139  (
140  cellDisplacement_,
141  pointDisplacement_
142  );
143 
144  if (pointLocation_.valid())
145  {
146  if (debug)
147  {
148  Info<< "displacementComponentLaplacianFvMotionSolver : applying "
149  << " boundary conditions on " << pointLocation_().name()
150  << " to new point location."
151  << endl;
152  }
153 
154  // Apply pointLocation_ b.c. to mesh points.
155 
156  pointLocation_().primitiveFieldRef() = fvMesh_.points();
157 
158  pointLocation_().primitiveFieldRef().replace
159  (
160  cmpt_,
161  points0_ + pointDisplacement_.primitiveField()
162  );
163 
164  pointLocation_().correctBoundaryConditions();
165 
166  // Implement frozen points
167  if (frozenPointsZone_ != -1)
168  {
169  const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
170 
171  forAll(pz, i)
172  {
173  label pointi = pz[i];
174 
175  pointLocation_()[pointi][cmpt_] = points0_[pointi];
176  }
177  }
178 
179  twoDCorrectPoints(pointLocation_().primitiveFieldRef());
180 
181  return tmp<pointField>(pointLocation_().primitiveField());
182  }
183  else
184  {
185  tmp<pointField> tcurPoints(new pointField(fvMesh_.points()));
186  pointField& curPoints = tcurPoints.ref();
187 
188  curPoints.replace
189  (
190  cmpt_,
191  points0_ + pointDisplacement_.primitiveField()
192  );
193 
194  // Implement frozen points
195  if (frozenPointsZone_ != -1)
196  {
197  const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
198 
199  forAll(pz, i)
200  {
201  label pointi = pz[i];
202 
203  curPoints[pointi][cmpt_] = points0_[pointi];
204  }
205  }
206 
207  twoDCorrectPoints(curPoints);
208 
209  return tcurPoints;
210  }
211 }
212 
213 
215 {
216  // The points have moved so before interpolation update
217  // the motionSolver accordingly
218  movePoints(fvMesh_.points());
219 
220  diffusivityPtr_->correct();
221  pointDisplacement_.boundaryFieldRef().updateCoeffs();
222 
224  (
226  (
227  diffusivityPtr_->operator()(),
228  cellDisplacement_,
229  "laplacian(diffusivity,cellDisplacement)"
230  )
231  );
232 }
233 
234 
236 (
237  const mapPolyMesh& mpm
238 )
239 {
241 
242  // Update diffusivity. Note two stage to make sure old one is de-registered
243  // before creating/registering new one.
244  diffusivityPtr_.reset(nullptr);
245  diffusivityPtr_ = motionDiffusivity::New
246  (
247  fvMesh_,
248  coeffDict().lookup("diffusivity")
249  );
250 }
251 
252 
253 // ************************************************************************* //
tmp< GeometricField< Type, pointPatchField, pointMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &) const
Interpolate volField using inverse distance weighting.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
virtual void updateMesh(const mapPolyMesh &)
Update local data for topology changes.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
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/any.
Definition: Switch.H:60
Calculate the matrix for the laplacian of the field.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
void replace(const direction, const UList< cmptType > &)
Replace a component field of the field.
Definition: Field.C:472
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
conserve primitiveFieldRef()+
Macros for easy insertion into run-time selection tables.
static const pointMesh & New(const polyMesh &mesh)
Definition: MeshObject.C:44
Base class for fvMesh based motionSolvers.
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
stressControl lookup("compactNormalStress") >> compactNormalStress
displacementComponentLaplacianFvMotionSolver(const polyMesh &, const dictionary &)
Construct from polyMesh and dictionary.
dynamicFvMesh & mesh
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
static const char nl
Definition: Ostream.H:265
const Time & time() const
Return time.
defineTypeNameAndDebug(combustionModel, 0)
A subset of mesh points. The labels of points in the zone can be obtained from the addressing() list...
Definition: pointZone.H:62
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
messageStream Info
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A class for managing temporary objects.
Definition: PtrList.H:53
static autoPtr< motionDiffusivity > New(const fvMesh &mesh, Istream &mdData)
Select null constructed.
Virtual base class for displacement motion solver.
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
bool found
Namespace for OpenFOAM.