engineValve.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 | Copyright (C) 2011-2013 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 "engineValve.H"
27 #include "engineTime.H"
28 #include "polyMesh.H"
29 #include "interpolateXY.H"
30 
31 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
32 
33 Foam::scalar Foam::engineValve::adjustCrankAngle(const scalar theta) const
34 {
35  if (theta < liftProfileStart_)
36  {
37  scalar adjustedTheta = theta;
38 
39  while (adjustedTheta < liftProfileStart_)
40  {
41  adjustedTheta += liftProfileEnd_ - liftProfileStart_;
42  }
43 
44  return adjustedTheta;
45  }
46  else if (theta > liftProfileEnd_)
47  {
48  scalar adjustedTheta = theta;
49 
50  while (adjustedTheta > liftProfileEnd_)
51  {
52  adjustedTheta -= liftProfileEnd_ - liftProfileStart_;
53  }
54 
55  return adjustedTheta;
56  }
57  else
58  {
59  return theta;
60  }
61 }
62 
63 
64 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
65 
66 // Construct from components
67 Foam::engineValve::engineValve
68 (
69  const word& name,
70  const polyMesh& mesh,
71  const autoPtr<coordinateSystem>& valveCS,
72  const word& bottomPatchName,
73  const word& poppetPatchName,
74  const word& stemPatchName,
75  const word& curtainInPortPatchName,
76  const word& curtainInCylinderPatchName,
77  const word& detachInCylinderPatchName,
78  const word& detachInPortPatchName,
79  const labelList& detachFaces,
80  const graph& liftProfile,
81  const scalar minLift,
82  const scalar minTopLayer,
83  const scalar maxTopLayer,
84  const scalar minBottomLayer,
85  const scalar maxBottomLayer,
86  const scalar diameter
87 )
88 :
89  name_(name),
90  mesh_(mesh),
91  engineDB_(refCast<const engineTime>(mesh.time())),
92  csPtr_(valveCS),
93  bottomPatch_(bottomPatchName, mesh.boundaryMesh()),
94  poppetPatch_(poppetPatchName, mesh.boundaryMesh()),
95  stemPatch_(stemPatchName, mesh.boundaryMesh()),
96  curtainInPortPatch_(curtainInPortPatchName, mesh.boundaryMesh()),
97  curtainInCylinderPatch_(curtainInCylinderPatchName, mesh.boundaryMesh()),
98  detachInCylinderPatch_(detachInCylinderPatchName, mesh.boundaryMesh()),
99  detachInPortPatch_(detachInPortPatchName, mesh.boundaryMesh()),
100  detachFaces_(detachFaces),
101  liftProfile_(liftProfile),
102  liftProfileStart_(min(liftProfile_.x())),
103  liftProfileEnd_(max(liftProfile_.x())),
104  minLift_(minLift),
105  minTopLayer_(minTopLayer),
106  maxTopLayer_(maxTopLayer),
107  minBottomLayer_(minBottomLayer),
108  maxBottomLayer_(maxBottomLayer),
109  diameter_(diameter)
110 {}
111 
112 
113 // Construct from dictionary
114 Foam::engineValve::engineValve
115 (
116  const word& name,
117  const polyMesh& mesh,
118  const dictionary& dict
119 )
120 :
121  name_(name),
122  mesh_(mesh),
123  engineDB_(refCast<const engineTime>(mesh_.time())),
124  csPtr_
125  (
127  (
128  mesh_,
129  dict.subDict("coordinateSystem")
130  )
131  ),
132  bottomPatch_(dict.lookup("bottomPatch"), mesh.boundaryMesh()),
133  poppetPatch_(dict.lookup("poppetPatch"), mesh.boundaryMesh()),
134  stemPatch_(dict.lookup("stemPatch"), mesh.boundaryMesh()),
135  curtainInPortPatch_
136  (
137  dict.lookup("curtainInPortPatch"),
138  mesh.boundaryMesh()
139  ),
140  curtainInCylinderPatch_
141  (
142  dict.lookup("curtainInCylinderPatch"),
143  mesh.boundaryMesh()
144  ),
145  detachInCylinderPatch_
146  (
147  dict.lookup("detachInCylinderPatch"),
148  mesh.boundaryMesh()
149  ),
150  detachInPortPatch_
151  (
152  dict.lookup("detachInPortPatch"),
153  mesh.boundaryMesh()
154  ),
155  detachFaces_(dict.lookup("detachFaces")),
156  liftProfile_("theta", "lift", name_, dict.lookup("liftProfile")),
157  liftProfileStart_(min(liftProfile_.x())),
158  liftProfileEnd_(max(liftProfile_.x())),
159  minLift_(readScalar(dict.lookup("minLift"))),
160  minTopLayer_(readScalar(dict.lookup("minTopLayer"))),
161  maxTopLayer_(readScalar(dict.lookup("maxTopLayer"))),
162  minBottomLayer_(readScalar(dict.lookup("minBottomLayer"))),
163  maxBottomLayer_(readScalar(dict.lookup("maxBottomLayer"))),
164  diameter_(readScalar(dict.lookup("diameter")))
165 {}
166 
167 
168 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
169 
170 
171 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
172 
173 Foam::scalar Foam::engineValve::lift(const scalar theta) const
174 {
175  return interpolateXY
176  (
177  adjustCrankAngle(theta),
178  liftProfile_.x(),
179  liftProfile_.y()
180  );
181 }
182 
183 
185 {
186  return lift(engineDB_.theta()) >= minLift_;
187 }
188 
189 
190 Foam::scalar Foam::engineValve::curLift() const
191 {
192  return max
193  (
194  lift(engineDB_.theta()),
195  minLift_
196  );
197 }
198 
199 
200 Foam::scalar Foam::engineValve::curVelocity() const
201 {
202  return
203  -(
204  curLift()
205  - max
206  (
207  lift(engineDB_.theta() - engineDB_.deltaTheta()),
208  minLift_
209  )
210  )/(engineDB_.deltaTValue() + VSMALL);
211 }
212 
213 
215 {
216  labelList mpIDs(2);
217  label nMpIDs = 0;
218 
219  if (bottomPatch_.active())
220  {
221  mpIDs[nMpIDs] = bottomPatch_.index();
222  nMpIDs++;
223  }
224 
225  if (poppetPatch_.active())
226  {
227  mpIDs[nMpIDs] = poppetPatch_.index();
228  nMpIDs++;
229  }
230 
231  mpIDs.setSize(nMpIDs);
232 
233  return mpIDs;
234 }
235 
236 
238 {
239  os << nl << name() << nl << token::BEGIN_BLOCK;
240 
241  cs().writeDict(os);
242 
243  os << "bottomPatch " << bottomPatch_.name() << token::END_STATEMENT << nl
244  << "poppetPatch " << poppetPatch_.name() << token::END_STATEMENT << nl
245  << "stemPatch " << stemPatch_.name() << token::END_STATEMENT << nl
246  << "curtainInPortPatch " << curtainInPortPatch_.name()
248  << "curtainInCylinderPatch " << curtainInCylinderPatch_.name()
250  << "detachInCylinderPatch " << detachInCylinderPatch_.name()
252  << "detachInPortPatch " << detachInPortPatch_.name()
254  << "detachFaces " << detachFaces_ << token::END_STATEMENT << nl
255  << "liftProfile " << nl << token::BEGIN_BLOCK
256  << liftProfile_ << token::END_BLOCK << token::END_STATEMENT << nl
257  << "minLift " << minLift_ << token::END_STATEMENT << nl
258  << "minTopLayer " << minTopLayer_ << token::END_STATEMENT << nl
259  << "maxTopLayer " << maxTopLayer_ << token::END_STATEMENT << nl
260  << "minBottomLayer " << minBottomLayer_ << token::END_STATEMENT << nl
261  << "maxBottomLayer " << maxBottomLayer_ << token::END_STATEMENT << nl
262  << "diameter " << diameter_ << token::END_STATEMENT << nl
263  << token::END_BLOCK << endl;
264 }
265 
266 
267 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:424
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
bool isOpen() const
Is the valve open?
Definition: engineValve.C:184
const scalarField & y() const
Definition: graph.C:137
const keyType & name() const
Return name.
Definition: DynamicID.H:95
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
void writeDict(Ostream &) const
Write dictionary.
Definition: engineValve.C:237
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:692
Class to create, store and output qgraph files.
Definition: graph.H:58
scalar theta() const
Return current crank-angle.
Definition: engineTime.C:139
scalar curVelocity() const
Return valve velocity for current time-step.
Definition: engineValve.C:200
A class for handling words, derived from string.
Definition: word.H:59
Interpolates y values from one curve to another with a different x distribution.
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:41
Field< Type > interpolateXY(const scalarField &xNew, const scalarField &xOld, const Field< Type > &yOld)
Definition: interpolateXY.C:38
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
label index() const
Return index of first matching zone.
Definition: DynamicID.H:107
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:262
const Time & time() const
Return time.
scalar curLift() const
Return current lift.
Definition: engineValve.C:190
scalar lift(const scalar theta) const
Return valve lift given crank angle in degrees.
Definition: engineValve.C:173
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
void setSize(const label)
Reset size of List.
Definition: List.C:281
static autoPtr< coordinateSystem > New(const objectRegistry &obr, const dictionary &dict)
Select constructed from dictionary and objectRegistry.
bool active() const
Has the zone been found.
Definition: DynamicID.H:113
const word & name() const
Return name.
Definition: engineValve.H:191
scalar deltaTheta() const
Return crank-angle increment.
Definition: engineTime.C:165
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
const coordinateSystem & cs() const
Return coordinate system.
Definition: engineValve.H:197
labelList movingPatchIDs() const
Return list of active patch labels for the valve head.
Definition: engineValve.C:214
const scalarField & x() const
Definition: graph.H:165
void writeDict(Ostream &, bool subDict=true) const
Write dictionary.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576