Evaluation Hessian of a Control paramater + Eigendecomposition

[casanovamaxime]

Hello everyone,
I have a question regarding the estimation of the hessian of a scalar cost function and its eigendecomposition,
To be more clear, lets define:

• The scalar cost function we want minimize J =0.5*(true_data - firedrake_data)**2
• The control parameter we optimize to fit the true data f
  So far i used the adjoint package to optimize the control f

Now i would like the compute the hessian of J w.r.t. f (dJ^2/df^2)

But since J do not depend explicitely of f i cant use the automatic differentiation, so is there a way to estimate such hessian ?

best,

[stephankramer]

What does your control space look like? For the general case where your control is a Function you generally don't want to construct the entire Hessian, as it would be too costly. pyadjoint does provide the functionality to efficiently calculate the Hessian in a certain direction, i.e. the Hessian matrix multiplied by a given vector: http://www.dolfin-adjoint.org/en/release/documentation/pyadjoint_api.html#pyadjoint.compute_hessian
This would be enough to iteratively approximate the largest eigenvalues, e.g. using slepc, although I don't ahve experience with that. For a small control space, say a small number of Constants, you could also use this to reconstruct the entire Hessian, by simply "applying the Hessian" in all standard basis directions.

[casanovamaxime]

Indeed the control is a Function, and it was exactly what i wanted to do. The fact is i do not see how to combine compute_hessian to slepc in order to solve H v = lambda v there are not a lot of documentation in Firedrake regarding that..

[dham]

The operation you are trying to do is above the level of Firedrake. Firedrake (and pyadjoint) can give you the hessian action, as described above. You need to wrap this in a PETSc matrix object in order to then have slepc do a matrix-free Eigen problem. You can take inspiration for the construction of the matrix from the way that we construct matrix-free petsc Mats using a form action. This is described in the documentation here: https://www.firedrakeproject.org/matrix-free.html#matrix-free-operators and that links to the code for ImplicitMatrix and ImplicitMatrix context. Note that you can't use that code directly because it's for the Firedrake case where the matrix action is given by assembling the action of a form. However the form of the code will be similar.

Similarly, you can take inspiration from Firedrakes LinearEigenSolver package to see how to put together a slepc eigensolver at Python level: https://www.firedrakeproject.org/_modules/firedrake/eigensolver.html#LinearEigensolver

It would be possible to wrap up all this into a higher level interface that took in a ReducedFunctional and solves the associated eigenproblem, but nobody has done that yet. If you'd like to contribute that to Firedrake we're more than happy to help.