Differentiating solution at a point with respect to parameters in PDE

[YanniPapandreou]

Hi, firstly thanks for maintaining this package, I have been finding it very useful!

I have been trying to use Autodiff to compute the sensitivities of the solution of a pde at a point with respect to parameters in the pde. For example say I have this elliptic pde:

$$ \begin{align} -\nabla(\kappa\nabla u) &= f , \quad x\in D\\ u &= 0, \quad x \in \partial D \end{align} $$

and I want to compute the derivative of the solution at some point with respect to the parameter $\kappa$, i.e. I want:

$$ \frac{\partial u}{\partial \kappa}(x_{0}) $$

for fixed $x_0$. Is it possible to do this with auto-diff in Gridap? I have implemented the following code to get the value of $u$ at an arbitrary fixed point:

using Gridap, Plots, Zygote
import Gridap: ∇

f(x) = 1.0

function U(k)
    order = 1
    n = 20
    domain = (0,1)
    partition = (n)
    model = CartesianDiscreteModel(domain, partition)

    reffe = ReferenceFE(lagrangian,Float64,order)
    V0 = TestFESpace(model,reffe,conformity=:H1,dirichlet_tags="boundary")
    U = TrialFESpace(V0,x->0)

    degree = order + 1

    Ω = Triangulation(model)
    dΩ = Measure(Ω,degree)
    
    a(u,v) = ∫( k * ∇(v)⊙∇(u) )*dΩ
    b(v) = ∫( v*f )*dΩ
    
    op = AffineFEOperator(a,b,U,V0)
    
    uh = solve(op)

    return uh(Point(0.5))
end

and I have tried using both Zygote.jl and ForwardDiff.jl to compute the derivative wrt k. Zygote just hangs and seems to not terminate and ForwardDiff produces the following error on running ForwardDiff.derivative(U,1.0):

MethodError: no method matching Float64(::ForwardDiff.Dual{ForwardDiff.Tag{typeof(U), Float64}, Float64, 1})

I am not sure if I have missed something in the docs. Any help would be much appreciated :)