thermophysicalPredictor.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 "film.H"
27 #include "fvmDdt.H"
28 #include "fvmDiv.H"
29 
30 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
31 
33 {
35 
36  fvScalarMatrix heEqn
37  (
39  - fvm::Sp(contErr(), he)
41  ==
42  fvModels().source(alpha, rho, he)
43  );
44 
45  heEqn.relax();
46 
47  fvConstraints().constrain(heEqn);
48 
49  heEqn.solve();
50 
52 
53  thermo_.correct();
54 }
55 
56 
57 // ************************************************************************* //
Generic GeometricField class.
virtual void correct()=0
Update properties.
virtual const volScalarField & he() const =0
Enthalpy/Internal energy [J/kg].
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
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:58
Finite volume models.
Definition: fvModels.H:65
Foam::fvConstraints & fvConstraints() const
Return the fvConstraints that are created on demand.
Definition: solver.C:107
virtual void thermophysicalPredictor()
Construct and solve the energy equation,.
autoPtr< filmThermophysicalTransportModel > thermophysicalTransport
Definition: film.H:70
rhoFluidThermo & thermo_
Reference to the fluid thermophysical properties.
const volScalarField & alpha
Film volume fraction in the cell layer.
tmp< volScalarField::Internal > contErr
Continuity error.
const surfaceScalarField & alphaRhoPhi
Reference to the film mass-flux field.
const volScalarField & rho
Reference to the thermodynamic density field.
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the divergence of the given field and flux.
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
thermo he()