working default Tuple type

Author: Brad Carman

src/parameters.jl

@@ -62,28 +62,50 @@ macro parameters(xs...)
         toparam) |> esc
 end
 
-function split_parameters_by_type(ps)
-    by = let set = Dict{Any, Int}(), counter = Ref(1)
+function find_types(array)
+    by = let set = Dict{Any, Int}(), counter = Ref(0)
         x -> begin
-            t = typeof(x)
+            # t = typeof(x)
+            
             get!(set, typeof(x)) do
-                if t == Float64
-                    1
-                else
-                    counter[] += 1
-                end
+                # if t == Float64
+                #     1
+                # else
+                counter[] += 1
+                # end
             end
         end
     end
-    idxs = by.(ps)
-    split_idxs = [Int[]]
-    for (i, idx) in enumerate(idxs)
-        if idx > length(split_idxs)
-            push!(split_idxs, Int[])
+    return by.(array)
+end
+
+
+function split_parameters_by_type(ps)
+
+    if ps === SciMLBase.NullParameters()
+        return Float64[],[] #use Float64 to avoid Any type warning
+    else
+        by = let set = Dict{Any, Int}(), counter = Ref(0)
+            x -> begin     
+                get!(set, typeof(x)) do
+                    counter[] += 1                    
+                end
+            end
+        end
+        idxs = by.(ps)
+        split_idxs = [Int[]]
+        for (i, idx) in enumerate(idxs)
+            if idx > length(split_idxs)
+                push!(split_idxs, Int[])
+            end
+            push!(split_idxs[idx], i)
+        end
+        tighten_types = x -> identity.(x)
+        split_ps = tighten_types.(Base.Fix1(getindex, ps).(split_idxs))
+        if length(split_ps) == 1  #Tuple not needed, only 1 type
+            return split_ps[1], split_idxs 
+        else
+            return (split_ps...,), split_idxs
         end
-        push!(split_idxs[idx], i)
     end
-    tighten_types = x -> identity.(x)
-    split_ps = tighten_types.(Base.Fix1(getindex, ps).(split_idxs))
-    (split_ps...,), split_idxs
 end

src/systems/diffeqs/abstractodesystem.jl

@@ -122,7 +122,6 @@ function generate_function(sys::AbstractODESystem, dvs = states(sys), ps = param
            nothing,
     isdde = false,
     has_difference = false,
-    split_parameters = false,
     kwargs...)
     if isdde
         eqs = delay_to_function(sys)
@@ -152,12 +151,14 @@ function generate_function(sys::AbstractODESystem, dvs = states(sys), ps = param
             build_function(rhss, ddvs, u, p, t; postprocess_fbody = pre,
                 states = sol_states,
                 kwargs...)
-        elseif split_parameters
-            build_function(rhss, u, p..., t; postprocess_fbody = pre, states = sol_states,
-                kwargs...)
         else
-            build_function(rhss, u, p, t; postprocess_fbody = pre, states = sol_states,
+            if p isa Tuple
+                build_function(rhss, u, p..., t; postprocess_fbody = pre, states = sol_states,
+                    kwargs...)
+            else
+                build_function(rhss, u, p, t; postprocess_fbody = pre, states = sol_states,
                 kwargs...)
+            end
         end
     end
 end
@@ -329,15 +330,14 @@ function DiffEqBase.ODEFunction{iip, specialize}(sys::AbstractODESystem, dvs = s
     checkbounds = false,
     sparsity = false,
     analytic = nothing,
-    split_parameters = false,
     kwargs...) where {iip, specialize}
     f_gen = generate_function(sys, dvs, ps; expression = Val{eval_expression},
-        expression_module = eval_module, checkbounds = checkbounds, split_parameters,
+        expression_module = eval_module, checkbounds = checkbounds, 
         kwargs...)
     f_oop, f_iip = eval_expression ?
                    (drop_expr(@RuntimeGeneratedFunction(eval_module, ex)) for ex in f_gen) :
                    f_gen
-    if split_parameters
+    if p isa Tuple
         g(u, p, t) = f_oop(u, p..., t)
         g(du, u, p, t) = f_iip(du, u, p..., t)
         f = g
@@ -696,9 +696,8 @@ Take dictionaries with initial conditions and parameters and convert them to num
 function get_u0_p(sys,
     u0map,
     parammap;
-    use_union = false,
-    tofloat = !use_union,
-    split_parameters = false,
+    use_union = true,
+    tofloat = true,
     symbolic_u0 = false)
     eqs = equations(sys)
     dvs = states(sys)
@@ -711,7 +710,7 @@ function get_u0_p(sys,
     if symbolic_u0
         u0 = varmap_to_vars(u0map, dvs; defaults = defs, tofloat = false, use_union = false)
     else
-        u0 = varmap_to_vars(u0map, dvs; defaults = defs, tofloat = !split_parameters)
+        u0 = varmap_to_vars(u0map, dvs; defaults = defs, tofloat = true)
     end
     p = varmap_to_vars(parammap, ps; defaults = defs, tofloat, use_union)
     p = p === nothing ? SciMLBase.NullParameters() : p
@@ -726,10 +725,10 @@ function process_DEProblem(constructor, sys::AbstractODESystem, u0map, parammap;
     simplify = false,
     linenumbers = true, parallel = SerialForm(),
     eval_expression = true,
-    use_union = false,
-    tofloat = !use_union,
+    use_union = true,
+    tofloat = true,
     symbolic_u0 = false,
-    split_parameters = false,
+    # split_parameters = true,
     kwargs...)
     eqs = equations(sys)
     dvs = states(sys)
@@ -741,12 +740,15 @@ function process_DEProblem(constructor, sys::AbstractODESystem, u0map, parammap;
         parammap;
         tofloat,
         use_union,
-        symbolic_u0,
-        split_parameters)
-    if split_parameters
+        symbolic_u0)
+
+    # if split_parameters
         p, split_idxs = split_parameters_by_type(p)
-        ps = Base.Fix1(getindex, parameters(sys)).(split_idxs)
-    end
+        if p isa Tuple
+            ps = Base.Fix1(getindex, parameters(sys)).(split_idxs)
+            ps = (ps...,) #if p is Tuple, ps should be Tuple
+        end
+    # end
 
     if implicit_dae && du0map !== nothing
         ddvs = map(Differential(iv), dvs)
@@ -763,7 +765,7 @@ function process_DEProblem(constructor, sys::AbstractODESystem, u0map, parammap;
     f = constructor(sys, dvs, ps, u0; ddvs = ddvs, tgrad = tgrad, jac = jac,
         checkbounds = checkbounds, p = p,
         linenumbers = linenumbers, parallel = parallel, simplify = simplify,
-        sparse = sparse, eval_expression = eval_expression, split_parameters,
+        sparse = sparse, eval_expression = eval_expression, 
         kwargs...)
     implicit_dae ? (f, du0, u0, p) : (f, u0, p)
 end

src/variables.jl

@@ -58,7 +58,7 @@ applicable.
 """
 function varmap_to_vars(varmap, varlist; defaults = Dict(), check = true,
     toterm = default_toterm, promotetoconcrete = nothing,
-    tofloat = true, use_union = false)
+    tofloat = true, use_union = true)
     varlist = collect(map(unwrap, varlist))
 
     # Edge cases where one of the arguments is effectively empty.
@@ -75,9 +75,10 @@ function varmap_to_vars(varmap, varlist; defaults = Dict(), check = true,
         end
     end
 
-    T = typeof(varmap)
-    # We respect the input type
-    container_type = T <: Dict ? Array : T
+    # T = typeof(varmap)
+    # We respect the input type (feature removed, not needed with Tuple support)
+    # container_type = T <: Union{Dict,Tuple} ? Array : T
+    container_type = Array
 
     vals = if eltype(varmap) <: Pair # `varmap` is a dict or an array of pairs
         varmap = todict(varmap)

test/odesystem.jl

@@ -734,18 +734,23 @@ let
     u0map = [A => 1.0]
     pmap = (k1 => 1.0, k2 => 1)
     tspan = (0.0, 1.0)
-    prob = ODEProblem(sys, u0map, tspan, pmap)
-    @test prob.p === Tuple([(Dict(pmap))[k] for k in values(parameters(sys))])
+    prob = ODEProblem(sys, u0map, tspan, pmap; tofloat=false)
+
+    @test prob.p == ([1], [1.0]) #Tuple([(Dict(pmap))[k] for k in values(parameters(sys))])
 
     pmap = [k1 => 1, k2 => 1]
     tspan = (0.0, 1.0)
     prob = ODEProblem(sys, u0map, tspan, pmap)
     @test eltype(prob.p) === Float64
-
-    pmap = Pair{Any, Union{Int, Float64}}[k1 => 1, k2 => 1.0]
-    tspan = (0.0, 1.0)
-    prob = ODEProblem(sys, u0map, tspan, pmap, use_union = true)
-    @test eltype(prob.p) === Union{Float64, Int}
+    
+    prob = ODEProblem(sys, u0map, tspan, pmap; tofloat=false)
+    @test eltype(prob.p) === Int
+
+    # No longer supported, Tuple used instead
+    # pmap = Pair{Any, Union{Int, Float64}}[k1 => 1, k2 => 1.0]
+    # tspan = (0.0, 1.0)
+    # prob = ODEProblem(sys, u0map, tspan, pmap, use_union = true)
+    # @test eltype(prob.p) === Union{Float64, Int}
 end
 
 let

test/split_parameters.jl

@@ -2,6 +2,10 @@ using ModelingToolkit, Test
 using ModelingToolkitStandardLibrary.Blocks
 using OrdinaryDiffEq
 
+
+
+# ------------------------ Mixed Single Values and Vector
+
 dt = 4e-4
 t_end = 10.0
 time = 0:dt:t_end
@@ -10,8 +14,31 @@ x = @. time^2 + 1.0
 @parameters t
 D = Differential(t)
 
+get_value(data, t, dt) = data[round(Int, t/dt+1)]
+@register_symbolic get_value(data, t, dt)
+
+
+function Sampled(; name, data=Float64[], dt=0.0)
+    pars = @parameters begin
+        data = data
+        dt = dt
+    end
+
+    vars = []
+    systems = @named begin
+        output = RealOutput()
+    end
+
+    eqs = [
+        output.u ~ get_value(data, t, dt)
+    ]
+
+    return ODESystem(eqs, t, vars, pars; name, systems,
+        defaults = [output.u => data[1]])
+end
+
 vars = @variables y(t)=1 dy(t)=0 ddy(t)=0
-@named src = SampledData(; data = Float64[], dt)
+@named src = Sampled(; data = Float64[], dt)
 @named int = Integrator()
 
 eqs = [y ~ src.output.u
@@ -22,10 +49,39 @@ eqs = [y ~ src.output.u
 @named sys = ODESystem(eqs, t, vars, []; systems = [int, src])
 s = complete(sys)
 sys = structural_simplify(sys)
-
-prob = ODEProblem(sys, [], (0.0, t_end), [s.src.data => x]; split_parameters = true)
+prob = ODEProblem(sys, [], (0.0, t_end), [s.src.data => x])
 @test prob.p isa Tuple{Vector{Float64}, Vector{Int}, Vector{Vector{Float64}}}
-@time sol = solve(prob, ImplicitEuler());
-prob2 = ODEProblem(sys, [], (0.0, t_end), [s.src.data => x]; to_float = false)
-@test prob2.p isa Vector{Union{Float64, Int64, Vector{Float64}}}
-@time sol2 = solve(prob2, ImplicitEuler());
+sol = solve(prob, ImplicitEuler());
+@test sol.retcode == ReturnCode.Success
+
+
+# ------------------------ Mixed Type Converted to float (default behavior)
+
+vars = @variables y(t)=1 dy(t)=0 ddy(t)=0
+pars = @parameters a=1.0 b=2.0 c=3
+eqs = [
+    D(y) ~ dy*a
+    D(dy) ~ ddy*b
+    ddy ~ sin(t)*c]
+
+@named sys = ODESystem(eqs, t, vars, pars)
+sys = structural_simplify(sys)
+
+tspan = (0.0, t_end)
+prob = ODEProblem(sys, [], tspan, [])
+
+@test prob.p isa Vector{Float64}
+sol = solve(prob, ImplicitEuler());
+@test sol.retcode == ReturnCode.Success
+
+
+# ------------------------ Mixed Type Conserved
+
+prob = ODEProblem(sys, [], tspan, []; tofloat=false)
+
+@test prob.p isa Tuple{Vector{Float64}, Vector{Int64}}
+sol = solve(prob, ImplicitEuler());
+@test sol.retcode == ReturnCode.Success
+
+
+