State Dependent Mass Matrix

Hello,

I’ve looked at the documentation for Handling Mass Matrices for DifferentialEquations.jl. In the example was a singular matrix. I would like to know how to implement a state dependent mass matrix, where the matrix M in equation
Mu’ =f(u,p,t)
has the dependent variables in the matrix itself.

The image below is the equation I’m trying to replicate:

I’m trying to integrate from Matlab to Julia for this type of problem. I thought replacing the mass_matrix option with a function would work but it produces the method error MethodError: no method matching issingular(::typeof(function goes here))

An example or helpful tip would be appreciated. The link relating to mass matrices is https://docs.sciml.ai/latest/tutorials/advanced_ode_example/#Handling-Mass-Matrices-1 Thank you!

I think this is not implemented right now ( but will be, I am sure ).

I planned on doing this after a discussion with @ChrisRackauckas and @YingboMa on Slack, but never found the time.

If I encounter such problems right now, I tend to use something along the lines of

function dudt(u, p, t)
    m = M(u,p,t)
    f = F(u,p,t)
    inv(m)*f
end

or solve it manually, if possible.

It’s implemented. Examples:

We just need to clean up our documentation of DiffEqOperators to make it more clear how to define dependent DiffEqOperators, since that’s all it is (the mass matrix can be any DiffEqOperator)

Ah, I did not know that :+1:

Thank you for that. I’m still new to Julia, could you link me to the equation your code snippet is solving to help visualize it for me? I can try to figure out the syntax from there.

I am really really slow with making matrices in LaTeX, so does this simplification help?

function update_func(A,u,p,t)
    A[1,1] = cos(t)
    A[2,1] = sin(t)*u[1]
    A[3,1] = t^2
    A[1,2] = cos(t)*sin(t)
    A[2,2] = cos(t)^2 + u[3]
    A[3,2] = sin(t)*u[2]
    A[1,3] = sin(t)
    A[2,3] = t^2
    A[3,3] = t*cos(t) + 1
end
M = DiffEqArrayOperator(ones(3,3),update_func=update_func)
u0 = ones(3)
tspan = (0.0, 1.0)
f(du,u,p,t) = (@. du = u + t; nothing)
prob = ODEProblem(ODEFunction(f, mass_matrix=M), u0, tspan)
sol = solve(prob,RadauIIA5())

The mass matrix is just the A defined in the update_func. I’ll trade you for LaTeX and together we can get this in the docs :slight_smile:

2 Likes

Has the solution to this problem updated? DifferentialEquations.jl documentation says that state dependent mass matrices are “not directly supported” and there is no documentation in DiffEqOperators.jl either. Would the current implementation support Complex variables? @ChrisRackauckas @Julius_Martensen

State-dependent mass matrices are inherently pretty unstable, so I would recommend a DAE formulation instead.

1 Like