Calculating Jacobian from two matrices

Specific Domains Modelling & Simulations
, ,
1

Hi everyone,

I would like to compute a Jacobian for one of my problem for a PDE, but I don’t have a function linking the unknown to my function. So I only have 2 matrices, 1 with the values of the function and 1 with the values of the unknowns.

From my shallow understanding of the problem, I cannot use automatic differentiation, because I would need the expression linking the function with the unknowns. Should I then go for finite differences to obtain a pure numerical jacobian?

I have this simple flux function, in 1D to simplify:

f(C, x, t) = v (C, x, t) C(x, t)

and I need to compute:

\frac{\partial f_i}{\partial C_i} = v_i + C_i \frac{\partial v_i}{\partial C_i}

which is a simple diagonal matrix, but I don’t know the relationship betweem v and C, so also \frac{\partial v}{\partial C}.

Is there an easy way to compute it directly without this relationship? For example with FiniteDiff.jl?

2

I’m not sure I understand your problem perfectly, but it sounds like you can use automatic differentiation.
It is not necessary to work with a “function” in the mathematical sense, that is, an explicit formula linking input to output. All you need is a “function” in the programming sense, that is, an algorithm computing the output from the input. It can use for loops, if/else statements and pretty much anything else you need.
However, it is unclear from your question if C is a function, an array or both?

3

Thx for your answer. I will try to clarify a bit my question with an example:

My problem is that I don’t know the relationship between the terms in the function I want to calculate the Jacobian from.

So, for example, for a function x * sin(x), let’s try to calculate \frac{\partial (x * sin(x))}{\partial x} = sin(x) + cos(x) * x

using ForwardDiff
using LinearAlgebra
using Plots

const nx = 20
const Lx = π
const Δx = Lx / (nx-1)
x = range(0, length=nx, stop= Lx)
v = zeros(nx)

v .= sin.(x)

function flux_false(x, v=v)
    x .* v
end

function flux_true(x)
    x .* sin.(x)
end

flux_value = flux_false(x)

Jacobian_wrong = ForwardDiff.jacobian(flux_false, x)
Jacobian_true = ForwardDiff.jacobian(flux_true, x)
deriv= sin.(x) .+ cos.(x) .* x

plot(flux_value, label="value function")
plot!(diag(Jacobian_wrong),label="jacobian wrong")
plot!(diag(Jacobian_true),label="jacobian true")
plot!(deriv, label="true derivatives", linestyle=:dash)

with the plot:

plot_29

So here, I can’t use automatic differentiation because it doesn’t know the relationship between v and x, and it is just considering that v is a constant, so \frac{\partial (x * v)}{\partial x} = v = sin(x). This is the case for my real case problem.

Am I missing something?

4

The actual problem is that you do not include the definition of v in the flux_false function. If instead you would write

function first_step(x)
    return sin.(x)
end

function second_step(x, v)
    return x .* v
end

function flux(x)
    v = first_step(x)
    y = second_step(x, v)
    return y
end

then everything would work smoothly. But your current flux_false treats v as an outside constant, with respect to which it doesn’t try to differentiate.

2 Likes
5

Thanks, I think I should be able to solve my problem with that idea!

1 Like