fvFieldDecomposer.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
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 Class
25  Foam::fvFieldDecomposer
26 
27 Description
28  Finite Volume volume and surface field decomposer.
29 
30 SourceFiles
31  fvFieldDecomposer.C
32  fvFieldDecomposerDecomposeFields.C
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #ifndef fvFieldDecomposer_H
37 #define fvFieldDecomposer_H
38 
39 #include "fvMesh.H"
41 #include "surfaceFields.H"
42 
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 
45 namespace Foam
46 {
47 
48 class IOobjectList;
49 
50 /*---------------------------------------------------------------------------*\
51  Class fvFieldDecomposer Declaration
52 \*---------------------------------------------------------------------------*/
53 
55 {
56 public:
57 
58  // Public Classes
59 
60  //- Patch field decomposer class
62  :
63  public labelList,
65  {
66  public:
67 
68  // Constructors
69 
70  //- Construct given addressing
72  };
73 
74 
75 private:
76 
77  // Private Data
78 
79  //- Reference to complete mesh
80  const fvMesh& completeMesh_;
81 
82  //- Reference to processor mesh
83  const fvMesh& procMesh_;
84 
85  //- Reference to face addressing
86  const labelList& faceAddressing_;
87 
88  //- Reference to cell addressing
89  const labelList& cellAddressing_;
90 
91  //- Reference to face addressing boundary field
92  const surfaceLabelField::Boundary& faceAddressingBf_;
93 
94  //- List of patch field decomposers
95  PtrList<patchFieldDecomposer> patchFieldDecomposers_;
96 
97 
98  // Private Member Functions
99 
100  //- Convert a processor patch to the corresponding complete patch index
101  label completePatchID(const label procPatchi) const;
102 
103  //- Map cell values to faces
104  template<class Type>
105  static tmp<Field<Type>> mapCellToFace
106  (
107  const labelUList& owner,
108  const labelUList& neighbour,
109  const Field<Type>& field,
110  const labelUList& addressing
111  );
112 
113  //- Map face values to faces
114  template<class Type>
115  static tmp<Field<Type>> mapFaceToFace
116  (
117  const Field<Type>& field,
118  const labelUList& addressing,
119  const bool isFlux
120  );
121 
122 
123 public:
124 
125  // Constructors
126 
127  //- Construct from components
129  (
130  const fvMesh& completeMesh,
131  const fvMesh& procMesh,
132  const labelList& faceAddressing,
133  const labelList& cellAddressing,
134  const surfaceLabelField::Boundary& faceAddressingBf
135  );
136 
137  //- Disallow default bitwise copy construction
138  fvFieldDecomposer(const fvFieldDecomposer&) = delete;
139 
140 
141  //- Destructor
143 
144 
145  // Member Functions
146 
147  //- Decompose volume field
148  template<class Type>
151  (
153  const bool allowUnknownPatchFields = false
154  ) const;
155 
156  //- Decompose surface field
157  template<class Type>
160  (
162  ) const;
163 
164  //- Decompose a list of fields
165  template<class GeoField>
166  void decomposeFields(const PtrList<GeoField>& fields) const;
167 
168 
169  // Member Operators
170 
171  //- Disallow default bitwise assignment
172  void operator=(const fvFieldDecomposer&) = delete;
173 };
174 
175 
176 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
177 
178 } // End namespace Foam
179 
180 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
181 
182 #ifdef NoRepository
184 #endif
185 
186 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
187 
188 #endif
189 
190 // ************************************************************************* //
Foam::surfaceFields.
Generic GeometricField class.
tmp< GeometricField< Type, fvPatchField, volMesh > > decomposeField(const GeometricField< Type, fvPatchField, volMesh > &field, const bool allowUnknownPatchFields=false) const
Decompose volume field.
fvFieldDecomposer(const fvMesh &completeMesh, const fvMesh &procMesh, const labelList &faceAddressing, const labelList &cellAddressing, const surfaceLabelField::Boundary &faceAddressingBf)
Construct from components.
Pre-declare SubField and related Field type.
Definition: Field.H:56
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
virtual const labelUList & addressing() const
Access to the direct map addressing.
Finite Volume volume and surface field decomposer.
patchFieldDecomposer(const labelUList &addressing)
Construct given addressing.
Generic GeometricBoundaryField class.
Definition: fvMesh.H:80
bool isFlux(const DimensionedField< Type, surfaceMesh > &df)
Check if surfaceField is a flux.
Definition: surfaceFields.H:55
void operator=(const UList< label > &)
Assignment to UList operator. Takes linear time.
Definition: List.C:376
Info<< "Reading field p_rgh\"<< endl;volScalarField p_rgh(IOobject("p_rgh", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);pressureReference pressureReference(p, p_rgh, pimple.dict(), thermo.incompressible());mesh.schemes().setFluxRequired(p_rgh.name());hydrostaticInitialisation(p_rgh, p, rho, U, gh, ghf, pRef, thermo, pimple.dict());Info<< "Creating field dpdt\"<< endl;volScalarField dpdt(IOobject("dpdt", runTime.timeName(), mesh), mesh, dimensionedScalar(p.dimensions()/dimTime, 0));Info<< "Creating field kinetic energy K\"<< endl;volScalarField K("K", 0.5 *magSqr(U));dimensionedScalar initialMass=fvc::domainIntegrate(rho);multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:131
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
~fvFieldDecomposer()
Destructor.
rDeltaTY field()
A class for managing temporary objects.
Definition: PtrList.H:53
void decomposeFields(const PtrList< GeoField > &fields) const
Decompose a list of fields.
Namespace for OpenFOAM.