fvMeshMoversEngine.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) 2021 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 "fvMeshMoversEngine.H"
27 #include "engineTime.H"
28 #include "unitConversion.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace fvMeshMovers
35 {
37 }
38 }
39 
40 
41 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42 
44 :
45  fvMeshMover(mesh),
46  meshCoeffs_(dict()),
47  rpm_
48  (
49  refCast<const userTimes::engine>(mesh.time().userTime()).rpm()
50  ),
51  conRodLength_("conRodLength", dimLength, meshCoeffs_),
52  bore_("bore", dimLength, meshCoeffs_),
53  stroke_("stroke", dimLength, meshCoeffs_),
54  clearance_("clearance", dimLength, meshCoeffs_),
55  pistonIndex_(-1),
56  linerIndex_(-1),
57  cylinderHeadIndex_(-1),
58  deckHeight_("deckHeight", dimLength, great),
59  pistonPosition_("pistonPosition", dimLength, -great)
60 {
61  bool foundPiston = false;
62  bool foundLiner = false;
63  bool foundCylinderHead = false;
64 
65  forAll(mesh.boundary(), i)
66  {
67  if (mesh.boundary()[i].name() == "piston")
68  {
69  pistonIndex_ = i;
70  foundPiston = true;
71  }
72  else if (mesh.boundary()[i].name() == "liner")
73  {
74  linerIndex_ = i;
75  foundLiner = true;
76  }
77  else if (mesh.boundary()[i].name() == "cylinderHead")
78  {
80  foundCylinderHead = true;
81  }
82  }
83 
84  reduce(foundPiston, orOp<bool>());
85  reduce(foundLiner, orOp<bool>());
86  reduce(foundCylinderHead, orOp<bool>());
87 
88  if (!foundPiston)
89  {
91  << "cannot find piston patch"
92  << exit(FatalError);
93  }
94 
95  if (!foundLiner)
96  {
98  << "cannot find liner patch"
99  << exit(FatalError);
100  }
101 
102  if (!foundCylinderHead)
103  {
105  << "cannot find cylinderHead patch"
106  << exit(FatalError);
107  }
108 
109  {
110  if (pistonIndex_ != -1)
111  {
112  pistonPosition_.value() = -great;
113  if (mesh.boundary()[pistonIndex_].patch().localPoints().size())
114  {
116  max(mesh.boundary()[pistonIndex_].patch().localPoints())
117  .z();
118  }
119  }
121 
122  if (cylinderHeadIndex_ != -1)
123  {
124  deckHeight_.value() = great;
125  if
126  (
127  mesh.boundary()[cylinderHeadIndex_].patch().localPoints().size()
128  )
129  {
130  deckHeight_.value() = min
131  (
132  mesh.boundary()[cylinderHeadIndex_].patch().localPoints()
133  ).z();
134  }
135  }
137 
138  Info<< "deckHeight: " << deckHeight_.value() << nl
139  << "piston position: " << pistonPosition_.value() << endl;
140  }
141 }
142 
143 
144 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
145 
147 {}
148 
149 
150 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
151 
153 {
154  return mesh().time().userTimeValue();
155 }
156 
157 
159 {
160  return mesh().time().timeToUserTime(mesh().time().deltaTValue());
161 }
162 
163 
165 (
166  const scalar theta
167 ) const
168 {
169  return
170  (
171  conRodLength_.value()
172  + stroke_.value()/2.0
173  + clearance_.value()
174  )
175  - (
176  stroke_.value()*::cos(degToRad(theta))/2.0
177  + ::sqrt
178  (
179  sqr(conRodLength_.value())
180  - sqr(stroke_.value()*::sin(degToRad(theta))/2.0)
181  )
182  );
183 }
184 
185 
187 {
188  return dimensionedScalar
189  (
190  "pistonPosition",
191  dimLength,
192  pistonPosition(theta())
193  );
194 }
195 
196 
198 {
199  return dimensionedScalar
200  (
201  "pistonDisplacement",
202  dimLength,
203  pistonPosition(theta() - deltaTheta()) - pistonPosition().value()
204  );
205 }
206 
207 
209 {
210  return dimensionedScalar
211  (
212  "pistonSpeed",
213  dimVelocity,
214  pistonDisplacement().value()/(mesh().time().deltaTValue() + vSmall)
215  );
216 }
217 
218 
219 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
const Type & value() const
Return const reference to value.
Abstract base class for fvMesh movers.
Definition: fvMeshMover.H:53
fvMesh & mesh()
Return the fvMesh.
Definition: fvMeshMover.H:132
Basic mesh motion specifically for engines.
virtual scalar theta() const
Return current crank-angle.
virtual ~engine()
Destructor.
dimensionedScalar pistonDisplacement() const
Return piston displacement for current time step.
virtual scalar deltaTheta() const
Return crank-angle increment.
dimensionedScalar pistonPosition() const
Return current piston position.
dimensionedScalar pistonPosition_
dimensionedScalar pistonSpeed() const
Return piston speed for current time step.
engine(fvMesh &mesh)
Construct from fvMesh.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:101
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:895
An abstract class for the user time description.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
defineTypeNameAndDebug(none, 0)
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
To & refCast(From &r)
Reference type cast template function.
Definition: typeInfo.H:111
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
dimensionedScalar sin(const dimensionedScalar &ds)
messageStream Info
const dimensionSet dimLength
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
dimensionedScalar sqrt(const dimensionedScalar &ds)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
const dimensionSet dimVelocity
error FatalError
static const char nl
Definition: Ostream.H:260
dimensionedScalar cos(const dimensionedScalar &ds)
scalar degToRad(const scalar deg)
Conversion from degrees to radians.
dictionary dict
Unit conversion functions.