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

alphaContactAngleFvPatchScalarField.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 "alphaContactAngleFvPatchScalarField.H"
00027 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00028 #include <finiteVolume/fvPatchFieldMapper.H>
00029 #include <finiteVolume/volMesh.H>
00030 
00031 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00032 
00033 namespace Foam
00034 {
00035     defineTypeNameAndDebug(alphaContactAngleFvPatchScalarField, 0);
00036 }
00037 
00038 template<>
00039 const char* Foam::NamedEnum
00040 <
00041     Foam::alphaContactAngleFvPatchScalarField::limitControls,
00042     4
00043 >::names[] =
00044 {
00045     "none",
00046     "gradient",
00047     "zeroGradient",
00048     "alpha"
00049 };
00050 
00051 const Foam::NamedEnum
00052 <
00053     Foam::alphaContactAngleFvPatchScalarField::limitControls,
00054     4
00055 > Foam::alphaContactAngleFvPatchScalarField::limitControlNames_;
00056 
00057 
00058 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00059 
00060 Foam::alphaContactAngleFvPatchScalarField::alphaContactAngleFvPatchScalarField
00061 (
00062     const fvPatch& p,
00063     const DimensionedField<scalar, volMesh>& iF
00064 )
00065 :
00066     fixedGradientFvPatchScalarField(p, iF),
00067     limit_(lcZeroGradient)
00068 {}
00069 
00070 
00071 Foam::alphaContactAngleFvPatchScalarField::alphaContactAngleFvPatchScalarField
00072 (
00073     const alphaContactAngleFvPatchScalarField& acpsf,
00074     const fvPatch& p,
00075     const DimensionedField<scalar, volMesh>& iF,
00076     const fvPatchFieldMapper& mapper
00077 )
00078 :
00079     fixedGradientFvPatchScalarField(acpsf, p, iF, mapper),
00080     limit_(acpsf.limit_)
00081 {}
00082 
00083 
00084 Foam::alphaContactAngleFvPatchScalarField::alphaContactAngleFvPatchScalarField
00085 (
00086     const fvPatch& p,
00087     const DimensionedField<scalar, volMesh>& iF,
00088     const dictionary& dict
00089 )
00090 :
00091     fixedGradientFvPatchScalarField(p, iF),
00092     limit_(limitControlNames_.read(dict.lookup("limit")))
00093 {
00094     if (dict.found("gradient"))
00095     {
00096         gradient() = scalarField("gradient", dict, p.size());
00097         fixedGradientFvPatchScalarField::updateCoeffs();
00098         fixedGradientFvPatchScalarField::evaluate();
00099     }
00100     else
00101     {
00102         fvPatchField<scalar>::operator=(patchInternalField());
00103         gradient() = 0.0;
00104     }
00105 }
00106 
00107 
00108 Foam::alphaContactAngleFvPatchScalarField::alphaContactAngleFvPatchScalarField
00109 (
00110     const alphaContactAngleFvPatchScalarField& acpsf
00111 )
00112 :
00113     fixedGradientFvPatchScalarField(acpsf),
00114     limit_(acpsf.limit_)
00115 {}
00116 
00117 
00118 Foam::alphaContactAngleFvPatchScalarField::alphaContactAngleFvPatchScalarField
00119 (
00120     const alphaContactAngleFvPatchScalarField& acpsf,
00121     const DimensionedField<scalar, volMesh>& iF
00122 )
00123 :
00124     fixedGradientFvPatchScalarField(acpsf, iF),
00125     limit_(acpsf.limit_)
00126 {}
00127 
00128 
00129 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00130 
00131 void Foam::alphaContactAngleFvPatchScalarField::evaluate
00132 (
00133     const Pstream::commsTypes
00134 )
00135 {
00136     if (limit_ == lcGradient)
00137     {
00138         gradient() =
00139         patch().deltaCoeffs()
00140        *(
00141            max(min
00142            (
00143                *this + gradient()/patch().deltaCoeffs(),
00144                scalar(1)), scalar(0)
00145            ) - *this
00146        );
00147     }
00148     else if (limit_ == lcZeroGradient)
00149     {
00150         gradient() = 0.0;
00151     }
00152 
00153     fixedGradientFvPatchScalarField::evaluate();
00154 
00155     if (limit_ == lcAlpha)
00156     {
00157         scalarField::operator=(max(min(*this, scalar(1)), scalar(0)));
00158     }
00159 }
00160 
00161 
00162 void Foam::alphaContactAngleFvPatchScalarField::write
00163 (
00164     Ostream& os
00165 ) const
00166 {
00167     fixedGradientFvPatchScalarField::write(os);
00168     os.writeKeyword("limit")
00169         << limitControlNames_[limit_] << token::END_STATEMENT << nl;
00170 }
00171 
00172 
00173 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines