FreeFOAM The Cross-Platform CFD Toolkit
resErrorLaplacian.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 
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace resError
36 {
37 
38 template<class Type>
39 tmp<errorEstimate<Type> >
41 (
43 )
44 {
45  surfaceScalarField Gamma
46  (
47  IOobject
48  (
49  "gamma",
50  vf.time().constant(),
51  vf.db(),
53  ),
54  vf.mesh(),
55  dimensionedScalar("1", dimless, 1.0)
56  );
57 
58  return resError::laplacian(Gamma, vf);
59 }
60 
61 
62 template<class Type>
65 (
66  const dimensionedScalar& gamma,
68 )
69 {
70  surfaceScalarField Gamma
71  (
72  IOobject
73  (
74  gamma.name(),
75  vf.time().timeName(),
76  vf.db(),
78  ),
79  vf.mesh(),
80  gamma
81  );
82 
83  return resError::laplacian(Gamma, vf);
84 }
85 
86 
87 template<class Type>
90 (
91  const volScalarField& gamma,
93 )
94 {
95  return resError::laplacian(fvc::interpolate(gamma), vf);
96 }
97 
98 
99 template<class Type>
101 laplacian
102 (
103  const tmp<volScalarField>& tgamma,
105 )
106 {
107  tmp<errorEstimate<Type> > Laplacian(resError::laplacian(tgamma(), vf));
108  tgamma.clear();
109  return Laplacian;
110 }
111 
112 
113 template<class Type>
115 laplacian
116 (
117  const surfaceScalarField& gamma,
119 )
120 {
121  const fvMesh& mesh = vf.mesh();
122 
123  const scalarField& vols = mesh.V();
124  const surfaceVectorField& Sf = mesh.Sf();
125  const surfaceScalarField magSf = mesh.magSf();
126  const fvPatchList& patches = mesh.boundary();
127  const unallocLabelList& owner = mesh.owner();
128  const unallocLabelList& neighbour = mesh.neighbour();
129 
130  const surfaceScalarField& delta =
131  mesh.surfaceInterpolation::deltaCoeffs();
132 
133  Field<Type> res(vols.size(), pTraits<Type>::zero);
134  scalarField aNorm(vols.size(), 0.0);
135 
136  // Calculate gradient of the solution
138  <
140  > gradVf = fvc::grad(vf);
141 
142  // Internal faces
143  forAll (owner, faceI)
144  {
145  // Owner
146 
147  // Subtract diffusion
148  res[owner[faceI]] -=
149  gamma[faceI]*(Sf[faceI] & gradVf[owner[faceI]]);
150 
151  aNorm[owner[faceI]] += delta[faceI]*gamma[faceI]*magSf[faceI];
152 
153  // Neighbour
154 
155  // Subtract diffusion
156  res[neighbour[faceI]] +=
157  gamma[faceI]*(Sf[faceI] & gradVf[neighbour[faceI]]);
158 
159  aNorm[neighbour[faceI]] += delta[faceI]*gamma[faceI]*magSf[faceI];
160 
161  }
162 
163  forAll (patches, patchI)
164  {
165  const vectorField& patchSf = Sf.boundaryField()[patchI];
166  const scalarField& patchMagSf = magSf.boundaryField()[patchI];
167  const scalarField& patchGamma = gamma.boundaryField()[patchI];
168  const scalarField& patchDelta = delta.boundaryField()[patchI];
169 
170  const labelList& fCells = patches[patchI].faceCells();
171 
172  forAll (fCells, faceI)
173  {
174  // Subtract diffusion
175  res[fCells[faceI]] -=
176  patchGamma[faceI]*
177  (
178  patchSf[faceI] & gradVf[fCells[faceI]]
179  );
180 
181  aNorm[fCells[faceI]] +=
182  patchDelta[faceI]*patchGamma[faceI]*patchMagSf[faceI];
183  }
184  }
185 
186  res /= vols;
187  aNorm /= vols;
188 
189  return tmp<errorEstimate<Type> >
190  (
191  new errorEstimate<Type>
192  (
193  vf,
194  delta.dimensions()*gamma.dimensions()*magSf.dimensions()
195  *vf.dimensions(),
196  res,
197  aNorm
198  )
199  );
200 }
201 
202 template<class Type>
203 tmp<errorEstimate<Type> >
204 laplacian
205 (
206  const tmp<surfaceScalarField>& tgamma,
207  const GeometricField<Type, fvPatchField, volMesh>& vf
208 )
209 {
210  tmp<errorEstimate<Type> > tresError(resError::laplacian(tgamma(), vf));
211  tgamma.clear();
212  return tresError;
213 }
214 
215 
216 template<class Type>
217 tmp<errorEstimate<Type> >
218 laplacian
219 (
220  const volTensorField& gamma,
221  const GeometricField<Type, fvPatchField, volMesh>& vf
222 )
223 {
224  const fvMesh& mesh = vf.mesh();
225 
226  return resError::laplacian
227  (
228  (mesh.Sf() & fvc::interpolate(gamma) & mesh.Sf())
229  /sqr(mesh.magSf()),
230  vf
231  );
232 }
233 
234 template<class Type>
235 tmp<errorEstimate<Type> >
236 laplacian
237 (
238  const tmp<volTensorField>& tgamma,
239  const GeometricField<Type, fvPatchField, volMesh>& vf
240 )
241 {
242  tmp<errorEstimate<Type> > Laplacian = resError::laplacian(tgamma(), vf);
243  tgamma.clear();
244  return Laplacian;
245 }
246 
247 
248 template<class Type>
249 tmp<errorEstimate<Type> >
250 laplacian
251 (
252  const surfaceTensorField& gamma,
253  const GeometricField<Type, fvPatchField, volMesh>& vf
254 )
255 {
256  const fvMesh& mesh = vf.mesh();
257 
258  return resError::laplacian
259  (
260  (mesh.Sf() & gamma & mesh.Sf())/sqr(mesh.magSf()),
261  vf
262  );
263 }
264 
265 template<class Type>
266 tmp<errorEstimate<Type> >
267 laplacian
268 (
269  const tmp<surfaceTensorField>& tgamma,
270  const GeometricField<Type, fvPatchField, volMesh>& vf
271 )
272 {
273  tmp<errorEstimate<Type> > Laplacian = resError::laplacian(tgamma(), vf);
274  tgamma.clear();
275  return Laplacian;
276 }
277 
278 
279 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
280 
281 } // End namespace resError
282 
283 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
284 
285 } // End namespace Foam
286 
287 // ************************ vim: set sw=4 sts=4 et: ************************ //