Go to the documentation of this file.
57 if (incoVars_.useSolverNameForFields())
60 <<
"useSolverNameForFields is set to true for primalSolver "
61 << solverName() <<
nl <<
tab
62 <<
"Appending variable names with the solver name" <<
nl <<
tab
63 <<
"Please adjust the necessary entries in fvSchemes and fvSolution"
75 mag(
contErr)().weightedAverage(mesh_.V()).value();
83 <<
", cumulative = " << cumulativeContErr_
93 const word& managerType,
99 incoVars_(allocateVars()),
101 cumulativeContErr_(
Zero),
112 solverControl_().dict(),
113 solverControl_().pRefCell(),
114 solverControl_().pRefValue()
141 incoVars_.turbulence();
142 label&
pRefCell = solverControl_().pRefCell();
143 scalar&
pRefValue = solverControl_().pRefValue();
148 MRF_.correctBoundaryVelocity(
U);
160 addOptimisationTypeSource(
UEqn);
164 fvOptions_().constrain(
UEqn);
170 fvOptions_().correct(
U);
184 if (solverControl_().consistent())
198 while (solverControl_().correctNonOrthogonal())
209 if (solverControl_().finalNonOrthogonalIter())
222 U.correctBoundaryConditions();
223 fvOptions_().correct(
U);
226 incoVars_.laminarTransport().correct();
229 solverControl_().write();
235 Info<< obj->objectiveName() <<
" : " << obj->J() <<
endl;
236 obj->accumulateJMean(solverControl_());
237 obj->writeInstantaneousValue();
241 incoVars_.computeMeanFields();
254 if (objectives_.empty())
256 objectives_ = getObjectiveFunctions();
261 incoVars_.resetMeanFields();
264 incoVars_.turbulence()->validate();
266 while (solverControl_().loop())
279 return solverControl_().loop();
285 incoVars_.restoreInitValues();
291 os.
writeEntry(
"averageIter", solverControl_().averageIter());
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
void constrainPressure(volScalarField &p, const RhoType &rho, const volVectorField &U, const surfaceScalarField &phiHbyA, const RAUType &rhorAU, const MRFType &MRF)
tmp< GeometricField< Type, fvPatchField, volMesh > > H() const
Return the H operation source.
tmp< volVectorField > constrainHbyA(const tmp< volVectorField > &tHbyA, const volVectorField &U, const volScalarField &p)
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
volScalarField rAtU(1.0/(1.0/rAU - UEqn.H1()))
volScalarField::Internal contErr((fvc::ddt(rho)+fvc::div(rhoPhi) -(fvOptions(alpha1, mixture.thermo1().rho())&rho1) -(fvOptions(alpha2, mixture.thermo2().rho())&rho2))())
bool adjustPhi(surfaceScalarField &phi, const volVectorField &U, volScalarField &p)
Adjust the balance of fluxes to obey continuity.
const dictionary & simple
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
const bool momentumPredictor
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
tmp< volScalarField > rAU
A class for handling words, derived from Foam::string.
autoPtr< SIMPLEControl > solverControl_
Solver control.
Info<< "Reading field U\n"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);volScalarField rho(IOobject("rho", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), thermo.rho());volVectorField rhoU(IOobject("rhoU", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *U);volScalarField rhoE(IOobject("rhoE", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *(e+0.5 *magSqr(U)));surfaceScalarField pos(IOobject("pos", runTime.timeName(), mesh), mesh, dimensionedScalar("pos", dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar("neg", dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\n"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
Find the reference cell nearest (in index) to the given cell but which is not on a cyclic,...
const Time & time() const
Return time.
fvMesh & mesh_
Reference to the mesh database.
static word timeName(const scalar t, const int precision=precision_)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Ostream & endl(Ostream &os)
Add newline and flush stream.
void setReference(const label celli, const Type &value, const bool forceReference=false)
Set reference level for solution.
For cases which do no have a pressure boundary adjust the balance of fluxes to obey continuity....
void continuityErrors()
Compute continuity errors.
virtual void solveIter()
Execute one iteration of the solution algorithm.
Ostream & printExecutionTime(OSstream &os) const
Print the elapsed ExecutionTime (cpu-time), ClockTime.
tmp< fvVectorMatrix > tUEqn(fvm::ddt(rho, U)+fvm::div(phi, U)+MRF.DDt(rho, U)+turbulence->divDevRhoReff(U)==fvOptions(rho, U))
virtual bool writeData(Ostream &os) const
Write average iteration.
messageStream Info
Information stream (uses stdout - output is on the master only)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< volScalarField > H1() const
Return H(1)
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
virtual bool readDict(const dictionary &dict)
Read dict if updated.
bool setRefCell(const volScalarField &field, const volScalarField &fieldRef, const dictionary &dict, label &refCelli, scalar &refValue, const bool forceReference=false)
If the field fieldRef needs referencing find the reference cell nearest.
SolverPerformance< Type > solve(faMatrix< Type > &, Istream &)
Solve returning the solution statistics given convergence tolerance.
virtual void restoreInitValues()
Restore initial field values if necessary.
Base class for primal incompressible solvers.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
tmp< volScalarField > A() const
Return the central coefficient.
virtual bool loop()
Looper (advances iters, time step)
Mesh data needed to do the Finite Volume discretisation.
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
dictionary subOrEmptyDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX, const bool mandatory=false) const
static autoPtr< SIMPLEControl > New(fvMesh &mesh, const word &managerType, const solver &solver)
Return a reference to the selected turbulence model.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
autoPtr< variablesSet > vars_
Base variableSet pointer.
dimensioned< Type > weightedAverage(const DimensionedField< scalar, GeoMesh > &weightField) const
Calculate and return weighted average.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
incompressibleVars & allocateVars()
Protected Member Functions.
void relax(const scalar alpha)
Relax field (for steady-state solution).
A special matrix type and solver, designed for finite volume solutions of scalar equations....
SolverPerformance< Type > solve(const dictionary &)
Solve returning the solution statistics.
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
virtual bool readDict(const dictionary &dict)
Read dict if updated.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
const incompressibleVars & getIncoVars() const
Access to the incompressible variables set.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Abstract base class for objective functions. No point in making this runTime selectable since its chi...
List of finite volume options.
defineTypeNameAndDebug(combustionModel, 0)
#define WarningInFunction
Report a warning using Foam::Warning.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux() const
Return the face-flux field from the matrix.
Base class for solution control classes.
virtual void solve()
Main control loop.