multivariateScheme.C
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 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 \*---------------------------------------------------------------------------*/
25 
26 #include "volFields.H"
27 #include "surfaceFields.H"
28 #include "upwind.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 template<class Type, class Scheme>
34 (
35  const fvMesh& mesh,
37  fieldTable& fields,
38  const surfaceScalarField& faceFlux,
39  Istream& schemeData
40 )
41 :
43  (
44  mesh,
45  fields,
46  faceFlux,
47  schemeData
48  ),
49  Scheme::LimiterType(schemeData),
50  faceFlux_(faceFlux),
51  weights_
52  (
53  IOobject
54  (
55  "multivariateWeights",
56  mesh.time().timeName(),
57  mesh
58  ),
59  mesh,
60  dimless
61  )
62 {
64  fieldTable::const_iterator iter = this->fields().begin();
65 
67  (
68  Scheme(mesh, faceFlux_, *this).limiter(*iter())
69  );
70 
71  for (++iter; iter != this->fields().end(); ++iter)
72  {
73  limiter = min
74  (
75  limiter,
76  Scheme(mesh, faceFlux_, *this).limiter(*iter())
77  );
78  }
79 
80  weights_ =
81  limiter*mesh.surfaceInterpolation::weights()
82  + (scalar(1) - limiter)*upwind<Type>(mesh, faceFlux_).weights();
83 }
84 
85 
86 // ************************************************************************* //
Foam::surfaceFields.
tmp< surfaceScalarField > weights() const
Return the interpolation weighting factors.
Definition: upwind.H:131
void limiter(scalarField &allLambda, const RdeltaTType &rDeltaT, const RhoType &rho, const volScalarField &psi, const surfaceScalarField &phiBD, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su, const scalar psiMax, const scalar psiMin)
Generic multi-variate discretisation scheme class which may be instantiated for any of the NVD...
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:104
dynamicFvMesh & mesh
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
Upwind differencing scheme class.
Definition: upwind.H:51
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Abstract base class for multi-variate surface interpolation schemes.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92