correctAlpha.C
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) 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 \*---------------------------------------------------------------------------*/
25 
26 #include "isothermalFilm.H"
27 #include "surfaceTensionModel.H"
28 #include "constrainHbyA.H"
29 #include "fvcFlux.H"
30 #include "fvcSnGrad.H"
31 #include "fvcReconstruct.H"
32 #include "fvmDiv.H"
33 #include "fvmLaplacian.H"
34 
35 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
36 
37 void Foam::solvers::isothermalFilm::correctAlpha()
38 {
40 
41  fvVectorMatrix& UEqn = tUEqn.ref();
42 
44 
45  const surfaceScalarField pbByAlphaf(this->pbByAlphaf());
46  const surfaceScalarField pbByAlphaGradRhof
47  (
48  constrainedField(this->pbByAlphaGradRhof()*mesh.magSf())
49  );
50 
51  const surfaceScalarField phip
52  (
53  constrainedField
54  (
55  fvc::snGrad(pe() + pc(surfaceTension->sigma()), "snGrad(p)")
56  *mesh.magSf()
57  - rhof*(g & mesh.Sf())
58  )
59  );
60 
61  while (pimple.correct())
62  {
63  const volScalarField rAU(1/UEqn.A());
65 
66  const surfaceScalarField alphaf
67  (
68  constrainedField(fvc::interpolate(alpha))
69  );
70 
71  const surfaceScalarField alpharAUf
72  (
73  constrainedField(fvc::interpolate(alpha*rAU))
74  );
75 
76  const surfaceScalarField phig("phig", phip + pbByAlphaGradRhof*alphaf);
77 
78  phi_ = constrainedField(fvc::flux(HbyA) - alpharAUf*phig);
79 
80  const surfaceScalarField phid("phid", rhof*phi);
81 
82  const surfaceScalarField rAUf
83  (
84  "rAUf",
85  alphaf*rhof*alpharAUf*pbByAlphaf
86  );
87 
88  fvScalarMatrix alphadEqn
89  (
91  + fvm::div(phid, alpha)
92  ==
93  fvModels().source(rho, alpha)
94  );
95 
97  {
98  fvScalarMatrix alphaEqn(alphadEqn - fvm::laplacian(rAUf, alpha));
99 
100  alphaEqn.solve();
101 
103  {
104  alphaRhoPhi_ = alphaEqn.flux();
105  }
106  }
107 
108  const surfaceScalarField phiGradAlpha
109  (
110  constrainedField(pbByAlphaf*fvc::snGrad(alpha)*mesh.magSf())
111  );
112 
113  phi_ -= alpharAUf*phiGradAlpha;
114 
115  U_ = HbyA - rAU*fvc::reconstruct(alphaf*(phig + phiGradAlpha));
116 
117  // Remove film-normal component of velocity
118  U_ -= nHat*(nHat & U);
119 
121 
123 
125 
126  continuityErrors();
127  }
128 
129  // Update film thickness
130  correctDelta();
131 
132  tUEqn.clear();
133 }
134 
135 
136 // ************************************************************************* //
void correctBoundaryConditions()
Correct boundary field.
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
const surfaceVectorField & Sf() const
Return cell face area vectors.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
bool finalNonOrthogonalIter() const
Flag to indicate the last non-orthogonal iteration.
bool correct()
Piso loop within outer loop.
bool correctNonOrthogonal()
Non-orthogonal corrector loop.
Foam::fvModels & fvModels() const
Return the fvModels that are created on demand.
Definition: solver.C:96
pimpleNoLoopControl pimple
PIMPLE inner-loop controls.
Definition: solver.H:107
Foam::fvConstraints & fvConstraints() const
Return the fvConstraints that are created on demand.
Definition: solver.C:107
const fvMesh & mesh
Region mesh.
Definition: solver.H:101
const surfaceScalarField & phi
Reference to the film volumetric-flux field.
volVectorField U_
Film velocity field.
tmp< fvVectorMatrix > tUEqn
Cached momentum matrix.
const volScalarField & alpha
Film volume fraction in the cell layer.
const volVectorField & U
Reference to the film velocity field.
autoPtr< surfaceTensionModel > surfaceTension
Pointer to the surface tension coefficient model.
const uniformDimensionedVectorField g
Acceleration due to gravity.
surfaceScalarField alphaRhoPhi_
Film mass-flux field.
surfaceScalarField phi_
Film volumetric-flux field.
const volScalarField & rho
Reference to the thermodynamic density field.
const volVectorField & nHat
Film wall normal.
volScalarField alpha_
Film volume fraction in the cell layer.
fvVectorMatrix & UEqn
Definition: UEqn.H:11
Calculate the face-flux of the given field.
Reconstruct volField from a face flux field.
Calculate the snGrad of the given volField.
Calculate the matrix for the divergence of the given field and flux.
Calculate the matrix for the laplacian of the field.
volScalarField rAU(1.0/UEqn.A())
volVectorField & HbyA
Definition: pEqn.H:13
tmp< SurfaceField< typename innerProduct< vector, Type >::type > > flux(const VolField< Type > &vf)
Return the face-flux field obtained from the given volVectorField.
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< VolField< typename outerProduct< vector, Type >::type > > reconstruct(const SurfaceField< Type > &ssf)
tmp< SurfaceField< Type > > snGrad(const VolField< Type > &vf, const word &name)
Definition: fvcSnGrad.C:45
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const VolField< Type > &vf, const word &name)
Definition: fvmDiv.C:48
tmp< fvMatrix< Type > > ddt(const VolField< Type > &vf)
Definition: fvmDdt.C:46
VolField< vector > volVectorField
Definition: volFieldsFwd.H:65
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
SurfaceField< scalar > surfaceScalarField
void constrainHbyA(volVectorField &HbyA, const volVectorField &U, const volScalarField &p)
Definition: constrainHbyA.C:34
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
fvMatrix< vector > fvVectorMatrix
Definition: fvMatricesFwd.H:45