Is there a specific reason why Float64 was enforced in PoincareMap?
BigFloat support in PoincareMap
Right now, PoincareMap is hard-wired to Float64 in several fields:
```
mutable struct PoincareMap{I<:ContinuousTimeDynamicalSystem, F, P, R, V} <: DiscreteTimeDynamicalSystem
ds::I
plane_distance::F
planecrossing::P
Tmax::Float64 # here
rootkw::R
state_on_plane::V
tcross::Float64 # here
t::Int
# These two fields are for setting the state of the pmap from the plane
# (i.e., given a D-1 dimensional state, create the full D-dimensional state)
dummy::Vector{Float64} # here
diffidxs::Vector{Int}
state_initial::V
end
```
and initializer is
```
function PoincareMap(
ds::DS, plane;
Tmax = 1e3,
direction = -1, u0 = nothing,
rootkw = (xrtol = 1e-6, atol = 1e-8)
) where {DS<:ContinuousTimeDynamicalSystem}
reinit!(ds, u0)
D = dimension(ds)
check_hyperplane_match(plane, D)
planecrossing = PlaneCrossing(plane, direction > 0)
plane_distance = (t) -> planecrossing(ds(t))
v = recursivecopy(current_state(ds))
dummy = zeros(D)
diffidxs = _indices_on_poincare_plane(plane, D)
pmap = PoincareMap(
ds, plane_distance, planecrossing, Tmax, rootkw,
v, current_time(ds), 0, dummy, diffidxs, recursivecopy(v), # problem occurs here!
)
step!(pmap) # this ensures that the state is on the hyperplane
pmap.state_initial = recursivecopy(current_state(pmap))
pmap.t = 0 # first step is 0
return pmap
end
```
However, current_time(ds) adheres to the type of the initial state, u0.
So if u0 is BigFloat, the integrator’s time is also BigFloat, which the struct doesn't accept
**Proposal:**
Make PoincaréMap parametric over the time type and state element type, e.g.
```
mutable struct PoincareMap{I<:ContinuousTimeDynamicalSystem,Tv<:AbstractFloat}
ds::I
...
Tmax::Tv
tcross::Tv
dummy::Vector{Tv}
...
end
```
and convert literals like 1e3 into Tv(1e3) inside the constructor.
Is there a reason Float64 was fixed originally, or would a PR making this change be welcome?
1 条评论