cutCellIso.H
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) 2016-2017 DHI
9 Copyright (C) 2016-2017 OpenCFD Ltd.
10 Copyright (C) 2018 Johan Roenby
11 Copyright (c) 2017-2019, German Aerospace Center (DLR)
12-------------------------------------------------------------------------------
13License
14 This file is part of OpenFOAM.
15
16 OpenFOAM is free software: you can redistribute it and/or modify it
17 under the terms of the GNU General Public License as published by
18 the Free Software Foundation, either version 3 of the License, or
19 (at your option) any later version.
20
21 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
22 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
23 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
24 for more details.
25
26 You should have received a copy of the GNU General Public License
27 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
28
29Class
30 Foam::cutCellIso
31
32Description
33 Class for cutting a cell, celli, of an fvMesh, mesh_, at its intersection
34 with an isosurface defined by the mesh point values f_ and the cutValue,
35 isoValue_.
36
37 Reference:
38 \verbatim
39 Roenby, J., Bredmose, H. and Jasak, H. (2016).
40 A computational method for sharp interface advection
41 Royal Society Open Science, 3
42 doi 10.1098/rsos.160405
43
44 Henning Scheufler, Johan Roenby,
45 Accurate and efficient surface reconstruction from volume
46 fraction data on general meshes,
47 Journal of Computational Physics, 2019,
48 doi 10.1016/j.jcp.2019.01.009
49
50 \endverbatim
51
52 Original code supplied by Johan Roenby, DHI (2016)
53
54SourceFiles
55 cutCellIso.C
56
57\*---------------------------------------------------------------------------*/
58
59#ifndef cutCellIso_H
60#define cutCellIso_H
61
62#include "cutCell.H"
63#include "cutFaceIso.H"
64#include "fvMesh.H"
65#include "surfaceFields.H"
66#include "volFields.H"
67
68// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
69
70namespace Foam
71{
72
73/*---------------------------------------------------------------------------*\
74 Class cutCellIso Declaration
75\*---------------------------------------------------------------------------*/
77class cutCellIso
78:
79 public cutCell
80{
81 // Private Data
82
83 //- Mesh whose cells and faces to cut at their intersection with an
84 // isosurface.
85 const fvMesh& mesh_;
86
87 //- Cell to cut
88 label cellI_;
89
90 //- Isofunction values at mesh points. f_size() = mesh_.nPoints().
91 scalarField& f_;
92
93 //- Cutvalue used to cut cell
94 scalar cutValue_;
95
96 //- An cutFaceIso object to get access to its face cutting functionality
97 cutFaceIso cutFace_;
98
99 //- List of face centres for CutFaces
100 DynamicList<point> cutFaceCentres_;
101
102 //- List of face area vectors for isoCutFaces
103 DynamicList<vector> cutFaceAreas_;
104
105 //- Storage for subFace edges belonging to isoFace
106 DynamicList<DynamicList<point>> isoFaceEdges_;
107
108 //- Points constituting the cell-isosurface intersection (isoface)
109 DynamicList<point> facePoints_;
110
111 //- Face centre of the cutFace
112 point faceCentre_;
113
114 //- Face normal of the isoface by convention pointing from high to low
115 // values (i.e. opposite of the gradient vector).
116 vector faceArea_;
117
118 //- Cell centre of the subcell of celli which is "fully submerged", i.e.
119 // where the function value is higher than the isoValue_
120 point subCellCentre_;
121
122 //- Volume of fully submerged subcell
123 scalar subCellVolume_;
124
125 //- Volume of Fluid for cellI (subCellVolume_/mesh_.V()[cellI])
126 scalar VOF_;
127
128 //- A cell status label taking one of the values:
129 //
130 // -1: cell is fully below the isosurface
131 // 0: cell is cut
132 // +1: cell is fully above the isosurface
133 label cellStatus_;
134
135
136 public:
137
138 // Constructors
139
140 //- Construct from fvMesh and a scalarField
141 // Length of scalarField should equal number of mesh points
143
144
145 // Member Functions
146
147 //- Sets internal values and returns face status
148 label calcSubCell(const label cellI, const scalar cutValue);
149
150 //- Returns subCellCentre
151 const point& subCellCentre() const noexcept
152 {
153 return subCellCentre_;
154 }
155
156 //- Returns subCellVolume
157 scalar subCellVolume() const noexcept
158 {
159 return subCellVolume_;
160 }
161
162 //- Returns the points of the cutting PLICface
164
165 //- Returns the centre of the cutting PLICface
166 const point& faceCentre() const noexcept
167 {
168 return faceCentre_;
169 }
170
171 //- Returns the area normal vector of the cutting PLICface
172 const vector& faceArea() const noexcept
173 {
174 return faceArea_;
175 }
176
177 //- Returns cellStatus
178 label cellStatus() const noexcept
179 {
180 return cellStatus_;
181 }
182
183 //- Returns volume of fluid value
184 scalar VolumeOfFluid() const noexcept
185 {
186 return VOF_;
187 }
188
189 //- Returns cutValue
190 scalar cutValue() const noexcept
191 {
192 return cutValue_;
193 }
194
195 //- Resets internal values
196 void clearStorage();
197};
198
199// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
200
201} // End namespace Foam
202
203// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
204
205#endif
206
207// ************************************************************************* //
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
Definition: DynamicList.H:72
Class for cutting a cell, celli, of an fvMesh, mesh_, at its intersection with an isosurface defined ...
Definition: cutCellIso.H:79
const point & faceCentre() const noexcept
Returns the centre of the cutting PLICface.
Definition: cutCellIso.H:165
const vector & faceArea() const noexcept
Returns the area normal vector of the cutting PLICface.
Definition: cutCellIso.H:171
scalar cutValue() const noexcept
Returns cutValue.
Definition: cutCellIso.H:189
label calcSubCell(const label cellI, const scalar cutValue)
Sets internal values and returns face status.
Definition: cutCellIso.C:60
const DynamicList< point > & facePoints()
Returns the points of the cutting PLICface.
Definition: cutCellIso.C:173
void clearStorage()
Resets internal values.
Definition: cutCellIso.C:191
scalar VolumeOfFluid() const noexcept
Returns volume of fluid value.
Definition: cutCellIso.H:183
label cellStatus() const noexcept
Returns cellStatus.
Definition: cutCellIso.H:177
scalar subCellVolume() const noexcept
Returns subCellVolume.
Definition: cutCellIso.H:156
const point & subCellCentre() const noexcept
Returns subCellCentre.
Definition: cutCellIso.H:150
Service routines for cutting a cell, celli, of an fvMesh, mesh_, at its intersection with a surface.
Definition: cutCell.H:60
Class for cutting a face, faceI, of an fvMesh, mesh_, at its intersection with an isosurface defined ...
Definition: cutFaceIso.H:71
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
dynamicFvMesh & mesh
Namespace for OpenFOAM.
const direction noexcept
Definition: Scalar.H:223
labelList f(nPoints)
Foam::surfaceFields.