blended.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Class
25  Foam::blended
26 
27 Description
28  linear/upwind blended differencing scheme.
29 
30 SourceFiles
31  blended.C
32 
33 \*---------------------------------------------------------------------------*/
34 
35 #ifndef blended_H
36 #define blended_H
37 
39 #include "blendedSchemeBase.H"
40 
41 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42 
43 namespace Foam
44 {
45 
46 /*---------------------------------------------------------------------------*\
47  Class blended Declaration
48 \*---------------------------------------------------------------------------*/
49 
50 template<class Type>
51 class blended
52 :
54  public blendedSchemeBase<Type>
55 {
56  // Private data
57 
58  const scalar blendingFactor_;
59 
60 
61  // Private Member Functions
62 
63  //- Disallow default bitwise copy construct
64  blended(const blended&);
65 
66  //- Disallow default bitwise assignment
67  void operator=(const blended&);
68 
69 
70 public:
71 
72  //- Runtime type information
73  TypeName("blended");
74 
75 
76  // Constructors
77 
78  //- Construct from mesh, faceFlux and blendingFactor
80  (
81  const fvMesh& mesh,
82  const surfaceScalarField& faceFlux,
83  const scalar blendingFactor
84  )
85  :
86  limitedSurfaceInterpolationScheme<Type>(mesh, faceFlux),
87  blendingFactor_(blendingFactor)
88  {}
89 
90  //- Construct from mesh and Istream.
91  // The name of the flux field is read from the Istream and looked-up
92  // from the mesh objectRegistry
94  (
95  const fvMesh& mesh,
96  Istream& is
97  )
98  :
99  limitedSurfaceInterpolationScheme<Type>(mesh, is),
100  blendingFactor_(readScalar(is))
101  {}
102 
103  //- Construct from mesh, faceFlux and Istream
105  (
106  const fvMesh& mesh,
107  const surfaceScalarField& faceFlux,
108  Istream& is
109  )
110  :
111  limitedSurfaceInterpolationScheme<Type>(mesh, faceFlux),
112  blendingFactor_(readScalar(is))
113  {}
114 
115 
116  //- Destructor
117  virtual ~blended()
118  {}
119 
120 
121  // Member Functions
122 
123  //- Return the face-based blending factor
125  (
127  ) const
128  {
130  (
132  (
133  IOobject
134  (
135  vf.name() + "BlendingFactor",
136  this->mesh().time().timeName(),
137  this->mesh(),
140  false
141  ),
142  this->mesh(),
144  (
145  "blendingFactor",
146  dimless,
147  blendingFactor_
148  )
149  )
150  );
151  }
152 
153  //- Return the interpolation limiter
155  (
157  ) const
158  {
160  (
162  (
163  IOobject
164  (
165  "blendedLimiter",
166  this->mesh().time().timeName(),
167  this->mesh(),
170  false
171  ),
172  this->mesh(),
174  (
175  "blendedLimiter",
176  dimless,
177  1 - blendingFactor_
178  )
179  )
180  );
181  }
182 
183  //- Return the interpolation weighting factors
185  (
187  ) const
188  {
189  return
190  blendingFactor_*this->mesh().surfaceInterpolation::weights()
191  + (1 - blendingFactor_)*pos0(this->faceFlux_);
192  }
193 };
194 
195 
196 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
197 
198 } // End namespace Foam
199 
200 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
201 
202 #endif
203 
204 // ************************************************************************* //
virtual ~blended()
Destructor.
Definition: blended.H:116
const word & name() const
Return name.
Definition: IOobject.H:291
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
virtual tmp< surfaceScalarField > limiter(const GeometricField< Type, fvPatchField, volMesh > &) const
Return the interpolation limiter.
Definition: blended.H:154
TypeName("blended")
Runtime type information.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:644
const fvMesh & mesh() const
Return mesh reference.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
virtual tmp< surfaceScalarField > blendingFactor(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the face-based blending factor.
Definition: blended.H:124
Base class for blended schemes to provide access to the blending factor surface field.
word timeName
Definition: getTimeIndex.H:3
linear/upwind blended differencing scheme.
Definition: blended.H:50
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
dimensionedScalar pos0(const dimensionedScalar &ds)
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Abstract base class for limited surface interpolation schemes.
A class for managing temporary objects.
Definition: PtrList.H:53
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
virtual tmp< surfaceScalarField > weights(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the interpolation weighting factors.
Definition: blended.H:184
Namespace for OpenFOAM.