FreeFOAM The Cross-Platform CFD Toolkit
GuldersEGR.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 "GuldersEGR.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace laminarFlameSpeedModels
34 {
35  defineTypeNameAndDebug(GuldersEGR, 0);
36 
38  (
39  laminarFlameSpeed,
40  GuldersEGR,
41  dictionary
42  );
43 }
44 }
45 
46 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 
48 Foam::laminarFlameSpeedModels::GuldersEGR::GuldersEGR
49 (
50  const dictionary& dict,
51  const hhuCombustionThermo& ct
52 )
53 :
54  laminarFlameSpeed(dict, ct),
55 
56  coeffsDict_(dict.subDict(typeName + "Coeffs").subDict(fuel_)),
57  W_(readScalar(coeffsDict_.lookup("W"))),
58  eta_(readScalar(coeffsDict_.lookup("eta"))),
59  xi_(readScalar(coeffsDict_.lookup("xi"))),
60  f_(readScalar(coeffsDict_.lookup("f"))),
61  alpha_(readScalar(coeffsDict_.lookup("alpha"))),
62  beta_(readScalar(coeffsDict_.lookup("beta")))
63 {}
64 
65 
66 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
67 
69 {}
70 
71 
72 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
73 
74 inline Foam::scalar Foam::laminarFlameSpeedModels::GuldersEGR::SuRef
75 (
76  scalar phi
77 ) const
78 {
79  if (phi > SMALL)
80  {
81  return W_*pow(phi, eta_)*exp(-xi_*sqr(phi - 1.075));
82  }
83  else
84  {
85  return 0.0;
86  }
87 }
88 
89 
90 inline Foam::scalar Foam::laminarFlameSpeedModels::GuldersEGR::Su0pTphi
91 (
92  scalar p,
93  scalar Tu,
94  scalar phi,
95  scalar Yres
96 ) const
97 {
98  static const scalar Tref = 300.0;
99  static const scalar pRef = 1.013e5;
100 
101  return SuRef(phi)*pow((Tu/Tref), alpha_)*pow((p/pRef), beta_)*(1 - f_*Yres);
102 }
103 
104 
106 Foam::laminarFlameSpeedModels::GuldersEGR::Su0pTphi
107 (
108  const volScalarField& p,
109  const volScalarField& Tu,
110  scalar phi
111 ) const
112 {
113  tmp<volScalarField> tSu0
114  (
115  new volScalarField
116  (
117  IOobject
118  (
119  "Su0",
120  p.time().timeName(),
121  p.db(),
124  ),
125  p.mesh(),
126  dimensionedScalar("Su0", dimVelocity, 0.0)
127  )
128  );
129 
130  volScalarField& Su0 = tSu0();
131 
132  forAll(Su0, celli)
133  {
134  Su0[celli] = Su0pTphi(p[celli], Tu[celli], phi, 0.0);
135  }
136 
137  forAll(Su0.boundaryField(), patchi)
138  {
139  forAll(Su0.boundaryField()[patchi], facei)
140  {
141  Su0.boundaryField()[patchi][facei] =
142  Su0pTphi
143  (
144  p.boundaryField()[patchi][facei],
145  Tu.boundaryField()[patchi][facei],
146  phi,
147  0.0
148  );
149  }
150  }
151 
152  return tSu0;
153 }
154 
155 
157 Foam::laminarFlameSpeedModels::GuldersEGR::Su0pTphi
158 (
159  const volScalarField& p,
160  const volScalarField& Tu,
161  const volScalarField& phi,
162  const volScalarField& egr
163 ) const
164 {
165  tmp<volScalarField> tSu0
166  (
167  new volScalarField
168  (
169  IOobject
170  (
171  "Su0",
172  p.time().timeName(),
173  p.db(),
176  ),
177  p.mesh(),
178  dimensionedScalar("Su0", dimVelocity, 0.0)
179  )
180  );
181 
182  volScalarField& Su0 = tSu0();
183 
184  forAll(Su0, celli)
185  {
186  Su0[celli] = Su0pTphi(p[celli], Tu[celli], phi[celli], egr[celli]);
187  }
188 
189  forAll(Su0.boundaryField(), patchi)
190  {
191  forAll(Su0.boundaryField()[patchi], facei)
192  {
193  Su0.boundaryField()[patchi][facei] =
194  Su0pTphi
195  (
196  p.boundaryField()[patchi][facei],
197  Tu.boundaryField()[patchi][facei],
198  phi.boundaryField()[patchi][facei],
199  egr.boundaryField()[patchi][facei]
200  );
201  }
202  }
203 
204  return tSu0;
205 }
206 
207 
210 {
211  if
212  (
213  hhuCombustionThermo_.composition().contains("ft")
214  && hhuCombustionThermo_.composition().contains("egr")
215  )
216  {
217  return Su0pTphi
218  (
219  hhuCombustionThermo_.p(),
220  hhuCombustionThermo_.Tu(),
222  (
223  hhuCombustionThermo_.lookup("stoichiometricAirFuelMassRatio")
224  )/
225  (
226  scalar(1)/hhuCombustionThermo_.composition().Y("ft")
227  - scalar(1)
228  ),
229  hhuCombustionThermo_.composition().Y("egr")
230  );
231  }
232  else
233  {
234  return Su0pTphi
235  (
236  hhuCombustionThermo_.p(),
237  hhuCombustionThermo_.Tu(),
238  equivalenceRatio_
239  );
240  }
241 }
242 
243 
244 // ************************ vim: set sw=4 sts=4 et: ************************ //