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-2023 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 "polyTopoChangeMap.H"
31 #include "volPointInterpolation.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
38 
40  (
44  );
45 }
46 
47 
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49 
52 (
53  const word& name,
54  const polyMesh& mesh,
55  const dictionary& dict
56 )
57 :
59  fvMotionSolver(mesh),
60  cellDisplacement_
61  (
62  IOobject
63  (
64  "cellDisplacement" + cmptName_,
65  mesh.time().name(),
66  mesh,
67  IOobject::READ_IF_PRESENT,
68  IOobject::AUTO_WRITE
69  ),
70  fvMesh_,
71  dimensionedScalar(pointDisplacement_.dimensions(), 0),
72  cellMotionBoundaryTypes<scalar>(pointDisplacement_.boundaryField())
73  ),
74  pointLocation_(nullptr),
75  diffusivityPtr_
76  (
77  motionDiffusivity::New(fvMesh_, coeffDict().lookup("diffusivity"))
78  ),
79  frozenPointsZone_
80  (
81  coeffDict().found("frozenPointsZone")
82  ? fvMesh_.pointZones().findIndex(coeffDict().lookup("frozenPointsZone"))
83  : -1
84  )
85 {
87  (
88  "pointLocation",
89  fvMesh_.time().name(),
90  fvMesh_,
93  );
94 
95  if (debug)
96  {
97  Info<< "displacementComponentLaplacianFvMotionSolver:" << nl
98  << " diffusivity : " << diffusivityPtr_().type() << nl
99  << " frozenPoints zone : " << frozenPointsZone_ << endl;
100  }
101 
102  if (io.headerOk())
103  {
104  pointLocation_.reset
105  (
106  new pointVectorField
107  (
108  IOobject
109  (
110  "pointLocation",
111  fvMesh_.time().name(),
112  fvMesh_,
115  ),
117  )
118  );
119 
120  if (debug)
121  {
122  Info<< "displacementComponentLaplacianFvMotionSolver :"
123  << " Read pointVectorField "
124  << pointLocation_().name()
125  << " to be used for boundary conditions on points."
126  << nl
127  << "Boundary conditions:"
128  << pointLocation_().boundaryField().types() << endl;
129  }
130  }
131 }
132 
133 
134 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
135 
138 {}
139 
140 
141 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
142 
145 {
147  (
148  cellDisplacement_,
149  pointDisplacement_
150  );
151 
152  if (pointLocation_.valid())
153  {
154  if (debug)
155  {
156  Info<< "displacementComponentLaplacianFvMotionSolver : applying "
157  << " boundary conditions on " << pointLocation_().name()
158  << " to new point location."
159  << endl;
160  }
161 
162  // Apply pointLocation_ b.c. to mesh points.
163 
164  pointLocation_().primitiveFieldRef() = fvMesh_.points();
165 
166  pointLocation_().primitiveFieldRef().replace
167  (
168  cmpt_,
169  points0_ + pointDisplacement_.primitiveField()
170  );
171 
172  pointLocation_().correctBoundaryConditions();
173 
174  // Implement frozen points
175  if (frozenPointsZone_ != -1)
176  {
177  const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
178 
179  forAll(pz, i)
180  {
181  label pointi = pz[i];
182 
183  pointLocation_()[pointi][cmpt_] = points0_[pointi];
184  }
185  }
186 
187  twoDCorrectPoints(pointLocation_().primitiveFieldRef());
188 
189  return tmp<pointField>(pointLocation_().primitiveField());
190  }
191  else
192  {
193  tmp<pointField> tcurPoints(new pointField(fvMesh_.points()));
194  pointField& curPoints = tcurPoints.ref();
195 
196  curPoints.replace
197  (
198  cmpt_,
199  points0_ + pointDisplacement_.primitiveField()
200  );
201 
202  // Implement frozen points
203  if (frozenPointsZone_ != -1)
204  {
205  const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
206 
207  forAll(pz, i)
208  {
209  label pointi = pz[i];
210 
211  curPoints[pointi][cmpt_] = points0_[pointi];
212  }
213  }
214 
215  twoDCorrectPoints(curPoints);
216 
217  return tcurPoints;
218  }
219 }
220 
221 
223 {
224  // The points have moved so before interpolation update
225  // the motionSolver accordingly
226  movePoints(fvMesh_.points());
227 
228  diffusivityPtr_->correct();
229  pointDisplacement_.boundaryFieldRef().updateCoeffs();
230 
232  (
234  (
235  diffusivityPtr_->operator()(),
236  cellDisplacement_,
237  "laplacian(diffusivity,cellDisplacement)"
238  )
239  );
240 }
241 
242 
244 (
245  const polyTopoChangeMap& map
246 )
247 {
249 
250  // Update diffusivity. Note two stage to make sure old one is de-registered
251  // before creating/registering new one.
252  diffusivityPtr_.reset(nullptr);
253  diffusivityPtr_ = motionDiffusivity::New
254  (
255  fvMesh_,
256  coeffDict().lookup("diffusivity")
257  );
258 }
259 
260 
262 (
263  const polyMeshMap& map
264 )
265 {
267 
268  // Update diffusivity. Note two stage to make sure old one is de-registered
269  // before creating/registering new one.
270  diffusivityPtr_.reset(nullptr);
271  diffusivityPtr_ = motionDiffusivity::New
272  (
273  fvMesh_,
274  coeffDict().lookup("diffusivity")
275  );
276 }
277 
278 
279 // ************************************************************************* //
bool found
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
static pointMesh & New(const word &name, const polyMesh &mesh)
Construct and return the named DemandDrivenMeshObject.
void replace(const direction, const UList< cmptType > &)
Replace a component field of the field.
Definition: Field.C:491
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
Virtual base class for displacement motion solver.
virtual void topoChange(const polyTopoChangeMap &)
Update local data for topology changes.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
const word & name() const
Return const reference to name.
Mesh motion solver for an fvMesh. Based on solving the cell-centre Laplacian for the given component ...
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
virtual void topoChange(const polyTopoChangeMap &)
Update corresponding to the given map.
displacementComponentLaplacianFvMotionSolver(const word &name, const polyMesh &, const dictionary &)
Construct from polyMesh and dictionary.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:418
Base class for fvMesh based motionSolvers.
const fvMesh & fvMesh_
The fvMesh to be moved.
Abstract base class for cell-centre mesh motion diffusivity.
static autoPtr< motionDiffusivity > New(const fvMesh &mesh, Istream &mdData)
Select null constructed.
Virtual base class for mesh motion solver.
Definition: motionSolver.H:57
A subset of mesh points. The labels of points in the zone can be obtained from the addressing() list.
Definition: pointZone.H:59
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:51
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
Templated form of IOobject providing type information for file reading and header type checking.
Definition: IOobject.H:531
bool headerOk()
Read header (uses typeGlobalFile to find file) and check.
tmp< PointField< Type > > interpolate(const VolField< Type > &) const
Interpolate volField using inverse distance weighting.
A class for handling words, derived from string.
Definition: word.H:62
Calculate the matrix for the laplacian of the field.
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
Namespace for OpenFOAM.
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
addToRunTimeSelectionTable(polyPatch, mergedCyclicPolyPatch, word)
const HashTable< dimensionSet > & dimensions()
Get the table of dimension sets.
Definition: dimensionSets.C:96
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
messageStream Info
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
defineTypeNameAndDebug(combustionModel, 0)
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
static const char nl
Definition: Ostream.H:266
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.
dictionary dict
conserve primitiveFieldRef()+