fvMeshSubset.H 10.8 KB
Newer Older
1
2
3
4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
5
    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
6
     \\/     M anipulation  | Copyright (C) 2016 OpenCFD Ltd.
7
8
9
10
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

11
12
13
14
    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.
15
16
17
18
19
20
21

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    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
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43

Class
    Foam::fvMeshSubset

Description
    Post-processing mesh subset tool.  Given the original mesh and the
    list of selected cells, it creates the mesh consisting only of the
    desired cells, with the mapping list for points, faces, and cells.

    Puts all exposed internal faces into either
    - a user supplied patch
    - a newly created patch "oldInternalFaces"

    - setCellSubset is for small subsets. Uses Maps to minimize memory.
    - setLargeCellSubset is for largish subsets (>10% of mesh).
      Uses labelLists instead.

    - setLargeCellSubset does coupled patch subsetting as well. If it detects
      a face on a coupled patch 'losing' its neighbour it will move the
      face into the oldInternalFaces patch.

44
45
46
47
48
49
    - if a user supplied patch is used it is up to the destination
      patchField to handle exposed internal faces (mapping from face -1).
      If not provided the default is to assign the internalField. All the
      basic patch field types (e.g. fixedValue) will give a warning and
      preferably derived patch field types should be used that know how to
      handle exposed faces (e.g. use uniformFixedValue instead of fixedValue)
50
51
52
53
54
55
56
57
58
59

SourceFiles
    fvMeshSubset.C

\*---------------------------------------------------------------------------*/

#ifndef fvMeshSubset_H
#define fvMeshSubset_H

#include "fvMesh.H"
Mattijs Janssens's avatar
Mattijs Janssens committed
60
#include "pointMesh.H"
61
#include "GeometricField.H"
62
#include "HashSet.H"
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
#include "surfaceMesh.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

namespace Foam
{

/*---------------------------------------------------------------------------*\
                        Class fvMeshSubset Declaration
\*---------------------------------------------------------------------------*/

class fvMeshSubset
{

private:

    // Private data

        //- Mesh to subset from
        const fvMesh& baseMesh_;

        //- Subset mesh pointer
Mattijs Janssens's avatar
Mattijs Janssens committed
85
86
        autoPtr<fvMesh> fvMeshSubsetPtr_;

87
88
89
90
91
92
93
94
95
96
97
98
        //- Point mapping array
        labelList pointMap_;

        //- Face mapping array
        labelList faceMap_;

        //- Cell mapping array
        labelList cellMap_;

        //- Patch mapping array
        labelList patchMap_;

99
        //- Optional face mapping array with flip encoded
100
101
        mutable autoPtr<labelList> faceFlipMapPtr_;

102
103
104
105
106
107
108

    // Private Member Functions

        //- Check if subset has been performed
        bool checkCellSubset() const;

        //- Mark points in Map
109
        static void markPoints(const labelList&, Map<label>&);
110
111

        //- Mark points (with 0) in labelList
112
        static void markPoints(const labelList&, labelList&);
113
114
115
116
117

        //- Adapt nCellsUsingFace for coupled faces becoming 'uncoupled'.
        void doCoupledPatches
        (
            const bool syncPar,
118
            Map<label>& facesToSubset,
119
120
121
122
123
124
125
126
127
128
129
130
131
132
            labelList& nCellsUsingFace
        ) const;

        //- Subset of subset
        static labelList subset
        (
            const label nElems,
            const labelList& selectedElements,
            const labelList& subsetMap
        );

        //- Create zones for submesh
        void subsetZones();

133
134
135
136
137
138
139
        //- Helper: extract cells-to-remove from cells-to-keep
        labelList getCellsToRemove
        (
            const labelList& region,
            const label currentRegion
        ) const;

140
141
        //- No copy construct
        fvMeshSubset(const fvMeshSubset&) = delete;
142

143
144
        //- No copy assignment
        void operator=(const fvMeshSubset&) = delete;
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164

public:

    // Constructors

        //- Construct given a mesh to subset
        explicit fvMeshSubset(const fvMesh&);


    // Member Functions

        // Edit

            //- Set the subset. Create "oldInternalFaces" patch for exposed
            //  internal faces (patchID==-1) or use supplied patch.
            //  Does not handle coupled patches correctly if only one side
            //  gets deleted.
            void setCellSubset
            (
                const labelHashSet& globalCellMap,
165
166
                const label patchID = -1,
                const bool syncPar = true
167
168
169
170
171
            );

