velocityLaplacianFvMotionSolver.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-2022 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 "volPointInterpolation.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
37 
39  (
43  );
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
50 (
51  const word& name,
52  const polyMesh& mesh,
53  const dictionary& dict
54 )
55 :
56  velocityMotionSolver(name, mesh, dict, typeName),
57  fvMotionSolver(mesh),
58  cellMotionU_
59  (
60  IOobject
61  (
62  "cellMotionU",
63  mesh.time().name(),
64  mesh,
65  IOobject::READ_IF_PRESENT,
66  IOobject::AUTO_WRITE
67  ),
68  fvMesh_,
70  (
71  "cellMotionU",
72  pointMotionU_.dimensions(),
73  Zero
74  ),
75  cellMotionBoundaryTypes<vector>(pointMotionU_.boundaryField())
76  ),
77  diffusivityPtr_
78  (
79  motionDiffusivity::New(fvMesh_, coeffDict().lookup("diffusivity"))
80  )
81 {}
82 
83 
84 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
85 
87 {}
88 
89 
90 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
91 
94 {
96  (
97  cellMotionU_,
98  pointMotionU_
99  );
100 
101  tmp<pointField> tcurPoints
102  (
103  fvMesh_.points()
104  + fvMesh_.time().deltaTValue()*pointMotionU_.primitiveField()
105  );
106 
107  twoDCorrectPoints(tcurPoints.ref());
108 
109  return tcurPoints;
110 }
111 
112 
114 {
115  // The points have moved so before interpolation update
116  // the fvMotionSolver accordingly
117  movePoints(fvMesh_.points());
118 
119  diffusivityPtr_->correct();
120  pointMotionU_.boundaryFieldRef().updateCoeffs();
121 
123  (
125  (
126  diffusivityPtr_->operator()(),
127  cellMotionU_,
128  "laplacian(diffusivity,cellMotionU)"
129  )
130  );
131 }
132 
133 
134 //void Foam::velocityLaplacianFvMotionSolver::movePoints(const pointField& p)
135 //{
136 // // Movement of pointMesh and volPointInterpolation already
137 // // done by polyMesh,fvMesh
138 //}
139 
140 
142 (
143  const polyTopoChangeMap& map
144 )
145 {
147 
148  // Update diffusivity. Note two stage to make sure old one is de-registered
149  // before creating/registering new one.
150  diffusivityPtr_.reset(nullptr);
151  diffusivityPtr_ = motionDiffusivity::New
152  (
153  fvMesh_,
154  coeffDict().lookup("diffusivity")
155  );
156 }
157 
158 
160 (
161  const polyMeshMap& map
162 )
163 {
165 
166  // Update diffusivity. Note two stage to make sure old one is de-registered
167  // before creating/registering new one.
168  diffusivityPtr_.reset(nullptr);
169  diffusivityPtr_ = motionDiffusivity::New
170  (
171  fvMesh_,
172  coeffDict().lookup("diffusivity")
173  );
174 }
175 
176 
177 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
static volPointInterpolation & New(const word &name, const fvMesh &mesh)
Construct and return the named DemandDrivenMeshObject.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Base class for fvMesh based motionSolvers.
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
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
Mesh motion solver for an fvMesh. Based on solving the cell-centre Laplacian for the motion velocity.
velocityLaplacianFvMotionSolver(const word &name, const polyMesh &, const dictionary &)
Construct from polyMesh and dictionary.
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.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
Virtual base class for velocity motion solver.
virtual void topoChange(const polyTopoChangeMap &)
Update corresponding to the given map.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
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.
static const zero Zero
Definition: zero.H:97
addToRunTimeSelectionTable(polyPatch, mergedCyclicPolyPatch, word)
const HashTable< dimensionSet > & dimensions()
Get the table of dimension sets.
Definition: dimensionSets.C:96
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
defineTypeNameAndDebug(combustionModel, 0)
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
dictionary dict