# Integral equations in ModelingToolkit

**URL:** <https://discourse.julialang.org/t/integral-equations-in-modelingtoolkit/74346>\
**Category:** Modelling & Simulations\
**Tags:** diffeq, sciml\
**Created:** [January 10, 2022, 5:27pm UTC](https://discourse.julialang.org/t/integral-equations-in-modelingtoolkit/74346 "2022-01-10T17:27:11Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![sdwfrost](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sdwfrost/32/2831_2.png) [@sdwfrost](https://discourse.julialang.org/u/sdwfrost)\
**Post date:** [January 10, 2022, 5:27pm UTC](https://discourse.julialang.org/t/integral-equations-in-modelingtoolkit/74346/1 "2022-01-10T17:27:11Z")

</div>

Hi Everyone,

I’d like to know how/if the following integral equation can be solved in ModelingToolkit:

```julia
using ModelingToolkit
using DomainSets

@parameters t β γ
@variables S(..) I(..) R(..)
It = Integral(t in DomainSets.ClosedInterval(0,t))

eqs = [S(t) ~ S(0.) - It(β*S(t)*I(t)),
       I(t) ~ I(0.) + It(β*S(t)*I(t)) - It(γ*I(t)),
       R(t) ~ R(0.) + It(γ*I(t))];

```

I know it’s a trivial example that can be rewritten as ODEs, but knowing how to do this would open up a lot.

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [January 10, 2022, 7:23pm UTC](https://discourse.julialang.org/t/integral-equations-in-modelingtoolkit/74346/2 "2022-01-10T19:23:55Z")

</div>

NeuralPDE.jl should be able to do it now? [https://github.com/SciML/NeuralPDE.jl/pull/409](https://github.com/SciML/NeuralPDE.jl/pull/409) was what was required. Give it a try.

---

<div class="post-metadata">

**Author:** ![sdwfrost](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sdwfrost/32/2831_2.png) [@sdwfrost](https://discourse.julialang.org/u/sdwfrost)\
**Post date:** [January 23, 2023, 9:36pm UTC](https://discourse.julialang.org/t/integral-equations-in-modelingtoolkit/74346/3 "2023-01-23T21:36:22Z")

</div>

I gave this a try in NeuralPDE.jl, but like all my other attempts, it didn’t work well at all - can you spot any issues?

```julia
using ModelingToolkit
using DomainSets
using NeuralPDE
using Flux
using Optimization
using OptimizationOptimisers
using OptimizationOptimJL
using Plots

@parameters t
@variables S(..) I(..) R(..)
It = Integral(t in DomainSets.ClosedInterval(0,t))

eqs = [S(t) ~ S(0.) - It(0.5*S(t)*I(t)),
       I(t) ~ I(0.) + It(0.5*S(t)*I(t)) - It(0.25*I(t)),
       R(t) ~ R(0.) + It(0.25*I(t))];

bcs = [S(0) ~ 0.99, I(0) ~ 0.01, R(0) ~ 0.0];

domains = [t ∈ (0.0, 40.0)];

@named pde_system = PDESystem(eqs, bcs, domains, [t], [S(t), I(t), R(t)]);

numhid = 4
chain = [Chain(Dense(1, numhid, Flux.σ),
               Dense(numhid, numhid, Flux.σ),
               Dense(numhid, 1, abs)) for i in 1:3] |> f64;

strategy_ = GridTraining(0.1);
discretization = PhysicsInformedNN(chain, strategy_);
prob = NeuralPDE.discretize(pde_system, discretization);

callback = function (p,l)
    println("Current loss is: $l")
    return false
end;
res = Optimization.solve(prob, OptimizationOptimisers.Adam(0.1); callback = callback, maxiters=20);
prob = remake(prob, u0 = res.minimizer)
res = Optimization.solve(prob, OptimizationOptimJL.BFGS(); callback = callback, maxiters=20);

ts = [infimum(d.domain):0.1:supremum(d.domain) for d in domains][1]
phi = discretization.phi

np = length(res.u) ÷ 3
acum = [0;accumulate(+, [np,np,np])]
sep = [acum[i]+1 : acum[i+1] for i in 1:length(acum)-1] 

Spred = hcat([phi[1]([t],res.u[sep[1]]) for t in ts]...)'
Ipred = hcat([phi[2]([t],res.u[sep[2]]) for t in ts]...)'
Rpred = hcat([phi[3]([t],res.u[sep[3]]) for t in ts]...)'

plot(ts, Spred, label="S")
plot!(ts, Ipred, label="I")
plot!(ts, Rpred, label="R")

```

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [January 30, 2023, 1:36pm UTC](https://discourse.julialang.org/t/integral-equations-in-modelingtoolkit/74346/4 "2023-01-30T13:36:26Z")

</div>

Open an issue. Part of the problem here is that some tuning of layer choices and such may be necessary.
