crankConRod.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) 2017-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 "crankConRod.H"
27 #include "unitConversion.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34  defineTypeNameAndDebug(crankConRod, 0);
35  addToRunTimeSelectionTable(engineTime, crankConRod, dictionary);
36 }
37 
38 
39 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
40 
41 void Foam::crankConRod::timeAdjustment()
42 {
43  deltaT_ = degToTime(deltaT_);
44  endTime_ = degToTime(endTime_);
45 
46  if
47  (
48  writeControl_ == wcRunTime
49  || writeControl_ == wcAdjustableRunTime
50  )
51  {
52  writeInterval_ = degToTime(writeInterval_);
53  }
54 }
55 
56 
57 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
58 
59 Foam::crankConRod::crankConRod
60 (
61  const word& name,
62  const fileName& rootPath,
63  const fileName& caseName,
64  const fileName& systemName,
65  const fileName& constantName,
66  const fileName& dictName
67 )
68 :
70  (
71  name,
72  rootPath,
73  caseName,
74  systemName,
75  constantName
76  ),
77  rpm_(dict_.lookup("rpm")),
78  conRodLength_(dimensionedScalar("conRodLength", dimLength, 0)),
79  bore_(dimensionedScalar("bore", dimLength, 0)),
80  stroke_(dimensionedScalar("stroke", dimLength, 0)),
81  clearance_(dimensionedScalar("clearance", dimLength, 0))
82 {
83  // geometric parameters are not strictly required for Time
84  dict_.readIfPresent("conRodLength", conRodLength_);
85  dict_.readIfPresent("bore", bore_);
86  dict_.readIfPresent("stroke", stroke_);
87  dict_.readIfPresent("clearance", clearance_);
88 
89  timeAdjustment();
90 
91  startTime_ = degToTime(startTime_);
92  value() = degToTime(value());
93 
94  deltaTSave_ = deltaT_;
95  deltaT0_ = deltaT_;
96 }
97 
98 
99 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
100 
102 {
103  Time::readDict();
104  timeAdjustment();
105 }
106 
107 
109 {
110  if (Time::read())
111  {
112  timeAdjustment();
113  return true;
114  }
115  else
116  {
117  return false;
118  }
119 }
120 
121 
122 Foam::scalar Foam::crankConRod::degToTime(const scalar theta) const
123 {
124  // 6 * rpm => deg/s
125  return theta/(6.0*rpm_.value());
126 }
127 
128 
129 Foam::scalar Foam::crankConRod::timeToDeg(const scalar t) const
130 {
131  // 6 * rpm => deg/s
132  return t*(6.0*rpm_.value());
133 }
134 
135 
136 Foam::scalar Foam::crankConRod::theta() const
137 {
138  return timeToDeg(value());
139 }
140 
141 
143 {
144  return " CAD";
145 }
146 
147 
149 {
150  scalar t = theta();
151 
152  while (t > 180.0)
153  {
154  t -= 360.0;
155  }
156 
157  while (t < -180.0)
158  {
159  t += 360.0;
160  }
161 
162  return t;
163 }
164 
165 
166 Foam::scalar Foam::crankConRod::deltaTheta() const
167 {
168  return timeToDeg(deltaTValue());
169 }
170 
171 
172 Foam::scalar Foam::crankConRod::pistonPosition(const scalar theta) const
173 {
174  return
175  (
176  conRodLength_.value()
177  + stroke_.value()/2.0
178  + clearance_.value()
179  )
180  - (
181  stroke_.value()*::cos(degToRad(theta))/2.0
182  + ::sqrt
183  (
184  sqr(conRodLength_.value())
185  - sqr(stroke_.value()*::sin(degToRad(theta))/2.0)
186  )
187  );
188 }
189 
190 
191 
192 Foam::scalar Foam::crankConRod::userTimeToTime(const scalar theta) const
193 {
194  return degToTime(theta);
195 }
196 
197 
198 Foam::scalar Foam::crankConRod::timeToUserTime(const scalar t) const
199 {
200  return timeToDeg(t);
201 }
202 
203 
204 // ************************************************************************* //
virtual void readDict()
Read the control dictionary and set the write controls etc.
Definition: crankConRod.C:101
A class for handling file names.
Definition: fileName.H:69
virtual scalar theta() const
Return current crank-angle.
Definition: crankConRod.C:136
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Unit conversion functions.
dimensionedScalar sqrt(const dimensionedScalar &ds)
bool readIfPresent(const dictionary &)
Update the value of dimensioned<Type> if found in the dictionary.
dimensionedScalar pistonPosition() const
Return current piston position.
Definition: engineTime.C:93
An abstract class for the time description of the piston motion.
Definition: engineTime.H:51
Macros for easy insertion into run-time selection tables.
scalar degToRad(const scalar deg)
Conversion from degrees to radians.
scalar timeToDeg(const scalar t) const
Convert seconds to degrees (for given engine speed in RPM)
Definition: crankConRod.C:129
virtual scalar deltaTheta() const
Return crank-angle increment.
Definition: crankConRod.C:166
virtual bool read()
Read the controlDict and set all the parameters.
Definition: crankConRod.C:108
dimensionedScalar cos(const dimensionedScalar &ds)
A class for handling words, derived from string.
Definition: word.H:59
virtual scalar timeToUserTime(const scalar t) const
Convert the real-time (s) into user-time (CA deg)
Definition: crankConRod.C:198
dimensionedScalar sin(const dimensionedScalar &ds)
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
virtual void readDict()
Read the control dictionary and set the write controls etc.
Definition: TimeIO.C:35
scalar degToTime(const scalar theta) const
Convert degrees to seconds (for given engine speed in RPM)
Definition: crankConRod.C:122
defineTypeNameAndDebug(combustionModel, 0)
virtual scalar userTimeToTime(const scalar theta) const
Convert the user-time (CA deg) to real-time (s).
Definition: crankConRod.C:192
virtual bool read()
Read control dictionary, update controls and time.
Definition: TimeIO.C:431
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual word unit() const
Return time unit.
Definition: crankConRod.C:142
scalar thetaRevolution() const
Return current crank-angle translated to a single revolution.
Definition: crankConRod.C:148
Namespace for OpenFOAM.