momentumPredictor.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 
27 #include "fvcDdt.H"
28 #include "fvcSnGrad.H"
29 #include "fvcReconstruct.H"
30 #include "fvmDiv.H"
31 #include "fvmDdt.H"
32 
33 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
34 
36 {
38 
39  tUcEqn =
40  (
43  + momentumTransport->divDevSigma(Uc)
44  ==
46  );
47  fvVectorMatrix& UcEqn = tUcEqn.ref();
48 
49  UcEqn.relax();
50 
51  fvConstraints().constrain(UcEqn);
52 
54  {
55  // Face buoyancy force
56  const surfaceScalarField Fgf(g & mesh.Sf());
57 
58  solve
59  (
60  UcEqn
61  ==
63  (
64  Fgf - fvc::snGrad(p)*mesh.magSf()
65  - Dcf()*(phic - phid())
66  )
67  + Dc()*fvc::reconstruct(phic - phid())
68  + Fd() - fvm::Sp(Dc(), Uc)
69  );
70 
72  }
73 }
74 
75 
76 // ************************************************************************* //
Generic GeometricField class.
bool momentumPredictor() const
Flag to indicate to solve for momentum.
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:603
const surfaceVectorField & Sf() const
Return cell face area vectors.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
tmp< fvMatrix< Type > > source(const VolField< Type > &field) const
Return source for an equation.
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
autoPtr< surfaceScalarField > phid
Effective volumetric flux of the dispersed phase.
autoPtr< phaseIncompressible::momentumTransportModel > momentumTransport
Pointer to the momentum transport model.
uniformDimensionedVectorField g
Acceleration due to gravity.
const volVectorField & Uc
Reference to the continuous phase velocity field.
const volScalarField & alphac
Reference continuous phase-fraction.
volVectorField Uc_
Continuous phase velocity field.
surfaceScalarField alphaPhic
Continuous phase volumetric-flux field.
tmp< fvVectorMatrix > tUcEqn
Cached momentum matrix.
autoPtr< volVectorField > Fd
Dispersed phase drag force.
virtual void momentumPredictor()
Construct and optionally solve the momentum equation.
autoPtr< surfaceScalarField > Dcf
Face continuous-dispersed phase drag coefficient.
const surfaceScalarField & phic
Reference to the continuous phase volumetric-flux field.
const volScalarField & p
Reference to the pressure field.
Calculate the first temporal derivative.
Reconstruct volField from a face flux field.
Calculate the snGrad of the given volField.
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the divergence of the given field and flux.
tmp< VolField< Type > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
tmp< VolField< typename outerProduct< vector, Type >::type > > reconstruct(const SurfaceField< Type > &ssf)
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
Definition: fvcDiv.C:47
tmp< SurfaceField< Type > > snGrad(const VolField< Type > &vf, const word &name)
Definition: fvcSnGrad.C:45
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
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.