Go to the documentation of this file.
43 Foam::optimisationManager::optimisationManager(
fvMesh&
mesh)
58 time_(const_cast<
Time&>(
mesh.time())),
60 adjointSolverManagers_(),
61 managerType_(get<
word>(
"optimisationManager")),
65 const wordList& primalSolverNames = primalSolversDict.
toc();
78 primalSolversDict.
subDict(primalSolverNames[solveri])
85 const wordList& adjointManagerNames = adjointManagersDict.
toc();
88 label nAdjointSolvers(0);
98 adjointManagersDict.
subDict(adjointManagerNames[manageri])
109 if (!solveri.useSolverNameForFields())
112 <<
"Multiple primal solvers are present but "
113 <<
"useSolverNameForFields is set to false in "
114 <<
"primal solver " << solveri.solverName() <<
nl
115 <<
"This is considered fatal."
121 if (nAdjointSolvers > 1)
128 if (!asI.useSolverNameForFields())
131 <<
"Multiple adjoint solvers are present but "
132 <<
"useSolverNameForFields is set to false in "
133 <<
"adjoint solver " << asI.solverName() <<
nl
134 <<
"This is considered fatal."
165 Info<<
"optimisationManager type : " << modelType <<
endl;
167 auto cstrIter = dictionaryConstructorTablePtr_->cfind(modelType);
169 if (!cstrIter.found())
174 "optimisationManager",
176 *dictionaryConstructorTablePtr_
191 const dictionary& primalSolversDict = subDict(
"primalSolvers");
194 sol.readDict(primalSolversDict.
subDict(sol.solverName()));
197 const dictionary& adjointManagersDict = subDict(
"adjointManagers");
200 man.readDict(adjointManagersDict.
subDict(man.managerName()));
212 return primalSolvers_;
219 return adjointSolverManagers_;
226 forAll(primalSolvers_, psI)
228 primalSolvers_[psI].solve();
236 forAll(adjointSolverManagers_, amI)
238 adjointSolverManagers_[amI].solveAdjointEquations();
246 forAll(adjointSolverManagers_, amI)
248 adjointSolverManagers_[amI].computeAllSensitivities();
255 forAll(adjointSolverManagers_, amI)
258 adjointSolverManagers_[amI].adjointSolvers();
260 forAll(adjointSolvers, asI)
262 adjointSolvers[asI].updatePrimalBasedQuantities();
Base class for adjoint solvers.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
PtrList< adjointSolverManager > adjointSolverManagers_
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
virtual void updatePrimalBasedQuantities()
Solve all primal equations.
int system(const std::string &command, const bool bg=false)
Execute the specified command via the shell.
virtual bool read()
Read object.
static autoPtr< primalSolver > New(fvMesh &mesh, const word &managerType, const dictionary &dict)
Return a reference to the selected primal solver.
Ostream & endl(Ostream &os)
Add newline and flush stream.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
#define forAll(list, i)
Loop across all elements in list.
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
virtual bool read()
Read object.
PtrList< primalSolver > primalSolvers_
messageStream Info
Information stream (uses stdout - output is on the master only)
virtual void solvePrimalEquations()
Solve all primal equations.
virtual PtrList< primalSolver > & primalSolvers()
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
static autoPtr< optimisationManager > New(fvMesh &mesh)
Return a reference to the selected turbulence model.
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.
virtual void solveAdjointEquations()
Solve all adjoint equations.
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)
virtual void computeSensitivities()
Compute all adjoint sensitivities.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
const word & system() const
Return system name.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Base class for primal solvers.
const Time & time() const
Return the top-level database.
wordList toc() const
Return the table of contents.
virtual PtrList< adjointSolverManager > & adjointSolverManagers()
defineTypeNameAndDebug(combustionModel, 0)