searchableSurfaceSelection.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) 2012-2016 OpenFOAM Foundation
9 Copyright (C) 2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
31#include "syncTools.H"
32#include "searchableSurface.H"
33#include "fvMesh.H"
34#include "Time.H"
35
36// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
37
38namespace Foam
39{
40namespace faceSelections
41{
42 defineTypeNameAndDebug(searchableSurfaceSelection, 0);
44 (
45 faceSelection,
46 searchableSurfaceSelection,
47 dictionary
48 );
49}
50}
51
52
53// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
54
56(
57 const word& name,
58 const fvMesh& mesh,
59 const dictionary& dict
60)
61:
62 faceSelection(name, mesh, dict),
63 surfacePtr_
64 (
65 searchableSurface::New
66 (
67 dict.get<word>("surface"),
68 IOobject
69 (
70 dict.getOrDefault("name", mesh.objectRegistry::db().name()),
71 mesh.time().constant(),
72 "triSurface",
73 mesh.objectRegistry::db(),
74 IOobject::MUST_READ,
75 IOobject::NO_WRITE
76 ),
77 dict
78 )
79 )
80{}
81
82
83// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
84
86{}
87
88
89// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
90
92(
93 const label zoneID,
94 labelList& faceToZoneID,
95 boolList& faceToFlip
96) const
97{
98 // Get cell-cell centre vectors
99
100 pointField start(mesh_.nFaces());
101 pointField end(mesh_.nFaces());
102
103 // Internal faces
104 for (label facei = 0; facei < mesh_.nInternalFaces(); facei++)
105 {
106 start[facei] = mesh_.cellCentres()[mesh_.faceOwner()[facei]];
107 end[facei] = mesh_.cellCentres()[mesh_.faceNeighbour()[facei]];
108 }
109
110 // Boundary faces
111 vectorField neighbourCellCentres;
113 (
114 mesh_,
115 mesh_.cellCentres(),
116 neighbourCellCentres
117 );
118
119 const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
120
121 forAll(pbm, patchi)
122 {
123 const polyPatch& pp = pbm[patchi];
124
125 if (pp.coupled())
126 {
127 forAll(pp, i)
128 {
129 label facei = pp.start()+i;
130 start[facei] = mesh_.cellCentres()[mesh_.faceOwner()[facei]];
131 end[facei] = neighbourCellCentres[facei-mesh_.nInternalFaces()];
132 }
133 }
134 else
135 {
136 forAll(pp, i)
137 {
138 label facei = pp.start()+i;
139 start[facei] = mesh_.cellCentres()[mesh_.faceOwner()[facei]];
140 end[facei] = mesh_.faceCentres()[facei];
141 }
142 }
143 }
144
145 List<pointIndexHit> hits;
146 surfacePtr_().findLine(start, end, hits);
147 pointField normals;
148 surfacePtr_().getNormal(hits, normals);
149
150 //- Note: do not select boundary faces.
151
152 for (label facei = 0; facei < mesh_.nInternalFaces(); facei++)
153 {
154 if (hits[facei].hit())
155 {
156 faceToZoneID[facei] = zoneID;
157 vector d = end[facei]-start[facei];
158 faceToFlip[facei] = ((normals[facei] & d) < 0);
159 }
160 }
161 forAll(pbm, patchi)
162 {
163 const polyPatch& pp = pbm[patchi];
164
165 if (pp.coupled())
166 {
167 forAll(pp, i)
168 {
169 label facei = pp.start()+i;
170 if (hits[facei].hit())
171 {
172 faceToZoneID[facei] = zoneID;
173 vector d = end[facei]-start[facei];
174 faceToFlip[facei] = ((normals[facei] & d) < 0);
175 }
176 }
177 }
178 }
179
180 faceSelection::select(zoneID, faceToZoneID, faceToFlip);
181}
182
183
184// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
virtual void select(const label, labelList &, boolList &) const =0
Selects all (internal or coupled) faces intersecting the searchableSurface.
virtual void select(const label zoneID, labelList &, boolList &) const
constant condensation/saturation model.
static void swapBoundaryCellPositions(const polyMesh &mesh, const UList< point > &cellData, List< point > &neighbourCellData)
Swap to obtain neighbour cell positions for all boundary faces.
Definition: syncTools.C:34
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
dynamicFvMesh & mesh
const labelIOList & zoneID
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
Namespace for OpenFOAM.
List< label > labelList
A List of labels.
Definition: List.H:66
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:44
Field< vector > vectorField
Specialisation of Field<T> for vector.
List< bool > boolList
A List of bools.
Definition: List.H:64
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
Definition: stdFoam.H:158
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333