FreeFOAM The Cross-Platform CFD Toolkit
gaussLaplacianSchemes.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) 1991-2010 OpenCFD Ltd.
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 "gaussLaplacianScheme.H"
27 #include <finiteVolume/fvMesh.H>
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace fv
34 {
35  makeFvLaplacianScheme(gaussLaplacianScheme)
36 }
37 }
38 
39 #define declareFvmLaplacianScalarGamma(Type) \
40  \
41 template<> \
42 Foam::tmp<Foam::fvMatrix<Foam::Type> > \
43 Foam::fv::gaussLaplacianScheme<Foam::Type, Foam::scalar>::fvmLaplacian \
44 ( \
45  const GeometricField<scalar, fvsPatchField, surfaceMesh>& gamma, \
46  GeometricField<Type, fvPatchField, volMesh>& vf \
47 ) \
48 { \
49  const fvMesh& mesh = this->mesh(); \
50  \
51  GeometricField<scalar, fvsPatchField, surfaceMesh> gammaMagSf \
52  ( \
53  gamma*mesh.magSf() \
54  ); \
55  \
56  tmp<fvMatrix<Type> > tfvm = fvmLaplacianUncorrected(gammaMagSf, vf); \
57  fvMatrix<Type>& fvm = tfvm(); \
58  \
59  if (this->tsnGradScheme_().corrected()) \
60  { \
61  if (mesh.fluxRequired(vf.name())) \
62  { \
63  fvm.faceFluxCorrectionPtr() = new \
64  GeometricField<Type, fvsPatchField, surfaceMesh> \
65  ( \
66  gammaMagSf*this->tsnGradScheme_().correction(vf) \
67  ); \
68  \
69  fvm.source() -= \
70  mesh.V()* \
71  fvc::div \
72  ( \
73  *fvm.faceFluxCorrectionPtr() \
74  )().internalField(); \
75  } \
76  else \
77  { \
78  fvm.source() -= \
79  mesh.V()* \
80  fvc::div \
81  ( \
82  gammaMagSf*this->tsnGradScheme_().correction(vf) \
83  )().internalField(); \
84  } \
85  } \
86  \
87  return tfvm; \
88 } \
89  \
90  \
91 template<> \
92 Foam::tmp<Foam::GeometricField<Foam::Type, Foam::fvPatchField, Foam::volMesh> >\
93 Foam::fv::gaussLaplacianScheme<Foam::Type, Foam::scalar>::fvcLaplacian \
94 ( \
95  const GeometricField<scalar, fvsPatchField, surfaceMesh>& gamma, \
96  const GeometricField<Type, fvPatchField, volMesh>& vf \
97 ) \
98 { \
99  const fvMesh& mesh = this->mesh(); \
100  \
101  tmp<GeometricField<Type, fvPatchField, volMesh> > tLaplacian \
102  ( \
103  fvc::div(gamma*this->tsnGradScheme_().snGrad(vf)*mesh.magSf()) \
104  ); \
105  \
106  tLaplacian().rename("laplacian(" + gamma.name() + ',' + vf.name() + ')');\
107  \
108  return tLaplacian; \
109 }
110 
111 
117 
118 
119 // ************************ vim: set sw=4 sts=4 et: ************************ //