All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
compressibilityEqns.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) 2022-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 "multiphaseEuler.H"
27 #include "fvcDdt.H"
28 #include "fvcDiv.H"
29 #include "fvcSup.H"
30 #include "fvmDdt.H"
31 #include "fvmDiv.H"
32 #include "fvmSup.H"
33 
34 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
35 
37 Foam::solvers::multiphaseEuler::compressibilityEqns
38 (
39  const PtrList<volScalarField>& dmdts,
40  const PtrList<volScalarField>& d2mdtdps
41 ) const
42 {
43  PtrList<fvScalarMatrix> pEqnComps(phases.size());
44 
45  forAll(phases_, phasei)
46  {
47  phaseModel& phase = phases_[phasei];
48  const volScalarField& alpha = phase;
49  volScalarField& rho = phase.rho();
50 
51  pEqnComps.set(phasei, new fvScalarMatrix(p_rgh, dimVolume/dimTime));
52  fvScalarMatrix& pEqnComp = pEqnComps[phasei];
53 
54  // Density variation
55  if (!phase.isochoric() || !phase.pure())
56  {
57  pEqnComp +=
58  (
59  fvc::ddt(alpha, rho) + fvc::div(phase.alphaRhoPhi())
60  - fvc::Sp(fvc::ddt(alpha) + fvc::div(phase.alphaPhi()), rho)
61  )/rho;
62  }
63 
64  // Mesh dilatation correction
65  if (mesh.moving())
66  {
67  pEqnComp += fvc::div(mesh.phi())*alpha;
68  }
69 
70  // Compressibility
71  if (!phase.incompressible())
72  {
73  if (pimple.transonic())
74  {
75  const surfaceScalarField phid
76  (
77  IOobject::groupName("phid", phase.name()),
78  fvc::interpolate(phase.fluidThermo().psi())*phase.phi()
79  );
80 
81  pEqnComp +=
83  (
84  (alpha/rho)*
85  (
86  phase.fluidThermo().psi()*fvm::ddt(p_rgh)
87  + fvm::div(phid, p_rgh)
88  - fvm::Sp(fvc::div(phid), p_rgh)
89  )
90  );
91 
92  pEqnComps[phasei].relax();
93  }
94  else
95  {
96  pEqnComp +=
97  (alpha*phase.fluidThermo().psi()/rho)
99  }
100  }
101 
102  // Option sources
103  if (fvModels().addsSupToField(rho.name()))
104  {
105  pEqnComp -= fvModels().sourceProxy(alpha, rho, p_rgh)/rho;
106  }
107 
108  // Mass transfer
109  if (dmdts.set(phasei))
110  {
111  pEqnComp -= dmdts[phasei]/rho;
112  }
113  if (d2mdtdps.set(phasei))
114  {
115  pEqnComp -= correction(fvm::Sp(d2mdtdps[phasei]/rho, p_rgh));
116  }
117  }
118 
119  return pEqnComps;
120 }
121 
122 
123 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
static word groupName(Name name, const word &group)
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
bool transonic() const
Flag to indicate to solve using transonic algorithm.
const surfaceScalarField & phi() const
Return cell face motion fluxes.
tmp< fvMatrix< Type > > sourceProxy(const VolField< Type > &eqnField) const
Return source for an equation.
virtual bool addsSupToField(const word &fieldName) const
Return true if an fvModel adds a source term to the given.
Definition: fvModels.C:229
bool moving() const
Is mesh moving.
Definition: polyMesh.H:476
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
const fvMesh & mesh
Region mesh.
Definition: solver.H:101
volScalarField & p_rgh
Reference to the buoyant pressure for buoyant cases.
phaseSystem::phaseModelList & phases_
const phaseSystem::phaseModelList & phases
Reference to the phases.
Calculate the first temporal derivative.
Calculate the divergence of the given field.
Calculate the field for explicit evaluation of implicit and explicit sources.
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the divergence of the given field and flux.
Calculate the matrix for implicit and explicit sources.
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
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< Type > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
tmp< VolField< Type > > Sp(const volScalarField &sp, const VolField< Type > &vf)
Definition: fvcSup.C:67
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
Definition: fvcDiv.C:47
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const VolField< Type > &vf, const word &name)
Definition: fvmDiv.C:48
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const VolField< Type > &)
tmp< fvMatrix< Type > > ddt(const VolField< Type > &vf)
Definition: fvmDdt.C:46
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
SurfaceField< scalar > surfaceScalarField
const dimensionSet dimTime
const dimensionSet dimVolume
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.