Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

ODE that depends nonlinearly on a jump process #885

Closed
cumberworth opened this issue Jul 26, 2022 · 7 comments
Closed

ODE that depends nonlinearly on a jump process #885

cumberworth opened this issue Jul 26, 2022 · 7 comments

Comments

@cumberworth
Copy link

cumberworth commented Jul 26, 2022

I need to solve a differential equation that has a nonlinear dependence on a jump process, whose rate is itself dependent on the solution of the ODE. This does seem to be possible with this package; it seems only possible to add a jump process to an ODE. Is there any plan to allow for a more general dependence of an ODE on a jump process than what is currently available?

@cumberworth cumberworth changed the title ODE that depends on a jump process ODE that depends nonlinearly on a jump process Jul 26, 2022
@ChrisRackauckas
Copy link
Member

Have you seen https://diffeq.sciml.ai/stable/tutorials/discrete_stochastic_example/ and https://diffeq.sciml.ai/stable/tutorials/jump_diffusion/? The rate is just a function of (u,p,t) so you can make the jump rate use u, which is from the combined system.

@cumberworth
Copy link
Author

I have seen that, but it doesn't allow for the ODE to use the jump process. What I want is to solve something of the form

$$ \frac{\mathrm{d} u}{\mathrm{d} t} = f(u, p, N, t) $$

where $p$ are parameters, and $N$ is a Poisson process with rate $\lambda(u, p, N, t)$. From the documentation, it is only possible to handle problems of the form

$$ \frac{\mathrm{d} u}{\mathrm{d} t} = f(u, p, t) + c(u, p, t) \mathrm{d} N $$

where $c(u, p, t)$ sets the jump size and the Poisson process has rate $\lambda(u, p, t)$.

@ChrisRackauckas
Copy link
Member

Oh yes, this would require using the random ordinary differential equation (RODE) solvers with a Poisson process.

@cumberworth
Copy link
Author

cumberworth commented Jul 27, 2022

I thought about this as well, but according to the documentation, the RODE problems seem to be restricted to the form

$$ \frac{\mathrm{d} u}{\mathrm{d} t} ​= f(u, p, t, W(t)) $$

where $W(t)$ is the noise/stochastic process. Reading through the documentation for noise processes, I cannot see that it is possible to define a noise process that has a distribution that depends on the solution to the differential equation that it is being used in, or even on it's own solution. Again, I would like to have a Poisson process who's rate depends on both the current value of the solution to the ODE and on the current solution of the Poisson process.

@isaacsas
Copy link
Member

Apologies if I've misunderstood, but if you just want to have a counting process as an input in the ODE, why not make it an extra state-variable within u? i.e. below u[1] is the ODE component and u[2] is the counting process

using DifferentialEquations
function f!(du,u,p,t)
    du[1] = 10 - u[2]*u[1]
    du[2] = 0
    nothing
end
# rate the counting process fires at 
rate(u,p,t) = u[1]*u[2] + 1
function affect!(integrator) 
    # when a jump occurs we increment that value of the counting process by 1
    integrator.u[2] += 1   
    nothing
end
vrj = VariableRateJump(rate, affect!)
oprob = ODEProblem(f!, [100.0, 0.0], (0.0, 10.0))
jprob = JumpProblem(oprob, Direct(), vrj)
sol = solve(jprob, Tsit5())

@isaacsas
Copy link
Member

Also, doesn't your general mathematical problem still fit within the jump-ODE form you wrote down above, i.e. just take your combined equation to be

$$ du(t) = \begin{pmatrix} du_1(t) \\ du_2(t) \end{pmatrix} = \begin{pmatrix} f(u,p,t) \\ 0 \end{pmatrix} + \begin{pmatrix} 0 \\ 1 \end{pmatrix} dN(t) $$

choosing $N(t)$ to have whatever intensity function you want via the rate(u,p,t) function.

@cumberworth
Copy link
Author

Ah, yes, you are right. Thank you both for helping to clarify this.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants