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

turbulentIntensityKineticEnergyInletFvPatchScalarField.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) 2006-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 "turbulentIntensityKineticEnergyInletFvPatchScalarField.H"
00027 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00028 #include <finiteVolume/fvPatchFieldMapper.H>
00029 #include <finiteVolume/surfaceFields.H>
00030 #include <finiteVolume/volFields.H>
00031 
00032 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00033 
00034 Foam::turbulentIntensityKineticEnergyInletFvPatchScalarField::
00035 turbulentIntensityKineticEnergyInletFvPatchScalarField
00036 (
00037     const fvPatch& p,
00038     const DimensionedField<scalar, volMesh>& iF
00039 )
00040 :
00041     inletOutletFvPatchScalarField(p, iF),
00042     intensity_(0.0),
00043     UName_("undefined-U"),
00044     phiName_("undefined-phi")
00045 {
00046     this->refValue() = 0.0;
00047     this->refGrad() = 0.0;
00048     this->valueFraction() = 0.0;
00049 }
00050 
00051 Foam::turbulentIntensityKineticEnergyInletFvPatchScalarField::
00052 turbulentIntensityKineticEnergyInletFvPatchScalarField
00053 (
00054     const turbulentIntensityKineticEnergyInletFvPatchScalarField& ptf,
00055     const fvPatch& p,
00056     const DimensionedField<scalar, volMesh>& iF,
00057     const fvPatchFieldMapper& mapper
00058 )
00059 :
00060     inletOutletFvPatchScalarField(ptf, p, iF, mapper),
00061     intensity_(ptf.intensity_),
00062     UName_(ptf.UName_),
00063     phiName_(ptf.phiName_)
00064 {}
00065 
00066 Foam::turbulentIntensityKineticEnergyInletFvPatchScalarField::
00067 turbulentIntensityKineticEnergyInletFvPatchScalarField
00068 (
00069     const fvPatch& p,
00070     const DimensionedField<scalar, volMesh>& iF,
00071     const dictionary& dict
00072 )
00073 :
00074     inletOutletFvPatchScalarField(p, iF),
00075     intensity_(readScalar(dict.lookup("intensity"))),
00076     UName_(dict.lookupOrDefault<word>("U", "U")),
00077     phiName_(dict.lookupOrDefault<word>("phi", "phi"))
00078 {
00079     if (intensity_ < 0 || intensity_ > 1)
00080     {
00081         FatalErrorIn
00082         (
00083             "turbulentIntensityKineticEnergyInletFvPatchScalarField::"
00084             "turbulentIntensityKineticEnergyInletFvPatchScalarField"
00085             "(const fvPatch& p, const DimensionedField<scalar, volMesh>& iF, "
00086             "const dictionary& dict)"
00087         )   << "Turbulence intensity should be specified as a fraction 0-1 "
00088                "of the mean velocity\n"
00089                "    value given is " << intensity_
00090             << "\n    on patch " << this->patch().name()
00091             << " of field " << this->dimensionedInternalField().name()
00092             << " in file " << this->dimensionedInternalField().objectPath()
00093             << exit(FatalError);
00094     }
00095 
00096     fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
00097 
00098     this->refValue() = 0.0;
00099     this->refGrad() = 0.0;
00100     this->valueFraction() = 0.0;
00101 }
00102 
00103 Foam::turbulentIntensityKineticEnergyInletFvPatchScalarField::
00104 turbulentIntensityKineticEnergyInletFvPatchScalarField
00105 (
00106     const turbulentIntensityKineticEnergyInletFvPatchScalarField& ptf
00107 )
00108 :
00109     inletOutletFvPatchScalarField(ptf),
00110     intensity_(ptf.intensity_),
00111     UName_(ptf.UName_),
00112     phiName_(ptf.phiName_)
00113 {}
00114 
00115 
00116 Foam::turbulentIntensityKineticEnergyInletFvPatchScalarField::
00117 turbulentIntensityKineticEnergyInletFvPatchScalarField
00118 (
00119     const turbulentIntensityKineticEnergyInletFvPatchScalarField& ptf,
00120     const DimensionedField<scalar, volMesh>& iF
00121 )
00122 :
00123     inletOutletFvPatchScalarField(ptf, iF),
00124     intensity_(ptf.intensity_),
00125     UName_(ptf.UName_),
00126     phiName_(ptf.phiName_)
00127 {}
00128 
00129 
00130 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00131 
00132 void Foam::turbulentIntensityKineticEnergyInletFvPatchScalarField::
00133 updateCoeffs()
00134 {
00135     if (updated())
00136     {
00137         return;
00138     }
00139 
00140     const fvPatchVectorField& Up =
00141         patch().lookupPatchField<volVectorField, vector>(UName_);
00142 
00143     const fvsPatchScalarField& phip =
00144         patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
00145 
00146     this->refValue() = 1.5*sqr(intensity_)*magSqr(Up);
00147     this->valueFraction() = 1.0 - pos(phip);
00148 
00149     inletOutletFvPatchScalarField::updateCoeffs();
00150 }
00151 
00152 
00153 void Foam::turbulentIntensityKineticEnergyInletFvPatchScalarField::write
00154 (
00155     Ostream& os
00156 ) const
00157 {
00158     fvPatchScalarField::write(os);
00159     os.writeKeyword("intensity") << intensity_ << token::END_STATEMENT << nl;
00160     os.writeKeyword("U") << UName_ << token::END_STATEMENT << nl;
00161     os.writeKeyword("phi") << phiName_ << token::END_STATEMENT << nl;
00162     writeEntry("value", os);
00163 }
00164 
00165 
00166 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00167 
00168 namespace Foam
00169 {
00170     makePatchTypeField
00171     (
00172         fvPatchScalarField,
00173         turbulentIntensityKineticEnergyInletFvPatchScalarField
00174     );
00175 }
00176 
00177 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines