Integration within ModelingToolkit Models

Specific Domains Modelling & Simulations
1

Hello,

Sorry I am asking so many questions on here. I would like to be able to integrate part of an equation within a MTK model so that each time step the integration is performed with the current variables and parameters. However, when building the model, the integration wants to perform at that moment therefore it is unsolvable as it tries to perform symbolic rather than numeric integration. This is compounded by the fact that a component within the integral is an inteprolation function which can’t be integrated any way other than numerically. Here is a MWE:

function func()
    @parameters a b
    @variables x
    return a*x + b
end

function theintegral(int)
    @variables x
    return SymbolicNumericIntegration.integrate(x*int(x)*func(),(x,-10,10),symbolic=false)
end

function OurODE(;name)
    @variables a(t) placeholder(t)
    
    eqs = D(a) ~ placeholder
end

function main()
    int = get_interpolate(stuff)
    @named ODE = OurODE
    integral = theintegral(int)
    connections = [ODE.placeholder ~ integral]
    connected = compose(ODESystem(connections,t,name=:connected,ODE)
    connected_simp = structural_simplify(connected)
    u0=[ODE.a => 0.0]
    p[ODE.b => 1]
    prob=ODEProblem(connected_simp,u0,(0.0,500),p
    sol=solve(prob,AutoVern7(Rodas5());abstol=1e-10,reltol=1e-10)
end

This should have enough information that it can solve everything numerically at the time I want it solved but it tries to solve the integral at the call to theintegral(int) which makes sense but I want it to return almost a placeholder num that says it will integrate numerically when it gets there. Is this too much to ask/is there a different way of doing this
Any help would be appreciated thank you.

2

I had a thought that it might be possible with some callback that tells the system to update the placeholder parameter every time step with the solution to the equation but I’m not sure if that would work or is possible

3

Why not just set it up with Integrals.jl and register the numerical approximation?

4

How would I go about doing this?

As the function call to Integrals.jl in the MWE wouldn’t have numerical values for the parameters a and b yet so it would only be able to evaluate the integral once the ODESystem has been built and solution begun

5

You need to write a function that calls solve and then @register it

6

Hi @ChrisRackauckas,

I’ve tried to do this is my real code except it is now passing nums to the integral which I don’t think it knows how to handle. Is there a keyword somewhere that tells register_symbolic to unpack Nums into floats when sending to these functions.
The documentation isn’t very clear and I couldn’t find a tutorial on how to use the macro.

This is the code for the integral component. kB,Tel and μ are parameters and variables while DOS is an interpolation that I pass to the function

function dFDdT(kB,Tel,μ,E)
    e = exp((E-μ)/kB*Tel)
    return (E-μ)*e/(kB*Tel^2*(e+1)^2)
end

function HeatCapacity(kB,Tel,μ,DOS)
    int(u,p)=dFDdT(kB,Tel,μ,u)*DOS(u)*u
    prob=IntegralProblem(int,0,20)
    sol = solve(prob,HCubatureJL(initdiv=2);reltol=1e-5,abstol=1e-5)
    return sol.u
end
@register_symbolic HeatCapacity(kB,Tel,μ,DOS)
7

DOS is a function and not symbolic? Then you need:

@register_symbolic HeatCapacity(kB,Tel,μ,DOS::Function)