FreeFOAM The Cross-Platform CFD Toolkit
engineTime.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) 1991-2010 OpenCFD Ltd.
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 "engineTime.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 void Foam::engineTime::timeAdjustment()
32 {
35 
36  if
37  (
40  )
41  {
43  }
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
49 //- Construct from objectRegistry arguments
50 Foam::engineTime::engineTime
51 (
52  const word& name,
53  const fileName& rootPath,
54  const fileName& caseName,
55  const fileName& systemName,
56  const fileName& constantName,
57  const fileName& dictName
58 )
59 :
60  Time
61  (
62  name,
63  rootPath,
64  caseName,
65  systemName,
66  constantName
67  ),
68  dict_
69  (
70  IOobject
71  (
72  "engineGeometry",
73  constant(),
74  *this,
77  false
78  )
79  ),
80  rpm_(dict_.lookup("rpm")),
81  conRodLength_(dimensionedScalar("conRodLength", dimLength, 0)),
82  bore_(dimensionedScalar("bore", dimLength, 0)),
83  stroke_(dimensionedScalar("stroke", dimLength, 0)),
84  clearance_(dimensionedScalar("clearance", dimLength, 0))
85 {
86  // geometric parameters are not strictly required for Time
87  dict_.readIfPresent("conRodLength", conRodLength_);
88  dict_.readIfPresent("bore", bore_);
89  dict_.readIfPresent("stroke", stroke_);
90  dict_.readIfPresent("clearance", clearance_);
91 
92  timeAdjustment();
93 
94  startTime_ = degToTime(startTime_);
95  value() = degToTime(value());
96  deltaT0_ = deltaT_;
97 }
98 
99 
100 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
101 
102 // Read the controlDict and set all the parameters
104 {
105  Time::readDict();
106  timeAdjustment();
107 }
108 
109 
110 // Read the controlDict and set all the parameters
112 {
113  if (Time::read())
114  {
115  timeAdjustment();
116  return true;
117  }
118  else
119  {
120  return false;
121  }
122 }
123 
124 
125 Foam::scalar Foam::engineTime::degToRad(const scalar deg) const
126 {
127  return mathematicalConstant::pi*deg/180.0;
128 }
129 
130 
131 Foam::scalar Foam::engineTime::degToTime(const scalar theta) const
132 {
133  // 6 * rpm => deg/s
134  return theta/(6.0*rpm_.value());
135 }
136 
137 
138 Foam::scalar Foam::engineTime::timeToDeg(const scalar t) const
139 {
140  // 6 * rpm => deg/s
141  return t*(6.0*rpm_.value());
142 }
143 
144 
145 Foam::scalar Foam::engineTime::theta() const
146 {
147  return timeToDeg(value());
148 }
149 
150 
151 // Return current crank-angle translated to a single revolution
152 // (value between -180 and 180 with 0 = top dead centre)
154 {
155  scalar t = theta();
156 
157  while (t > 180.0)
158  {
159  t -= 360.0;
160  }
161 
162  while (t < -180.0)
163  {
164  t += 360.0;
165  }
166 
167  return t;
168 }
169 
170 
171 Foam::scalar Foam::engineTime::deltaTheta() const
172 {
173  return timeToDeg(deltaT().value());
174 }
175 
176 
177 Foam::scalar Foam::engineTime::pistonPosition(const scalar theta) const
178 {
179  return
180  (
181  conRodLength_.value()
182  + stroke_.value()/2.0
183  + clearance_.value()
184  )
185  - (
186  stroke_.value()*::cos(degToRad(theta))/2.0
187  + ::sqrt
188  (
189  sqr(conRodLength_.value())
190  - sqr(stroke_.value()*::sin(degToRad(theta))/2.0)
191  )
192  );
193 }
194 
195 
197 {
198  return dimensionedScalar
199  (
200  "pistonPosition",
201  dimLength,
202  pistonPosition(theta())
203  );
204 }
205 
206 
208 {
209  return dimensionedScalar
210  (
211  "pistonDisplacement",
212  dimLength,
213  pistonPosition(theta() - deltaTheta()) - pistonPosition().value()
214  );
215 }
216 
217 
219 {
220  return dimensionedScalar
221  (
222  "pistonSpeed",
223  dimVelocity,
224  pistonDisplacement().value()/(deltaT().value() + VSMALL)
225  );
226 }
227 
228 
229 Foam::scalar Foam::engineTime::userTimeToTime(const scalar theta) const
230 {
231  return degToTime(theta);
232 }
233 
234 
235 Foam::scalar Foam::engineTime::timeToUserTime(const scalar t) const
236 {
237  return timeToDeg(t);
238 }
239 
240 
241 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
242 
243 // ************************ vim: set sw=4 sts=4 et: ************************ //