FreeFOAM The Cross-Platform CFD Toolkit
extendedLeastSquaresGrad.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 
28 #include <finiteVolume/gaussGrad.H>
29 #include <finiteVolume/fvMesh.H>
30 #include <finiteVolume/volMesh.H>
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 
40 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41 
42 namespace fv
43 {
44 
45 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
46 
47 template<class Type>
48 tmp
49 <
50  GeometricField
51  <
52  typename outerProduct<vector, Type>::type, fvPatchField, volMesh
53  >
54 >
56 (
57  const GeometricField<Type, fvPatchField, volMesh>& vsf
58 ) const
59 {
60  typedef typename outerProduct<vector, Type>::type GradType;
61 
62  const fvMesh& mesh = vsf.mesh();
63 
64  tmp<GeometricField<GradType, fvPatchField, volMesh> > tlsGrad
65  (
66  new GeometricField<GradType, fvPatchField, volMesh>
67  (
68  IOobject
69  (
70  "grad("+vsf.name()+')',
71  vsf.instance(),
72  mesh,
75  ),
76  mesh,
77  dimensioned<GradType>
78  (
79  "zero",
80  vsf.dimensions()/dimLength,
81  pTraits<GradType>::zero
82  ),
83  zeroGradientFvPatchField<GradType>::typeName
84  )
85  );
86  GeometricField<GradType, fvPatchField, volMesh>& lsGrad = tlsGrad();
87 
88  // Get reference to least square vectors
89  const extendedLeastSquaresVectors& lsv = extendedLeastSquaresVectors::New
90  (
91  mesh,
92  minDet_
93  );
94 
95  const surfaceVectorField& ownLs = lsv.pVectors();
96  const surfaceVectorField& neiLs = lsv.nVectors();
97 
98  const unallocLabelList& owner = mesh.owner();
99  const unallocLabelList& neighbour = mesh.neighbour();
100 
101  forAll(owner, facei)
102  {
103  label own = owner[facei];
104  label nei = neighbour[facei];
105 
106  Type deltaVsf = vsf[nei] - vsf[own];
107 
108  lsGrad[own] += ownLs[facei]*deltaVsf;
109  lsGrad[nei] -= neiLs[facei]*deltaVsf;
110  }
111 
112  // Boundary faces
113  forAll(vsf.boundaryField(), patchi)
114  {
115  const fvsPatchVectorField& patchOwnLs = ownLs.boundaryField()[patchi];
116 
117  const unallocLabelList& faceCells =
118  lsGrad.boundaryField()[patchi].patch().faceCells();
119 
120  if (vsf.boundaryField()[patchi].coupled())
121  {
122  Field<Type> neiVsf =
123  vsf.boundaryField()[patchi].patchNeighbourField();
124 
125  forAll(neiVsf, patchFaceI)
126  {
127  lsGrad[faceCells[patchFaceI]] +=
128  patchOwnLs[patchFaceI]
129  *(neiVsf[patchFaceI] - vsf[faceCells[patchFaceI]]);
130  }
131  }
132  else
133  {
134  const fvPatchField<Type>& patchVsf = vsf.boundaryField()[patchi];
135 
136  forAll(patchVsf, patchFaceI)
137  {
138  lsGrad[faceCells[patchFaceI]] +=
139  patchOwnLs[patchFaceI]
140  *(patchVsf[patchFaceI] - vsf[faceCells[patchFaceI]]);
141  }
142  }
143  }
144 
145 
146  const List<labelPair>& additionalCells = lsv.additionalCells();
147  const vectorField& additionalVectors = lsv.additionalVectors();
148 
149  forAll(additionalCells, i)
150  {
151  lsGrad[additionalCells[i][0]] +=
152  additionalVectors[i]
153  *(vsf[additionalCells[i][1]] - vsf[additionalCells[i][0]]);
154  }
155 
156 
157  lsGrad.correctBoundaryConditions();
159 
160  return tlsGrad;
161 }
162 
163 
164 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
165 
166 } // End namespace fv
167 
168 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
169 
170 } // End namespace Foam
171 
172 // ************************ vim: set sw=4 sts=4 et: ************************ //