MRFZoneTemplates.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2018 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "MRFZone.H"
29#include "fvMesh.H"
30#include "volFields.H"
31#include "surfaceFields.H"
32#include "fvMatrices.H"
33
34// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
35
36template<class RhoFieldType>
37void Foam::MRFZone::makeRelativeRhoFlux
38(
39 const RhoFieldType& rho,
40 surfaceScalarField& phi
41) const
42{
43 if (!active_)
44 {
45 return;
46 }
47
48 const surfaceVectorField& Cf = mesh_.Cf();
49 const surfaceVectorField& Sf = mesh_.Sf();
50
51 const vector Omega = omega_->value(mesh_.time().timeOutputValue())*axis_;
52
53 const vectorField& Cfi = Cf;
54 const vectorField& Sfi = Sf;
56
57 // Internal faces
58 forAll(internalFaces_, i)
59 {
60 label facei = internalFaces_[i];
61 phii[facei] -= rho[facei]*(Omega ^ (Cfi[facei] - origin_)) & Sfi[facei];
62 }
63
64 makeRelativeRhoFlux(rho.boundaryField(), phi.boundaryFieldRef());
65}
66
67
68template<class RhoFieldType>
69void Foam::MRFZone::makeRelativeRhoFlux
70(
71 const RhoFieldType& rho,
72 FieldField<fvsPatchField, scalar>& phi
73) const
74{
75 if (!active_)
76 {
77 return;
78 }
79
80 const surfaceVectorField& Cf = mesh_.Cf();
81 const surfaceVectorField& Sf = mesh_.Sf();
82
83 const vector Omega = omega_->value(mesh_.time().timeOutputValue())*axis_;
84
85 // Included patches
86 forAll(includedFaces_, patchi)
87 {
88 forAll(includedFaces_[patchi], i)
89 {
90 label patchFacei = includedFaces_[patchi][i];
91
92 phi[patchi][patchFacei] = 0.0;
93 }
94 }
95
96 // Excluded patches
97 forAll(excludedFaces_, patchi)
98 {
99 forAll(excludedFaces_[patchi], i)
100 {
101 label patchFacei = excludedFaces_[patchi][i];
102
103 phi[patchi][patchFacei] -=
104 rho[patchi][patchFacei]
105 * (Omega ^ (Cf.boundaryField()[patchi][patchFacei] - origin_))
106 & Sf.boundaryField()[patchi][patchFacei];
107 }
108 }
109}
110
111
112template<class RhoFieldType>
113void Foam::MRFZone::makeRelativeRhoFlux
114(
115 const RhoFieldType& rho,
116 Field<scalar>& phi,
117 const label patchi
118) const
119{
120 if (!active_)
121 {
122 return;
123 }
124
125 const surfaceVectorField& Cf = mesh_.Cf();
126 const surfaceVectorField& Sf = mesh_.Sf();
127
128 const vector Omega = omega_->value(mesh_.time().timeOutputValue())*axis_;
129
130 // Included patches
131 forAll(includedFaces_[patchi], i)
132 {
133 label patchFacei = includedFaces_[patchi][i];
134
135 phi[patchFacei] = 0.0;
136 }
137
138 // Excluded patches
139 forAll(excludedFaces_[patchi], i)
140 {
141 label patchFacei = excludedFaces_[patchi][i];
142
143 phi[patchFacei] -=
144 rho[patchFacei]
145 * (Omega ^ (Cf.boundaryField()[patchi][patchFacei] - origin_))
146 & Sf.boundaryField()[patchi][patchFacei];
147 }
148}
149
150
151template<class RhoFieldType>
152void Foam::MRFZone::makeAbsoluteRhoFlux
153(
154 const RhoFieldType& rho,
156) const
157{
158 if (!active_)
159 {
160 return;
161 }
162
163 const surfaceVectorField& Cf = mesh_.Cf();
164 const surfaceVectorField& Sf = mesh_.Sf();
165
166 const vector Omega = omega_->value(mesh_.time().timeOutputValue())*axis_;
167
168 const vectorField& Cfi = Cf;
169 const vectorField& Sfi = Sf;
171
172 // Internal faces
173 forAll(internalFaces_, i)
174 {
175 label facei = internalFaces_[i];
176 phii[facei] += rho[facei]*(Omega ^ (Cfi[facei] - origin_)) & Sfi[facei];
177 }
178
180
181
182 // Included patches
183 forAll(includedFaces_, patchi)
184 {
185 forAll(includedFaces_[patchi], i)
186 {
187 label patchFacei = includedFaces_[patchi][i];
188
189 phibf[patchi][patchFacei] +=
190 rho.boundaryField()[patchi][patchFacei]
191 * (Omega ^ (Cf.boundaryField()[patchi][patchFacei] - origin_))
192 & Sf.boundaryField()[patchi][patchFacei];
193 }
194 }
195
196 // Excluded patches
197 forAll(excludedFaces_, patchi)
198 {
199 forAll(excludedFaces_[patchi], i)
200 {
201 label patchFacei = excludedFaces_[patchi][i];
202
203 phibf[patchi][patchFacei] +=
204 rho.boundaryField()[patchi][patchFacei]
205 * (Omega ^ (Cf.boundaryField()[patchi][patchFacei] - origin_))
206 & Sf.boundaryField()[patchi][patchFacei];
207 }
208 }
209}
210
211
212template<class Type>
214(
216) const
217{
218 if (!active_)
219 {
220 return;
221 }
222
224
225 forAll(internalFaces_, i)
226 {
227 phii[internalFaces_[i]] = Zero;
228 }
229
230 auto& phibf = phi.boundaryFieldRef();
231
232 forAll(includedFaces_, patchi)
233 {
234 forAll(includedFaces_[patchi], i)
235 {
236 phibf[patchi][includedFaces_[patchi][i]] = Zero;
237 }
238 }
239
240 forAll(excludedFaces_, patchi)
241 {
242 forAll(excludedFaces_[patchi], i)
243 {
244 phibf[patchi][excludedFaces_[patchi][i]] = Zero;
245 }
246 }
247}
248
249
250// ************************************************************************* //
surfaceScalarField & phi
Generic templated field type.
Definition: Field.H:82
Generic GeometricField class.
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
GeometricBoundaryField< scalar, fvsPatchField, surfaceMesh > Boundary
Type of boundary fields.
vector Omega() const
Return the current Omega vector.
Definition: MRFZone.C:264
scalar timeOutputValue() const
Return current time value.
Definition: TimeStateI.H:31
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
const surfaceVectorField & Cf() const
Return face centres as surfaceVectorField.
const surfaceVectorField & Sf() const
Return cell face area vectors.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
A special matrix type and solver, designed for finite volume solutions of scalar equations.
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Field< vector > vectorField
Specialisation of Field<T> for vector.
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
A non-counting (dummy) refCount.
Definition: refCount.H:59
Foam::surfaceFields.