Go to the documentation of this file.
43 incompressibleAdjointSolver,
70 if (adjointVars_.useSolverNameForFields())
73 <<
"useSolverNameForFields is set to true for adjointSolver "
74 << solverName() <<
nl <<
tab
75 <<
"Appending variable names with the solver name" <<
nl <<
tab
76 <<
"Please adjust the necessary entries in fvSchemes and fvSolution"
88 mag(
contErr)().weightedAverage(mesh_.V()).value();
96 <<
", cumulative = " << cumulativeContErr_
103 Foam::adjointSimple::adjointSimple
106 const word& managerType,
108 const word& primalSolverName
113 adjointVars_(allocateVars()),
114 cumulativeContErr_(
Zero),
115 adjointSensitivity_(
nullptr)
131 adjointVars_.paInst(),
132 solverControl_().dict(),
133 solverControl_().pRefCell(),
134 solverControl_().pRefValue()
137 if (computeSensitivities_)
142 adjointSensitivity_.reset
150 objectiveManagerPtr_(),
164 if (adjointSensitivity_.valid())
169 adjointSensitivity_().readDict
194 adjointVars_.adjointTurbulence();
195 const label& paRefCell = solverControl_().pRefCell();
196 const scalar& paRefValue = solverControl_().pRefValue();
207 fvOptionsAdjoint_(Ua)
215 objectiveManagerPtr_().addUaEqnSource(UaEqn);
218 ATCModel_->addATC(UaEqn);
221 addOptimisationTypeSource(UaEqn);
225 fvOptionsAdjoint_.constrain(UaEqn);
231 fvOptionsAdjoint_.correct(Ua);
246 if (solverControl_().consistent())
248 rAtUa = 1.0/(1.0/rAUa - UaEqn.H1());
260 while (solverControl_().correctNonOrthogonal())
273 if (solverControl_().finalNonOrthogonalIter())
275 phia = phiaHbyA - paEqn.
flux();
287 fvOptionsAdjoint_.correct(Ua);
293 if (solverControl_().printMaxMags())
297 Info<<
"Max mag of adjoint velocity = " << maxUa.
value() <<
endl;
298 Info<<
"Max mag of adjoint pressure = " << maxpa.
value() <<
endl;
301 solverControl_().
write();
304 adjointVars_.computeMeanFields();
316 adjointVars_.resetMeanFields();
319 while (solverControl_().loop())
329 return solverControl_().loop();
335 if (computeSensitivities_)
337 adjointSensitivity_->accumulateIntegrand(scalar(1));
338 const scalarField& sens = adjointSensitivity_->calculateSensitivities();
339 if (sensitivities_.empty())
343 sensitivities_.ref() = sens;
354 if (!sensitivities_.valid())
356 computeObjectiveSensitivities();
359 return sensitivities_();
365 if (computeSensitivities_)
367 adjointSensitivity_->clearSensitivities();
375 if (!adjointSensitivity_.valid())
378 <<
"Sensitivity object not allocated" <<
nl
379 <<
"Turn computeSensitivities on in "
385 return adjointSensitivity_();
394 objectiveManagerPtr_->updateAndWrite();
400 os.
writeEntry(
"averageIter", solverControl_().averageIter());
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
void continuityErrors()
Compute continuity errors.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
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)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
void boundaryManipulate(typename GeometricField< Type, fvPatchField, volMesh >::Boundary &values)
Manipulate based on a boundary field.
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.
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.
virtual void clearSensitivities()
Base class for incompressibleAdjoint solvers.
A class for handling words, derived from Foam::string.
virtual void solveIter()
Execute one iteration of the solution algorithm.
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,...
virtual void computeObjectiveSensitivities()
Compute sensitivities of the underlaying objectives.
virtual const incompressibleAdjointVars & getAdjointVars() const
Access to the incompressible adjoint variables set.
const Time & time() const
Return time.
fvMesh & mesh_
Reference to the mesh database.
static word timeName(const scalar t, const int precision=precision_)
virtual void correct()
Solve the adjoint turbulence equations.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Ostream & endl(Ostream &os)
Add newline and flush stream.
const Type & value() const
Return const reference to value.
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....
Class including all adjoint fields for incompressible flows.
autoPtr< objectiveManager > objectiveManagerPtr_
Object to manage objective functions.
Ostream & printExecutionTime(OSstream &os) const
Print the elapsed ExecutionTime (cpu-time), ClockTime.
messageStream Info
Information stream (uses stdout - output is on the master only)
virtual sensitivity & getSensitivityBase()
Return the base sensitivity object.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
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.
static autoPtr< adjointSensitivity > New(const fvMesh &mesh, const dictionary &dict, incompressibleVars &primalVars, incompressibleAdjointVars &adjointVars, objectiveManager &objectiveManager, fv::optionAdjointList &fvOptionsAdjoint)
Return a reference to the selected turbulence model.
SolverPerformance< Type > solve(faMatrix< Type > &, Istream &)
Solve returning the solution statistics given convergence tolerance.
const Type & lookupObject(const word &name, const bool recursive=false) const
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
virtual void updatePrimalBasedQuantities()
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual void clearSensitivities()
Clears the sensitivity field known by the adjoint solver.
Macros for easy insertion into run-time selection tables.
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
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.
virtual void solve()
Main control loop.
static autoPtr< SIMPLEControl > New(fvMesh &mesh, const word &managerType, const solver &solver)
Return a reference to the selected turbulence model.
errorManipArg< error, int > exit(error &err, const int errNo=1)
void correctBoundaryConditions()
Correct boundary field.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
virtual bool readDict(const dictionary &dict)
Read dict if updated.
autoPtr< variablesSet > vars_
Base variableSet pointer.
dimensioned< Type > weightedAverage(const DimensionedField< scalar, GeoMesh > &weightField) const
Calculate and return weighted average.
virtual void updatePrimalBasedQuantities()
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
virtual tmp< volVectorField > adjointMeanFlowSource()=0
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
incompressibleVars & primalVars_
Primal variable set.
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 writeData(Ostream &os) const
Write average iteration.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
virtual const scalarField & getObjectiveSensitivities()
Grab a reference to the computed sensitivities.
Abstract base class for adjoint sensitivities.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
virtual bool readDict(const dictionary &dict)
Read dict if updated.
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.
virtual bool loop()
Looper (advances iters, time step)
autoPtr< SIMPLEControl > solverControl_
Solver control.
static autoPtr< ATCModel > New(const fvMesh &mesh, const incompressibleVars &primalVars, const incompressibleAdjointVars &adjointVars, const dictionary &dict)
Return a reference to the selected turbulence model.
incompressibleAdjointVars & allocateVars()
virtual tmp< fvVectorMatrix > divDevReff(volVectorField &U) const =0
Return the diffusion term for the momentum equation.