displacementLaplacianFvMotionSolver.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2015-2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
30#include "motionInterpolation.H"
31#include "motionDiffusivity.H"
32#include "fvmLaplacian.H"
34#include "OFstream.H"
35#include "meshTools.H"
36#include "mapPolyMesh.H"
37#include "fvOptions.H"
38
39// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40
41namespace Foam
42{
44
46 (
50 );
51
53 (
56 displacement
57 );
58}
59
60
61// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
62
64(
65 const polyMesh& mesh,
66 const IOdictionary& dict
67)
68:
71 cellDisplacement_
72 (
74 (
75 "cellDisplacement",
76 mesh.time().timeName(),
77 mesh,
78 IOobject::READ_IF_PRESENT,
79 IOobject::AUTO_WRITE
80 ),
81 fvMesh_,
82 dimensionedVector(pointDisplacement_.dimensions(), Zero),
83 cellMotionBoundaryTypes<vector>(pointDisplacement_.boundaryField())
84 ),
85 pointLocation_(nullptr),
86 interpolationPtr_
87 (
88 coeffDict().found("interpolation")
89 ? motionInterpolation::New(fvMesh_, coeffDict().lookup("interpolation"))
90 : motionInterpolation::New(fvMesh_)
91 ),
92 diffusivityPtr_
93 (
94 motionDiffusivity::New(fvMesh_, coeffDict().lookup("diffusivity"))
95 ),
96 frozenPointsZone_
97 (
98 coeffDict().found("frozenPointsZone")
99 ? fvMesh_.pointZones().findZoneID
100 (
101 coeffDict().get<word>("frozenPointsZone")
102 )
103 : -1
104 )
105{
107 (
108 "pointLocation",
110 fvMesh_,
113 );
114
115 if (debug)
116 {
117 Info<< "displacementLaplacianFvMotionSolver:" << nl
118 << " diffusivity : " << diffusivityPtr_().type() << nl
119 << " frozenPoints zone : " << frozenPointsZone_ << endl;
120 }
121
122
124 {
125 pointLocation_.reset
126 (
128 (
129 io,
131 )
132 );
133
134 if (debug)
135 {
136 Info<< "displacementLaplacianFvMotionSolver :"
137 << " Read pointVectorField "
138 << io.name()
139 << " to be used for boundary conditions on points."
140 << nl
141 << "Boundary conditions:"
142 << pointLocation_().boundaryField().types() << endl;
143 }
144 }
145}
146
147
148Foam::displacementLaplacianFvMotionSolver::
149displacementLaplacianFvMotionSolver
150(
151 const polyMesh& mesh,
152 const IOdictionary& dict,
153 const pointVectorField& pointDisplacement,
154 const pointIOField& points0
155)
156:
157 displacementMotionSolver(mesh, dict, pointDisplacement, points0, typeName),
159 cellDisplacement_
160 (
162 (
163 "cellDisplacement",
164 mesh.time().timeName(),
165 mesh,
166 IOobject::READ_IF_PRESENT,
167 IOobject::AUTO_WRITE
168 ),
169 fvMesh_,
170 dimensionedVector(pointDisplacement_.dimensions(), Zero),
171 cellMotionBoundaryTypes<vector>(pointDisplacement_.boundaryField())
172 ),
173 pointLocation_(nullptr),
174 interpolationPtr_
175 (
176 coeffDict().found("interpolation")
177 ? motionInterpolation::New(fvMesh_, coeffDict().lookup("interpolation"))
178 : motionInterpolation::New(fvMesh_)
179 ),
180 diffusivityPtr_
181 (
182 motionDiffusivity::New(fvMesh_, coeffDict().lookup("diffusivity"))
183 ),
184 frozenPointsZone_
185 (
186 coeffDict().found("frozenPointsZone")
187 ? fvMesh_.pointZones().findZoneID
188 (
189 coeffDict().get<word>("frozenPointsZone")
190 )
191 : -1
192 )
193{
195 (
196 "pointLocation",
198 fvMesh_,
201 );
202
203 if (debug)
204 {
205 Info<< "displacementLaplacianFvMotionSolver:" << nl
206 << " diffusivity : " << diffusivityPtr_().type() << nl
207 << " frozenPoints zone : " << frozenPointsZone_ << endl;
208 }
209
210
212 {
213 pointLocation_.reset
214 (
216 (
217 io,
219 )
220 );
221
222 if (debug)
223 {
224 Info<< "displacementLaplacianFvMotionSolver :"
225 << " Read pointVectorField "
226 << io.name()
227 << " to be used for boundary conditions on points."
228 << nl
229 << "Boundary conditions:"
230 << pointLocation_().boundaryField().types() << endl;
231 }
232 }
233}
234
235
236// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
237
240{}
241
242
243// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
244
247{
248 if (!diffusivityPtr_)
249 {
250 diffusivityPtr_ = motionDiffusivity::New
251 (
252 fvMesh_,
253 coeffDict().lookup("diffusivity")
254 );
255 }
256
257 return *diffusivityPtr_;
258}
259
260
263{
264 interpolationPtr_->interpolate
265 (
266 cellDisplacement_,
267 pointDisplacement_
268 );
269
270 if (pointLocation_)
271 {
272 if (debug)
273 {
274 Info<< "displacementLaplacianFvMotionSolver : applying "
275 << " boundary conditions on " << pointLocation_().name()
276 << " to new point location."
277 << endl;
278 }
279
280 pointLocation_().primitiveFieldRef() =
281 points0()
282 + pointDisplacement_.primitiveField();
283
284 pointLocation_().correctBoundaryConditions();
285
286 // Implement frozen points
287 if (frozenPointsZone_ != -1)
288 {
289 const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
290
291 forAll(pz, i)
292 {
293 pointLocation_()[pz[i]] = points0()[pz[i]];
294 }
295 }
296
297 twoDCorrectPoints(pointLocation_().primitiveFieldRef());
298
299 return tmp<pointField>(pointLocation_().primitiveField());
300 }
301 else
302 {
303 tmp<pointField> tcurPoints
304 (
305 points0() + pointDisplacement_.primitiveField()
306 );
307 pointField& curPoints = tcurPoints.ref();
308
309 // Implement frozen points
310 if (frozenPointsZone_ != -1)
311 {
312 const pointZone& pz = fvMesh_.pointZones()[frozenPointsZone_];
313
314 forAll(pz, i)
315 {
316 curPoints[pz[i]] = points0()[pz[i]];
317 }
318 }
319
320 twoDCorrectPoints(curPoints);
321
322 return tcurPoints;
323 }
324}
325
326
328{
329 // The points have moved so before interpolation update
330 // the motionSolver accordingly
331 movePoints(fvMesh_.points());
332
333 diffusivity().correct();
334 pointDisplacement_.boundaryFieldRef().updateCoeffs();
335
337
338 // We explicitly do NOT want to interpolate the motion inbetween
339 // different regions so bypass all the matrix manipulation.
341 (
343 (
344 dimensionedScalar("viscosity", dimViscosity, 1.0)
345 *diffusivity().operator()(),
346 cellDisplacement_,
347 "laplacian(diffusivity,cellDisplacement)"
348 )
349 ==
350 fvOptions(cellDisplacement_)
351 );
352
355 fvOptions.correct(cellDisplacement_);
356}
357
358
360(
361 const mapPolyMesh& mpm
362)
363{
365
366 // Update diffusivity. Note two stage to make sure old one is de-registered
367 // before creating/registering new one.
368 diffusivityPtr_.clear();
369}
370
371
372// ************************************************************************* //
bool found
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
fv::options & fvOptions
A primitive field of type <T> with automated input and output.
Definition: IOField.H:58
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:65
bool typeHeaderOk(const bool checkType=true, const bool search=true, const bool verbose=true)
Read header (uses typeFilePath to find file) and check its info.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Mesh motion solver for an fvMesh. Based on solving the cell-centre Laplacian for the motion displacem...
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
motionDiffusivity & diffusivity()
Return reference to the diffusivity field.
Virtual base class for displacement motion solver.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
SolverPerformance< Type > solveSegregatedOrCoupled(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:62
const dictionary & solverDict() const
Return the solver dictionary taking into account finalIteration.
Definition: fvMatrix.C:1558
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
Base class for fvMesh based motionSolvers.
const fvMesh & fvMesh_
The fvMesh to be moved.
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
Finite-volume options.
Definition: fvOptions.H:59
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:162
Abstract base class for cell-centre mesh motion diffusivity.
Base class for interpolation of cell displacement fields, generated by fvMotionSolvers,...
void updateMesh()
Update for new mesh topology.
Virtual base class for mesh motion solver.
Definition: motionSolver.H:61
A subset of mesh points.
Definition: pointZone.H:68
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
Lookup type of boundary radiation properties.
Definition: lookup.H:66
A class for managing temporary objects.
Definition: tmp.H:65
T & ref() const
Definition: tmpI.H:227
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
dynamicFvMesh & mesh
Calculate the matrix for the laplacian of the field.
fvScalarMatrix TEqn(fvm::ddt(T)+fvm::div(phi, T) - fvm::laplacian(alphaEff, T)==radiation->ST(rhoCpRef, T)+fvOptions(T))
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, false)
word timeName
Definition: getTimeIndex.H:3
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:48
Namespace for OpenFOAM.
const dimensionSet dimViscosity
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
pointField points0(pointIOField(IOobject("points", mesh.time().constant(), polyMesh::meshSubDir, mesh, IOobject::MUST_READ, IOobject::NO_WRITE, false)))
conserve primitiveFieldRef()+