            //- Set the subset from all cells with region == currentRegion.
            //  Create "oldInternalFaces" patch for exposed
            //  internal faces (patchID==-1) or use supplied patch.
andy's avatar
andy committed
172
            //  Handles coupled patches by if necessary making coupled patch
173
174
175
176
177
178
179
180
181
            //  face part of patchID (so uncoupled)
            void setLargeCellSubset
            (
                const labelList& region,
                const label currentRegion,
                const label patchID = -1,
                const bool syncCouples = true
            );

182
183
184
            //- setLargeCellSubset but only marking certain cells
            void setLargeCellSubset
            (
185
                const labelUList& globalCellMap,
186
187
188
189
                const label patchID = -1,
                const bool syncPar = true
            );

190
191
192
193
194
195
196
197
198
            //- setLargeCellSubset but with labelHashSet.
            void setLargeCellSubset
            (
                const labelHashSet& globalCellMap,
                const label patchID = -1,
                const bool syncPar = true
            );


199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
            //- Two step subsetting

                //- Get labels of exposed faces.
                //  These are
                //  - internal faces that become boundary faces
                //  - coupled faces that become uncoupled (since one of the
                //    sides gets deleted)
                labelList getExposedFaces
                (
                    const labelList& region,
                    const label currentRegion,
                    const bool syncCouples = true
                ) const;

                //- For every exposed face (from above getExposedFaces)
                //  used supplied (existing!) patch
                void setLargeCellSubset
                (
                    const labelList& region,
                    const label currentRegion,
                    const labelList& exposedFaces,
                    const labelList& patchIDs,
                    const bool syncCouples = true
                );


225
226
227
228
229
230
231
232
        // Access

            //- Original mesh
            const fvMesh& baseMesh() const
            {
                return baseMesh_;
            }

mattijs's avatar
mattijs committed
233
234
235
            //- Have subMesh?
            bool hasSubMesh() const;

236
237
            //- Return reference to subset mesh
            const fvMesh& subMesh() const;
Mattijs Janssens's avatar
Mattijs Janssens committed
238

239
240
241
242
243
244
245
246
            fvMesh& subMesh();

            //- Return point map
            const labelList& pointMap() const;

            //- Return face map
            const labelList& faceMap() const;

247
248
249
            //- Return face map with sign to encode flipped faces
            const labelList& faceFlipMap() const;

250
251
252
253
254
255
256
257
258
259
260
            //- Return cell map
            const labelList& cellMap() const;

            //- Return patch map
            const labelList& patchMap() const;


        // Field mapping

            //- Map volume field
            template<class Type>
261
            static tmp<GeometricField<Type, fvPatchField, volMesh>>
262
263
264
265
266
267
268
269
270
271
            interpolate
            (
                const GeometricField<Type, fvPatchField, volMesh>&,
                const fvMesh& sMesh,
                const labelList& patchMap,
                const labelList& cellMap,
                const labelList& faceMap
            );

            template<class Type>
272
            tmp<GeometricField<Type, fvPatchField, volMesh>>
273
274
275
276
277
            interpolate
            (
                const GeometricField<Type, fvPatchField, volMesh>&
            ) const;

278
279
            //- Map surface field. Optionally negates value if flipping
            //  a face (from exposing an internal face)
280
            template<class Type>
281
            static tmp<GeometricField<Type, fvsPatchField, surfaceMesh>>
Mattijs Janssens's avatar
Mattijs Janssens committed
282
            interpolate
283
284
285
286
            (
                const GeometricField<Type, fvsPatchField, surfaceMesh>&,
                const fvMesh& sMesh,
                const labelList& patchMap,
287
                const labelList& cellMap,
288
                const labelList& faceMap
289
290
291
            );

            template<class Type>
292
            tmp<GeometricField<Type, fvsPatchField, surfaceMesh>>
293
294
            interpolate
            (
295
                const GeometricField<Type, fvsPatchField, surfaceMesh>&
296
            ) const;
Mattijs Janssens's avatar
Mattijs Janssens committed
297
298
299

            //- Map point field
            template<class Type>
300
            static tmp<GeometricField<Type, pointPatchField, pointMesh>>
Mattijs Janssens's avatar
Mattijs Janssens committed
301
302
303
304
305
306
307
308
309
            interpolate
            (
                const GeometricField<Type, pointPatchField, pointMesh>&,
                const pointMesh& sMesh,
                const labelList& patchMap,
                const labelList& pointMap
            );

            template<class Type>
310
            tmp<GeometricField<Type, pointPatchField, pointMesh>>
Mattijs Janssens's avatar
Mattijs Janssens committed
311
312
313
314
            interpolate
            (
                const GeometricField<Type, pointPatchField, pointMesh>&
            ) const;
315

316
            //- Map dimensioned field
317
            template<class Type>
318
            static tmp<DimensionedField<Type, volMesh>>
319
320
321
322
323
324
325
326
            interpolate
            (
                const DimensionedField<Type, volMesh>&,
                const fvMesh& sMesh,
                const labelList& cellMap
            );

            template<class Type>
327
            tmp<DimensionedField<Type, volMesh>>
328
            interpolate(const DimensionedField<Type, volMesh>&) const;
329
330
331
332
333
334
335
336
337
338
};


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

} // End namespace Foam

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#ifdef NoRepository
339
    #include "fvMeshSubsetInterpolate.C"
340
341
342
343
344
345
346
#endif

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#endif

// ************************************************************************* //