Differentiating solution at a point with respect to parameters in PDE
- Dominant language
- Julia
- Stars
- 879
- Forks
- 119
- Avg merge
- 6d 3h
- Merged PRs (30d)
- 4
Description
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:
```julia
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 :)
Contributor guide
Research direction
Start with the supplied U(k) example, especially the Float64 conversion reported by ForwardDiff.derivative(U, 1.0), and compare it with the Zygote behavior. Determine whether Gridap can support differentiating the point evaluation with respect to k; done means the example returns the requested derivative without hanging or the reported MethodError.
Written by the indexing model from the issue text.
Assessment
- Tech stack
- julia
- Domain
- backend
- Issue type
- Feature
- Difficulty
- 5/5
- Estimated time
- Over a week
- Activity status
- Stale
- Clarity
- Needs clarification
- Newbie friendliness
- 25/100