layeredEngineMesh.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-2018 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 
26 #include "layeredEngineMesh.H"
28 #include "fvcMeshPhi.H"
29 #include "surfaceInterpolate.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  defineTypeNameAndDebug(layeredEngineMesh, 0);
36  addToRunTimeSelectionTable(engineMesh, layeredEngineMesh, IOobject);
37 }
38 
39 
40 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41 
43 :
44  engineMesh(io),
45  pistonLayers_("pistonLayers", dimLength, 0.0)
46 {
47  engineDB_.engineDict().readIfPresent("pistonLayers", pistonLayers_);
48 }
49 
50 
51 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
52 
54 {}
55 
56 
57 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
58 
60 {
61  scalar deltaZ = engineDB_.pistonDisplacement().value();
62  Info<< "deltaZ = " << deltaZ << endl;
63 
64  // Position of the top of the static mesh layers above the piston
65  scalar pistonPlusLayers = pistonPosition_.value() + pistonLayers_.value();
66 
67  pointField newPoints(points());
68 
69  forAll(newPoints, pointi)
70  {
71  point& p = newPoints[pointi];
72 
73  if (p.z() < pistonPlusLayers) // In piston bowl
74  {
75  p.z() += deltaZ;
76  }
77  else if (p.z() < deckHeight_.value()) // In liner region
78  {
79  p.z() +=
80  deltaZ
81  *(deckHeight_.value() - p.z())
82  /(deckHeight_.value() - pistonPlusLayers);
83  }
84  }
85 
87  {
90 
91  const volScalarField& rho =
93 
94  const volVectorField& U =
96 
97  bool absolutePhi = false;
98  if (moving())
99  {
100  phi += fvc::interpolate(rho)*fvc::meshPhi(rho, U);
101  absolutePhi = true;
102  }
103 
104  movePoints(newPoints);
105 
106  if (absolutePhi)
107  {
108  phi -= fvc::interpolate(rho)*fvc::meshPhi(rho, U);
109  }
110  }
111  else
112  {
113  movePoints(newPoints);
114  }
115 
116  pistonPosition_.value() += deltaZ;
117  scalar pistonSpeed = deltaZ/engineDB_.deltaTValue();
118 
119  Info<< "clearance: " << deckHeight_.value() - pistonPosition_.value() << nl
120  << "Piston speed = " << pistonSpeed << " m/s" << endl;
121 }
122 
123 
124 // ************************************************************************* //
~layeredEngineMesh()
Destructor.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
bool moving() const
Is mesh moving.
Definition: polyMesh.H:512
Type & lookupObjectRef(const word &name) const
Lookup and return the object reference of the given Type.
const surfaceScalarField & phi() const
Return cell face motion fluxes.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
bool movePoints()
Do what is necessary if the mesh has moved.
bool foundObject(const word &name) const
Is the named Type found?
const Cmpt & z() const
Definition: VectorI.H:87
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
const IOdictionary & engineDict() const
Return the engine geometry dictionary.
Definition: engineTime.H:129
Macros for easy insertion into run-time selection tables.
const dimensionSet dimLength
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1131
dimensionedScalar deckHeight_
Definition: engineMesh.H:65
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:41
const Type & value() const
Return const reference to value.
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
const engineTime & engineDB_
Definition: engineMesh.H:59
layeredEngineMesh(const IOobject &io)
Construct from IOobject.
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
static const char nl
Definition: Ostream.H:260
Foam::engineMesh.
Definition: engineMesh.H:51
defineTypeNameAndDebug(combustionModel, 0)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
U
Definition: pEqn.H:72
tmp< surfaceScalarField > meshPhi(const volVectorField &U)
Definition: fvcMeshPhi.C:34
messageStream Info
volScalarField & p
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
dimensionedScalar pistonPosition_
Definition: engineMesh.H:66
Namespace for OpenFOAM.
dimensionedScalar pistonDisplacement() const
Return piston displacement for current time step.
Definition: engineTime.C:104