gaussConvectionScheme.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 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 "gaussConvectionScheme.H"
27 #include "fvcSurfaceIntegrate.H"
28 #include "fvMatrices.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 namespace fv
38 {
39 
40 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41 
42 template<class Type>
43 const surfaceInterpolationScheme<Type>&
45 {
46  return tinterpScheme_();
47 }
48 
49 
50 template<class Type>
53 (
54  const surfaceScalarField&,
56 ) const
57 {
58  return tinterpScheme_().interpolate(vf);
59 }
60 
61 
62 template<class Type>
65 (
66  const surfaceScalarField& faceFlux,
68 ) const
69 {
70  return faceFlux*interpolate(faceFlux, vf);
71 }
72 
73 
74 template<class Type>
77 (
78  const surfaceScalarField& faceFlux,
80 ) const
81 {
82  tmp<surfaceScalarField> tweights = tinterpScheme_().weights(vf);
83  const surfaceScalarField& weights = tweights();
84 
85  tmp<fvMatrix<Type>> tfvm
86  (
87  new fvMatrix<Type>
88  (
89  vf,
90  faceFlux.dimensions()*vf.dimensions()
91  )
92  );
93  fvMatrix<Type>& fvm = tfvm.ref();
94 
95  fvm.lower() = -weights.primitiveField()*faceFlux.primitiveField();
96  fvm.upper() = fvm.lower() + faceFlux.primitiveField();
97  fvm.negSumDiag();
98 
100  {
101  const fvPatchField<Type>& psf = vf.boundaryField()[patchi];
102  const fvsPatchScalarField& patchFlux = faceFlux.boundaryField()[patchi];
103  const fvsPatchScalarField& pw = weights.boundaryField()[patchi];
104 
105  fvm.internalCoeffs()[patchi] = patchFlux*psf.valueInternalCoeffs(pw);
106  fvm.boundaryCoeffs()[patchi] = -patchFlux*psf.valueBoundaryCoeffs(pw);
107  }
108 
109  if (tinterpScheme_().corrected())
110  {
111  fvm += fvc::surfaceIntegrate(faceFlux*tinterpScheme_().correction(vf));
112  }
113 
114  return tfvm;
115 }
116 
117 
118 template<class Type>
121 (
122  const surfaceScalarField& faceFlux,
124 ) const
125 {
127  (
128  fvc::surfaceIntegrate(flux(faceFlux, vf))
129  );
130 
131  tConvection.ref().rename
132  (
133  "convection(" + faceFlux.name() + ',' + vf.name() + ')'
134  );
135 
136  return tConvection;
137 }
138 
139 
140 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
141 
142 } // End namespace fv
143 
144 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
145 
146 } // End namespace Foam
147 
148 // ************************************************************************* //
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
Surface integrate surfaceField creating a volField. Surface sum a surfaceField creating a volField...
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
const word & name() const
Return name.
Definition: IOobject.H:315
void surfaceIntegrate(Field< Type > &ivf, const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
FieldField< Field, Type > & internalCoeffs()
fvBoundary scalar field containing pseudo-matrix coeffs
Definition: fvMatrix.H:312
const Boundary & boundaryField() const
Return const-reference to the boundary field.
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
virtual tmp< Field< Type > > valueBoundaryCoeffs(const tmp< Field< scalar >> &) const
Return the matrix source coefficients corresponding to the.
Definition: fvPatchField.H:472
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
scalarField & upper()
Definition: lduMatrix.C:197
const surfaceInterpolationScheme< Type > & interpScheme() const
const dimensionSet & dimensions() const
Return dimensions.
labelList fv(nPoints)
FieldField< Field, Type > & boundaryCoeffs()
fvBoundary scalar field containing pseudo-matrix coeffs
Definition: fvMatrix.H:319
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDiv(const surfaceScalarField &, const GeometricField< Type, fvPatchField, volMesh > &) const
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
virtual tmp< Field< Type > > valueInternalCoeffs(const tmp< Field< scalar >> &) const
Return the matrix diagonal coefficients corresponding to the.
Definition: fvPatchField.H:461
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const surfaceScalarField &, const GeometricField< Type, fvPatchField, volMesh > &) const
label patchi
tmp< fvMatrix< Type > > fvmDiv(const surfaceScalarField &, const GeometricField< Type, fvPatchField, volMesh > &) const
scalarField & lower()
Definition: lduMatrix.C:168
A special matrix type and solver, designed for finite volume solutions of scalar equations.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux(const surfaceScalarField &, const GeometricField< Type, fvPatchField, volMesh > &) const
A class for managing temporary objects.
Definition: PtrList.H:53
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
Definition: fvcFlux.C:32
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
Namespace for OpenFOAM.