FreeFOAM The Cross-Platform CFD Toolkit
PhiScheme.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 <finiteVolume/volFields.H>
28 #include <finiteVolume/fvcGrad.H>
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
39 template<class Type, class PhiLimiter>
40 tmp<surfaceScalarField> PhiScheme<Type, PhiLimiter>::limiter
41 (
43 ) const
44 {
45  const fvMesh& mesh = this->mesh();
46 
48  (
50  (
51  IOobject
52  (
53  "PhiLimiter",
54  mesh.time().timeName(),
55  mesh
56  ),
57  mesh,
58  dimless
59  )
60  );
61  surfaceScalarField& Limiter = tLimiter();
62 
63  const surfaceScalarField& CDweights = mesh.surfaceInterpolation::weights();
64 
65  const surfaceVectorField& Sf = mesh.Sf();
66  const surfaceScalarField& magSf = mesh.magSf();
67 
68  const unallocLabelList& owner = mesh.owner();
69  const unallocLabelList& neighbour = mesh.neighbour();
70 
71  tmp<surfaceScalarField> tUflux = this->faceFlux_;
72 
73  if (this->faceFlux_.dimensions() == dimDensity*dimVelocity*dimArea)
74  {
75  const volScalarField& rho =
76  phi.db().objectRegistry::lookupObject<volScalarField>("rho");
77  tUflux = this->faceFlux_/fvc::interpolate(rho);
78  }
79  else if (this->faceFlux_.dimensions() != dimVelocity*dimArea)
80  {
82  (
83  "PhiScheme<PhiLimiter>::limiter"
84  "(const GeometricField<Type, fvPatchField, volMesh>& phi)"
85  ) << "dimensions of faceFlux are not correct"
86  << exit(FatalError);
87  }
88 
89  const surfaceScalarField& Uflux = tUflux();
90 
91  scalarField& pLimiter = Limiter.internalField();
92 
93  forAll(pLimiter, face)
94  {
95  pLimiter[face] = PhiLimiter::limiter
96  (
97  CDweights[face],
98  Uflux[face],
99  phi[owner[face]],
100  phi[neighbour[face]],
101  Sf[face],
102  magSf[face]
103  );
104  }
105 
106 
107  surfaceScalarField::GeometricBoundaryField& bLimiter =
108  Limiter.boundaryField();
109 
110  forAll(bLimiter, patchI)
111  {
112  scalarField& pLimiter = bLimiter[patchI];
113 
114  if (bLimiter[patchI].coupled())
115  {
116  const scalarField& pCDweights = CDweights.boundaryField()[patchI];
117  const vectorField& pSf = Sf.boundaryField()[patchI];
118  const scalarField& pmagSf = magSf.boundaryField()[patchI];
119  const scalarField& pFaceFlux = Uflux.boundaryField()[patchI];
120  Field<Type> pphiP =
121  phi.boundaryField()[patchI].patchInternalField();
122  Field<Type> pphiN =
123  phi.boundaryField()[patchI].patchNeighbourField();
124 
125  forAll(pLimiter, face)
126  {
127  pLimiter[face] = PhiLimiter::limiter
128  (
129  pCDweights[face],
130  pFaceFlux[face],
131  pphiP[face],
132  pphiN[face],
133  pSf[face],
134  pmagSf[face]
135  );
136  }
137  }
138  else
139  {
140  pLimiter = 1.0;
141  }
142  }
143 
144  return tLimiter;
145 }
146 
147 
148 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
149 
150 } // End namespace Foam
151 
152 // ************************ vim: set sw=4 sts=4 et: ************************ //