# Rational Type support for DifferentialEquations

**URL:** <https://discourse.julialang.org/t/rational-type-support-for-differentialequations/90445>\
**Category:** General Usage\
**Tags:** question, package, diffeq\
**Created:** [November 18, 2022, 12:28pm UTC](https://discourse.julialang.org/t/rational-type-support-for-differentialequations/90445 "2022-11-18T12:28:17Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![svretina](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/svretina/32/206157_2.png) [@svretina](https://discourse.julialang.org/u/svretina)\
**Post date:** [November 18, 2022, 12:28pm UTC](https://discourse.julialang.org/t/rational-type-support-for-differentialequations/90445/1 "2022-11-18T12:28:17Z")

</div>

Hello,

I am trying to use the DifferentialEquations functions with Rational{BigInt} as data type.

for simpler cases it works as it’s shown in [https://tutorials.sciml.ai/html/type\_handling/01-number\_types.html](https://tutorials.sciml.ai/html/type_handling/01-number_types.html)

Even though I turn off the adaptive option in the solve command, i get an

```julia
InexactError: BigInt(1//1000000)

```

and after searching around I found out that DiffEq is trying to set a tolerance (1//1000000) and convert it to BigInt. In order to reproduce the issue:

```julia
function f(u, p, t)
       return 3//1*u
end

u0 = ones(Rational{BigInt}, 10, 3)
prob = ODEProblem(f, u0, (0//1,1//1))
sol = solve(prob, RK4(), dt=1//4, adaptive=false)

```

This will result in the following error message:

```julia
ERROR: InexactError: BigInt(1//1000000)
Stacktrace:
 [1] Integer
   @ ./rational.jl:110 [inlined]
 [2] convert
   @ ./number.jl:7 [inlined]
 [3] __init(prob::ODEProblem{Matrix{Rational{BigInt}}, Tuple{Rational{Int64}, Rational{Int64}}, false, SciMLBase.NullParameters, ODEFunction{false, SciMLBase.AutoSpecialize, typeof(d), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, alg::RK4, timeseries_init::Tuple{}, ts_init::Tuple{}, ks_init::Tuple{}, recompile::Type{Val{true}}; saveat::Tuple{}, tstops::Tuple{}, d_discontinuities::Tuple{}, save_idxs::Nothing, save_everystep::Bool, save_on::Bool, save_start::Bool, save_end::Nothing, callback::Nothing, dense::Bool, calck::Bool, dt::Rational{Int64}, dtmin::Nothing, dtmax::Rational{Int64}, force_dtmin::Bool, adaptive::Bool, gamma::Rational{Int64}, abstol::Nothing, reltol::Nothing, qmin::Rational{Int64}, qmax::Int64, qsteady_min::Int64, qsteady_max::Int64, beta1::Nothing, beta2::Nothing, qoldinit::Rational{Int64}, controller::Nothing, fullnormalize::Bool, failfactor::Int64, maxiters::Int64, internalnorm::typeof(DiffEqBase.ODE_DEFAULT_NORM), internalopnorm::typeof(LinearAlgebra.opnorm), isoutofdomain::typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN), unstable_check::typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK), verbose::Bool, timeseries_errors::Bool, dense_errors::Bool, advance_to_tstop::Bool, stop_at_next_tstop::Bool, initialize_save::Bool, progress::Bool, progress_steps::Int64, progress_name::String, progress_message::typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE), userdata::Nothing, allow_extrapolation::Bool, initialize_integrator::Bool, alias_u0::Bool, alias_du0::Bool, initializealg::OrdinaryDiffEq.DefaultInit, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
   @ OrdinaryDiffEq ~/.julia/packages/OrdinaryDiffEq/8rtul/src/solve.jl:166
 [4] __solve(::ODEProblem{Matrix{Rational{BigInt}}, Tuple{Rational{Int64}, Rational{Int64}}, false, SciMLBase.NullParameters, ODEFunction{false, SciMLBase.AutoSpecialize, typeof(d), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, ::RK4; kwargs::Base.Pairs{Symbol, Real, Tuple{Symbol, Symbol}, NamedTuple{(:dt, :adaptive), Tuple{Rational{Int64}, Bool}}})
   @ OrdinaryDiffEq ~/.julia/packages/OrdinaryDiffEq/8rtul/src/solve.jl:5
 [5] #solve_call#26
   @ ~/.julia/packages/DiffEqBase/Lq1gG/src/solve.jl:473 [inlined]
 [6] solve_up(prob::ODEProblem{Matrix{Rational{BigInt}}, Tuple{Rational{Int64}, Rational{Int64}}, false, SciMLBase.NullParameters, ODEFunction{false, SciMLBase.AutoSpecialize, typeof(d), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, sensealg::Nothing, u0::Matrix{Rational{BigInt}}, p::SciMLBase.NullParameters, args::RK4; kwargs::Base.Pairs{Symbol, Real, Tuple{Symbol, Symbol}, NamedTuple{(:dt, :adaptive), Tuple{Rational{Int64}, Bool}}})
   @ DiffEqBase ~/.julia/packages/DiffEqBase/Lq1gG/src/solve.jl:835
 [7] #solve#31
   @ ~/.julia/packages/DiffEqBase/Lq1gG/src/solve.jl:802 [inlined]
 [8] top-level scope
   @ REPL[68]:1

```

is it possible to overcome this and have something working purely with Rationals, even the solution?

Thanks in advance

---

<div class="post-metadata">

**Author:** ![albheim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/albheim/32/34660_2.png) [@albheim](https://discourse.julialang.org/u/albheim)\
**Post date:** [November 18, 2022, 12:47pm UTC](https://discourse.julialang.org/t/rational-type-support-for-differentialequations/90445/2 "2022-11-18T12:47:30Z")

</div>

It doesn’t error for me using Julia v1.8.2 and DifferentialEquations v7.6.0

```julia
julia> function f(u, p, t)
              return 3//1*u
       end
f (generic function with 1 method)

julia> u0 = ones(Rational{BigInt}, 10, 3)
10×3 Matrix{Rational{BigInt}}:
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1

julia> prob = ODEProblem(f, u0, (0//1,1//1))
ODEProblem with uType Matrix{Rational{BigInt}} and tType Rational{Int64}. In-place: false
timespan: (0//1, 1//1)
u0: 10×3 Matrix{Rational{BigInt}}:
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1
 1//1 1//1 1//1

julia> sol = solve(prob, RK4(), dt=1//4, adaptive=false)
retcode: Success
Interpolation: 3rd order Hermite
t: 5-element Vector{Rational{Int64}}:
 0//1
 1//4
 1//2
 3//4
 1//1
u: 5-element Vector{Matrix{Rational{BigInt}}}:
 [1//1 1//1 1//1; 1//1 1//1 1//1; … ; 1//1 1//1 1//1; 1//1 1//1 1//1]
 [4331//2048 4331//2048 4331//2048; 4331//2048 4331//2048 4331//2048; … ; 4331//2048 4331//2048 4331//2048; 4331//2048 4331//2048 4331//2048]
 [18757561//4194304 18757561//4194304 18757561//4194304; 18757561//4194304 18757561//4194304 18757561//4194304; … ; 18757561//4194304 18757561//4194304 18757561//4194304; 18757561//4194304 18757561//4194304 18757561//4194304]
 [81238996691//8589934592 81238996691//8589934592 81238996691//8589934592; 81238996691//8589934592 81238996691//8589934592 81238996691//8589934592; … ; 81238996691//8589934592 81238996691//8589934592 81238996691//8589934592; 81238996691//8589934592 81238996691//8589934592 81238996691//8589934592]
 [351846094668721//17592186044416 351846094668721//17592186044416 351846094668721//17592186044416; 351846094668721//17592186044416 351846094668721//17592186044416 351846094668721//17592186044416; … ; 351846094668721//17592186044416 351846094668721//17592186044416 351846094668721//17592186044416; 351846094668721//17592186044416 351846094668721//17592186044416 351846094668721//17592186044416]

```

---

<div class="post-metadata">

**Author:** ![John\_Gibson](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/john_gibson/32/5321_2.png) [@John\_Gibson](https://discourse.julialang.org/u/John_Gibson)\
**Post date:** [November 18, 2022, 1:27pm UTC](https://discourse.julialang.org/t/rational-type-support-for-differentialequations/90445/3 "2022-11-18T13:27:31Z")

</div>

What’s the motivation to simulate ODEs using rationals? Those numerators and denominators are going to get large very fast, even when you cancel common factors.

Is it just to investigate or illustrate the time-stepping scheme exactly for a few steps?

---

<div class="post-metadata">

**Author:** ![svretina](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/svretina/32/206157_2.png) [@svretina](https://discourse.julialang.org/u/svretina)\
**Post date:** [November 18, 2022, 1:36pm UTC](https://discourse.julialang.org/t/rational-type-support-for-differentialequations/90445/4 "2022-11-18T13:36:35Z")

</div>

I have some operators that need to satisfy a relation exactly. This would imply exact energy conservation for my simulation. I tried it with floats and the relation hold up to floating point error but the simulation gets unstable. I want to show exact energy conservation even for some time steps.

Now lower orders is fine since I can write a RK4 easily, but for an 8th order operator I need something of order 9 or higher to exclude errors from the time stepper.

It seems to be a problem with my version of the package, I am trying to resolve it and will come back to you.

Thanks for the help!

---

<div class="post-metadata">

**Author:** ![Oscar\_Smith](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/oscar_smith/32/25343_2.png) [@Oscar\_Smith](https://discourse.julialang.org/u/Oscar_Smith)\
**Post date:** [November 18, 2022, 2:41pm UTC](https://discourse.julialang.org/t/rational-type-support-for-differentialequations/90445/5 "2022-11-18T14:41:11Z")

</div>

This sounds like you might want to try using floats and a symplectic solver.

---

<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:** [November 19, 2022, 7:14pm UTC](https://discourse.julialang.org/t/rational-type-support-for-differentialequations/90445/6 "2022-11-19T19:14:20Z")

</div>

Even with 0 floating point error, methods like RK4 will still exhibit lack of energy conservation. It has nothing to do with numerical error if there’s discretization error.