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-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 "volPointInterpolation.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  defineTypeNameAndDebug(velocityLaplacianFvMotionSolver, 0);
37 
39  (
40  motionSolver,
41  velocityLaplacianFvMotionSolver,
42  dictionary
43  );
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
50 (
51  const polyMesh& mesh,
52  const dictionary& dict
53 )
54 :
55  velocityMotionSolver(mesh, dict, typeName),
56  fvMotionSolver(mesh),
57  cellMotionU_
58  (
59  IOobject
60  (
61  "cellMotionU",
62  mesh.time().timeName(),
63  mesh,
66  ),
67  fvMesh_,
69  (
70  "cellMotionU",
71  pointMotionU_.dimensions(),
72  Zero
73  ),
74  cellMotionBoundaryTypes<vector>(pointMotionU_.boundaryField())
75  ),
76  diffusivityPtr_
77  (
78  motionDiffusivity::New(fvMesh_, coeffDict().lookup("diffusivity"))
79  )
80 {}
81 
82 
83 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
84 
86 {}
87 
88 
89 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
90 
93 {
95  (
96  cellMotionU_,
97  pointMotionU_
98  );
99 
100  tmp<pointField> tcurPoints
101  (
102  fvMesh_.points()
103  + fvMesh_.time().deltaTValue()*pointMotionU_.primitiveField()
104  );
105 
106  twoDCorrectPoints(tcurPoints.ref());
107 
108  return tcurPoints;
109 }
110 
111 
113 {
114  // The points have moved so before interpolation update
115  // the fvMotionSolver accordingly
116  movePoints(fvMesh_.points());
117 
118  diffusivityPtr_->correct();
119  pointMotionU_.boundaryFieldRef().updateCoeffs();
120 
122  (
124  (
125  diffusivityPtr_->operator()(),
126  cellMotionU_,
127  "laplacian(diffusivity,cellMotionU)"
128  )
129  );
130 }
131 
132 
133 //void Foam::velocityLaplacianFvMotionSolver::movePoints(const pointField& p)
134 //{
135 // // Movement of pointMesh and volPointInterpolation already
136 // // done by polyMesh,fvMesh
137 //}
138 
139 
141 (
142  const mapPolyMesh& mpm
143 )
144 {
146 
147  // Update diffusivity. Note two stage to make sure old one is de-registered
148  // before creating/registering new one.
149  diffusivityPtr_.reset(nullptr);
150  diffusivityPtr_ = motionDiffusivity::New
151  (
152  fvMesh_,
153  coeffDict().lookup("diffusivity")
154  );
155 }
156 
157 
158 // ************************************************************************* //
tmp< GeometricField< Type, pointPatchField, pointMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &) const
Interpolate volField using inverse distance weighting.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
Virtual base class for velocity motion solver.
velocityLaplacianFvMotionSolver(const polyMesh &, const dictionary &)
Construct from polyMesh and dictionary.
Calculate the matrix for the laplacian of the field.
virtual void updateMesh(const mapPolyMesh &)
Update local data for topology changes.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
Macros for easy insertion into run-time selection tables.
static const volPointInterpolation & New(const fvMesh &mesh)
Definition: MeshObject.C:44
Base class for fvMesh based motionSolvers.
stressControl lookup("compactNormalStress") >> compactNormalStress
dynamicFvMesh & mesh
static const zero Zero
Definition: zero.H:97
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
const Time & time() const
Return time.
defineTypeNameAndDebug(combustionModel, 0)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
virtual void updateMesh(const mapPolyMesh &)
Update topology.
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.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Namespace for OpenFOAM.