diff --git a/lib/ODEProblemLibrary/Project.toml b/lib/ODEProblemLibrary/Project.toml index babc503..fb1f4ec 100644 --- a/lib/ODEProblemLibrary/Project.toml +++ b/lib/ODEProblemLibrary/Project.toml @@ -4,23 +4,17 @@ version = "0.1.12" [deps] DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e" -Latexify = "23fbe1c1-3f47-55db-b15f-69d7ec21a316" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Markdown = "d6f4376e-aef5-505a-96c1-9c027394607a" -ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" -RuntimeGeneratedFunctions = "7e49a35a-f44a-4d26-94aa-eba1b4ca6b47" [compat] Aqua = "0.5" DiffEqBase = "6" -Latexify = "0.16" -ModelingToolkit = "10" -RuntimeGeneratedFunctions = "0.5" julia = "1.10" [extras] Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595" [targets] -test = ["Aqua"] +test = ["Aqua"] \ No newline at end of file diff --git a/lib/ODEProblemLibrary/src/ODEProblemLibrary.jl b/lib/ODEProblemLibrary/src/ODEProblemLibrary.jl index ba3be4a..9dc24e4 100644 --- a/lib/ODEProblemLibrary/src/ODEProblemLibrary.jl +++ b/lib/ODEProblemLibrary/src/ODEProblemLibrary.jl @@ -1,17 +1,10 @@ module ODEProblemLibrary using DiffEqBase -using Latexify -using ModelingToolkit -using ModelingToolkit: t_nounits as t, D_nounits as D - using LinearAlgebra using Markdown using Random -import RuntimeGeneratedFunctions -RuntimeGeneratedFunctions.init(@__MODULE__) - Random.seed!(100) #ODE Example Problems @@ -30,19 +23,13 @@ export prob_ode_linear, prob_ode_bigfloatlinear, prob_ode_2Dlinear, prob_ode_chen, prob_ode_rossler, prob_ode_rabinovich_fabrikant, prob_ode_sprott, prob_ode_hindmarsh_rose -# For MTKv9 compatibility -@static if !isdefined(ModelingToolkit, :mtkcompile) - function mtkcompile(args...; kwargs...) - structural_simplify(args...; kwargs...) - end -end - +# Include all problem definitions include("ode_linear_prob.jl") +include("nonlinchem.jl") include("ode_simple_nonlinear_prob.jl") include("brusselator_prob.jl") include("pollution_prob.jl") include("filament_prob.jl") -include("nonlinchem.jl") include("strange_attractors.jl") -end # module +end # module \ No newline at end of file diff --git a/lib/ODEProblemLibrary/src/ode_simple_nonlinear_prob.jl b/lib/ODEProblemLibrary/src/ode_simple_nonlinear_prob.jl index 7f45788..c0be973 100644 --- a/lib/ODEProblemLibrary/src/ode_simple_nonlinear_prob.jl +++ b/lib/ODEProblemLibrary/src/ode_simple_nonlinear_prob.jl @@ -48,13 +48,15 @@ with initial condition ``v=w=1`` prob_ode_fitzhughnagumo = ODEProblem(fitz, [1.0; 1.0], (0.0, 1.0), (0.7, 0.8, 1 / 12.5, 0.5)) -#Van der Pol Equations -@parameters μ -@variables x(t) y(t) +## Van der Pol Equations -eqs = [D(y) ~ μ * ((1 - x^2) * y - x), - D(x) ~ y] -van = System(eqs, t; name = :van_der_pol) |> mtkcompile +function vanderpol(du, u, p, t) + x = u[1] + y = u[2] + μ = p[1] + du[1] = y + du[2] = μ * ((1 - x^2) * y - x) +end """ Van der Pol Equations @@ -66,12 +68,11 @@ Van der Pol Equations \\frac{dy}{dt} = μ((1-x^2)y -x) ``` -with ``μ=1.0`` and ``u_0=[x => \\sqrt{3}, y => 0]`` +with ``μ=1.0`` and ``u_0=[\\sqrt{3}, 0]`` (where ``u[1] = x``, ``u[2] = y``) Non-stiff parameters. """ -prob_ode_vanderpol = ODEProblem(van, [y => 0, x => sqrt(3), μ => 1.0], (0.0, 1.0), - jac = true, eval_module = @__MODULE__) +prob_ode_vanderpol = ODEProblem(vanderpol, [sqrt(3), 0.0], (0.0, 1.0), [1.0]) """ Van der Pol Equations @@ -83,21 +84,25 @@ Van der Pol Equations \\frac{dy}{dt} = μ((1-x^2)y -x) ``` -with ``μ=10^6`` and ``u_0=[x => \\sqrt{3}, y => 0]`` +with ``μ=10^6`` and ``u_0=[\\sqrt{3}, 0]`` (where ``u[1] = x``, ``u[2] = y``) Stiff parameters. """ -prob_ode_vanderpol_stiff = ODEProblem(van, [y => 0, x => sqrt(3), μ => 1e6], (0.0, 1.0), - jac = true, eval_module = @__MODULE__) - -# ROBER -@parameters k₁ k₂ k₃ -@variables y₁(t) y₂(t) y₃(t) - -eqs = [D(y₁) ~ -k₁ * y₁ + k₃ * y₂ * y₃, - D(y₂) ~ k₁ * y₁ - k₂ * y₂^2 - k₃ * y₂ * y₃, - D(y₃) ~ k₂ * y₂^2] -rober = System(eqs, t; name = :rober) |> mtkcompile +prob_ode_vanderpol_stiff = ODEProblem(vanderpol, [sqrt(3), 0.0], (0.0, 1.0), [1e6]) + +## ROBER + +function rober(du, u, p, t) + y₁ = u[1] + y₂ = u[2] + y₃ = u[3] + k₁ = p[1] + k₂ = p[2] + k₃ = p[3] + du[1] = -k₁ * y₁ + k₃ * y₂ * y₃ + du[2] = k₁ * y₁ - k₂ * y₂^2 - k₃ * y₂ * y₃ + du[3] = k₂ * y₂^2 +end """ The Robertson biochemical reactions: (Stiff) @@ -118,9 +123,7 @@ Hairer Norsett Wanner Solving Ordinary Differential Equations I - Nonstiff Probl Usually solved on ``[0,1e11]`` """ -prob_ode_rober = ODEProblem( - rober, [[y₁, y₂, y₃] .=> [1.0; 0.0; 0.0]; [k₁, k₂, k₃] .=> (0.04, 3e7, 1e4)], - (0.0, 1e11), jac = true, eval_module = @__MODULE__) +prob_ode_rober = ODEProblem(rober, [1.0, 0.0, 0.0], (0.0, 1e11), [0.04, 3e7, 1e4]) # Three Body const threebody_μ = big(0.012277471); @@ -174,15 +177,19 @@ prob_ode_threebody = ODEProblem(threebody, [0.994, 0.0, 0.0, big(-2.00158510637908252240537862224)], (big(0.0), big(17.0652165601579625588917206249))) -# Rigid Body Equations - -@parameters I₁ I₂ I₃ -@variables y₁(t) y₂(t) y₃(t) - -eqs = [D(y₁) ~ I₁ * y₂ * y₃, - D(y₂) ~ I₂ * y₁ * y₃, - D(y₃) ~ I₃ * y₁ * y₂] -rigid = System(eqs, t; name = :rigid_body) |> mtkcompile +## Rigid Body Equations + +function rigidbody(du, u, p, t) + y₁ = u[1] + y₂ = u[2] + y₃ = u[3] + I₁ = p[1] + I₂ = p[2] + I₃ = p[3] + du[1] = I₁ * y₂ * y₃ + du[2] = I₂ * y₁ * y₃ + du[3] = I₃ * y₁ * y₂ +end """ Rigid Body Equations (Non-stiff) @@ -207,10 +214,7 @@ or Hairer Norsett Wanner Solving Ordinary Differential Equations I - Nonstiff Pr Usually solved from 0 to 20. """ -prob_ode_rigidbody = ODEProblem( - rigid, [[y₁, y₂, y₃] .=> [1.0, 0.0, 0.9]; [I₁, I₂, I₃] .=> (-2.0, 1.25, -0.5)], - (0.0, 20.0), - jac = true, eval_module = @__MODULE__) +prob_ode_rigidbody = ODEProblem(rigidbody, [1.0, 0.0, 0.9], (0.0, 20.0), [-2.0, 1.25, -0.5]) # Pleiades Problem @@ -356,28 +360,27 @@ mm_f = ODEFunction(mm_linear; analytic = (u0, p, t) -> exp(inv(MM_linear) * mm_A mass_matrix = MM_linear) prob_ode_mm_linear = ODEProblem(mm_f, rand(4), (0.0, 1.0)) -@parameters p1 p2 p3 p4 p5 p6 p7 p8 p9 p10 p11 p12 -@variables y1(t) y2(t) y3(t) y4(t) y5(t) y6(t) y7(t) y8(t) - -eqs = [D(y1) ~ -p1 * y1 + p2 * y2 + p3 * y3 + p4, - D(y2) ~ p1 * y1 - p5 * y2, - D(y3) ~ -p6 * y3 + p2 * y4 + p7 * y5, - D(y4) ~ p3 * y2 + p1 * y3 - p8 * y4, - D(y5) ~ -p9 * y5 + p2 * y6 + p2 * y7, - D(y6) ~ -p10 * y6 * y8 + p11 * y4 + p1 * y5 - - p2 * y6 + p11 * y7, - D(y7) ~ p10 * y6 * y8 - p12 * y7, - D(y8) ~ -p10 * y6 * y8 + p12 * y7] -hires = System(eqs, t; name = :hires) |> mtkcompile +## Hires Problem + +function hires(du, u, p, t) + y1, y2, y3, y4, y5, y6, y7, y8 = u + p1, p2, p3, p4, p5, p6, p7, p8, p9, p10, p11, p12 = p + + du[1] = -p1 * y1 + p2 * y2 + p3 * y3 + p4 + du[2] = p1 * y1 - p5 * y2 + du[3] = -p6 * y3 + p2 * y4 + p7 * y5 + du[4] = p3 * y2 + p1 * y3 - p8 * y4 + du[5] = -p9 * y5 + p2 * y6 + p2 * y7 + du[6] = -p10 * y6 * y8 + p11 * y4 + p1 * y5 - p2 * y6 + p11 * y7 + du[7] = p10 * y6 * y8 - p12 * y7 + du[8] = -p10 * y6 * y8 + p12 * y7 +end u0 = zeros(8) u0[1] = 1 u0[8] = 0.0057 -u0 = [y1, y2, y3, y4, y5, y6, y7, y8] .=> u0 -p = [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) +p = (1.71, 0.43, 8.32, 0.0007, 8.75, 10.03, 0.035, 1.12, 1.745, 280.0, 0.69, 1.81) """ Hires Problem (Stiff) @@ -401,16 +404,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; p], (0.0, 321.8122), jac = true, eval_module = @__MODULE__) +prob_ode_hires = ODEProblem(hires, u0, (0.0, 321.8122), p) -@parameters p1 p2 p3 -@variables y1(t) y2(t) y3(t) +## Orego Problem -eqs = [D(y1) ~ p1 * (y2 + y1 * (1 - p2 * y1 - y2)), - D(y2) ~ (y3 - (1 + y1) * y2) / p1, - D(y3) ~ p3 * (y1 - y3)] -orego = System(eqs, t; name = :orego) |> mtkcompile +function orego(du, u, p, t) + y1, y2, y3 = u + p1, p2, p3 = p + + du[1] = p1 * (y2 + y1 * (1 - p2 * y1 - y2)) + du[2] = (y3 - (1 + y1) * y2) / p1 + du[3] = p3 * (y1 - y3) +end """ Orego Problem (Stiff) @@ -430,6 +435,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, [[y1, y2, y3] .=> [1.0, 2.0, 3.0]; [p1, p2, p3] .=> [77.27, 8.375e-6, 0.161]], - (0.0, 30.0), jac = true, eval_module = @__MODULE__) +prob_ode_orego = ODEProblem(orego, [1.0, 2.0, 3.0], (0.0, 30.0), [77.27, 8.375e-6, 0.161]) \ No newline at end of file diff --git a/lib/ODEProblemLibrary/src/strange_attractors.jl b/lib/ODEProblemLibrary/src/strange_attractors.jl index a47ec01..a3341fb 100644 --- a/lib/ODEProblemLibrary/src/strange_attractors.jl +++ b/lib/ODEProblemLibrary/src/strange_attractors.jl @@ -2,206 +2,259 @@ # https://en.wikipedia.org/wiki/List_of_chaotic_maps # Opted for the equations as reported in papers -# Thomas -@parameters b = 0.208186 -@variables x(t)=1 y(t)=0 z(t)=0 +## Thomas -eqs = [D(x) ~ sin(y) - b * x, - D(y) ~ sin(z) - b * y, - D(z) ~ sin(x) - b * z] - -@mtkcompile thomas = System(eqs, t) +function thomas_eqs(du, u, p, t) + x, y, z = u + b = p[1] + du[1] = sin(y) - b * x + du[2] = sin(z) - b * y + du[3] = sin(x) - b * z +end """ Thomas' cyclically symmetric attractor equations ```math -$(latexify(thomas; mult_symbol="")) +\\begin{align} +\\frac{dx}{dt} &= \\sin(y) - bx \\\\ +\\frac{dy}{dt} &= \\sin(z) - by \\\\ +\\frac{dz}{dt} &= \\sin(x) - bz +\\end{align} ``` +with parameter ``b = 0.208186`` and initial conditions ``x(0)=1, y(0)=0, z(0)=0`` + [Reference](https://www.worldscientific.com/doi/abs/10.1142/S0218127499001383) [Wikipedia](https://en.wikipedia.org/wiki/Thomas%27_cyclically_symmetric_attractor) """ -prob_ode_thomas = ODEProblem(thomas, [], (0.0, 1.0), jac = true, eval_module = @__MODULE__) +thomas = ODEFunction(thomas_eqs) +prob_ode_thomas = ODEProblem(thomas, [1.0, 0.0, 0.0], (0.0, 1.0), [0.208186]) -# Lorenz -@parameters σ=10 ρ=28 β=8/3 -@variables x(t)=1 y(t)=0 z(t)=0 +## Lorenz -eqs = [D(x) ~ σ * (y - x), - D(y) ~ x * (ρ - z) - y, - D(z) ~ x * y - β * z] - -@mtkcompile lorenz = System(eqs, t) +function lorenz_eqs(du, u, p, t) + x, y, z = u + σ, ρ, β = p + du[1] = σ * (y - x) + du[2] = x * (ρ - z) - y + du[3] = x * y - β * z +end """ Lorenz equations ```math -$(latexify(lorenz; mult_symbol="")) +\\begin{align} +\\frac{dx}{dt} &= σ(y - x) \\\\ +\\frac{dy}{dt} &= x(ρ - z) - y \\\\ +\\frac{dz}{dt} &= xy - βz +\\end{align} ``` +with parameters ``σ=10, ρ=28, β=8/3`` and initial conditions ``x(0)=1, y(0)=0, z(0)=0`` + [Reference](https://journals.ametsoc.org/view/journals/atsc/20/2/1520-0469_1963_020_0130_dnf_2_0_co_2.xml) [Wikipedia](https://en.wikipedia.org/wiki/Lorenz_system) """ -prob_ode_lorenz = ODEProblem(lorenz, [], (0.0, 1.0), jac = true, eval_module = @__MODULE__) +lorenz = ODEFunction(lorenz_eqs) +prob_ode_lorenz = ODEProblem(lorenz, [1.0, 0.0, 0.0], (0.0, 1.0), [10.0, 28.0, 8/3]) -# Aizawa -@parameters a=0.95 b=0.7 c=0.6 d=3.5 e=0.25 f=0.1 -@variables x(t)=1 y(t)=0 z(t)=0 +## Aizawa -eqs = [D(x) ~ (z - b) * x - d * y, - D(y) ~ d * x + (z - b) * y, - D(z) ~ c + a * z - z^3 / 3 - (x^2 + y^2) * (1 + e * z) + f * z * x^3] - -@mtkcompile aizawa = System(eqs, t) +function aizawa_eqs(du, u, p, t) + x, y, z = u + a, b, c, d, e, f = p + du[1] = (z - b) * x - d * y + du[2] = d * x + (z - b) * y + du[3] = c + a * z - z^3 / 3 - (x^2 + y^2) * (1 + e * z) + f * z * x^3 +end """ Aizawa equations ```math -$(latexify(aizawa; mult_symbol="")) +\\begin{align} +\\frac{dx}{dt} &= (z - b)x - dy \\\\ +\\frac{dy}{dt} &= dx + (z - b)y \\\\ +\\frac{dz}{dt} &= c + az - \\frac{z^3}{3} - (x^2 + y^2)(1 + ez) + fzx^3 +\\end{align} ``` -[Reference](https://journals.ametsoc.org/view/journals/atsc/20/2/1520-0469_1963_020_0130_dnf_2_0_co_2.xml) +with parameters ``a=0.95, b=0.7, c=0.6, d=3.5, e=0.25, f=0.1`` and initial conditions ``x(0)=1, y(0)=0, z(0)=0`` +[Reference](https://journals.ametsoc.org/view/journals/atsc/20/2/1520-0469_1963_020_0130_dnf_2_0_co_2.xml) """ -prob_ode_aizawa = ODEProblem(aizawa, [], (0.0, 1.0), jac = true, eval_module = @__MODULE__) +aizawa = ODEFunction(aizawa_eqs) +prob_ode_aizawa = ODEProblem(aizawa, [1.0, 0.0, 0.0], (0.0, 1.0), [0.95, 0.7, 0.6, 3.5, 0.25, 0.1]) -# Dadras -@parameters a=3 b=2.7 c=1.7 d=2 e=9 -@variables x(t)=1 y(t)=0 z(t)=0 +## Dadras -eqs = [D(x) ~ y - a * x + b * y * z, - D(y) ~ c * y - x * z + z, - D(z) ~ d * x * y - e * z] - -@mtkcompile dadras = System(eqs, t) +function dadras_eqs(du, u, p, t) + x, y, z = u + a, b, c, d, e = p + du[1] = y - a * x + b * y * z + du[2] = c * y - x * z + z + du[3] = d * x * y - e * z +end """ Dadras equations ```math -$(latexify(dadras; mult_symbol="")) +\\begin{align} +\\frac{dx}{dt} &= y - ax + byz \\\\ +\\frac{dy}{dt} &= cy - xz + z \\\\ +\\frac{dz}{dt} &= dxy - ez +\\end{align} ``` -[Reference](https://www.sciencedirect.com/science/article/abs/pii/S0375960109009591) +with parameters ``a=3, b=2.7, c=1.7, d=2, e=9`` and initial conditions ``x(0)=1, y(0)=0, z(0)=0`` +[Reference](https://www.sciencedirect.com/science/article/abs/pii/S0375960109009591) """ -prob_ode_dadras = ODEProblem(dadras, [], (0.0, 1.0), jac = true, eval_module = @__MODULE__) +dadras = ODEFunction(dadras_eqs) +prob_ode_dadras = ODEProblem(dadras, [1.0, 0.0, 0.0], (0.0, 1.0), [3.0, 2.7, 1.7, 2.0, 9.0]) -# chen -@parameters a=35 b=3 c=28 -@variables x(t)=1 y(t)=0 z(t)=0 +## Chen -eqs = [D(x) ~ a * (y - x), - D(y) ~ (c - a) * x - x * z + c * y, - D(z) ~ x * y - b * z] - -@mtkcompile chen = System(eqs, t) +function chen_eqs(du, u, p, t) + x, y, z = u + a, b, c = p + du[1] = a * (y - x) + du[2] = (c - a) * x - x * z + c * y + du[3] = x * y - b * z +end """ -chen equations +Chen equations ```math -$(latexify(chen; mult_symbol="")) +\\begin{align} +\\frac{dx}{dt} &= a(y - x) \\\\ +\\frac{dy}{dt} &= (c - a)x - xz + cy \\\\ +\\frac{dz}{dt} &= xy - bz +\\end{align} ``` -[Reference](https://www.worldscientific.com/doi/abs/10.1142/S0218127499001024) +with parameters ``a=35, b=3, c=28`` and initial conditions ``x(0)=1, y(0)=0, z(0)=0`` +[Reference](https://www.worldscientific.com/doi/abs/10.1142/S0218127499001024) """ -prob_ode_chen = ODEProblem(chen, [], (0.0, 1.0), jac = true, eval_module = @__MODULE__) +chen = ODEFunction(chen_eqs) +prob_ode_chen = ODEProblem(chen, [1.0, 0.0, 0.0], (0.0, 1.0), [35.0, 3.0, 28.0]) -# rossler -@parameters a=0.2 b=0.2 c=5.7 -@variables x(t)=1 y(t)=0 z(t)=0 +## Rössler -eqs = [D(x) ~ -(y + z), - D(y) ~ x + a * y, - D(z) ~ b + z * (x - c)] - -@mtkcompile rossler = System(eqs, t) +function rossler_eqs(du, u, p, t) + x, y, z = u + a, b, c = p + du[1] = -(y + z) + du[2] = x + a * y + du[3] = b + z * (x - c) +end """ -rossler equations +Rössler equations ```math -$(latexify(rossler; mult_symbol="")) +\\begin{align} +\\frac{dx}{dt} &= -(y + z) \\\\ +\\frac{dy}{dt} &= x + ay \\\\ +\\frac{dz}{dt} &= b + z(x - c) +\\end{align} ``` +with parameters ``a=0.2, b=0.2, c=5.7`` and initial conditions ``x(0)=1, y(0)=0, z(0)=0`` + [Reference](https://www.sciencedirect.com/science/article/abs/pii/0375960176901018) [Wikipedia](https://en.wikipedia.org/wiki/R%C3%B6ssler_attractor) - """ -prob_ode_rossler = ODEProblem( - rossler, [], (0.0, 1.0), jac = true, eval_module = @__MODULE__) +rossler = ODEFunction(rossler_eqs) +prob_ode_rossler = ODEProblem(rossler, [1.0, 0.0, 0.0], (0.0, 1.0), [0.2, 0.2, 5.7]) -# rabinovich_fabrikant -@parameters a=0.14 b=0.10 -@variables x(t)=1 y(t)=0 z(t)=0 +## Rabinovich-Fabrikant -eqs = [D(x) ~ y * (z - 1 + x^2) + b * x, - D(y) ~ x * (3 * z + 1 - x^2) + b * y, - D(z) ~ -2 * z * (a + x * y)] - -@mtkcompile rabinovich_fabrikant = System(eqs, t) +function rabinovich_fabrikant_eqs(du, u, p, t) + x, y, z = u + a, b = p + du[1] = y * (z - 1 + x^2) + b * x + du[2] = x * (3 * z + 1 - x^2) + b * y + du[3] = -2 * z * (a + x * y) +end """ -rabinovich_fabrikant equations +Rabinovich-Fabrikant equations ```math -$(latexify(rabinovich_fabrikant; mult_symbol="")) +\\begin{align} +\\frac{dx}{dt} &= y(z - 1 + x^2) + bx \\\\ +\\frac{dy}{dt} &= x(3z + 1 - x^2) + by \\\\ +\\frac{dz}{dt} &= -2z(a + xy) +\\end{align} ``` -[Reference](https://en.wikipedia.org/wiki/Rabinovich%E2%80%93Fabrikant_equations) +with parameters ``a=0.14, b=0.10`` and initial conditions ``x(0)=1, y(0)=0, z(0)=0`` +[Reference](https://en.wikipedia.org/wiki/Rabinovich%E2%80%93Fabrikant_equations) """ -prob_ode_rabinovich_fabrikant = ODEProblem( - rabinovich_fabrikant, [], (0.0, 1.0), jac = true, eval_module = @__MODULE__) +rabinovich_fabrikant = ODEFunction(rabinovich_fabrikant_eqs) +prob_ode_rabinovich_fabrikant = ODEProblem(rabinovich_fabrikant, [1.0, 0.0, 0.0], (0.0, 1.0), [0.14, 0.10]) -# sprott -@parameters a=2.07 b=1.79 -@variables x(t)=1 y(t)=0 z(t)=0 +## Sprott -eqs = [D(x) ~ y + a * x * y + x * z, - D(y) ~ 1 - b * x^2 + y * z, - D(z) ~ x - x^2 - y^2] - -@mtkcompile sprott = System(eqs, t) +function sprott_eqs(du, u, p, t) + x, y, z = u + a, b = p + du[1] = y + a * x * y + x * z + du[2] = 1 - b * x^2 + y * z + du[3] = x - x^2 - y^2 +end """ -sprott equations +Sprott equations ```math -$(latexify(sprott; mult_symbol="")) +\\begin{align} +\\frac{dx}{dt} &= y + axy + xz \\\\ +\\frac{dy}{dt} &= 1 - bx^2 + yz \\\\ +\\frac{dz}{dt} &= x - x^2 - y^2 +\\end{align} ``` -[Reference](https://sprott.physics.wisc.edu/pubs/paper423.pdf) +with parameters ``a=2.07, b=1.79`` and initial conditions ``x(0)=1, y(0)=0, z(0)=0`` +[Reference](https://sprott.physics.wisc.edu/pubs/paper423.pdf) """ -prob_ode_sprott = ODEProblem(sprott, [], (0.0, 1.0), jac = true, eval_module = @__MODULE__) +sprott = ODEFunction(sprott_eqs) +prob_ode_sprott = ODEProblem(sprott, [1.0, 0.0, 0.0], (0.0, 1.0), [2.07, 1.79]) -# hindmarsh_rose -@parameters a=1 b=3 c=1 d=5 r=1e-2 s=4 xr=-8/5 i=5 -@variables x(t)=1 y(t)=0 z(t)=0 +## Hindmarsh-Rose -eqs = [D(x) ~ y - a * x^3 + b * x^2 - z + i, - D(y) ~ c - d * x^2 - y, - D(z) ~ r * (s * (x - xr) - z)] - -@mtkcompile hindmarsh_rose = System(eqs, t) +function hindmarsh_rose_eqs(du, u, p, t) + x, y, z = u + a, b, c, d, r, s, xr, i = p + du[1] = y - a * x^3 + b * x^2 - z + i + du[2] = c - d * x^2 - y + du[3] = r * (s * (x - xr) - z) +end """ -hindmarsh_rose equations +Hindmarsh-Rose equations ```math -$(latexify(hindmarsh_rose; mult_symbol="")) +\\begin{align} +\\frac{dx}{dt} &= y - ax^3 + bx^2 - z + i \\\\ +\\frac{dy}{dt} &= c - dx^2 - y \\\\ +\\frac{dz}{dt} &= r(s(x - x_r) - z) +\\end{align} ``` -[Reference](https://en.wikipedia.org/wiki/Hindmarsh%E2%80%93Rose_model) +with parameters ``a=1, b=3, c=1, d=5, r=1e-2, s=4, x_r=-8/5, i=5`` and initial conditions ``x(0)=1, y(0)=0, z(0)=0`` +[Reference](https://en.wikipedia.org/wiki/Hindmarsh%E2%80%93Rose_model) """ -prob_ode_hindmarsh_rose = ODEProblem( - hindmarsh_rose, [], (0.0, 1.0), jac = true, eval_module = @__MODULE__) +hindmarsh_rose = ODEFunction(hindmarsh_rose_eqs) +prob_ode_hindmarsh_rose = ODEProblem(hindmarsh_rose, [1.0, 0.0, 0.0], (0.0, 1.0), [1.0, 3.0, 1.0, 5.0, 1e-2, 4.0, -8/5, 5.0]) \ No newline at end of file