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)
138 createZeroFieldPtr<vector>
157 createZeroFieldPtr<scalar>
176 createZeroFieldPtr<scalar>
195 createZeroFieldPtr<scalar>
198 (
"dJdTMvar1_"+
type()),
214 createZeroFieldPtr<scalar>
217 (
"dJdTMvar2_"+
type()),
231 if (bdJdvPtr_.empty())
233 bdJdvPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
235 return bdJdvPtr_()[patchI];
244 if (bdJdvnPtr_.empty())
246 bdJdvnPtr_.reset(createZeroBoundaryPtr<scalar>(mesh_));
248 return bdJdvnPtr_()[patchI];
257 if (bdJdvtPtr_.empty())
259 bdJdvtPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
261 return bdJdvtPtr_()[patchI];
270 if (bdJdpPtr_.empty())
272 bdJdpPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
274 return bdJdpPtr_()[patchI];
283 if (bdJdTPtr_.empty())
285 bdJdTPtr_.reset(createZeroBoundaryPtr<scalar>(mesh_));
287 return bdJdTPtr_()[patchI];
296 if (bdJdTMvar1Ptr_.empty())
298 bdJdTMvar1Ptr_.reset(createZeroBoundaryPtr<scalar>(mesh_));
300 return bdJdTMvar1Ptr_()[patchI];
309 if (bdJdTMvar2Ptr_.empty())
311 bdJdTMvar2Ptr_.reset(createZeroBoundaryPtr<scalar>(mesh_));
313 return bdJdTMvar2Ptr_()[patchI];
autoPtr< boundaryScalarField > bdJdTPtr_
Adjoint outlet temperature.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
virtual void nullify()
Nullify adjoint contributions.
virtual void write() const
Write objective function history.
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.
bool hasBoundarydJdp() const
virtual void update_boundarydJdb()
Update objective function derivative term.
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
virtual void update_boundarydJdTMvar2()
autoPtr< boundaryScalarField > bdJdTMvar2Ptr_
Adjoint outlet turbulence model var 2.
virtual void update_dJdp()
virtual void update_gradDxDbMultiplier()
Update grad( dx/db multiplier). Volume-based sensitivity term.
autoPtr< boundaryScalarField > bdJdvnPtr_
Adjoint outlet pressure.
virtual void update_dJdTMvar1()
autoPtr< boundaryVectorField > bdJdvPtr_
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
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()
Update d (x) / db multiplier. Surface-based sensitivity term.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
virtual void update_dSdbMultiplier()
Update d (normal dS) / db multiplier. Surface-based sensitivity term.
const boundaryVectorField & boundarydJdv()
Objective partial deriv wrt velocity for all patches.
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
virtual void update_dndbMultiplier()
Update d (normal) / db multiplier. Surface-based sensitivity term.
autoPtr< volScalarField > dJdpPtr_
messageStream Info
Information stream (uses stdout - output is on the master only)
virtual void update_boundaryEdgeContribution()
Update boundary edge contributions.
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_
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()
Update div( dx/db multiplier). Volume-based sensitivity term.
autoPtr< volScalarField > dJdTMvar2Ptr_
Second turbulence model variable.
static const Vector< scalar > zero
const volScalarField & dJdTMvar2()
Contribution to field adjoint turbulence model variable 2.
virtual void update_dJdStressMultiplier()
Update dJ/dStress field.
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 write() const
Write objective function history.