FreeFOAM The Cross-Platform CFD Toolkit
ConeInjectionMP.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) 2008-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 "ConeInjectionMP.H"
29 
30 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
31 
32 template<class CloudType>
34 (
35  const scalar time0,
36  const scalar time1
37 ) const
38 {
39  if ((time0 >= 0.0) && (time0 < duration_))
40  {
41  const scalar targetVolume = volumeFlowRate_().integrate(0, time1);
42 
43  const label targetParcels =
44  parcelsPerInjector_*targetVolume/this->volumeTotal_;
45 
46  const label nToInject = targetParcels - nInjected_;
47 
48  nInjected_ += nToInject;
49 
50  return positions_.size()*nToInject;
51  }
52  else
53  {
54  return 0;
55  }
56 }
57 
58 
59 template<class CloudType>
61 (
62  const scalar time0,
63  const scalar time1
64 ) const
65 {
66  if ((time0 >= 0.0) && (time0 < duration_))
67  {
68  return volumeFlowRate_().integrate(time0, time1);
69  }
70  else
71  {
72  return 0.0;
73  }
74 }
75 
76 
77 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
78 
79 template<class CloudType>
81 (
82  const dictionary& dict,
83  CloudType& owner
84 )
85 :
86  InjectionModel<CloudType>(dict, owner, typeName),
87  positionsFile_(this->coeffDict().lookup("positionsFile")),
88  positions_
89  (
90  IOobject
91  (
92  positionsFile_,
93  owner.db().time().constant(),
94  owner.mesh(),
95  IOobject::MUST_READ,
96  IOobject::NO_WRITE
97  )
98  ),
99  injectorCells_(positions_.size()),
100  axesFile_(this->coeffDict().lookup("axesFile")),
101  axes_
102  (
103  IOobject
104  (
105  axesFile_,
106  owner.db().time().constant(),
107  owner.mesh(),
108  IOobject::MUST_READ,
109  IOobject::NO_WRITE
110  )
111  ),
112  duration_(readScalar(this->coeffDict().lookup("duration"))),
113  parcelsPerInjector_
114  (
115  readScalar(this->coeffDict().lookup("parcelsPerInjector"))
116  ),
117  volumeFlowRate_
118  (
120  (
121  "volumeFlowRate",
122  this->coeffDict()
123  )
124  ),
125  Umag_
126  (
128  (
129  "Umag",
130  this->coeffDict()
131  )
132  ),
133  thetaInner_
134  (
136  (
137  "thetaInner",
138  this->coeffDict()
139  )
140  ),
141  thetaOuter_
142  (
144  (
145  "thetaOuter",
146  this->coeffDict()
147  )
148  ),
149  parcelPDF_
150  (
151  pdfs::pdf::New
152  (
153  this->coeffDict().subDict("parcelPDF"),
154  owner.rndGen()
155  )
156  ),
157  nInjected_(this->parcelsAddedTotal()),
158  tanVec1_(positions_.size()),
159  tanVec2_(positions_.size())
160 {
161  // Normalise direction vector and determine direction vectors
162  // tangential to direction
163  forAll(axes_, i)
164  {
165  axes_[i] /= mag(axes_[i]);
166 
167  vector tangent = vector::zero;
168  scalar magTangent = 0.0;
169 
170  while (magTangent < SMALL)
171  {
172  vector v = this->owner().rndGen().vector01();
173 
174  tangent = v - (v & axes_[i])*axes_[i];
175  magTangent = mag(tangent);
176  }
177 
178  tanVec1_[i] = tangent/magTangent;
179  tanVec2_[i] = axes_[i]^tanVec1_[i];
180  }
181 
182  // Set total volume to inject
183  this->volumeTotal_ = volumeFlowRate_().integrate(0.0, duration_);
184 
185  // Set/cache the injector cells
186  forAll(positions_, i)
187  {
188  this->findCellAtPosition
189  (
190  injectorCells_[i],
191  positions_[i]
192  );
193  }
194 }
195 
196 
197 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
198 
199 template<class CloudType>
201 {}
202 
203 
204 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
205 
206 template<class CloudType>
208 {
209  return true;
210 }
211 
212 
213 template<class CloudType>
215 {
216  return this->SOI_ + duration_;
217 }
218 
219 
220 template<class CloudType>
222 (
223  const label parcelI,
224  const label,
225  const scalar,
226  vector& position,
227  label& cellOwner
228 )
229 {
230  const label i = parcelI%positions_.size();
231 
232  position = positions_[i];
233  cellOwner = injectorCells_[i];
234 }
235 
236 
237 template<class CloudType>
239 (
240  const label parcelI,
241  const label,
242  const scalar time,
243  typename CloudType::parcelType& parcel
244 )
245 {
246  // set particle velocity
247  const label i = parcelI%positions_.size();
248 
249  const scalar deg2Rad = mathematicalConstant::pi/180.0;
250 
251  scalar t = time - this->SOI_;
252  scalar ti = thetaInner_().value(t);
253  scalar to = thetaOuter_().value(t);
254  scalar coneAngle = this->owner().rndGen().scalar01()*(to - ti) + ti;
255 
256  coneAngle *= deg2Rad;
257  scalar alpha = sin(coneAngle);
258  scalar dcorr = cos(coneAngle);
259  scalar beta = mathematicalConstant::twoPi*this->owner().rndGen().scalar01();
260 
261  vector normal = alpha*(tanVec1_[i]*cos(beta) + tanVec2_[i]*sin(beta));
262  vector dirVec = dcorr*axes_[i];
263  dirVec += normal;
264 
265  dirVec /= mag(dirVec);
266 
267  parcel.U() = Umag_().value(t)*dirVec;
268 
269  // set particle diameter
270  parcel.d() = parcelPDF_().sample();
271 }
272 
273 
274 template<class CloudType>
276 {
277  return false;
278 }
279 
280 
281 template<class CloudType>
283 {
284  return true;
285 }
286 
287 
288 // ************************ vim: set sw=4 sts=4 et: ************************ //