fvcCellReduce.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) 2013-2016 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 "fvcCellReduce.H"
29#include "fvMesh.H"
30#include "volFields.H"
31#include "surfaceFields.H"
33
34// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35
36template<class Type, class CombineOp>
39(
40 const GeometricField<Type, fvsPatchField, surfaceMesh>& ssf,
41 const CombineOp& cop,
42 const Type& nullValue
43)
44{
45 typedef GeometricField<Type, fvPatchField, volMesh> volFieldType;
46
47 const fvMesh& mesh = ssf.mesh();
48
49 tmp<volFieldType> tresult
50 (
51 new volFieldType
52 (
53 IOobject
54 (
55 "cellReduce(" + ssf.name() + ')',
56 ssf.instance(),
57 mesh,
58 IOobject::NO_READ,
59 IOobject::NO_WRITE
60 ),
61 mesh,
62 dimensioned<Type>("initialValue", ssf.dimensions(), nullValue),
63 extrapolatedCalculatedFvPatchField<Type>::typeName
64 )
65 );
66
67 volFieldType& result = tresult.ref();
68
69 const labelUList& own = mesh.owner();
70 const labelUList& nbr = mesh.neighbour();
71
72 forAll(own, i)
73 {
74 label celli = own[i];
75 cop(result[celli], ssf[i]);
76 }
77 forAll(nbr, i)
78 {
79 label celli = nbr[i];
80 cop(result[celli], ssf[i]);
81 }
82
83 result.correctBoundaryConditions();
84
85 return tresult;
86}
87
88
89template<class Type, class CombineOp>
92(
93 const tmp<GeometricField<Type, fvsPatchField, surfaceMesh>>& tssf,
94 const CombineOp& cop,
95 const Type& nullValue
96)
97{
98 tmp<GeometricField<Type, fvPatchField, volMesh>> tvf
99 (
100 cellReduce
101 (
102 tssf,
103 cop,
104 nullValue
105 )
106 );
107
108 tssf.clear();
109
110 return tvf;
111}
112
113
114// ************************************************************************* //
A class for managing temporary objects.
Definition: tmp.H:65
dynamicFvMesh & mesh
Construct a volume field from a surface field using a combine operator.
tmp< GeometricField< Type, fvPatchField, volMesh > > cellReduce(const GeometricField< Type, fvsPatchField, surfaceMesh > &, const CombineOp &cop, const Type &nullValue=pTraits< Type >::zero)
UList< label > labelUList
A UList of labels.
Definition: UList.H:85
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Foam::surfaceFields.