https://mooseframework.inl.gov
TransientAndAdjoint.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 // MOOSE includes
11 #include "TransientAndAdjoint.h"
12 #include "OptimizationAppTypes.h"
13 #include "FEProblemBase.h"
14 
15 registerMooseObject("OptimizationApp", TransientAndAdjoint);
16 
19 {
22  params.addClassDescription("Executioner for evaluating transient simulations and their adjoint.");
23 
24  // We need to full matrix for the adjoint solve, so set this to NEWTON
25  params.set<MooseEnum>("solve_type") = "newton";
26  params.suppressParameter<MooseEnum>("solve_type");
27 
28  // The adjoint system (second one) is solved by _adjoint_solve
29  // This is a parameter of the MultiSystemSolveObject, which we set from here, the executioner.
30  // We seek to prevent the MultiSystemSolveObject from solving both systems
31  // This is abusing input parameters, but SolveObjects do not have their own syntax
32  // and we need to send this parameter from the executioner to the default nested SolveObject
33  params.renameParam("system_names", "forward_system", "");
34 
35  return params;
36 }
37 
39  : Transient(parameters),
40  _adjoint_solve(*this),
41  _forward_times(declareRecoverableData<std::vector<Real>>("forward_times"))
42 {
43  // Error out on unsupported time integration schemes
44  switch (getTimeScheme())
45  {
47  break;
48  default:
49  paramError(
50  "scheme", getParam<MooseEnum>("scheme"), " is not supported for computing adjoint.");
51  break;
52  }
53 }
54 
55 void
57 {
59 
60  // Save the forward initial condition
62 
63  // Save the time of the initial condition
64  // The vector of times should be reset unless we are recovering
65  if (!_app.isRecovering())
67  else
68  _forward_times.push_back(_time);
69 }
70 
71 void
73 {
75 
76  // Save the converged forward solution and time
77  if (lastSolveConverged())
78  {
80  _forward_times.push_back(_time);
81  }
82 }
83 
84 void
86 {
88 
89  // If it is a half transient, then the app is meant to be run with recovery. Therefore, it doesn't
90  // make sense to run the adjoint calculation since we aren't getting to the final time required
91  // for a consistent adjoint.
93  return;
94 
95  // Looping backward through forward time steps
96  for (const auto n : make_range(_forward_times.size() - 1))
97  {
98  // Set important time information so the Jacobian is properly computed in the forward system
99  _t_step = _forward_times.size() - 1 - n;
102  _dt = _time - _time_old;
104 
105  // Set the forward solution to the time step that is currently being solved
107 
108  // Incase there are some residual objects that need this call
111 
112  // Solve the adjoint system
114  if (!lastSolveConverged())
115  break;
116 
117  // FIXME: This works well enough for console and CSV output, but exodus is a mess since I don't
118  // think it knows what to do with backward timestepping or duplicate outputs at a given time
120  }
121 }
void renameParam(const std::string &old_name, const std::string &new_name, const std::string &new_docstring)
virtual void preExecute() override
void timestepSetup() override
virtual void postExecute() override
virtual void preExecute() override
For storing the initial condition.
bool & _last_solve_converged
TransientAndAdjoint(const InputParameters &parameters)
T & set(const std::string &name, bool quiet_mode=false)
virtual void postStep() override
For storing the converged solution.
Real & _time
static InputParameters validParams()
registerMooseObject("OptimizationApp", TransientAndAdjoint)
void suppressParameter(const std::string &name)
FEProblemBase & _problem
std::vector< Real > & _forward_times
Cached forward time points so we can properly loop backward in time.
TI_IMPLICIT_EULER
virtual void postExecute() override
For solving the adjoint problem.
static InputParameters validParams()
Real & _time_old
Real & _dt
virtual void postStep()
bool testCheckpointHalfTransient() const
AdjointTransientSolve _adjoint_solve
The transient adjoint solve object responsible for storing forward solutions and solving the adjoint ...
void paramError(const std::string &param, Args... args) const
void insertForwardSolution(int tstep)
This function should be called after every converged forward time step.
const ExecFlagType EXEC_ADJOINT_TIMESTEP_END
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
MooseApp & _app
void setForwardSolution(int tstep)
Takes the previously saved forward solutions residing in the adjoint system and copies them to the av...
Real & _dt_old
virtual void onTimestepBegin() override
IntRange< T > make_range(T beg, T end)
void addClassDescription(const std::string &doc_string)
static InputParameters validParams()
Moose::TimeIntegratorType getTimeScheme() const
virtual bool lastSolveConverged() const override
bool isRecovering() const
virtual bool solve() override
Overriding parent class so the time-derivative residual is stored for the next time step...
virtual void outputStep(ExecFlagType type)
int & _t_step