FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

constInjector.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 \*---------------------------------------------------------------------------*/
00025 
00026 #include "constInjector.H"
00027 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00028 #include <OpenFOAM/mathematicalConstants.H>
00029 
00030 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00031 
00032 namespace Foam
00033 {
00034 
00035 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
00036 
00037 defineTypeNameAndDebug(constInjector, 0);
00038 
00039 addToRunTimeSelectionTable
00040 (
00041     injectorModel,
00042     constInjector,
00043     dictionary
00044 );
00045 
00046 
00047 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00048 
00049 // Construct from components
00050 constInjector::constInjector
00051 (
00052     const dictionary& dict,
00053     spray& sm
00054 )
00055 :
00056     injectorModel(dict, sm),
00057     specDict_(dict.subDict(typeName + "Coeffs")),
00058     dropletNozzleDiameterRatio_(specDict_.lookup("dropletNozzleDiameterRatio")),
00059     sprayAngle_(specDict_.lookup("sprayAngle"))
00060 {
00061     if (sm.injectors().size() != dropletNozzleDiameterRatio_.size())
00062     {
00063         FatalError << "constInjector::constInjector"
00064             << "(const dictionary& dict, spray& sm)\n"
00065             << "Wrong number of entries in dropletNozzleDiameterRatio"
00066             << abort(FatalError);
00067     }
00068 
00069     if (sm.injectors().size() != sprayAngle_.size())
00070     {
00071         FatalError << "constInjector::constInjector"
00072             << "(const dictionary& dict, spray& sm)\n"
00073             << "Wrong number of entries in sprayAngle"
00074             << abort(FatalError);
00075     }
00076 
00077     scalar referencePressure = sm.p().average().value();
00078 
00079     // correct velocity and pressure profiles
00080     forAll(sm.injectors(), i)
00081     {
00082         sm.injectors()[i].properties()->correctProfiles(sm.fuels(), referencePressure);
00083     }
00084 
00085 }
00086 
00087 
00088 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
00089 
00090 constInjector::~constInjector()
00091 {}
00092 
00093 
00094 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00095 
00096 scalar constInjector::d0
00097 (
00098     const label n, 
00099     const scalar
00100 ) const
00101 {
00102     return injectors_[n].properties()->d()*dropletNozzleDiameterRatio_[n];
00103 }
00104 
00105 
00106 vector constInjector::direction
00107 (
00108     const label n,
00109     const label hole,
00110     const scalar time,
00111     const scalar d
00112 ) const
00113 {
00114 
00115     /*
00116         randomly distribute parcels in a solid cone
00117         with angle = sprayAngle,
00118         alpha = radius of the two normal vectors,
00119         = maximum sin(sprayAngle/2)
00120         beta = angle in the normal plane
00121         
00122                         o                    / (beta)
00123                         |\                  /
00124                         | \                /)
00125                         |  \              o-----------> (x-axis)
00126                         |   \
00127                         v  (alpha)
00128         */
00129 
00130     scalar angle = rndGen_.scalar01()*sprayAngle_[n]*mathematicalConstant::pi/360.0;
00131     scalar alpha = sin(angle);
00132     scalar dcorr = cos(angle);
00133 
00134     scalar beta = 2.0*mathematicalConstant::pi*rndGen_.scalar01();
00135 
00136     // randomly distributed vector normal to the injection vector
00137     vector normal = vector::zero;
00138     
00139     if (sm_.twoD())
00140     {
00141         scalar reduce = 0.01;
00142         // correct beta if this is a 2D run
00143         // map it onto the 'angleOfWedge'
00144         beta *= (1.0-2.0*reduce)*0.5*sm_.angleOfWedge()/mathematicalConstant::pi;
00145         beta += reduce*sm_.angleOfWedge();
00146 
00147         normal = alpha*
00148         (
00149             sm_.axisOfWedge()*cos(beta) +
00150             sm_.axisOfWedgeNormal()*sin(beta)
00151         );
00152 
00153     }
00154     else
00155     {
00156         normal = alpha*
00157         (
00158             injectors_[n].properties()->tan1(hole)*cos(beta) +
00159             injectors_[n].properties()->tan2(hole)*sin(beta)
00160         );
00161     }
00162     
00163     // set the direction of injection by adding the normal vector
00164     vector dir = dcorr*injectors_[n].properties()->direction(n, time) + normal;
00165     dir /= mag(dir);
00166 
00167     return dir;
00168 }
00169 
00170 scalar constInjector::velocity
00171 (
00172     const label i,
00173     const scalar time
00174 ) const
00175 {
00176     const injectorType& it = sm_.injectors()[i].properties();
00177     if (it.pressureIndependentVelocity())
00178     {
00179         return it.getTableValue(it.velocityProfile(), time);
00180     }
00181     else
00182     {
00183         scalar Pref = sm_.ambientPressure();
00184         scalar Pinj = it.getTableValue(it.injectionPressureProfile(), time);
00185         scalar rho = sm_.fuels().rho(Pinj, it.T(time), it.X());
00186         scalar dp = max(0.0, Pinj - Pref);
00187         return sqrt(2.0*dp/rho);
00188     }
00189 }
00190 
00191 scalar constInjector::averageVelocity
00192 (
00193     const label i
00194 ) const
00195 {    
00196     const injectorType& it = sm_.injectors()[i].properties();
00197     scalar dt = it.teoi() - it.tsoi();
00198 
00199     return it.integrateTable(it.velocityProfile())/dt;
00200 }
00201 
00202 } // End namespace Foam
00203 
00204 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines