Properly use symbolic constructors

This was using an old form that didn't necessarily always get the order right. This PR fixes the constructors to ensure they are safe for any equation reordering, and thus fixes downstream.

Author: ChrisRackauckas

lib/ODEProblemLibrary/Project.toml

@@ -1,6 +1,6 @@
 name = "ODEProblemLibrary"
 uuid = "fdc4e326-1af4-4b90-96e7-779fcce2daa5"
-version = "0.1.8"
+version = "0.1.9"
 
 [deps]
 DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e"

lib/ODEProblemLibrary/src/ode_simple_nonlinear_prob.jl

@@ -54,8 +54,7 @@ prob_ode_fitzhughnagumo = ODEProblem(fitz, [1.0; 1.0], (0.0, 1.0),
 
 eqs = [D(y) ~ μ * ((1 - x^2) * y - x),
     D(x) ~ y]
-de = ODESystem(eqs, t; name = :van_der_pol) |> structural_simplify |> complete
-van = ODEFunction(de, [y, x], [μ], jac = true, eval_module = @__MODULE__)
+van = ODESystem(eqs, t; name = :van_der_pol) |> structural_simplify |> complete
 
 """
 Van der Pol Equations
@@ -67,11 +66,11 @@ Van der Pol Equations
 \\frac{dy}{dt} = μ((1-x^2)y -x)
 ```
 
-with ``μ=1.0`` and ``u_0=[0,\\sqrt{3}]``
+with ``μ=1.0`` and ``u_0=[x => \\sqrt{3}, y => 0]``
 
 Non-stiff parameters.
 """
-prob_ode_vanderpol = ODEProblem(van, [0; sqrt(3)], (0.0, 1.0), 1.0)
+prob_ode_vanderpol = ODEProblem(van, [y => 0, x => sqrt(3)], (0.0, 1.0), [μ => 1.0], jac=true, eval_module = @__MODULE__)
 
 """
 Van der Pol Equations
@@ -83,11 +82,11 @@ Van der Pol Equations
 \\frac{dy}{dt} = μ((1-x^2)y -x)
 ```
 
-with ``μ=10^6`` and ``u_0=[0,\\sqrt{3}]``
+with ``μ=10^6`` and ``u_0=[x => \\sqrt{3}, y => 0]``
 
 Stiff parameters.
 """
-prob_ode_vanderpol_stiff = ODEProblem(van, [0; sqrt(3)], (0.0, 1.0), 1e6)
+prob_ode_vanderpol_stiff = ODEProblem(van, [y => 0, x => sqrt(3)], (0.0, 1.0), [μ => 1e6], jac=true, eval_module = @__MODULE__)
 
 # ROBER
 @parameters k₁ k₂ k₃
@@ -96,8 +95,7 @@ prob_ode_vanderpol_stiff = ODEProblem(van, [0; sqrt(3)], (0.0, 1.0), 1e6)
 eqs = [D(y₁) ~ -k₁ * y₁ + k₃ * y₂ * y₃,
     D(y₂) ~ k₁ * y₁ - k₂ * y₂^2 - k₃ * y₂ * y₃,
     D(y₃) ~ k₂ * y₂^2]
-de = ODESystem(eqs, t; name = :rober) |> structural_simplify |> complete
-rober = ODEFunction(de, [y₁, y₂, y₃], [k₁, k₂, k₃], jac = true, eval_module = @__MODULE__)
+rober = ODESystem(eqs, t; name = :rober) |> structural_simplify |> complete
 
 """
 The Robertson biochemical reactions: (Stiff)
@@ -118,7 +116,7 @@ Hairer Norsett Wanner Solving Ordinary Differential Equations I - Nonstiff Probl
 
 Usually solved on ``[0,1e11]``
 """
-prob_ode_rober = ODEProblem(rober, [1.0; 0.0; 0.0], (0.0, 1e11), (0.04, 3e7, 1e4))
+prob_ode_rober = ODEProblem(rober, [y₁, y₂, y₃] .=> [1.0; 0.0; 0.0], (0.0, 1e11), [k₁, k₂, k₃] .=> (0.04, 3e7, 1e4), jac = true, eval_module = @__MODULE__)
 
 # Three Body
 const threebody_μ = big(0.012277471);
@@ -177,8 +175,7 @@ prob_ode_threebody = ODEProblem(threebody,
 eqs = [D(y₁) ~ I₁ * y₂ * y₃,
     D(y₂) ~ I₂ * y₁ * y₃,
     D(y₃) ~ I₃ * y₁ * y₂]
-de = ODESystem(eqs, t; name = :rigid_body) |> structural_simplify |> complete
-rigid = ODEFunction(de, [y₁, y₂, y₃], [I₁, I₂, I₃], jac = true, eval_module = @__MODULE__)
+rigid = ODESystem(eqs, t; name = :rigid_body) |> structural_simplify |> complete
 
 """
 Rigid Body Equations (Non-stiff)
@@ -203,7 +200,9 @@ or Hairer Norsett Wanner Solving Ordinary Differential Equations I - Nonstiff Pr
 
 Usually solved from 0 to 20.
 """
-prob_ode_rigidbody = ODEProblem(rigid, [1.0, 0.0, 0.9], (0.0, 20.0), (-2.0, 1.25, -0.5))
+prob_ode_rigidbody = ODEProblem(rigid, [y₁, y₂, y₃] .=> [1.0, 0.0, 0.9], (0.0, 20.0),
+                                       [I₁, I₂, I₃] .=> (-2.0, 1.25, -0.5), jac = true,
+                                       eval_module = @__MODULE__)
 
 # Pleiades Problem
 
@@ -360,10 +359,7 @@ eqs = [D(y1) ~ -p1 * y1 + p2 * y2 + p3 * y3 + p4,
             p2 * y6 + p11 * y7,
     D(y7) ~ p10 * y6 * y8 - p12 * y7,
     D(y8) ~ -p10 * y6 * y8 + p12 * y7]
-de = ODESystem(eqs, t; name = :hires) |> structural_simplify |> complete
-hires = ODEFunction(de, [y1, y2, y3, y4, y5, y6, y7, y8],
-    [p1, p2, p3, p4, p5, p6, p7, p8, p9, p10, p11, p12],
-    jac = true)
+hires = ODESystem(eqs, t; name = :hires) |> structural_simplify |> complete
 
 u0 = zeros(8)
 u0[1] = 1
@@ -391,20 +387,18 @@ where ``f`` is defined by
 Reference: [demohires.pdf](http://www.radford.edu/~thompson/vodef90web/problems/demosnodislin/Demos_Pitagora/DemoHires/demohires.pdf)
 Notebook: [Hires.ipynb](http://nbviewer.jupyter.org/github/JuliaDiffEq/DiffEqBenchmarks.jl/blob/master/StiffODE/Hires.ipynb)
 """
-prob_ode_hires = ODEProblem(hires, u0, (0.0, 321.8122),
-    (1.71, 0.43, 8.32, 0.0007, 8.75,
-        10.03, 0.035, 1.12, 1.745, 280.0,
-        0.69, 1.81))
+prob_ode_hires = ODEProblem(hires, [y1, y2, y3, y4, y5, y6, y7, y8] .=> u0, (0.0, 321.8122),
+                            [p1, p2, p3, p4, p5, p6, p7, p8, p9, p10, p11, p12] .=> (1.71, 0.43, 8.32, 0.0007, 8.75,
+                                    10.03, 0.035, 1.12, 1.745, 280.0,
+                                    0.69, 1.81), jac = true, eval_module = @__MODULE__)
 
 @parameters p1 p2 p3
 @variables y1(t) y2(t) y3(t)
 
 eqs = [D(y1) ~ p1 * (y2 + y1 * (1 - p2 * y1 - y2)),
     D(y2) ~ (y3 - (1 + y1) * y2) / p1,
     D(y3) ~ p3 * (y1 - y3)]
-de = ODESystem(eqs, t; name = :orego) |> structural_simplify |> complete
-jac = calculate_jacobian(de)
-orego = ODEFunction(de, [y1, y2, y3], [p1, p2, p3], jac = true, eval_module = @__MODULE__)
+orego = ODESystem(eqs, t; name = :orego) |> structural_simplify |> complete
 
 """
 Orego Problem (Stiff)
@@ -424,4 +418,4 @@ where ``s=77.27``, ``w=0.161`` and ``q=8.375⋅10^{-6}``.
 Reference: [demoorego.pdf](http://www.radford.edu/~thompson/vodef90web/problems/demosnodislin/Demos_Pitagora/DemoOrego/demoorego.pdf)
 Notebook: [Orego.ipynb](http://nbviewer.jupyter.org/github/JuliaDiffEq/DiffEqBenchmarks.jl/blob/master/StiffODE/Orego.ipynb)
 """
-prob_ode_orego = ODEProblem(orego, [1.0, 2.0, 3.0], (0.0, 30.0), [77.27, 8.375e-6, 0.161])
+prob_ode_orego = ODEProblem(orego, [y1, y2, y3] .=> [1.0, 2.0, 3.0], (0.0, 30.0), [p1, p2, p3] .=> [77.27, 8.375e-6, 0.161], jac = true, eval_module = @__MODULE__)