OpenFOAM logo
Open Source CFD Toolkit

blended.H

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-2005 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 the
00013     Free Software Foundation; either version 2 of the License, or (at your
00014     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, write to the Free Software Foundation,
00023     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
00024 
00025 Class
00026     blended
00027 
00028 Description
00029     linear/upwind blended differencing scheme.
00030 
00031 SourceFiles
00032     blended.C
00033 
00034 \*---------------------------------------------------------------------------*/
00035 
00036 #ifndef blended_H
00037 #define blended_H
00038 
00039 #include "limitedSurfaceInterpolationScheme.H"
00040 
00041 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00042 
00043 namespace Foam
00044 {
00045 
00046 /*---------------------------------------------------------------------------*\
00047                            Class blended Declaration
00048 \*---------------------------------------------------------------------------*/
00049 
00050 template<class Type>
00051 class blended
00052 :
00053     public limitedSurfaceInterpolationScheme<Type>
00054 {
00055     // Private data
00056 
00057         const scalar blendingFactor_;
00058 
00059 
00060     // Private Member Functions
00061 
00062         //- Disallow default bitwise copy construct
00063         blended(const blended&);
00064 
00065         //- Disallow default bitwise assignment
00066         void operator=(const blended&);
00067 
00068 
00069 public:
00070 
00071     //- Runtime type information
00072     TypeName("blended");
00073 
00074 
00075     // Constructors
00076 
00077         //- Construct from mesh, faceFlux and blendingFactor
00078         blended
00079         (
00080             const fvMesh& mesh,
00081             const surfaceScalarField& faceFlux,
00082             const scalar blendingFactor
00083         )
00084         :
00085             limitedSurfaceInterpolationScheme<Type>(mesh, faceFlux),
00086             blendingFactor_(blendingFactor)
00087         {}
00088 
00089         //- Construct from mesh and Istream. 
00090         //  The name of the flux field is read from the Istream and looked-up
00091         //  from the mesh objectRegistry
00092         blended
00093         (
00094             const fvMesh& mesh,
00095             Istream& is
00096         )
00097         :
00098             limitedSurfaceInterpolationScheme<Type>(mesh, is),
00099             blendingFactor_(readScalar(is))
00100         {}
00101 
00102         //- Construct from mesh, faceFlux and Istream
00103         blended
00104         (
00105             const fvMesh& mesh,
00106             const surfaceScalarField& faceFlux,
00107             Istream& is
00108         )
00109         :
00110             limitedSurfaceInterpolationScheme<Type>(mesh, faceFlux),
00111             blendingFactor_(readScalar(is))
00112         {}
00113 
00114 
00115     // Member Functions
00116 
00117         //- Return the interpolation limiter
00118         virtual tmp<surfaceScalarField> limiter
00119         (
00120             const GeometricField<Type, fvPatchField, volMesh>&
00121         ) const
00122         {
00123             return tmp<surfaceScalarField>
00124             (
00125                 new surfaceScalarField
00126                 (
00127                     IOobject
00128                     (
00129                         "limiter",
00130                         this->mesh().time().timeName(),
00131                         this->mesh()
00132                     ),
00133                     this->mesh(),
00134                     dimensionedScalar
00135                     (
00136                         "limiter",
00137                         dimless,
00138                         1 - blendingFactor_
00139                     )
00140                 )
00141             );
00142         }
00143 
00144         //- Return the interpolation weighting factors
00145         virtual tmp<surfaceScalarField> weights
00146         (
00147             const GeometricField<Type, fvPatchField, volMesh>& vf
00148         ) const
00149         {
00150             return
00151                 blendingFactor_*this->mesh().surfaceInterpolation::weights()
00152               + (1 - blendingFactor_)*pos(this->faceFlux_);
00153         }
00154 };
00155 
00156 
00157 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00158 
00159 } // End namespace Foam
00160 
00161 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00162 
00163 #endif
00164 
00165 // ************************************************************************* //
For further information go to www.openfoam.org