MPLIC.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) 2020-2023 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::MPLIC
26 
27 Description
28  Multicut Piecewise-Linear Interface Calculation (MPLIC) corrected scheme is
29  a surface interpolation scheme for flux calculation in advection of a
30  bounded variable, e.g. phase fraction and for interface capturing in the
31  volume of fluid (VoF) method.
32 
33  The interface is represented by multiple cuts which split each cell to match
34  the volume fraction of the phase in the cell. The cut planes are oriented
35  according to the point field of the local phase fraction. The phase
36  fraction at each cell face - the interpolated value - is then calculated
37  from the face area on either side of the cuts.
38 
39  Three progressively more complex algorithms are used to ensure the cell
40  volume fraction is accurately reproduced:
41  -# single cut: cuts all the cell faces regardless the order
42  -# multi cut: topological face-edge-face walk which can split cell into
43  multiple sub-volumes
44  -# tetrahedron cut: decomposes cell into tetrahedrons which are cut
45 
46  Example:
47  \verbatim
48  divSchemes
49  {
50  .
51  .
52  div(phi,alpha) Gauss MPLIC;
53  .
54  .
55  }
56  \endverbatim
57 
58 See also
59  Foam::MPLICU
60  Foam::PLIC
61  Foam::PLICU
62  Foam::interfaceCompression
63 
64 SourceFiles
65  MPLIC.C
66 
67 \*---------------------------------------------------------------------------*/
68 
69 #ifndef MPLIC_H
70 #define MPLIC_H
71 
73 #include "DynamicList.H"
74 
75 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
76 
77 namespace Foam
78 {
79 
80 /*---------------------------------------------------------------------------*\
81  Class MPLIC Declaration
82 \*---------------------------------------------------------------------------*/
83 
84 class MPLIC
85 :
86  public surfaceInterpolationScheme<scalar>
87 {
88 protected:
89 
90  // Private member data
91 
92  const surfaceScalarField& phi_;
93 
94 
95  // Protected member functions
96 
97  //- Set alphaPhi for the faces of the given cell
98  void setCellAlphaf
99  (
100  const label celli,
101  const scalarField& phi,
102  scalarField& alphaf,
103  boolList& correctedFaces,
104  const DynamicList<scalar>& cellAlphaf,
105  const fvMesh& mesh
106  ) const;
107 
108  //- Return alpha interpolation
110  (
111  const volScalarField& alpha,
112  const surfaceScalarField& phi,
113  scalarField& splicedTvff,
114  const bool unweighted,
115  const scalar tol,
116  const bool isMPLIC = true
117  ) const;
118 
119 
120 public:
121 
122  //- Runtime type information
123  TypeName("MPLIC");
124 
125 
126  // Constructors
127 
128  //- Construct from faceFlux and Istream
129  MPLIC
130  (
131  const fvMesh& mesh,
132  const surfaceScalarField& faceFlux,
133  Istream& is
134  )
135  :
137  phi_(faceFlux)
138  {}
139 
140 
141  // Member Functions
142 
143  //- Return the interpolation weighting factors
145  (
146  const VolField<scalar>&
147  ) const
148  {
150 
151  return tmp<surfaceScalarField>(nullptr);
152  }
153 
154  //- Return the face-interpolate of the given cell field
156  (
157  const VolField<scalar>& vf
158  ) const;
159 
160 
161  // Member Operators
162 
163  //- Disallow default bitwise assignment
164  void operator=(const MPLIC&) = delete;
165 };
166 
167 
168 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
169 
170 } // End namespace Foam
171 
172 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
173 
174 #endif
175 
176 // ************************************************************************* //
Generic GeometricField class.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
Multicut Piecewise-Linear Interface Calculation (MPLIC) corrected scheme is a surface interpolation s...
Definition: MPLIC.H:86
TypeName("MPLIC")
Runtime type information.
void operator=(const MPLIC &)=delete
Disallow default bitwise assignment.
const surfaceScalarField & phi_
Definition: MPLIC.H:91
virtual tmp< surfaceScalarField > interpolate(const VolField< scalar > &vf) const
Return the face-interpolate of the given cell field.
Definition: MPLIC.C:237
tmp< surfaceScalarField > surfaceAlpha(const volScalarField &alpha, const surfaceScalarField &phi, scalarField &splicedTvff, const bool unweighted, const scalar tol, const bool isMPLIC=true) const
Return alpha interpolation.
Definition: MPLIC.C:82
MPLIC(const fvMesh &mesh, const surfaceScalarField &faceFlux, Istream &is)
Construct from faceFlux and Istream.
Definition: MPLIC.H:129
virtual tmp< surfaceScalarField > weights(const VolField< scalar > &) const
Return the interpolation weighting factors.
Definition: MPLIC.H:144
void setCellAlphaf(const label celli, const scalarField &phi, scalarField &alphaf, boolList &correctedFaces, const DynamicList< scalar > &cellAlphaf, const fvMesh &mesh) const
Set alphaPhi for the faces of the given cell.
Definition: MPLIC.C:47
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Abstract base class for surface interpolation schemes.
const fvMesh & mesh() const
Return mesh reference.
A class for managing temporary objects.
Definition: tmp.H:55
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:381
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
Namespace for OpenFOAM.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59