Go to the documentation of this file.
47 objectiveIncompressible,
53 objectiveIncompressible::objectiveIncompressible
57 const word& adjointSolverName,
58 const word& primalSolverName
77 dJdTMvar1Ptr_(
nullptr),
78 dJdTMvar2Ptr_(
nullptr),
86 bdJdTMvar1Ptr_(
nullptr),
87 bdJdTMvar2Ptr_(
nullptr)
89 weight_ =
dict.
get<scalar>(
"weight");
90 computeMeanFields_ = vars_.computeMeanFields();
100 const word& adjointSolverName,
101 const word& primalSolverName
107 <<
" of type " << modelType <<
endl;
109 auto cstrIter = dictionaryConstructorTablePtr_->cfind(modelType);
111 if (!cstrIter.found())
116 "objectiveIncompressible",
118 *dictionaryConstructorTablePtr_
124 cstrIter()(
mesh,
dict, adjointSolverName, primalSolverName)
139 dJdvPtr_().primitiveFieldRef() *= oneOverNorm;
143 dJdpPtr_().primitiveFieldRef() *= oneOverNorm;
147 dJdTPtr_().primitiveFieldRef() *= oneOverNorm;
187 objective::doNormalization();
199 createZeroFieldPtr<vector>
218 createZeroFieldPtr<scalar>
237 createZeroFieldPtr<scalar>
256 createZeroFieldPtr<scalar>
259 (
"dJdTMvar1_"+
type()),
275 createZeroFieldPtr<scalar>
278 (
"dJdTMvar2_"+
type()),
292 if (bdJdvPtr_.empty())
294 bdJdvPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
296 return bdJdvPtr_()[patchI];
305 if (bdJdvnPtr_.empty())
307 bdJdvnPtr_.reset(createZeroBoundaryPtr<scalar>(mesh_));
309 return bdJdvnPtr_()[patchI];
318 if (bdJdvtPtr_.empty())
320 bdJdvtPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
322 return bdJdvtPtr_()[patchI];
331 if (bdJdpPtr_.empty())
333 bdJdpPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
335 return bdJdpPtr_()[patchI];
344 if (bdJdTPtr_.empty())
346 bdJdTPtr_.reset(createZeroBoundaryPtr<scalar>(mesh_));
348 return bdJdTPtr_()[patchI];
357 if (bdJdTMvar1Ptr_.empty())
359 bdJdTMvar1Ptr_.reset(createZeroBoundaryPtr<scalar>(mesh_));
361 return bdJdTMvar1Ptr_()[patchI];
370 if (bdJdTMvar2Ptr_.empty())
372 bdJdTMvar2Ptr_.reset(createZeroBoundaryPtr<scalar>(mesh_));
374 return bdJdTMvar2Ptr_()[patchI];
480 update_boundaryEdgeContribution();
481 update_dJdStressMultiplier();
544 objective::nullify();
autoPtr< boundaryScalarField > bdJdTPtr_
Adjoint outlet temperature.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
autoPtr< scalar > normFactor_
Normalization factor.
const volVectorField & dJdv()
Contribution to field adjoint momentum eqs.
virtual scalar J()=0
Return the objective function value.
A class for handling words, derived from Foam::string.
virtual void update_boundarydJdTMvar1()
bool hasBoundarydJdvn() const
static constexpr const zero Zero
Global zero (0)
bool hasBoundarydJdp() const
virtual void update_boundarydJdb()
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
virtual void update_boundarydJdTMvar2()
autoPtr< boundaryScalarField > bdJdTMvar2Ptr_
Adjoint outlet turbulence model var 2.
virtual void update_dJdp()
virtual void update_gradDxDbMultiplier()
autoPtr< boundaryScalarField > bdJdvnPtr_
Adjoint outlet pressure.
virtual void update_dJdTMvar1()
autoPtr< boundaryVectorField > bdJdvPtr_
bool valid() const noexcept
True if the managed pointer is non-null.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const volScalarField & dJdTMvar1()
Contribution to field adjoint turbulence model variable 1.
virtual void update_dJdb()
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
virtual bool write(const bool valid=true) const
Write objective function history.
Dimension set for the base types.
autoPtr< boundaryVectorField > bdJdvtPtr_
Adjoint outlet velocity.
bool hasBoundarydJdv() const
virtual void update_dJdv()
Update vol and boundary fields and derivatives.
virtual void update_dxdbMultiplier()
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
virtual void update_dSdbMultiplier()
const boundaryVectorField & boundarydJdv()
Objective partial deriv wrt velocity for all patches.
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
virtual void update_dndbMultiplier()
autoPtr< volScalarField > dJdpPtr_
messageStream Info
Information stream (uses stdout - output is on the master only)
virtual void update_boundarydJdv()
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual void update_dJdTMvar2()
const Type & lookupObject(const word &name, const bool recursive=false) const
static autoPtr< objectiveIncompressible > New(const fvMesh &mesh, const dictionary &dict, const word &adjointSolverName, const word &primalSolverName)
Return a reference to the selected turbulence model.
virtual void update_dxdbDirectMultiplier()
word dictName() const
The local dictionary name (final part of scoped name)
Base class for primal incompressible solvers.
virtual void update_boundarydJdT()
bool hasBoundarydJdT() const
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
autoPtr< boundaryScalarField > bdJdTMvar1Ptr_
Adjoint outlet turbulence model var 1.
Macros for easy insertion into run-time selection tables.
virtual void update_boundarydJdp()
const volScalarField & dJdp()
Contribution to field adjoint continuity eq.
const boundaryScalarField & boundarydJdT()
Objective partial deriv wrt temperature for all patches.
autoPtr< volScalarField > dJdTPtr_
Mesh data needed to do the Finite Volume discretisation.
virtual void update_boundarydJdvn()
autoPtr< volScalarField > dJdTMvar1Ptr_
First turbulence model variable.
errorManipArg< error, int > exit(error &err, const int errNo=1)
const volScalarField & dJdT()
Contribution to field adjoint energy eq.
const boundaryScalarField & boundarydJdvn()
Objective partial deriv wrt normal velocity for all patches.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
bool hasBoundarydJdTMVar2() const
autoPtr< volVectorField > dJdvPtr_
void write(Ostream &os, const bool subDict=true) const
Write dictionary, normally with sub-dictionary formatting.
virtual void nullify()
Update objective function derivatives.
virtual void update_meanValues()
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
bool hasdJdTMVar2() const
const boundaryScalarField & boundarydJdTMvar2()
Objective partial deriv wrt turbulence model var 2 for all patches.
virtual void update_divDxDbMultiplier()
autoPtr< volScalarField > dJdTMvar2Ptr_
Second turbulence model variable.
static const Vector< scalar > zero
const volScalarField & dJdTMvar2()
Contribution to field adjoint turbulence model variable 2.
const boundaryVectorField & boundarydJdp()
Objective partial deriv wrt pressure (times normal) for all patches.
bool hasBoundarydJdTMVar1() const
virtual void update_boundarydJdvt()
const boundaryScalarField & boundarydJdTMvar1()
Objective partial deriv wrt turbulence model var 1 for all patches.
virtual void update()
Update objective function derivatives.
bool hasBoundarydJdvt() const
autoPtr< boundaryVectorField > bdJdpPtr_
Adjoint (intlet,wall) velocity.
bool hasdJdv() const
Inline functions for checking whether pointers are set or not.
const boundaryVectorField & boundarydJdvt()
Objective partial deriv wrt tangent velocity for all patches.
Abstract base class for objective functions. No point in making this runTime selectable since its chi...
virtual void update_dJdT()
defineTypeNameAndDebug(combustionModel, 0)
bool hasdJdTMVar1() const
virtual void doNormalization()
Normalize all fields allocated by the objective.