Go to the documentation of this file.
39 namespace incompressible
49 optimisationType::optimisationType
68 dict_.subDict(
"updateMethod").subOrEmptyDict(
"lineSearch"),
75 label nConstraints(0);
78 nConstraints += adjManagerI.nConstraints();
86 && !isA<constrainedOptimisationMethod>(updateMethod_())
89 const auto& cnstrTable =
90 *(constrainedOptimisationMethod::dictionaryConstructorTablePtr_);
94 <<
"Found " << nConstraints <<
" adjoint solvers corresponding to "
95 <<
"constraints but the optimisation method ("
96 << updateMethod_().type()
97 <<
") is not a constrainedOptimisationMethod." <<
nl
98 <<
"Available constrainedOptimisationMethods:" <<
nl
99 << cnstrTable.sortedToc()
105 && isA<constrainedOptimisationMethod>(updateMethod_())
110 <<
"Did not find any adjoint solvers corresponding to "
111 <<
"constraints but the optimisation method ("
112 << updateMethod_().type()
113 <<
") is a constrainedOptimisationMethod." <<
nl <<
nl
114 <<
"This can cause some constraintOptimisationMethods to misbehave."
116 <<
"Either the isConstraint bool is not set in one of the adjoint "
117 <<
"solvers or you should consider using an updateMethod "
118 <<
"that is not a constrainedOptimisationMethod"
133 const word modelType(
dict.subDict(
"optimisationType").get<
word>(
"type"));
135 Info<<
"optimisationType type : " << modelType <<
endl;
137 auto* ctorPtr = dictionaryConstructorTable(modelType);
146 *dictionaryConstructorTablePtr_
194 scalar objectiveValue(
Zero);
228 scalar& objectiveValue,
234 const scalar opWeight = adjSolvManager.operatingPointWeight();
238 adjSolvManager.aggregateSensitivities();
240 if (objectiveSens.empty())
242 objectiveSens.setSize(tadjointSolverManagerSens().size(),
Zero);
245 objectiveSens += opWeight*tadjointSolverManagerSens();
246 objectiveValue += opWeight*adjSolvManager.objectiveValue();
250 adjSolvManager.constraintSensitivities();
254 if (constraintSens.empty())
256 constraintSens.
setSize(adjointSolverManagerConstSens.size());
257 forAll(constraintSens, cI)
264 adjointSolverManagerConstSens[cI].size(),
268 constraintValues.setSize(cValues().size());
269 constraintValues =
Zero;
273 forAll(constraintSens, cI)
275 constraintSens[cI] += opWeight*adjointSolverManagerConstSens[cI];
277 constraintValues += opWeight*cValues();
286 scalar objectiveValue(
Zero);
291 const scalar opWeight = adjSolvManager.operatingPointWeight();
293 objectiveValue += opWeight*adjSolvManager.objectiveValue();
296 if (constraintValues.empty())
298 constraintValues.setSize(cValues().size(),
Zero);
300 constraintValues += opWeight*cValues();
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
autoPtr< lineSearch > & getLineSearch()
Get a reference to the line search object.
A class for handling words, derived from Foam::string.
defineTypeNameAndDebug(adjointEikonalSolver, 0)
virtual scalar meritFunctionDirectionalDerivative()
Derivative of the merit function.
virtual void write()
Write useful quantities to files.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
const autoPtr< volScalarField > & sourcePtr()
Get source term.
autoPtr< volScalarField > sourcePtr_
Ostream & endl(Ostream &os)
Add newline and flush stream.
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
const T * set(const label i) const
Return const pointer to element (can be nullptr),.
virtual tmp< scalarField > computeDirection()
Compute update direction.
defineRunTimeSelectionTable(adjointSensitivity, dictionary)
#define forAll(list, i)
Loop across all elements in list.
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
virtual void updateDesignVariables(scalarField &correction)=0
Update the design variables given their correction.
messageStream Info
Information stream (stdout output on master, null elsewhere)
static autoPtr< lineSearch > New(const dictionary &dict, const Time &time)
Return a reference to the selected turbulence model.
static autoPtr< updateMethod > New(const fvMesh &mesh, const dictionary &dict)
Return a reference to the selected turbulence model.
void setSize(const label newLen)
Same as resize()
PtrList< adjointSolverManager > & adjointSolverManagers
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Class for managing adjoint solvers, which may be more than one per operating point.
PtrList< adjointSolverManager > & adjointSolvManagers_
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Mesh data needed to do the Finite Volume discretisation.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
autoPtr< updateMethod > updateMethod_
Macros to ease declaration of run-time selection tables.
virtual void computeEta(scalarField &correction)=0
Compute eta if not set in the first step.
virtual void updateGradientsAndValues(scalarField &objectiveSens, PtrList< scalarField > &constraintSens, scalar &objectiveValue, scalarField &constraintValues)
Compute cumulative objective and constraint gradients.
static autoPtr< optimisationType > New(fvMesh &mesh, const dictionary &dict, PtrList< adjointSolverManager > &adjointSolverManagers)
Return a reference to the selected turbulence model.
virtual void update()
Update design variables.
#define WarningInFunction
Report a warning using Foam::Warning.
virtual void updateOldCorrection(const scalarField &)
Update old correction. Needed for quasi-Newton Methods.
autoPtr< lineSearch > lineSearch_
virtual scalar computeMeritFunction()
Compute the merit function of the optimisation problem.