Merge commit 'd3b269b7c6ffa0cd68845adfecdfb849316dba71' into nextRelease
[foam-extend-3.2.git] / src / lagrangian / dieselSpray / spray / spray.C
blobd9bcfe2336e2619ebefc6913e8b04902af9277dd
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | foam-extend: Open Source CFD
4    \\    /   O peration     |
5     \\  /    A nd           | For copyright notice see file Copyright
6      \\/     M anipulation  |
7 -------------------------------------------------------------------------------
8 License
9     This file is part of foam-extend.
11     foam-extend is free software: you can redistribute it and/or modify it
12     under the terms of the GNU General Public License as published by the
13     Free Software Foundation, either version 3 of the License, or (at your
14     option) any later version.
16     foam-extend is distributed in the hope that it will be useful, but
17     WITHOUT ANY WARRANTY; without even the implied warranty of
18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
19     General Public License for more details.
21     You should have received a copy of the GNU General Public License
22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
24 \*---------------------------------------------------------------------------*/
26 #include "spray.H"
28 #include "atomizationModel.H"
29 #include "breakupModel.H"
30 #include "collisionModel.H"
31 #include "dispersionModel.H"
32 #include "dragModel.H"
33 #include "evaporationModel.H"
34 #include "heatTransferModel.H"
35 #include "injectorModel.H"
36 #include "wallModel.H"
38 #include "basicMultiComponentMixture.H"
40 #include "symmetryPolyPatch.H"
41 #include "wedgePolyPatch.H"
43 #include "mathematicalConstants.H"
45 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
47 defineTemplateTypeNameAndDebug(IOPtrList<injector>, 0);
49 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
51 // Construct from components
52 Foam::spray::spray
54     const volVectorField& U,
55     const volScalarField& rho,
56     const volScalarField& p,
57     const volScalarField& T,
58     const basicMultiComponentMixture& composition,
59     const PtrList<gasThermoPhysics>& gasProperties,
60     const dictionary&,
61     const dimensionedVector& g,
62     bool readFields
65     Cloud<parcel>(U.mesh(), false), // suppress className checking on positions
66     runTime_(U.time()),
67     time0_(runTime_.value()),
68     mesh_(U.mesh()),
69     rndGen_(label(0)),
70     g_(g.value()),
72     U_(U),
73     rho_(rho),
74     p_(p),
75     T_(T),
77     sprayProperties_
78     (
79         IOobject
80         (
81             "sprayProperties",
82             U.time().constant(),
83             U.db(),
84             IOobject::MUST_READ,
85             IOobject::NO_WRITE
86         )
87     ),
89     ambientPressure_(p_.average().value()),
90     ambientTemperature_(T_.average().value()),
92     injectors_
93     (
94         IOobject
95         (
96             "injectorProperties",
97             U.time().constant(),
98             U.db(),
99             IOobject::MUST_READ,
100             IOobject::NO_WRITE
101         ),
102         injector::iNew(U.time())
103     ),
104     atomization_
105     (
106         atomizationModel::New
107         (
108             sprayProperties_,
109             *this
110         )
111     ),
112     drag_
113     (
114         dragModel::New
115         (
116             sprayProperties_
117         )
118     ),
119     evaporation_
120     (
121         evaporationModel::New
122         (
123             sprayProperties_
124         )
125     ),
126     heatTransfer_
127     (
128         heatTransferModel::New
129         (
130             sprayProperties_
131         )
132     ),
133     wall_
134     (
135         wallModel::New
136         (
137             sprayProperties_,
138             U,
139             *this
140         )
141     ),
142     breakupModel_
143     (
144         breakupModel::New
145         (
146             sprayProperties_,
147             *this
148         )
149     ),
150     collisionModel_
151     (
152         collisionModel::New
153         (
154             sprayProperties_,
155             *this,
156             rndGen_
157         )
158     ),
159     dispersionModel_
160     (
161         dispersionModel::New
162         (
163             sprayProperties_,
164             *this
165         )
166     ),
168     fuels_
169     (
170         liquidMixture::New
171         (
172             mesh_.lookupObject<dictionary>("thermophysicalProperties")
173         )
174     ),
175     injectorModel_
176     (
177         injectorModel::New
178         (
179             sprayProperties_,
180             *this
181         )
182     ),
184     sprayIteration_(sprayProperties_.subDict("sprayIteration")),
185     sprayIterate_(readLabel(sprayIteration_.lookup("sprayIterate"))),
186     sprayRelaxFactor_(readScalar(sprayIteration_.lookup("sprayRelaxFactor"))),
187     minimumParcelMass_
188     (
189         readScalar(sprayIteration_.lookup("minimumParcelMass"))
190     ),
192     subCycles_(readLabel(sprayProperties_.lookup("subCycles"))),
194     gasProperties_(gasProperties),
195     composition_(composition),
197     liquidToGasIndex_(fuels_->components().size(), -1),
198     gasToLiquidIndex_(composition.Y().size(), -1),
199     isLiquidFuel_(composition.Y().size(), false),
201     twoD_(0),
202     axisOfSymmetry_(vector::zero),
203     axisOfWedge_(vector(0,0,0)),
204     axisOfWedgeNormal_(vector(0,0,0)),
205     angleOfWedge_(0.0),
207     interpolationSchemes_(sprayProperties_.subDict("interpolationSchemes")),
208     UInterpolator_(NULL),
209     rhoInterpolator_(NULL),
210     pInterpolator_(NULL),
211     TInterpolator_(NULL),
213     sms_(mesh_.nCells(), vector::zero),
214     shs_(mesh_.nCells(), 0.0),
215     srhos_(fuels_->components().size()),
217     totalInjectedLiquidMass_(0.0),
218     injectedLiquidKE_(0.0)
221     // create the evaporation source fields
222     forAll(srhos_, i)
223     {
224         srhos_.set(i, new scalarField(mesh_.nCells(), 0.0));
225     }
227     // Write some information about injection parameters
228     forAll(injectors_, i)
229     {
230         const injectorType& it = injectors_[i].properties();
232         scalar v = injection().averageVelocity(i);
234         scalar ip = it.integrateTable(it.injectionPressureProfile());
235         scalar dt = it.teoi() - it.tsoi();
236         Info<< "Average Velocity for injector " << i << ": " << v << " m/s"
237             << ", injection pressure = "
238             << 1.0e-5*ip/dt << " bar"
239             << endl;
240     }
242     // Check if the case is 2D wedge
243     const polyBoundaryMesh& bMesh = mesh().boundaryMesh();
244     bool symPlaneExist = false;
245     bool wedgeExist = false;
246     label patches[2];
247     label n=0;
249     // check for the type of boundary condition
250     forAll(bMesh, patchi)
251     {
252         if (isA<symmetryPolyPatch>(bMesh[patchi]))
253         {
254             symPlaneExist = true;
255         }
256         else if (isA<wedgePolyPatch>(bMesh[patchi]))
257         {
258             wedgeExist = true;
259             patches[n++] = patchi;
260         }
261     }
263     // if wedge exist we assume that this is a 2D run.
264     twoD_ = wedgeExist;
266     if (twoD_)
267     {
268         if (n<2)
269         {
270             FatalErrorIn
271             (
272                 "spray::spray(const volVectorField& U, "
273                 "const volScalarField& rho, const volScalarField& p, "
274                 "const volScalarField& T, const combustionMixture& composition,"
275                 "const PtrList<gasThermoPhsyics>& gaseousFuelProperties, "
276                 "const dictionary& thermophysicalProperties, "
277                 "const dimensionedScalar& g)"
278             )   << "spray::(...) only one wedgePolyPatch found. "
279                    "Please check you BC-setup."
280                 << abort(FatalError);
281         }
283         Info<< "Constructing two dimensional spray injection.";
285         vector v1 = bMesh[patches[0]].faceAreas()[0];
286         vector v2 = bMesh[patches[1]].faceAreas()[0];
287         v1 /= mag(v1);
288         v2 /= mag(v2);
289         axisOfSymmetry_ = v1 ^ v2;
290         axisOfSymmetry_ /= mag(axisOfSymmetry_);
292         // assuming that 'v2' is the 'front' face
293         axisOfWedge_ = axisOfSymmetry_ ^ v2;
294         axisOfWedge_ /= mag(axisOfWedge_);
296         axisOfWedgeNormal_ = axisOfSymmetry_ ^ axisOfWedge_;
297         axisOfWedgeNormal_ /= mag(axisOfWedgeNormal_);
299         scalar arcCos = (v1 & v2)/mag(v1);
300         angleOfWedge_ = mathematicalConstant::pi - acos(arcCos);
302         Info<< "Calculated angle of wedge is "
303             << angleOfWedge_*180/mathematicalConstant::pi << " deg."
304             << endl;
305     }
306     else
307     {
308         if (symPlaneExist)
309         {
310             angleOfWedge_ = mathematicalConstant::pi;
311             Info<< "Constructing 180 deg three dimensional spray injection."
312                 << endl;
313         }
314         else
315         {
316             Info<< "Constructing three dimensional spray injection." << endl;
317         }
319     }
321     // find index mapping between liquid indeces and gas indeces
322     label Ns = composition_.Y().size();
324     forAll(fuels_->components(), i)
325     {
326         word liquidName(fuels_->components()[i]);
328         for (label j=0; j<Ns; j++)
329         {
330             word specieName(composition_.Y()[j].name());
332             if (specieName == liquidName)
333             {
334                 liquidToGasIndex_[i] = j;
335                 gasToLiquidIndex_[j] = i;
336                 isLiquidFuel_[j] = true;
337             }
338         }
339         if (liquidToGasIndex_[i] == -1)
340         {
341             Info << "In composition:" << endl;
342             for (label k=0; k<Ns; k++)
343             {
344                 word specieName(composition_.Y()[k].name());
345                 Info << specieName << endl;
346             }
348             FatalError<<
349                 "The liquid component " << liquidName
350                 << " does not exist in the species composition.Y() list.\n"
351                 << "(Probably not defined in <chem.inp>)"
352                 << abort(FatalError);
353         }
354     }
356     if (readFields)
357     {
358         parcel::readFields(*this);
359     }
363 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
365 Foam::spray::~spray()
369 // ************************************************************************* //