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

smoluchowskiJumpTFvPatchScalarField.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 "smoluchowskiJumpTFvPatchScalarField.H"
00027 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00028 #include <finiteVolume/fvPatchFieldMapper.H>
00029 #include <finiteVolume/volFields.H>
00030 #include <OpenFOAM/mathematicalConstants.H>
00031 
00032 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00033 
00034 namespace Foam
00035 {
00036 
00037 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00038 
00039 smoluchowskiJumpTFvPatchScalarField::smoluchowskiJumpTFvPatchScalarField
00040 (
00041     const fvPatch& p,
00042     const DimensionedField<scalar, volMesh>& iF
00043 )
00044 :
00045     mixedFvPatchScalarField(p, iF),
00046     accommodationCoeff_(1.0),
00047     Twall_(p.size(), 0.0),
00048     gamma_(1.4)
00049 {
00050     refValue() = 0.0;
00051     refGrad() = 0.0;
00052     valueFraction() = 0.0;
00053 }
00054 
00055 smoluchowskiJumpTFvPatchScalarField::smoluchowskiJumpTFvPatchScalarField
00056 (
00057     const smoluchowskiJumpTFvPatchScalarField& ptf,
00058     const fvPatch& p,
00059     const DimensionedField<scalar, volMesh>& iF,
00060     const fvPatchFieldMapper& mapper
00061 )
00062 :
00063     mixedFvPatchScalarField(ptf, p, iF, mapper),
00064     accommodationCoeff_(ptf.accommodationCoeff_),
00065     Twall_(ptf.Twall_),
00066     gamma_(ptf.gamma_)
00067 {}
00068 
00069 
00070 smoluchowskiJumpTFvPatchScalarField::smoluchowskiJumpTFvPatchScalarField
00071 (
00072     const fvPatch& p,
00073     const DimensionedField<scalar, volMesh>& iF,
00074     const dictionary& dict
00075 )
00076 :
00077     mixedFvPatchScalarField(p, iF),
00078     accommodationCoeff_(readScalar(dict.lookup("accommodationCoeff"))),
00079     Twall_("Twall", dict, p.size()),
00080     gamma_(dict.lookupOrDefault<scalar>("gamma", 1.4))
00081 {
00082     if
00083     (
00084         mag(accommodationCoeff_) < SMALL
00085      || mag(accommodationCoeff_) > 2.0
00086     )
00087     {
00088         FatalIOErrorIn
00089         (
00090             "smoluchowskiJumpTFvPatchScalarField::"
00091             "smoluchowskiJumpTFvPatchScalarField"
00092             "("
00093             "    const fvPatch&,"
00094             "    const DimensionedField<scalar, volMesh>&,"
00095             "    const dictionary&"
00096             ")",
00097             dict
00098         )   << "unphysical accommodationCoeff specified"
00099             << "(0 < accommodationCoeff <= 1)" << endl
00100             << exit(FatalError);
00101     }
00102 
00103     if (dict.found("value"))
00104     {
00105         fvPatchField<scalar>::operator=
00106         (
00107             scalarField("value", dict, p.size())
00108         );
00109     }
00110     else
00111     {
00112         fvPatchField<scalar>::operator=(patchInternalField());
00113     }
00114 
00115     refValue() = *this;
00116     refGrad() = 0.0;
00117     valueFraction() = 0.0;
00118 }
00119 
00120 
00121 smoluchowskiJumpTFvPatchScalarField::smoluchowskiJumpTFvPatchScalarField
00122 (
00123     const smoluchowskiJumpTFvPatchScalarField& ptpsf,
00124     const DimensionedField<scalar, volMesh>& iF
00125 )
00126 :
00127     mixedFvPatchScalarField(ptpsf, iF),
00128     accommodationCoeff_(ptpsf.accommodationCoeff_),
00129     Twall_(ptpsf.Twall_),
00130     gamma_(ptpsf.gamma_)
00131 {}
00132 
00133 
00134 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00135 
00136 // Map from self
00137 void smoluchowskiJumpTFvPatchScalarField::autoMap
00138 (
00139     const fvPatchFieldMapper& m
00140 )
00141 {
00142     mixedFvPatchScalarField::autoMap(m);
00143 }
00144 
00145 
00146 // Reverse-map the given fvPatchField onto this fvPatchField
00147 void smoluchowskiJumpTFvPatchScalarField::rmap
00148 (
00149     const fvPatchField<scalar>& ptf,
00150     const labelList& addr
00151 )
00152 {
00153     mixedFvPatchField<scalar>::rmap(ptf, addr);
00154 }
00155 
00156 
00157 // Update the coefficients associated with the patch field
00158 void smoluchowskiJumpTFvPatchScalarField::updateCoeffs()
00159 {
00160     if (updated())
00161     {
00162         return;
00163     }
00164 
00165     const fvPatchScalarField& pmu =
00166         patch().lookupPatchField<volScalarField, scalar>("mu");
00167     const fvPatchScalarField& prho =
00168         patch().lookupPatchField<volScalarField, scalar>("rho");
00169     const fvPatchField<scalar>& ppsi =
00170         patch().lookupPatchField<volScalarField, scalar>("psi");
00171     const fvPatchVectorField& pU =
00172         patch().lookupPatchField<volVectorField, vector>("U");
00173 
00174     // Prandtl number reading consistent with rhoCentralFoam
00175     const dictionary& thermophysicalProperties =
00176         db().lookupObject<IOdictionary>("thermophysicalProperties");
00177     dimensionedScalar Pr = dimensionedScalar("Pr", dimless, 1.0);
00178     if (thermophysicalProperties.found("Pr"))
00179     {
00180         Pr = thermophysicalProperties.lookup("Pr");
00181     }
00182 
00183     Field<scalar> C2 = pmu/prho
00184         *sqrt(ppsi*mathematicalConstant::pi/2.0)
00185         *2.0*gamma_/Pr.value()/(gamma_ + 1.0)
00186         *(2.0 - accommodationCoeff_)/accommodationCoeff_;
00187 
00188     Field<scalar> aCoeff = prho.snGrad() - prho/C2;
00189     Field<scalar> KEbyRho = 0.5*magSqr(pU);
00190 
00191     valueFraction() = (1.0/(1.0 + patch().deltaCoeffs()*C2));
00192     refValue() = Twall_;
00193     refGrad() = 0.0;
00194 
00195     mixedFvPatchScalarField::updateCoeffs();
00196 }
00197 
00198 
00199 // Write
00200 void smoluchowskiJumpTFvPatchScalarField::write(Ostream& os) const
00201 {
00202     fvPatchScalarField::write(os);
00203     os.writeKeyword("accommodationCoeff")
00204         << accommodationCoeff_ << token::END_STATEMENT << nl;
00205     Twall_.writeEntry("Twall", os);
00206     os.writeKeyword("gamma")
00207         << gamma_ << token::END_STATEMENT << nl;
00208     writeEntry("value", os);
00209 }
00210 
00211 
00212 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00213 
00214 makePatchTypeField(fvPatchScalarField, smoluchowskiJumpTFvPatchScalarField);
00215 
00216 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00217 
00218 } // End namespace Foam
00219 
00220 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines