Line data    Source code

       1     1272178 : undefined(::Type{SVector{N, T}}) where {N, T<:Real} =
       2     2570303 :     SVector{N, T}(ntuple(i -> T(NaN), N))
       3             : 
       4           1 : undefined(x::SVector{N, T}) where {T<:Real, N} = undefined(SVector{N, T})
       5             : 
       6             : """
       7             :     val = undefined(x)
       8             :     val = undefined(tyepof(x))
       9             : 
      10             : Helper tthat ceates undefined value of `typeof(x)`.
      11             : 
      12             : See also [`Maybe`](@ref)
      13             : """ undefined
      14             : 
      15             : """
      16             : Encode value of a computation of the ODE solver or failed
      17             : execution [`RK43.State`](@ref).
      18             : 
      19             : See also [`maybe`](@ref), [`RK43.State`](@ref)
      20             : """
      21             : struct Maybe{V}
      22             :     value::V
      23             :     state::State
      24             : 
      25     8656852 :     Maybe{V}(value::V, state::State=Ok) where {V} = new(value, state)
      26             : 
      27      809182 :     function Maybe{V}(state::State) where {V}
      28           0 :         @massert !issuccess(state) "missing value"
      29      809182 :         new(undefined(V), state)
      30             :     end
      31             : 
      32     3888779 :     Maybe{V}(val::Maybe{V}) where {V}= new(val.value, val.state)
      33             : end
      34             : 
      35           2 : function Base.show(io::IO, v::Maybe{V}) where {V}
      36           2 :     issuccess(v) && return show(io, v[])
      37           0 :     show(io, state(v))
      38             : end
      39             : 
      40          22 : maybe(value::V, state::State=Ok,::Type{V}=V) where {V} = Maybe{V}(value, state)
      41           2 : maybe(state::State, ::Type{V}) where {V} = Maybe{V}(state)
      42             : # NOTE: make sure to not confuse the two methods above
      43           2 : maybe(state::State) = error("missing type")
      44             : 
      45     5163322 : state(value::Maybe) = value.state
      46             : 
      47           4 : isundefined(value::Maybe) = !issuccess(value.state)
      48     1284120 : issuccess(value::Maybe) = issuccess(value.state)
      49             : 
      50      879287 : function Base.getindex(value::Maybe{V}) where {V}
      51           0 :     @massert !isundefined(value) "missing value, state=$(state(value))"
      52      879287 :     value.value
      53             : end
      54             : 
      55           0 : requirefinite(val::Maybe{V}) where {V} = val
      56             : 
      57      404590 : function requirefinite(s::State)
      58           0 :     @massert s != Ok "missing value"
      59      404590 :     s
      60             : end
      61             : 
      62     3888779 : @inline function requirefinite(val::SVector{N, T}) where {N, T}
      63     3888779 :     isfinite(sum(val)) && return Maybe{SVector{N, T}}(val)
      64           0 :     Maybe{SVector{N, T}}(Failed)
      65             : end
      66             : 
      67             : """
      68             :     requirefinite(s::State) == s        # requires s != Ok
      69             :     requirefinite(v::Maybe{...}) == v   # no checks
      70             : 
      71             :     v = requirefinite(v::SVector{...}}) :: Maybe{...}
      72             :     @massert isfinite(v.value) || !issucces(v)
      73             : 
      74             : Helper that checks result evaluated as `SVector` in [`WrapDy`](@ref).
      75             : 
      76             : See also [`WrapDy`](@ref)
      77             : """ requirefinite
      78             : 
      79             : 
      80             : """
      81             :     wdy = wrap(dy, solver)
      82             :     wdy(t, y) :: Maybe{...}
      83             : 
      84             : Warp function `dy` in type. `dy(t, y)` returns either
      85             : - an `SVector` (will be checked by  [`requirefinite`](@ref))
      86             : - a `State` other than `Ok` or
      87             : - [`Maybe`](@ref)
      88             : 
      89             : 
      90             : See also [`GuardDy`](@ref), [`Maybe`](@ref), [`requirefinite`](@ref)
      91             : """
      92             : struct WrapDy{F, N, T}
      93             :     dy::F
      94          28 :     WrapDy{F, N, T}(dy) where {F, N, T} = new(dy)
      95             : end
      96             : 
      97     4697959 : (dy::WrapDy{F, N, T})(t::T, y::SVector{N, T}) where {F, N, T} =
      98             :     dy.dy(t, y) |> requirefinite |> Maybe{SVector{N, T}}
      99             : 
     100             : # TODO: requirefinite is slow?!
     101             : 
     102          28 : wrap(dy::F, ::Type{SVector{N, T}}) where {F, N, T} = WrapDy{F, N, T}(dy)
     103           0 : wrap(dy::WrapDy{F, N, T}, ::Type{SVector{N, T}}) where {F, N, T} = dy
     104             : 
     105             : """
     106             :     gdy = guard(wrap(dy, solver))
     107             :     gdy = guard(dy, solver)        # not evaluated if above call failed
     108             : 
     109             :     gdy(t, y) :: SVector{...}
     110             : 
     111             :     issuccess(gdy)
     112             :     s = state(gdy)
     113             : 
     114             : Closure for lazy evaluations of `dy` (within the ODE solver's
     115             :  stepper): Evaluate only while [`issuccess`](@ref).
     116             : 
     117             : See also [`WrapDy`](@ref), [`guard`](@ref)
     118             : """
     119             : mutable struct GuardDy{F, N ,T}
     120             :     const dy::WrapDy{F, N ,T}
     121             :     thestate::State
     122             : 
     123     2616794 :     GuardDy{F, N, T}(dy::WrapDy{F, N, T}, s::State=Ok) where {F, N, T} =
     124             :         new(dy, s)
     125             : end
     126             : 
     127             : """
     128             :     guardedevaluation = guard(warp(dy))
     129             : 
     130             : Construct a [`GuardDy`](@ref) instance for a [`wrap`](@ref)ped
     131             : function `dy(t, y)`.
     132             : 
     133             : See also [`GuardDy`](@ref), [`wrap`](@ref)
     134             : """
     135     1308397 : guard(dy::WrapDy{F, N ,T}) where {F, N, T} = GuardDy{F, N, T}(dy)
     136             : 
     137     4755432 : function (c::GuardDy{F, N, T})(t::T, y::SVector{N, T})::SVector{N, T} where {F, N, T}
     138     4755432 :     (c.thestate != Ok) && return undefined(SVector{N, T})
     139             : 
     140     4697959 :     result = c.dy(t, y)
     141     4293369 :     c.thestate = state(result)
     142     4293369 :     result.value
     143             : end
     144             : 
     145     2187682 : issuccess(guard::GuardDy) = issuccess(state(guard))
     146     2592272 : state(guard::GuardDy) = guard.thestate