Home
Downloads
Documentation
Installation
User Guide
man-pages
API Documentation
README
Release Notes
Changes
License
Support
SourceForge Project
Main Page
Related Pages
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Pages
src
finiteVolume
interpolation
surfaceInterpolation
limitedSchemes
LimitedScheme
LimitedScheme.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
>
27
#include <
finiteVolume/surfaceFields.H
>
28
#include <
finiteVolume/fvcGrad.H
>
29
#include <
finiteVolume/coupledFvPatchFields.H
>
30
31
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32
33
namespace
Foam
34
{
35
36
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
37
38
template
<
class
Type,
class
Limiter,
template
<
class
>
class
LimitFunc>
39
tmp<surfaceScalarField>
LimitedScheme<Type, Limiter, LimitFunc>::limiter
40
(
41
const
GeometricField<Type, fvPatchField, volMesh>
&
phi
42
)
const
43
{
44
const
fvMesh
&
mesh
= this->
mesh
();
45
46
tmp<surfaceScalarField>
tLimiter
47
(
48
new
surfaceScalarField
49
(
50
IOobject
51
(
52
type
() +
"Limiter("
+ phi.
name
() +
')'
,
53
mesh.
time
().
timeName
(),
54
mesh
55
),
56
mesh,
57
dimless
58
)
59
);
60
surfaceScalarField
& lim = tLimiter();
61
62
tmp<GeometricField<typename Limiter::phiType, fvPatchField, volMesh>
>
63
tlPhi = LimitFunc<Type>()(phi);
64
65
const
GeometricField<typename Limiter::phiType, fvPatchField, volMesh>
&
66
lPhi = tlPhi();
67
68
GeometricField<typename Limiter::gradPhiType, fvPatchField, volMesh>
69
gradc(
fvc::grad
(lPhi));
70
71
const
surfaceScalarField
& CDweights = mesh.surfaceInterpolation::weights();
72
73
const
unallocLabelList
& owner = mesh.
owner
();
74
const
unallocLabelList
& neighbour = mesh.
neighbour
();
75
76
const
vectorField
&
C
= mesh.
C
();
77
78
scalarField
& pLim = lim.internalField();
79
80
forAll
(pLim,
face
)
81
{
82
label own = owner[
face
];
83
label nei = neighbour[
face
];
84
85
pLim[
face
] =
Limiter::limiter
86
(
87
CDweights[
face
],
88
this->faceFlux_[face],
89
lPhi[own],
90
lPhi[nei],
91
gradc[own],
92
gradc[nei],
93
C[nei] - C[own]
94
);
95
}
96
97
surfaceScalarField::GeometricBoundaryField& bLim = lim.boundaryField();
98
99
forAll
(bLim,
patchi
)
100
{
101
scalarField
& pLim = bLim[
patchi
];
102
103
if
(bLim[
patchi
].coupled())
104
{
105
const
scalarField
& pCDweights = CDweights.
boundaryField
()[
patchi
];
106
const
scalarField
& pFaceFlux =
107
this->faceFlux_.boundaryField()[
patchi
];
108
Field<typename Limiter::phiType>
plPhiP =
109
lPhi.
boundaryField
()[
patchi
].patchInternalField();
110
Field<typename Limiter::phiType>
plPhiN =
111
lPhi.
boundaryField
()[
patchi
].patchNeighbourField();
112
Field<typename Limiter::gradPhiType>
pGradcP =
113
gradc.
boundaryField
()[
patchi
].patchInternalField();
114
Field<typename Limiter::gradPhiType>
pGradcN =
115
gradc.
boundaryField
()[
patchi
].patchNeighbourField();
116
117
// Build the d-vectors
118
vectorField
pd =
119
mesh.
Sf
().
boundaryField
()[
patchi
]
120
/(
121
mesh.
magSf
().
boundaryField
()[
patchi
]
122
*mesh.
deltaCoeffs
().
boundaryField
()[
patchi
]
123
);
124
125
if
(!mesh.
orthogonal
())
126
{
127
pd -= mesh.
correctionVectors
().
boundaryField
()[
patchi
]
128
/mesh.
deltaCoeffs
().
boundaryField
()[
patchi
];
129
}
130
131
forAll
(pLim,
face
)
132
{
133
pLim[
face
] =
Limiter::limiter
134
(
135
pCDweights[
face
],
136
pFaceFlux[face],
137
plPhiP[face],
138
plPhiN[face],
139
pGradcP[face],
140
pGradcN[face],
141
pd[face]
142
);
143
}
144
}
145
else
146
{
147
pLim = 1.0;
148
}
149
}
150
151
return
tLimiter;
152
}
153
154
155
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
156
157
}
// End namespace Foam
158
159
// ************************ vim: set sw=4 sts=4 et: ************************ //