FreeFOAM The Cross-Platform CFD Toolkit
surfaceInterpolationScheme.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 Description
25  Abstract base class for surface interpolation schemes.
26 
27 \*---------------------------------------------------------------------------*/
28 
30 #include <finiteVolume/volFields.H>
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 
39 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
40 
41 // Return weighting factors for scheme given by name in dictionary
42 template<class Type>
43 tmp<surfaceInterpolationScheme<Type> > surfaceInterpolationScheme<Type>::New
44 (
45  const fvMesh& mesh,
46  Istream& schemeData
47 )
48 {
49  if (schemeData.eof())
50  {
52  (
53  "surfaceInterpolationScheme<Type>::New(const fvMesh&, Istream&)",
54  schemeData
55  ) << "Discretisation scheme not specified"
56  << endl << endl
57  << "Valid schemes are :" << endl
58  << MeshConstructorTablePtr_->sortedToc()
59  << exit(FatalIOError);
60  }
61 
62  word schemeName(schemeData);
63 
65  {
66  Info<< "surfaceInterpolationScheme<Type>::New"
67  "(const fvMesh&, Istream&)"
68  " : discretisation scheme = "
69  << schemeName
70  << endl;
71  }
72 
73  typename MeshConstructorTable::iterator constructorIter =
74  MeshConstructorTablePtr_->find(schemeName);
75 
76  if (constructorIter == MeshConstructorTablePtr_->end())
77  {
79  (
80  "surfaceInterpolationScheme<Type>::New(const fvMesh&, Istream&)",
81  schemeData
82  ) << "Unknown discretisation scheme " << schemeName
83  << endl << endl
84  << "Valid schemes are :" << endl
85  << MeshConstructorTablePtr_->sortedToc()
86  << exit(FatalIOError);
87  }
88 
89  return constructorIter()(mesh, schemeData);
90 }
91 
92 
93 // Return weighting factors for scheme given by name in dictionary
94 template<class Type>
96 (
97  const fvMesh& mesh,
98  const surfaceScalarField& faceFlux,
99  Istream& schemeData
100 )
101 {
102  if (schemeData.eof())
103  {
105  (
106  "surfaceInterpolationScheme<Type>::New"
107  "(const fvMesh&, const surfaceScalarField&, Istream&)",
108  schemeData
109  ) << "Discretisation scheme not specified"
110  << endl << endl
111  << "Valid schemes are :" << endl
112  << MeshConstructorTablePtr_->sortedToc()
113  << exit(FatalIOError);
114  }
115 
116  word schemeName(schemeData);
117 
119  {
120  Info<< "surfaceInterpolationScheme<Type>::New"
121  "(const fvMesh&, const surfaceScalarField&, Istream&)"
122  " : discretisation scheme = "
123  << schemeName
124  << endl;
125  }
126 
127  typename MeshFluxConstructorTable::iterator constructorIter =
128  MeshFluxConstructorTablePtr_->find(schemeName);
129 
130  if (constructorIter == MeshFluxConstructorTablePtr_->end())
131  {
133  (
134  "surfaceInterpolationScheme<Type>::New"
135  "(const fvMesh&, const surfaceScalarField&, Istream&)",
136  schemeData
137  ) << "Unknown discretisation scheme " << schemeName
138  << endl << endl
139  << "Valid schemes are :" << endl
140  << MeshFluxConstructorTablePtr_->sortedToc()
141  << exit(FatalIOError);
142  }
143 
144  return constructorIter()(mesh, faceFlux, schemeData);
145 }
146 
147 
148 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
149 
150 template<class Type>
152 {}
153 
154 
155 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
156 
157 //- Return the face-interpolate of the given cell field
158 // with the given owner and neighbour weighting factors
159 template<class Type>
162 (
164  const tmp<surfaceScalarField>& tlambdas,
165  const tmp<surfaceScalarField>& tys
166 )
167 {
169  {
170  Info<< "surfaceInterpolationScheme<Type>::uncorrectedInterpolate"
171  "(const GeometricField<Type, fvPatchField, volMesh>&, "
172  "const tmp<surfaceScalarField>&, "
173  "const tmp<surfaceScalarField>&) : "
174  "interpolating volTypeField from cells to faces "
175  "without explicit correction"
176  << endl;
177  }
178 
179  const surfaceScalarField& lambdas = tlambdas();
180  const surfaceScalarField& ys = tys();
181 
182  const Field<Type>& vfi = vf.internalField();
183  const scalarField& lambda = lambdas.internalField();
184  const scalarField& y = ys.internalField();
185 
186  const fvMesh& mesh = vf.mesh();
187  const unallocLabelList& P = mesh.owner();
188  const unallocLabelList& N = mesh.neighbour();
189 
191  (
193  (
194  IOobject
195  (
196  "interpolate("+vf.name()+')',
197  vf.instance(),
198  vf.db()
199  ),
200  mesh,
201  vf.dimensions()
202  )
203  );
205 
206  Field<Type>& sfi = sf.internalField();
207 
208  for (register label fi=0; fi<P.size(); fi++)
209  {
210  sfi[fi] = lambda[fi]*vfi[P[fi]] + y[fi]*vfi[N[fi]];
211  }
212 
213 
214  // Interpolate across coupled patches using given lambdas and ys
215 
216  forAll (lambdas.boundaryField(), pi)
217  {
218  const fvsPatchScalarField& pLambda = lambdas.boundaryField()[pi];
219  const fvsPatchScalarField& pY = ys.boundaryField()[pi];
220 
221  if (vf.boundaryField()[pi].coupled())
222  {
223  sf.boundaryField()[pi] =
224  pLambda*vf.boundaryField()[pi].patchInternalField()
225  + pY*vf.boundaryField()[pi].patchNeighbourField();
226  }
227  else
228  {
229  sf.boundaryField()[pi] = vf.boundaryField()[pi];
230  }
231  }
232 
233  tlambdas.clear();
234  tys.clear();
235 
236  return tsf;
237 }
238 
239 
240 //- Return the face-interpolate of the given cell field
241 // with the given weigting factors
242 template<class Type>
245 (
247  const tmp<surfaceScalarField>& tlambdas
248 )
249 {
251  {
252  Info<< "surfaceInterpolationScheme<Type>::interpolate"
253  "(const GeometricField<Type, fvPatchField, volMesh>&, "
254  "const tmp<surfaceScalarField>&) : "
255  "interpolating volTypeField from cells to faces "
256  "without explicit correction"
257  << endl;
258  }
259 
260  const surfaceScalarField& lambdas = tlambdas();
261 
262  const Field<Type>& vfi = vf.internalField();
263  const scalarField& lambda = lambdas.internalField();
264 
265  const fvMesh& mesh = vf.mesh();
266  const unallocLabelList& P = mesh.owner();
267  const unallocLabelList& N = mesh.neighbour();
268 
270  (
272  (
273  IOobject
274  (
275  "interpolate("+vf.name()+')',
276  vf.instance(),
277  vf.db()
278  ),
279  mesh,
280  vf.dimensions()
281  )
282  );
284 
285  Field<Type>& sfi = sf.internalField();
286 
287  for (register label fi=0; fi<P.size(); fi++)
288  {
289  sfi[fi] = lambda[fi]*(vfi[P[fi]] - vfi[N[fi]]) + vfi[N[fi]];
290  }
291 
292  // Interpolate across coupled patches using given lambdas
293 
294  forAll (lambdas.boundaryField(), pi)
295  {
296  const fvsPatchScalarField& pLambda = lambdas.boundaryField()[pi];
297 
298  if (vf.boundaryField()[pi].coupled())
299  {
300  tsf().boundaryField()[pi] =
301  pLambda*vf.boundaryField()[pi].patchInternalField()
302  + (1.0 - pLambda)*vf.boundaryField()[pi].patchNeighbourField();
303  }
304  else
305  {
306  sf.boundaryField()[pi] = vf.boundaryField()[pi];
307  }
308  }
309 
310  tlambdas.clear();
311 
312  return tsf;
313 }
314 
315 
316 //- Return the face-interpolate of the given cell field
317 // with explicit correction
318 template<class Type>
321 (
323 ) const
324 {
326  {
327  Info<< "surfaceInterpolationScheme<Type>::interpolate"
328  "(const GeometricField<Type, fvPatchField, volMesh>&) : "
329  << "interpolating volTypeField from cells to faces"
330  << endl;
331  }
332 
334  = interpolate(vf, weights(vf));
335 
336  if (corrected())
337  {
338  tsf() += correction(vf);
339  }
340 
341  return tsf;
342 }
343 
344 
345 //- Return the face-interpolate of the given cell field
346 // with explicit correction
347 template<class Type>
350 (
352 ) const
353 {
355  = interpolate(tvf());
356  tvf.clear();
357  return tinterpVf;
358 }
359 
360 
361 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
362 
363 } // End namespace Foam
364 
365 // ************************ vim: set sw=4 sts=4 et: ************************ //