interfaceCompression.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-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::interfaceCompression
26 
27 Description
28  Interface compression corrected scheme, based on counter-gradient
29  transport, to maintain sharp interfaces during VoF simulations.
30 
31  The interface compression is applied to the face interpolated field from a
32  suitable 2nd-order shape-preserving NVD or TVD scheme, e.g. vanLeer or
33  vanAlbada. A coefficient is supplied to control the degree of compression,
34  with a value of 1 suitable for most VoF cases to ensure interface integrity.
35  A value larger than 1 can be used but the additional compression can bias
36  the interface to follow the mesh more closely while a value smaller than 1
37  can lead to interface smearing.
38 
39  Example:
40  \verbatim
41  divSchemes
42  {
43  .
44  .
45  div(phi,alpha) Gauss interfaceCompression vanLeer 1;
46  .
47  .
48  }
49  \endverbatim
50 
51 See also
52  Foam::vanLeer
53  Foam::vanAlbada
54  Foam::PLIC
55  Foam::PLICU
56  Foam::MPLIC
57  Foam::MPLICU
58 
59 SourceFiles
60  interfaceCompression.C
61 
62 \*---------------------------------------------------------------------------*/
63 
64 #ifndef interfaceCompression_H
65 #define interfaceCompression_H
66 
68 
69 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
70 
71 namespace Foam
72 {
73 
74 /*---------------------------------------------------------------------------*\
75  Class interfaceCompression Declaration
76 \*---------------------------------------------------------------------------*/
77 
79 :
80  public surfaceInterpolationScheme<scalar>
81 {
82  // Private member data
83 
84  const surfaceScalarField& phi_;
85 
86  //- Base scheme to which the compression is applied
88 
89  //- Compression factor
90  const scalar cAlpha_;
91 
92 
93 public:
94 
95  //- Runtime type information
96  TypeName("interfaceCompression");
97 
98 
99  // Constructors
100 
101  //- Construct from faceFlux and Istream
103  (
104  const fvMesh& mesh,
105  const surfaceScalarField& faceFlux,
106  Istream& is
107  )
108  :
110  phi_(faceFlux),
111  tScheme_
112  (
113  surfaceInterpolationScheme<scalar>::New(mesh, faceFlux, is)
114  ),
115  cAlpha_(readScalar(is))
116  {}
117 
118 
119  // Member Functions
120 
121  //- Return the interpolation weighting factors
123  (
124  const volScalarField&
125  ) const
126  {
128 
129  return tmp<surfaceScalarField>(nullptr);
130  }
131 
132  //- Return the face-interpolate of the given cell field
134  (
135  const volScalarField& vf
136  ) const;
137 
138 
139  // Member Operators
140 
141  //- Disallow default bitwise assignment
142  void operator=(const interfaceCompressionNew&) = delete;
143 };
144 
145 
146 extern const wordHashSet compressionSchemes;
147 
148 
149 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
150 
151 } // End namespace Foam
152 
153 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
154 
155 #endif
156 
157 // ************************************************************************* //
Generic GeometricField class.
A HashTable with keys but without contents.
Definition: HashSet.H:62
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
virtual tmp< surfaceScalarField > interpolate(const volScalarField &vf) const
Return the face-interpolate of the given cell field.
TypeName("interfaceCompression")
Runtime type information.
virtual tmp< surfaceScalarField > weights(const volScalarField &) const
Return the interpolation weighting factors.
interfaceCompressionNew(const fvMesh &mesh, const surfaceScalarField &faceFlux, Istream &is)
Construct from faceFlux and Istream.
void operator=(const interfaceCompressionNew &)=delete
Disallow default bitwise assignment.
Abstract base class for surface interpolation schemes.
const fvMesh & mesh() const
Return mesh reference.
static tmp< surfaceInterpolationScheme< scalar > > New(const fvMesh &mesh, Istream &schemeData)
Return new tmp interpolation scheme.
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
Namespace for OpenFOAM.
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:75
const wordHashSet compressionSchemes