FreeFOAM The Cross-Platform CFD Toolkit
gaussGrad.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 "gaussGrad.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 namespace fv
37 {
38 
39 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
40 
41 template<class Type>
42 tmp
43 <
44  GeometricField
45  <
46  typename outerProduct<vector, Type>::type, fvPatchField, volMesh
47  >
48 >
50 (
51  const GeometricField<Type, fvsPatchField, surfaceMesh>& ssf
52 )
53 {
54  typedef typename outerProduct<vector, Type>::type GradType;
55 
56  const fvMesh& mesh = ssf.mesh();
57 
58  tmp<GeometricField<GradType, fvPatchField, volMesh> > tgGrad
59  (
60  new GeometricField<GradType, fvPatchField, volMesh>
61  (
62  IOobject
63  (
64  "grad("+ssf.name()+')',
65  ssf.instance(),
66  mesh,
69  ),
70  mesh,
71  dimensioned<GradType>
72  (
73  "0",
74  ssf.dimensions()/dimLength,
75  pTraits<GradType>::zero
76  ),
77  zeroGradientFvPatchField<GradType>::typeName
78  )
79  );
80  GeometricField<GradType, fvPatchField, volMesh>& gGrad = tgGrad();
81 
82  const unallocLabelList& owner = mesh.owner();
83  const unallocLabelList& neighbour = mesh.neighbour();
84  const vectorField& Sf = mesh.Sf();
85 
86  Field<GradType>& igGrad = gGrad;
87  const Field<Type>& issf = ssf;
88 
89  forAll(owner, facei)
90  {
91  GradType Sfssf = Sf[facei]*issf[facei];
92 
93  igGrad[owner[facei]] += Sfssf;
94  igGrad[neighbour[facei]] -= Sfssf;
95  }
96 
97  forAll(mesh.boundary(), patchi)
98  {
99  const unallocLabelList& pFaceCells =
100  mesh.boundary()[patchi].faceCells();
101 
102  const vectorField& pSf = mesh.Sf().boundaryField()[patchi];
103 
104  const fvsPatchField<Type>& pssf = ssf.boundaryField()[patchi];
105 
106  forAll(mesh.boundary()[patchi], facei)
107  {
108  igGrad[pFaceCells[facei]] += pSf[facei]*pssf[facei];
109  }
110  }
111 
112  igGrad /= mesh.V();
113 
114  gGrad.correctBoundaryConditions();
115 
116  return tgGrad;
117 }
118 
119 
120 template<class Type>
121 tmp
122 <
123  GeometricField
124  <
125  typename outerProduct<vector, Type>::type, fvPatchField, volMesh
126  >
127 >
129 (
130  const GeometricField<Type, fvPatchField, volMesh>& vsf
131 ) const
132 {
133  typedef typename outerProduct<vector, Type>::type GradType;
134 
135  tmp<GeometricField<GradType, fvPatchField, volMesh> > tgGrad
136  (
137  grad(tinterpScheme_().interpolate(vsf))
138  );
139  GeometricField<GradType, fvPatchField, volMesh>& gGrad = tgGrad();
140 
141  gGrad.rename("grad(" + vsf.name() + ')');
142  correctBoundaryConditions(vsf, gGrad);
143 
144  return tgGrad;
145 }
146 
147 
148 template<class Type>
150 (
153  <
155  >& gGrad
156 )
157 {
158  forAll(vsf.boundaryField(), patchi)
159  {
160  if (!vsf.boundaryField()[patchi].coupled())
161  {
162  vectorField n =
163  vsf.mesh().Sf().boundaryField()[patchi]
164  /vsf.mesh().magSf().boundaryField()[patchi];
165 
166  gGrad.boundaryField()[patchi] += n *
167  (
168  vsf.boundaryField()[patchi].snGrad()
169  - (n & gGrad.boundaryField()[patchi])
170  );
171  }
172  }
173 }
174 
175 
176 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
177 
178 } // End namespace fv
179 
180 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
181 
182 } // End namespace Foam
183 
184 // ************************ vim: set sw=4 sts=4 et: ************************ //