fvMotionSolverEngineMesh.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 
28 #include "fvcMeshPhi.H"
29 #include "surfaceInterpolate.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  defineTypeNameAndDebug(fvMotionSolverEngineMesh, 0);
36  addToRunTimeSelectionTable(engineMesh, fvMotionSolverEngineMesh, IOobject);
37 }
38 
39 
40 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41 
43 :
44  engineMesh(io),
45  pistonLayers_("pistonLayers", dimLength, 0.0),
46  motionSolver_
47  (
48  *this,
49  engineDB_.engineDict()
50  )
51 {
52  engineDB_.engineDict().readIfPresent("pistonLayers", pistonLayers_);
53 }
54 
55 
56 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
57 
59 {}
60 
61 
62 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
63 
65 {
66  scalar deltaZ = engineDB_.pistonDisplacement().value();
67  Info<< "deltaZ = " << deltaZ << endl;
68 
69  // Position of the top of the static mesh layers above the piston
70  scalar pistonPlusLayers = pistonPosition_.value() + pistonLayers_.value();
71 
72  scalar pistonSpeed = deltaZ/engineDB_.deltaTValue();
73 
74  motionSolver_.pointMotionU().boundaryFieldRef()[pistonIndex_] ==
75  pistonSpeed;
76 
77  {
78  scalarField linerPoints
79  (
80  boundary()[linerIndex_].patch().localPoints().component(vector::Z)
81  );
82 
83  motionSolver_.pointMotionU().boundaryFieldRef()[linerIndex_] ==
84  pistonSpeed*pos0(deckHeight_.value() - linerPoints)
85  *(deckHeight_.value() - linerPoints)
86  /(deckHeight_.value() - pistonPlusLayers);
87  }
88 
89  motionSolver_.solve();
90 
92  {
95 
96  const volScalarField& rho =
98 
99  const volVectorField& U =
101 
102  bool absolutePhi = false;
103  if (moving())
104  {
105  phi += fvc::interpolate(rho)*fvc::meshPhi(rho, U);
106  absolutePhi = true;
107  }
108 
109  movePoints(motionSolver_.curPoints());
110 
111  if (absolutePhi)
112  {
113  phi -= fvc::interpolate(rho)*fvc::meshPhi(rho, U);
114  }
115  }
116  else
117  {
118  movePoints(motionSolver_.curPoints());
119  }
120 
121 
122  pistonPosition_.value() += deltaZ;
123 
124  Info<< "clearance: " << deckHeight_.value() - pistonPosition_.value() << nl
125  << "Piston speed = " << pistonSpeed << " m/s" << endl;
126 }
127 
128 
129 // ************************************************************************* //
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:256
bool movePoints()
Do what is necessary if the mesh has moved.
bool foundObject(const word &name) const
Is the named Type found?
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
fvMotionSolverEngineMesh(const IOobject &io)
Construct from IOobject.
const IOdictionary & engineDict() const
Return the engine geometry dictionary.
Definition: engineTime.H:128
Macros for easy insertion into run-time selection tables.
dimensionedScalar deckHeight_
Definition: engineMesh.H:65
label pistonIndex_
Definition: engineMesh.H:61
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
dimensionedScalar pos0(const dimensionedScalar &ds)
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
static const char nl
Definition: Ostream.H:265
Foam::engineMesh.
Definition: engineMesh.H:51
defineTypeNameAndDebug(combustionModel, 0)
pointScalarField & pointMotionU()
Non-const access to the pointMotionU in order to allow changes.
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
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
tmp< surfaceScalarField > meshPhi(const volVectorField &U)
Definition: fvcMeshPhi.C:34
messageStream Info
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
dimensionedScalar pistonPosition_
Definition: engineMesh.H:66
Namespace for OpenFOAM.
dimensionedScalar pistonDisplacement() const
Return piston displacement for current time step.
Definition: engineTime.C:104
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540