From 6b0f03c215aea4e938078dd122f4340674176465 Mon Sep 17 00:00:00 2001 From: fchen121 Date: Thu, 29 May 2025 12:36:31 -0400 Subject: [PATCH 01/19] Created change of variable for SDE --- src/ModelingToolkit.jl | 2 +- src/systems/diffeqs/basic_transformations.jl | 114 +++++++++++++++++++ test/changeofvariables.jl | 96 ++++++++++++++++ 3 files changed, 211 insertions(+), 1 deletion(-) create mode 100644 test/changeofvariables.jl diff --git a/src/ModelingToolkit.jl b/src/ModelingToolkit.jl index 63468917d1..84547e5698 100644 --- a/src/ModelingToolkit.jl +++ b/src/ModelingToolkit.jl @@ -296,7 +296,7 @@ export isinput, isoutput, getbounds, hasbounds, getguess, hasguess, isdisturbanc hasunit, getunit, hasconnect, getconnect, hasmisc, getmisc, state_priority export liouville_transform, change_independent_variable, substitute_component, - add_accumulations, noise_to_brownians + add_accumulations, noise_to_brownians, changeofvariables, change_of_variable_SDE export PDESystem export Differential, expand_derivatives, @derivatives export Equation, ConstrainedEquation diff --git a/src/systems/diffeqs/basic_transformations.jl b/src/systems/diffeqs/basic_transformations.jl index 6243ca2405..b19b6911e5 100644 --- a/src/systems/diffeqs/basic_transformations.jl +++ b/src/systems/diffeqs/basic_transformations.jl @@ -53,6 +53,120 @@ function liouville_transform(sys::System; kwargs...) ) end +""" +$(TYPEDSIGNATURES) + +Generates the set of ODEs after change of variables. + + +Example: + +```julia +using ModelingToolkit, OrdinaryDiffEq, Test + +# Change of variables: z = log(x) +# (this implies that x = exp(z) is automatically non-negative) + +@parameters t α +@variables x(t) +D = Differential(t) +eqs = [D(x) ~ α*x] + +tspan = (0., 1.) +u0 = [x => 1.0] +p = [α => -0.5] + +@named sys = ODESystem(eqs; defaults=u0) +prob = ODEProblem(sys, [], tspan, p) +sol = solve(prob, Tsit5()) + +@variables z(t) +forward_subs = [log(x) => z] +backward_subs = [x => exp(z)] + +@named new_sys = changeofvariables(sys, forward_subs, backward_subs) +@test equations(new_sys)[1] == (D(z) ~ α) + +new_prob = ODEProblem(new_sys, [], tspan, p) +new_sol = solve(new_prob, Tsit5()) + +@test isapprox(new_sol[x][end], sol[x][end], atol=1e-4) +``` + +""" +function changeofvariables(sys::System, forward_subs, backward_subs; simplify=false, t0=missing) + t = independent_variable(sys) + + old_vars = first.(backward_subs) + new_vars = last.(forward_subs) + kept_vars = setdiff(states(sys), old_vars) + rhs = [eq.rhs for eq in equations(sys)] + + # use: dz/dt = ∂z/∂x dx/dt + ∂z/∂t + dzdt = Symbolics.derivative( first.(forward_subs), t ) + new_eqs = Equation[] + for (new_var, ex) in zip(new_vars, dzdt) + for ode_eq in equations(sys) + ex = substitute(ex, ode_eq.lhs => ode_eq.rhs) + end + ex = substitute(ex, Dict(forward_subs)) + ex = substitute(ex, Dict(backward_subs)) + if simplify + ex = Symbolics.simplify(ex, expand=true) + end + push!(new_eqs, Differential(t)(new_var) ~ ex) + end + + defs = get_defaults(sys) + new_defs = Dict() + for f_sub in forward_subs + #TODO call value(...)? + ex = substitute(first(f_sub), defs) + if !ismissing(t0) + ex = substitute(ex, t => t0) + end + new_defs[last(f_sub)] = ex + end + return ODESystem(new_eqs; + defaults=new_defs, + observed=vcat(observed(sys),first.(backward_subs) .~ last.(backward_subs)) + ) +end + +function change_of_variable_SDE(sys::System, forward_subs, backward_subs, iv; simplify=false, t0=missing) + t = independent_variable(sys) + + old_vars = first.(backward_subs) + new_vars = last.(forward_subs) + + # use: f = Y(t, X) + # use: dY = (∂f/∂t + μ∂f/∂x + (1/2)*σ^2*∂2f/∂x2)dt + σ∂f/∂xdW + old_eqs = equations(sys) + old_noise = get_noiseeqs(sys) + ∂f∂t = Symbolics.derivative( first.(forward_subs), t ) + ∂f∂x = Symbolics.derivative( first.(forward_subs), old_vars ) + ∂2f∂x2 = Symbolics.derivative( ∂f∂x, old_vars ) + new_eqs = Equation[] + + for (new_var, eq, noise, first, second, third) in zip(new_vars, old_eqs, old_noise, ∂f∂t, ∂f∂x, ∂2f∂x2) + ex = first + eq.rhs * second + 1/2 * noise^2 * third + for eqs in old_eqs + ex = substitute(ex, eqs.lhs => eqs.rhs) + end + ex = substitute(ex, Dict(forward_subs)) + ex = substitute(ex, Dict(backward_subs)) + if simplify + ex = Symbolics.simplify(ex, expand=true) + end + push!(new_eqs, Differential(t)(new_var) ~ ex) + end + new_noise = [noise * div for (noise, div) in zip(old_noise, ∂f∂x)] + + return SDESystem(new_eqs, new_noise; + observed=vcat(observed(sys),first.(backward_subs) .~ last.(backward_subs)) + ) +end + """ change_independent_variable( sys::System, iv, eqs = []; diff --git a/test/changeofvariables.jl b/test/changeofvariables.jl new file mode 100644 index 0000000000..f9a4568fa0 --- /dev/null +++ b/test/changeofvariables.jl @@ -0,0 +1,96 @@ +using ModelingToolkit, OrdinaryDiffEq +using Test, LinearAlgebra + + +# Change of variables: z = log(x) +# (this implies that x = exp(z) is automatically non-negative) + +@parameters t α +@variables x(t) +D = Differential(t) +eqs = [D(x) ~ α*x] + +tspan = (0., 1.) +u0 = [x => 1.0] +p = [α => -0.5] + +sys = ODESystem(eqs; defaults=u0) +prob = ODEProblem(sys, [], tspan, p) +sol = solve(prob, Tsit5()) + +@variables z(t) +forward_subs = [log(x) => z] +backward_subs = [x => exp(z)] +new_sys = changeofvariables(sys, forward_subs, backward_subs) +@test equations(new_sys)[1] == (D(z) ~ α) + +new_prob = ODEProblem(new_sys, [], tspan, p) +new_sol = solve(new_prob, Tsit5()) + +@test isapprox(new_sol[x][end], sol[x][end], atol=1e-4) + + + +# Riccati equation +@parameters t α +@variables x(t) +D = Differential(t) +eqs = [D(x) ~ t^2 + α - x^2] +sys = ODESystem(eqs, defaults=[x=>1.]) + +@variables z(t) +forward_subs = [t + α/(x+t) => z ] +backward_subs = [ x => α/(z-t) - t] + +new_sys = changeofvariables(sys, forward_subs, backward_subs; simplify=true, t0=0.) +# output should be equivalent to +# t^2 + α - z^2 + 2 (but this simplification is not found automatically) + +tspan = (0., 1.) +p = [α => 1.] +prob = ODEProblem(sys,[],tspan,p) +new_prob = ODEProblem(new_sys,[],tspan,p) + +sol = solve(prob, Tsit5()) +new_sol = solve(new_prob, Tsit5()) + +@test isapprox(sol[x][end], new_sol[x][end], rtol=1e-4) + + +# Linear transformation to diagonal system +@parameters t +@variables x[1:3](t) +D = Differential(t) +A = [0. -1. 0.; -0.5 0.5 0.; 0. 0. -1.] +eqs = D.(x) .~ A*x + +tspan = (0., 10.) +u0 = x .=> [1.0, 2.0, -1.0] + +sys = ODESystem(eqs; defaults=u0) +prob = ODEProblem(sys,[],tspan) +sol = solve(prob, Tsit5()) + +T = eigen(A).vectors + +@variables z[1:3](t) +forward_subs = T \ x .=> z +backward_subs = x .=> T*z + +new_sys = changeofvariables(sys, forward_subs, backward_subs; simplify=true) + +new_prob = ODEProblem(new_sys, [], tspan, p) +new_sol = solve(new_prob, Tsit5()) + +# test RHS +new_rhs = [eq.rhs for eq in equations(new_sys)] +new_A = Symbolics.value.(Symbolics.jacobian(new_rhs, z)) +@test isapprox(diagm(eigen(A).values), new_A, rtol = 1e-10) +@test isapprox( new_sol[x[1],end], sol[x[1],end], rtol=1e-4) + +# Change of variables for sde +@Browian B +@parameters μ σ +@variables x(t) y(t) + + From e0bb5d20c4a7ff7615b62442465fd6093dfb9648 Mon Sep 17 00:00:00 2001 From: fchen121 Date: Tue, 3 Jun 2025 15:54:55 -0400 Subject: [PATCH 02/19] Implement and fix change of variable for ODE --- Project.toml | 5 + src/systems/diffeqs/basic_transformations.jl | 36 +++++--- test/changeofvariables.jl | 96 +++++++++++--------- 3 files changed, 82 insertions(+), 55 deletions(-) diff --git a/Project.toml b/Project.toml index 1f51599a51..a6620c2f28 100644 --- a/Project.toml +++ b/Project.toml @@ -44,7 +44,9 @@ NaNMath = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3" NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec" OffsetArrays = "6fe1bfb0-de20-5000-8ca7-80f57d26f881" OrderedCollections = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" +OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed" OrdinaryDiffEqCore = "bbf590c4-e513-4bbe-9b18-05decba2e5d8" +OrdinaryDiffEqDefault = "50262376-6c5a-4cf5-baba-aaf4f84d72d7" PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a" RecursiveArrayTools = "731186ca-8d62-57ce-b412-fbd966d074cd" Reexport = "189a3867-3050-52da-a836-e630ba90ab69" @@ -59,9 +61,11 @@ SimpleNonlinearSolve = "727e6d20-b764-4bd8-a329-72de5adea6c7" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" +StochasticDiffEq = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0" SymbolicIndexingInterface = "2efcf032-c050-4f8e-a9bb-153293bab1f5" SymbolicUtils = "d1185830-fcd6-423d-90d6-eec64667417b" Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7" +Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" URIs = "5c2747f8-b7ea-4ff2-ba2e-563bfd36b1d4" UnPack = "3a884ed6-31ef-47d7-9d2a-63182c4928ed" Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" @@ -159,6 +163,7 @@ StochasticDiffEq = "6.72.1" SymbolicIndexingInterface = "0.3.39" SymbolicUtils = "3.26.1" Symbolics = "6.40" +Test = "1.11.0" URIs = "1" UnPack = "0.1, 1.0" Unitful = "1.1" diff --git a/src/systems/diffeqs/basic_transformations.jl b/src/systems/diffeqs/basic_transformations.jl index dedc7bc7a7..a0102891a5 100644 --- a/src/systems/diffeqs/basic_transformations.jl +++ b/src/systems/diffeqs/basic_transformations.jl @@ -94,13 +94,13 @@ new_sol = solve(new_prob, Tsit5()) ``` """ -function changeofvariables(sys::System, forward_subs, backward_subs; simplify=false, t0=missing) - t = independent_variable(sys) +function changeofvariables(sys::System, iv, forward_subs, backward_subs; simplify=true, t0=missing) + t = iv old_vars = first.(backward_subs) new_vars = last.(forward_subs) - kept_vars = setdiff(states(sys), old_vars) - rhs = [eq.rhs for eq in equations(sys)] + # kept_vars = setdiff(states(sys), old_vars) + # rhs = [eq.rhs for eq in equations(sys)] # use: dz/dt = ∂z/∂x dx/dt + ∂z/∂t dzdt = Symbolics.derivative( first.(forward_subs), t ) @@ -120,21 +120,29 @@ function changeofvariables(sys::System, forward_subs, backward_subs; simplify=fa defs = get_defaults(sys) new_defs = Dict() for f_sub in forward_subs - #TODO call value(...)? ex = substitute(first(f_sub), defs) if !ismissing(t0) ex = substitute(ex, t => t0) end new_defs[last(f_sub)] = ex end - return ODESystem(new_eqs; + for para in parameters(sys) + if haskey(defs, para) + new_defs[para] = defs[para] + end + end + @named new_sys = System(new_eqs, t; defaults=new_defs, observed=vcat(observed(sys),first.(backward_subs) .~ last.(backward_subs)) ) + if simplify + return mtkcompile(new_sys) + end + return new_sys end -function change_of_variable_SDE(sys::System, forward_subs, backward_subs, iv; simplify=false, t0=missing) - t = independent_variable(sys) +function change_of_variable_SDE(sys::System, iv, nvs, forward_subs, backward_subs; simplify=false, t0=missing) + t = iv old_vars = first.(backward_subs) new_vars = last.(forward_subs) @@ -142,14 +150,14 @@ function change_of_variable_SDE(sys::System, forward_subs, backward_subs, iv; si # use: f = Y(t, X) # use: dY = (∂f/∂t + μ∂f/∂x + (1/2)*σ^2*∂2f/∂x2)dt + σ∂f/∂xdW old_eqs = equations(sys) - old_noise = get_noiseeqs(sys) + old_noise = ModelingToolkit.get_noise_eqs(sys) ∂f∂t = Symbolics.derivative( first.(forward_subs), t ) - ∂f∂x = Symbolics.derivative( first.(forward_subs), old_vars ) + ∂f∂x = [Symbolics.derivative( first(f_sub), old_var )] ∂2f∂x2 = Symbolics.derivative( ∂f∂x, old_vars ) new_eqs = Equation[] - for (new_var, eq, noise, first, second, third) in zip(new_vars, old_eqs, old_noise, ∂f∂t, ∂f∂x, ∂2f∂x2) - ex = first + eq.rhs * second + 1/2 * noise^2 * third + for (new_var, eq, noise, nv, first, second, third) in zip(new_vars, old_eqs, old_noise, nvs, ∂f∂t, ∂f∂x, ∂2f∂x2) + ex = first + eq.rhs * second + 1/2 * noise^2 * third + noise*second*nv for eqs in old_eqs ex = substitute(ex, eqs.lhs => eqs.rhs) end @@ -160,11 +168,11 @@ function change_of_variable_SDE(sys::System, forward_subs, backward_subs, iv; si end push!(new_eqs, Differential(t)(new_var) ~ ex) end - new_noise = [noise * div for (noise, div) in zip(old_noise, ∂f∂x)] - return SDESystem(new_eqs, new_noise; + @named new_sys = System(new_eqs; observed=vcat(observed(sys),first.(backward_subs) .~ last.(backward_subs)) ) + return new_sys end """ diff --git a/test/changeofvariables.jl b/test/changeofvariables.jl index f9a4568fa0..1ce05c0ec5 100644 --- a/test/changeofvariables.jl +++ b/test/changeofvariables.jl @@ -1,30 +1,35 @@ -using ModelingToolkit, OrdinaryDiffEq +using ModelingToolkit, OrdinaryDiffEq, StochasticDiffEq using Test, LinearAlgebra # Change of variables: z = log(x) # (this implies that x = exp(z) is automatically non-negative) - -@parameters t α +@independent_variables t +# @variables z(t)[1:2, 1:2] +# D = Differential(t) +# eqs = [D(D(z)) ~ ones(2, 2)] +# @mtkcompile sys = System(eqs, t) +# @test_nowarn ODEProblem(sys, [z => zeros(2, 2), D(z) => ones(2, 2)], (0.0, 10.0)) + +@parameters α @variables x(t) D = Differential(t) eqs = [D(x) ~ α*x] tspan = (0., 1.) -u0 = [x => 1.0] -p = [α => -0.5] +def = [x => 1.0, α => -0.5] -sys = ODESystem(eqs; defaults=u0) -prob = ODEProblem(sys, [], tspan, p) +@mtkcompile sys = System(eqs, t;defaults=def) +prob = ODEProblem(sys, [], tspan) sol = solve(prob, Tsit5()) @variables z(t) forward_subs = [log(x) => z] backward_subs = [x => exp(z)] -new_sys = changeofvariables(sys, forward_subs, backward_subs) +new_sys = changeofvariables(sys, t, forward_subs, backward_subs) @test equations(new_sys)[1] == (D(z) ~ α) -new_prob = ODEProblem(new_sys, [], tspan, p) +new_prob = ODEProblem(new_sys, [], tspan) new_sol = solve(new_prob, Tsit5()) @test isapprox(new_sol[x][end], sol[x][end], atol=1e-4) @@ -32,24 +37,24 @@ new_sol = solve(new_prob, Tsit5()) # Riccati equation -@parameters t α +@parameters α @variables x(t) D = Differential(t) eqs = [D(x) ~ t^2 + α - x^2] -sys = ODESystem(eqs, defaults=[x=>1.]) +def = [x=>1., α => 1.] +@mtkcompile sys = System(eqs, t; defaults=def) @variables z(t) forward_subs = [t + α/(x+t) => z ] backward_subs = [ x => α/(z-t) - t] -new_sys = changeofvariables(sys, forward_subs, backward_subs; simplify=true, t0=0.) +new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=true, t0=0.) # output should be equivalent to # t^2 + α - z^2 + 2 (but this simplification is not found automatically) tspan = (0., 1.) -p = [α => 1.] -prob = ODEProblem(sys,[],tspan,p) -new_prob = ODEProblem(new_sys,[],tspan,p) +prob = ODEProblem(sys,[],tspan) +new_prob = ODEProblem(new_sys,[],tspan) sol = solve(prob, Tsit5()) new_sol = solve(new_prob, Tsit5()) @@ -57,40 +62,49 @@ new_sol = solve(new_prob, Tsit5()) @test isapprox(sol[x][end], new_sol[x][end], rtol=1e-4) -# Linear transformation to diagonal system -@parameters t -@variables x[1:3](t) -D = Differential(t) -A = [0. -1. 0.; -0.5 0.5 0.; 0. 0. -1.] -eqs = D.(x) .~ A*x +# # Linear transformation to diagonal system +# @variables x(t)[1:3] +# D = Differential(t) +# A = [0. -1. 0.; -0.5 0.5 0.; 0. 0. -1.] +# right = A.*transpose(x) +# eqs = [D(x[1]) ~ sum(right[1, 1:3]), D(x[2]) ~ sum(right[2, 1:3]), D(x[3]) ~ sum(right[3, 1:3])] -tspan = (0., 10.) -u0 = x .=> [1.0, 2.0, -1.0] +# tspan = (0., 10.) +# u0 = [x[1] => 1.0, x[2] => 2.0, x[3] => -1.0] -sys = ODESystem(eqs; defaults=u0) -prob = ODEProblem(sys,[],tspan) -sol = solve(prob, Tsit5()) +# @mtkcompile sys = System(eqs, t; defaults=u0) +# prob = ODEProblem(sys,[],tspan) +# sol = solve(prob, Tsit5()) -T = eigen(A).vectors +# T = eigen(A).vectors -@variables z[1:3](t) -forward_subs = T \ x .=> z -backward_subs = x .=> T*z +# @variables z(t)[1:3] +# forward_subs = T \ x .=> z +# backward_subs = x .=> T*z -new_sys = changeofvariables(sys, forward_subs, backward_subs; simplify=true) +# new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=true) -new_prob = ODEProblem(new_sys, [], tspan, p) -new_sol = solve(new_prob, Tsit5()) +# new_prob = ODEProblem(new_sys, [], tspan) +# new_sol = solve(new_prob, Tsit5()) -# test RHS -new_rhs = [eq.rhs for eq in equations(new_sys)] -new_A = Symbolics.value.(Symbolics.jacobian(new_rhs, z)) -@test isapprox(diagm(eigen(A).values), new_A, rtol = 1e-10) -@test isapprox( new_sol[x[1],end], sol[x[1],end], rtol=1e-4) +# # test RHS +# new_rhs = [eq.rhs for eq in equations(new_sys)] +# new_A = Symbolics.value.(Symbolics.jacobian(new_rhs, z)) +# @test isapprox(diagm(eigen(A).values), new_A, rtol = 1e-10) +# @test isapprox( new_sol[x[1],end], sol[x[1],end], rtol=1e-4) # Change of variables for sde -@Browian B -@parameters μ σ -@variables x(t) y(t) +# @independent_variables t +# @brownian B +# @parameters μ σ +# @variables x(t) y(t) +# D = Differential(t) +# eqs = [D(x) ~ μ*x + σ*x*B] + +# def = [x=>0., μ => 2., σ=>1.] +# @mtkcompile sys = System(eqs, t; defaults=def) +# forward_subs = [log(x) => y] +# backward_subs = [x => exp(y)] +# new_sys = change_of_variable_SDE(sys, t, [B], forward_subs, backward_subs) From 22ab62432b85a20de0c56f1b6025d50488f3c6f8 Mon Sep 17 00:00:00 2001 From: fchen121 Date: Wed, 4 Jun 2025 17:44:24 -0400 Subject: [PATCH 03/19] Implement change of variables for sde --- src/systems/diffeqs/basic_transformations.jl | 37 ++++++++++++++++---- test/changeofvariables.jl | 34 +++++++++--------- 2 files changed, 49 insertions(+), 22 deletions(-) diff --git a/src/systems/diffeqs/basic_transformations.jl b/src/systems/diffeqs/basic_transformations.jl index a0102891a5..fbf03eacda 100644 --- a/src/systems/diffeqs/basic_transformations.jl +++ b/src/systems/diffeqs/basic_transformations.jl @@ -141,8 +141,10 @@ function changeofvariables(sys::System, iv, forward_subs, backward_subs; simplif return new_sys end -function change_of_variable_SDE(sys::System, iv, nvs, forward_subs, backward_subs; simplify=false, t0=missing) +function change_of_variable_SDE(sys::System, iv, forward_subs, backward_subs; simplify=true, t0=missing) + sys = mtkcompile(sys) t = iv + @brownian B old_vars = first.(backward_subs) new_vars = last.(forward_subs) @@ -151,13 +153,17 @@ function change_of_variable_SDE(sys::System, iv, nvs, forward_subs, backward_sub # use: dY = (∂f/∂t + μ∂f/∂x + (1/2)*σ^2*∂2f/∂x2)dt + σ∂f/∂xdW old_eqs = equations(sys) old_noise = ModelingToolkit.get_noise_eqs(sys) + + # Is there a function to find partial derivative? ∂f∂t = Symbolics.derivative( first.(forward_subs), t ) - ∂f∂x = [Symbolics.derivative( first(f_sub), old_var )] - ∂2f∂x2 = Symbolics.derivative( ∂f∂x, old_vars ) + ∂f∂t = [substitute(f_t, Differential(t)(old_var) => 0) for (f_t, old_var) in zip(∂f∂t, old_vars)] + + ∂f∂x = [Symbolics.derivative( first(f_sub), old_var ) for (f_sub, old_var) in zip(forward_subs, old_vars)] + ∂2f∂x2 = Symbolics.derivative.( ∂f∂x, old_vars ) new_eqs = Equation[] - for (new_var, eq, noise, nv, first, second, third) in zip(new_vars, old_eqs, old_noise, nvs, ∂f∂t, ∂f∂x, ∂2f∂x2) - ex = first + eq.rhs * second + 1/2 * noise^2 * third + noise*second*nv + for (new_var, eq, noise, first, second, third) in zip(new_vars, old_eqs, old_noise, ∂f∂t, ∂f∂x, ∂2f∂x2) + ex = first + eq.rhs * second + 1/2 * noise^2 * third + noise*second*B for eqs in old_eqs ex = substitute(ex, eqs.lhs => eqs.rhs) end @@ -169,9 +175,28 @@ function change_of_variable_SDE(sys::System, iv, nvs, forward_subs, backward_sub push!(new_eqs, Differential(t)(new_var) ~ ex) end - @named new_sys = System(new_eqs; + defs = get_defaults(sys) + new_defs = Dict() + for f_sub in forward_subs + ex = substitute(first(f_sub), defs) + if !ismissing(t0) + ex = substitute(ex, t => t0) + end + new_defs[last(f_sub)] = ex + end + for para in parameters(sys) + if haskey(defs, para) + new_defs[para] = defs[para] + end + end + + @named new_sys = System(new_eqs, t; + defaults=new_defs, observed=vcat(observed(sys),first.(backward_subs) .~ last.(backward_subs)) ) + if simplify + return mtkcompile(new_sys) + end return new_sys end diff --git a/test/changeofvariables.jl b/test/changeofvariables.jl index 1ce05c0ec5..2667e143f2 100644 --- a/test/changeofvariables.jl +++ b/test/changeofvariables.jl @@ -5,11 +5,11 @@ using Test, LinearAlgebra # Change of variables: z = log(x) # (this implies that x = exp(z) is automatically non-negative) @independent_variables t -# @variables z(t)[1:2, 1:2] -# D = Differential(t) -# eqs = [D(D(z)) ~ ones(2, 2)] -# @mtkcompile sys = System(eqs, t) -# @test_nowarn ODEProblem(sys, [z => zeros(2, 2), D(z) => ones(2, 2)], (0.0, 10.0)) +@variables z(t)[1:2, 1:2] +D = Differential(t) +eqs = [D(D(z)) ~ ones(2, 2)] +@mtkcompile sys = System(eqs, t) +@test_nowarn ODEProblem(sys, [z => zeros(2, 2), D(z) => ones(2, 2)], (0.0, 10.0)) @parameters α @variables x(t) @@ -94,17 +94,19 @@ new_sol = solve(new_prob, Tsit5()) # @test isapprox( new_sol[x[1],end], sol[x[1],end], rtol=1e-4) # Change of variables for sde -# @independent_variables t -# @brownian B -# @parameters μ σ -# @variables x(t) y(t) -# D = Differential(t) -# eqs = [D(x) ~ μ*x + σ*x*B] +@independent_variables t +@brownian B +@parameters μ σ +@variables x(t) y(t) +D = Differential(t) +eqs = [D(x) ~ μ*x + σ*x*B] -# def = [x=>0., μ => 2., σ=>1.] -# @mtkcompile sys = System(eqs, t; defaults=def) -# forward_subs = [log(x) => y] -# backward_subs = [x => exp(y)] -# new_sys = change_of_variable_SDE(sys, t, [B], forward_subs, backward_subs) +def = [x=>0., μ => 2., σ=>1.] +@mtkcompile sys = System(eqs, t; defaults=def) +forward_subs = [log(x) => y] +backward_subs = [x => exp(y)] +new_sys = change_of_variable_SDE(sys, t, forward_subs, backward_subs) +@test equations(new_sys)[1] == (D(y) ~ μ - 1/2*σ^2) +@test ModelingToolkit.get_noise_eqs(new_sys)[1] === ModelingToolkit.value(σ) From 50a4fc3d5b493584d7db42f746824ba2eae589b2 Mon Sep 17 00:00:00 2001 From: fchen121 Date: Thu, 5 Jun 2025 12:33:07 -0400 Subject: [PATCH 04/19] Fix Linear transformation to diagonal system test --- test/changeofvariables.jl | 53 ++++++++++++++++++++++----------------- 1 file changed, 30 insertions(+), 23 deletions(-) diff --git a/test/changeofvariables.jl b/test/changeofvariables.jl index 2667e143f2..69f79768e3 100644 --- a/test/changeofvariables.jl +++ b/test/changeofvariables.jl @@ -62,36 +62,43 @@ new_sol = solve(new_prob, Tsit5()) @test isapprox(sol[x][end], new_sol[x][end], rtol=1e-4) -# # Linear transformation to diagonal system -# @variables x(t)[1:3] -# D = Differential(t) -# A = [0. -1. 0.; -0.5 0.5 0.; 0. 0. -1.] -# right = A.*transpose(x) -# eqs = [D(x[1]) ~ sum(right[1, 1:3]), D(x[2]) ~ sum(right[2, 1:3]), D(x[3]) ~ sum(right[3, 1:3])] +# Linear transformation to diagonal system +@independent_variables t +@variables x(t)[1:3] +x = reshape(x, 3, 1) +D = Differential(t) +A = [0. -1. 0.; -0.5 0.5 0.; 0. 0. -1.] +right = A*x +eqs = vec(D.(x) .~ right) -# tspan = (0., 10.) -# u0 = [x[1] => 1.0, x[2] => 2.0, x[3] => -1.0] +tspan = (0., 10.) +u0 = [x[1] => 1.0, x[2] => 2.0, x[3] => -1.0] -# @mtkcompile sys = System(eqs, t; defaults=u0) -# prob = ODEProblem(sys,[],tspan) -# sol = solve(prob, Tsit5()) +@mtkcompile sys = System(eqs, t; defaults=u0) +prob = ODEProblem(sys,[],tspan) +sol = solve(prob, Tsit5()) -# T = eigen(A).vectors +T = eigen(A).vectors +T_inv = inv(T) -# @variables z(t)[1:3] -# forward_subs = T \ x .=> z -# backward_subs = x .=> T*z +@variables z(t)[1:3] +z = reshape(z, 3, 1) +forward_subs = vec(T_inv*x .=> z) +backward_subs = vec(x .=> T*z) -# new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=true) +new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=true) -# new_prob = ODEProblem(new_sys, [], tspan) -# new_sol = solve(new_prob, Tsit5()) +new_prob = ODEProblem(new_sys, [], tspan) +new_sol = solve(new_prob, Tsit5()) -# # test RHS -# new_rhs = [eq.rhs for eq in equations(new_sys)] -# new_A = Symbolics.value.(Symbolics.jacobian(new_rhs, z)) -# @test isapprox(diagm(eigen(A).values), new_A, rtol = 1e-10) -# @test isapprox( new_sol[x[1],end], sol[x[1],end], rtol=1e-4) +# test RHS +new_rhs = [eq.rhs for eq in equations(new_sys)] +new_A = Symbolics.value.(Symbolics.jacobian(new_rhs, z)) +A = diagm(eigen(A).values) +A = sortslices(A, dims=1) +new_A = sortslices(new_A, dims=1) +@test isapprox(A, new_A, rtol = 1e-10) +@test isapprox( new_sol[x[1],end], sol[x[1],end], rtol=1e-4) # Change of variables for sde @independent_variables t From a1e99446edd556f92923c1c353675adb6a36b1c2 Mon Sep 17 00:00:00 2001 From: fchen121 Date: Thu, 5 Jun 2025 12:36:16 -0400 Subject: [PATCH 05/19] Update Project.toml --- Project.toml | 2 ++ 1 file changed, 2 insertions(+) diff --git a/Project.toml b/Project.toml index 6632f626a2..525e14534b 100644 --- a/Project.toml +++ b/Project.toml @@ -65,6 +65,7 @@ StochasticDiffEq = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0" SymbolicIndexingInterface = "2efcf032-c050-4f8e-a9bb-153293bab1f5" SymbolicUtils = "d1185830-fcd6-423d-90d6-eec64667417b" Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7" +TermInterface = "8ea1fca8-c5ef-4a55-8b96-4e9afe9c9a3c" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" URIs = "5c2747f8-b7ea-4ff2-ba2e-563bfd36b1d4" UnPack = "3a884ed6-31ef-47d7-9d2a-63182c4928ed" @@ -163,6 +164,7 @@ StochasticDiffEq = "6.72.1" SymbolicIndexingInterface = "0.3.39" SymbolicUtils = "3.26.1" Symbolics = "6.40" +TermInterface = "2.0.0" Test = "1.11.0" URIs = "1" UnPack = "0.1, 1.0" From e33bcc1f6ac5baa94278ca1bf520f9224977e49b Mon Sep 17 00:00:00 2001 From: fchen121 Date: Thu, 5 Jun 2025 13:39:14 -0400 Subject: [PATCH 06/19] Change backward subs from observed to equations --- src/systems/diffeqs/basic_transformations.jl | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/systems/diffeqs/basic_transformations.jl b/src/systems/diffeqs/basic_transformations.jl index fbf03eacda..3061ae20ea 100644 --- a/src/systems/diffeqs/basic_transformations.jl +++ b/src/systems/diffeqs/basic_transformations.jl @@ -131,9 +131,9 @@ function changeofvariables(sys::System, iv, forward_subs, backward_subs; simplif new_defs[para] = defs[para] end end - @named new_sys = System(new_eqs, t; + @named new_sys = System(vcat(new_eqs, first.(backward_subs) .~ last.(backward_subs)), t; defaults=new_defs, - observed=vcat(observed(sys),first.(backward_subs) .~ last.(backward_subs)) + observed=observed(sys) ) if simplify return mtkcompile(new_sys) @@ -190,9 +190,9 @@ function change_of_variable_SDE(sys::System, iv, forward_subs, backward_subs; si end end - @named new_sys = System(new_eqs, t; + @named new_sys = System(vcat(new_eqs, first.(backward_subs) .~ last.(backward_subs)), t; defaults=new_defs, - observed=vcat(observed(sys),first.(backward_subs) .~ last.(backward_subs)) + observed=observed(sys) ) if simplify return mtkcompile(new_sys) From 192872ce19f3fc0cb6e2af2c99b9be2bef2ca939 Mon Sep 17 00:00:00 2001 From: fchen121 Date: Thu, 5 Jun 2025 17:54:27 -0400 Subject: [PATCH 07/19] Change of variables for multiple Brownian SDE --- src/systems/diffeqs/basic_transformations.jl | 24 +-- test/changeofvariables.jl | 157 ++++++++++--------- 2 files changed, 99 insertions(+), 82 deletions(-) diff --git a/src/systems/diffeqs/basic_transformations.jl b/src/systems/diffeqs/basic_transformations.jl index 3061ae20ea..26131554ed 100644 --- a/src/systems/diffeqs/basic_transformations.jl +++ b/src/systems/diffeqs/basic_transformations.jl @@ -144,7 +144,6 @@ end function change_of_variable_SDE(sys::System, iv, forward_subs, backward_subs; simplify=true, t0=missing) sys = mtkcompile(sys) t = iv - @brownian B old_vars = first.(backward_subs) new_vars = last.(forward_subs) @@ -152,20 +151,27 @@ function change_of_variable_SDE(sys::System, iv, forward_subs, backward_subs; si # use: f = Y(t, X) # use: dY = (∂f/∂t + μ∂f/∂x + (1/2)*σ^2*∂2f/∂x2)dt + σ∂f/∂xdW old_eqs = equations(sys) - old_noise = ModelingToolkit.get_noise_eqs(sys) + neqs = get_noise_eqs(sys) + neqs = [neqs[i,:] for i in 1:size(neqs,1)] - # Is there a function to find partial derivative? - ∂f∂t = Symbolics.derivative( first.(forward_subs), t ) - ∂f∂t = [substitute(f_t, Differential(t)(old_var) => 0) for (f_t, old_var) in zip(∂f∂t, old_vars)] + brownvars = map([Symbol(:B, :_, i) for i in 1:length(neqs[1])]) do name + unwrap(only(@brownian $name)) + end + # df/dt = ∂f/∂x dx/dt + ∂f/∂t + dfdt = Symbolics.derivative( first.(forward_subs), t ) ∂f∂x = [Symbolics.derivative( first(f_sub), old_var ) for (f_sub, old_var) in zip(forward_subs, old_vars)] ∂2f∂x2 = Symbolics.derivative.( ∂f∂x, old_vars ) new_eqs = Equation[] - for (new_var, eq, noise, first, second, third) in zip(new_vars, old_eqs, old_noise, ∂f∂t, ∂f∂x, ∂2f∂x2) - ex = first + eq.rhs * second + 1/2 * noise^2 * third + noise*second*B - for eqs in old_eqs - ex = substitute(ex, eqs.lhs => eqs.rhs) + for (new_var, ex, first, second) in zip(new_vars, dfdt, ∂f∂x, ∂2f∂x2) + for (eqs, neq) in zip(old_eqs, neqs) + if occursin(value(eqs.lhs), value(ex)) + ex = substitute(ex, eqs.lhs => eqs.rhs) + for (noise, B) in zip(neq, brownvars) + ex = ex + 1/2 * noise^2 * second + noise*first*B + end + end end ex = substitute(ex, Dict(forward_subs)) ex = substitute(ex, Dict(backward_subs)) diff --git a/test/changeofvariables.jl b/test/changeofvariables.jl index 69f79768e3..d5eb7d8dd4 100644 --- a/test/changeofvariables.jl +++ b/test/changeofvariables.jl @@ -4,101 +4,101 @@ using Test, LinearAlgebra # Change of variables: z = log(x) # (this implies that x = exp(z) is automatically non-negative) -@independent_variables t -@variables z(t)[1:2, 1:2] -D = Differential(t) -eqs = [D(D(z)) ~ ones(2, 2)] -@mtkcompile sys = System(eqs, t) -@test_nowarn ODEProblem(sys, [z => zeros(2, 2), D(z) => ones(2, 2)], (0.0, 10.0)) +# @independent_variables t +# @variables z(t)[1:2, 1:2] +# D = Differential(t) +# eqs = [D(D(z)) ~ ones(2, 2)] +# @mtkcompile sys = System(eqs, t) +# @test_nowarn ODEProblem(sys, [z => zeros(2, 2), D(z) => ones(2, 2)], (0.0, 10.0)) -@parameters α -@variables x(t) -D = Differential(t) -eqs = [D(x) ~ α*x] +# @parameters α +# @variables x(t) +# D = Differential(t) +# eqs = [D(x) ~ α*x] -tspan = (0., 1.) -def = [x => 1.0, α => -0.5] +# tspan = (0., 1.) +# def = [x => 1.0, α => -0.5] -@mtkcompile sys = System(eqs, t;defaults=def) -prob = ODEProblem(sys, [], tspan) -sol = solve(prob, Tsit5()) +# @mtkcompile sys = System(eqs, t;defaults=def) +# prob = ODEProblem(sys, [], tspan) +# sol = solve(prob, Tsit5()) -@variables z(t) -forward_subs = [log(x) => z] -backward_subs = [x => exp(z)] -new_sys = changeofvariables(sys, t, forward_subs, backward_subs) -@test equations(new_sys)[1] == (D(z) ~ α) +# @variables z(t) +# forward_subs = [log(x) => z] +# backward_subs = [x => exp(z)] +# new_sys = changeofvariables(sys, t, forward_subs, backward_subs) +# @test equations(new_sys)[1] == (D(z) ~ α) -new_prob = ODEProblem(new_sys, [], tspan) -new_sol = solve(new_prob, Tsit5()) +# new_prob = ODEProblem(new_sys, [], tspan) +# new_sol = solve(new_prob, Tsit5()) -@test isapprox(new_sol[x][end], sol[x][end], atol=1e-4) +# @test isapprox(new_sol[x][end], sol[x][end], atol=1e-4) -# Riccati equation -@parameters α -@variables x(t) -D = Differential(t) -eqs = [D(x) ~ t^2 + α - x^2] -def = [x=>1., α => 1.] -@mtkcompile sys = System(eqs, t; defaults=def) +# # Riccati equation +# @parameters α +# @variables x(t) +# D = Differential(t) +# eqs = [D(x) ~ t^2 + α - x^2] +# def = [x=>1., α => 1.] +# @mtkcompile sys = System(eqs, t; defaults=def) -@variables z(t) -forward_subs = [t + α/(x+t) => z ] -backward_subs = [ x => α/(z-t) - t] +# @variables z(t) +# forward_subs = [t + α/(x+t) => z ] +# backward_subs = [ x => α/(z-t) - t] -new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=true, t0=0.) -# output should be equivalent to -# t^2 + α - z^2 + 2 (but this simplification is not found automatically) +# new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=true, t0=0.) +# # output should be equivalent to +# # t^2 + α - z^2 + 2 (but this simplification is not found automatically) -tspan = (0., 1.) -prob = ODEProblem(sys,[],tspan) -new_prob = ODEProblem(new_sys,[],tspan) +# tspan = (0., 1.) +# prob = ODEProblem(sys,[],tspan) +# new_prob = ODEProblem(new_sys,[],tspan) -sol = solve(prob, Tsit5()) -new_sol = solve(new_prob, Tsit5()) +# sol = solve(prob, Tsit5()) +# new_sol = solve(new_prob, Tsit5()) -@test isapprox(sol[x][end], new_sol[x][end], rtol=1e-4) +# @test isapprox(sol[x][end], new_sol[x][end], rtol=1e-4) -# Linear transformation to diagonal system -@independent_variables t -@variables x(t)[1:3] -x = reshape(x, 3, 1) -D = Differential(t) -A = [0. -1. 0.; -0.5 0.5 0.; 0. 0. -1.] -right = A*x -eqs = vec(D.(x) .~ right) +# # Linear transformation to diagonal system +# @independent_variables t +# @variables x(t)[1:3] +# x = reshape(x, 3, 1) +# D = Differential(t) +# A = [0. -1. 0.; -0.5 0.5 0.; 0. 0. -1.] +# right = A*x +# eqs = vec(D.(x) .~ right) -tspan = (0., 10.) -u0 = [x[1] => 1.0, x[2] => 2.0, x[3] => -1.0] +# tspan = (0., 10.) +# u0 = [x[1] => 1.0, x[2] => 2.0, x[3] => -1.0] -@mtkcompile sys = System(eqs, t; defaults=u0) -prob = ODEProblem(sys,[],tspan) -sol = solve(prob, Tsit5()) +# @mtkcompile sys = System(eqs, t; defaults=u0) +# prob = ODEProblem(sys,[],tspan) +# sol = solve(prob, Tsit5()) -T = eigen(A).vectors -T_inv = inv(T) +# T = eigen(A).vectors +# T_inv = inv(T) -@variables z(t)[1:3] -z = reshape(z, 3, 1) -forward_subs = vec(T_inv*x .=> z) -backward_subs = vec(x .=> T*z) +# @variables z(t)[1:3] +# z = reshape(z, 3, 1) +# forward_subs = vec(T_inv*x .=> z) +# backward_subs = vec(x .=> T*z) -new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=true) +# new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=true) -new_prob = ODEProblem(new_sys, [], tspan) -new_sol = solve(new_prob, Tsit5()) +# new_prob = ODEProblem(new_sys, [], tspan) +# new_sol = solve(new_prob, Tsit5()) -# test RHS -new_rhs = [eq.rhs for eq in equations(new_sys)] -new_A = Symbolics.value.(Symbolics.jacobian(new_rhs, z)) -A = diagm(eigen(A).values) -A = sortslices(A, dims=1) -new_A = sortslices(new_A, dims=1) -@test isapprox(A, new_A, rtol = 1e-10) -@test isapprox( new_sol[x[1],end], sol[x[1],end], rtol=1e-4) +# # test RHS +# new_rhs = [eq.rhs for eq in equations(new_sys)] +# new_A = Symbolics.value.(Symbolics.jacobian(new_rhs, z)) +# A = diagm(eigen(A).values) +# A = sortslices(A, dims=1) +# new_A = sortslices(new_A, dims=1) +# @test isapprox(A, new_A, rtol = 1e-10) +# @test isapprox( new_sol[x[1],end], sol[x[1],end], rtol=1e-4) # Change of variables for sde @independent_variables t @@ -116,4 +116,15 @@ new_sys = change_of_variable_SDE(sys, t, forward_subs, backward_subs) @test equations(new_sys)[1] == (D(y) ~ μ - 1/2*σ^2) @test ModelingToolkit.get_noise_eqs(new_sys)[1] === ModelingToolkit.value(σ) - +#Multiple Brownian and equations +@independent_variables t +@brownian Bx By +@parameters μ σ α +@variables x(t) y(t) z(t) w(t) u(t) v(t) +D = Differential(t) +eqs = [D(x) ~ μ*x + σ*x*Bx, D(y) ~ α*By, D(u) ~ μ*u + σ*u*Bx + α*u*By] +def = [x=>0., y=> 0., u=>0., μ => 2., σ=>1., α=>3.] +@mtkcompile sys = System(eqs, t; defaults=def) +forward_subs = [log(x) => z, y^2 => w, log(u) => v] +backward_subs = [x => exp(z), y => w^.5, u => exp(v)] +new_sys = change_of_variable_SDE(sys, t, forward_subs, backward_subs) \ No newline at end of file From 9a5312189533952593cc29d5bd40a5d19480022a Mon Sep 17 00:00:00 2001 From: fchen121 Date: Thu, 5 Jun 2025 18:21:53 -0400 Subject: [PATCH 08/19] Improve change of variables test --- test/changeofvariables.jl | 16 ++++++++++++++-- 1 file changed, 14 insertions(+), 2 deletions(-) diff --git a/test/changeofvariables.jl b/test/changeofvariables.jl index d5eb7d8dd4..68a8a143b6 100644 --- a/test/changeofvariables.jl +++ b/test/changeofvariables.jl @@ -101,6 +101,9 @@ using Test, LinearAlgebra # @test isapprox( new_sol[x[1],end], sol[x[1],end], rtol=1e-4) # Change of variables for sde +noise_eqs = ModelingToolkit.get_noise_eqs +value = ModelingToolkit.value + @independent_variables t @brownian B @parameters μ σ @@ -114,7 +117,7 @@ forward_subs = [log(x) => y] backward_subs = [x => exp(y)] new_sys = change_of_variable_SDE(sys, t, forward_subs, backward_subs) @test equations(new_sys)[1] == (D(y) ~ μ - 1/2*σ^2) -@test ModelingToolkit.get_noise_eqs(new_sys)[1] === ModelingToolkit.value(σ) +@test noise_eqs(new_sys)[1] === value(σ) #Multiple Brownian and equations @independent_variables t @@ -127,4 +130,13 @@ def = [x=>0., y=> 0., u=>0., μ => 2., σ=>1., α=>3.] @mtkcompile sys = System(eqs, t; defaults=def) forward_subs = [log(x) => z, y^2 => w, log(u) => v] backward_subs = [x => exp(z), y => w^.5, u => exp(v)] -new_sys = change_of_variable_SDE(sys, t, forward_subs, backward_subs) \ No newline at end of file +new_sys = change_of_variable_SDE(sys, t, forward_subs, backward_subs) +@test equations(new_sys)[1] == (D(z) ~ μ - 1/2*σ^2) +@test equations(new_sys)[2] == (D(w) ~ α^2) +@test equations(new_sys)[3] == (D(v) ~ μ - 1/2*(α^2 + σ^2)) +@test noise_eqs(new_sys)[1,1] === value(σ) +@test noise_eqs(new_sys)[1,2] === value(0) +@test noise_eqs(new_sys)[2,1] === value(0) +@test noise_eqs(new_sys)[2,2] === value(substitute(2*α*y, backward_subs[2])) +@test noise_eqs(new_sys)[3,1] === value(σ) +@test noise_eqs(new_sys)[3,2] === value(α) \ No newline at end of file From 8140b6254c38219e2bd0ab1195c805818103c04d Mon Sep 17 00:00:00 2001 From: fchen121 Date: Fri, 6 Jun 2025 18:49:11 -0400 Subject: [PATCH 09/19] Combined change of variable function for ODE and SDE --- src/systems/diffeqs/basic_transformations.jl | 71 +++------ test/changeofvariables.jl | 148 +++++++++---------- 2 files changed, 92 insertions(+), 127 deletions(-) diff --git a/src/systems/diffeqs/basic_transformations.jl b/src/systems/diffeqs/basic_transformations.jl index 26131554ed..e28149819e 100644 --- a/src/systems/diffeqs/basic_transformations.jl +++ b/src/systems/diffeqs/basic_transformations.jl @@ -94,55 +94,13 @@ new_sol = solve(new_prob, Tsit5()) ``` """ -function changeofvariables(sys::System, iv, forward_subs, backward_subs; simplify=true, t0=missing) - t = iv - - old_vars = first.(backward_subs) - new_vars = last.(forward_subs) - # kept_vars = setdiff(states(sys), old_vars) - # rhs = [eq.rhs for eq in equations(sys)] - - # use: dz/dt = ∂z/∂x dx/dt + ∂z/∂t - dzdt = Symbolics.derivative( first.(forward_subs), t ) - new_eqs = Equation[] - for (new_var, ex) in zip(new_vars, dzdt) - for ode_eq in equations(sys) - ex = substitute(ex, ode_eq.lhs => ode_eq.rhs) - end - ex = substitute(ex, Dict(forward_subs)) - ex = substitute(ex, Dict(backward_subs)) - if simplify - ex = Symbolics.simplify(ex, expand=true) - end - push!(new_eqs, Differential(t)(new_var) ~ ex) - end - - defs = get_defaults(sys) - new_defs = Dict() - for f_sub in forward_subs - ex = substitute(first(f_sub), defs) - if !ismissing(t0) - ex = substitute(ex, t => t0) - end - new_defs[last(f_sub)] = ex - end - for para in parameters(sys) - if haskey(defs, para) - new_defs[para] = defs[para] - end - end - @named new_sys = System(vcat(new_eqs, first.(backward_subs) .~ last.(backward_subs)), t; - defaults=new_defs, - observed=observed(sys) - ) - if simplify - return mtkcompile(new_sys) +function changeofvariables( + sys::System, iv, forward_subs, backward_subs; + simplify=true, t0=missing, isSDE=false +) + if !iscomplete(sys) + sys = mtkcompile(sys) end - return new_sys -end - -function change_of_variable_SDE(sys::System, iv, forward_subs, backward_subs; simplify=true, t0=missing) - sys = mtkcompile(sys) t = iv old_vars = first.(backward_subs) @@ -152,10 +110,15 @@ function change_of_variable_SDE(sys::System, iv, forward_subs, backward_subs; si # use: dY = (∂f/∂t + μ∂f/∂x + (1/2)*σ^2*∂2f/∂x2)dt + σ∂f/∂xdW old_eqs = equations(sys) neqs = get_noise_eqs(sys) - neqs = [neqs[i,:] for i in 1:size(neqs,1)] + if neqs !== nothing + isSDE = true + neqs = [neqs[i,:] for i in 1:size(neqs,1)] - brownvars = map([Symbol(:B, :_, i) for i in 1:length(neqs[1])]) do name - unwrap(only(@brownian $name)) + brownvars = map([Symbol(:B, :_, i) for i in 1:length(neqs[1])]) do name + unwrap(only(@brownian $name)) + end + else + neqs = ones(1, length(old_eqs)) end # df/dt = ∂f/∂x dx/dt + ∂f/∂t @@ -168,8 +131,10 @@ function change_of_variable_SDE(sys::System, iv, forward_subs, backward_subs; si for (eqs, neq) in zip(old_eqs, neqs) if occursin(value(eqs.lhs), value(ex)) ex = substitute(ex, eqs.lhs => eqs.rhs) - for (noise, B) in zip(neq, brownvars) - ex = ex + 1/2 * noise^2 * second + noise*first*B + if isSDE + for (noise, B) in zip(neq, brownvars) + ex = ex + 1/2 * noise^2 * second + noise*first*B + end end end end diff --git a/test/changeofvariables.jl b/test/changeofvariables.jl index 68a8a143b6..1c0204e8fd 100644 --- a/test/changeofvariables.jl +++ b/test/changeofvariables.jl @@ -4,101 +4,101 @@ using Test, LinearAlgebra # Change of variables: z = log(x) # (this implies that x = exp(z) is automatically non-negative) -# @independent_variables t -# @variables z(t)[1:2, 1:2] -# D = Differential(t) -# eqs = [D(D(z)) ~ ones(2, 2)] -# @mtkcompile sys = System(eqs, t) -# @test_nowarn ODEProblem(sys, [z => zeros(2, 2), D(z) => ones(2, 2)], (0.0, 10.0)) +@independent_variables t +@variables z(t)[1:2, 1:2] +D = Differential(t) +eqs = [D(D(z)) ~ ones(2, 2)] +@mtkcompile sys = System(eqs, t) +@test_nowarn ODEProblem(sys, [z => zeros(2, 2), D(z) => ones(2, 2)], (0.0, 10.0)) -# @parameters α -# @variables x(t) -# D = Differential(t) -# eqs = [D(x) ~ α*x] +@parameters α +@variables x(t) +D = Differential(t) +eqs = [D(x) ~ α*x] -# tspan = (0., 1.) -# def = [x => 1.0, α => -0.5] +tspan = (0., 1.) +def = [x => 1.0, α => -0.5] -# @mtkcompile sys = System(eqs, t;defaults=def) -# prob = ODEProblem(sys, [], tspan) -# sol = solve(prob, Tsit5()) +@mtkcompile sys = System(eqs, t;defaults=def) +prob = ODEProblem(sys, [], tspan) +sol = solve(prob, Tsit5()) -# @variables z(t) -# forward_subs = [log(x) => z] -# backward_subs = [x => exp(z)] -# new_sys = changeofvariables(sys, t, forward_subs, backward_subs) -# @test equations(new_sys)[1] == (D(z) ~ α) +@variables z(t) +forward_subs = [log(x) => z] +backward_subs = [x => exp(z)] +new_sys = changeofvariables(sys, t, forward_subs, backward_subs) +@test equations(new_sys)[1] == (D(z) ~ α) -# new_prob = ODEProblem(new_sys, [], tspan) -# new_sol = solve(new_prob, Tsit5()) +new_prob = ODEProblem(new_sys, [], tspan) +new_sol = solve(new_prob, Tsit5()) -# @test isapprox(new_sol[x][end], sol[x][end], atol=1e-4) +@test isapprox(new_sol[x][end], sol[x][end], atol=1e-4) -# # Riccati equation -# @parameters α -# @variables x(t) -# D = Differential(t) -# eqs = [D(x) ~ t^2 + α - x^2] -# def = [x=>1., α => 1.] -# @mtkcompile sys = System(eqs, t; defaults=def) +# Riccati equation +@parameters α +@variables x(t) +D = Differential(t) +eqs = [D(x) ~ t^2 + α - x^2] +def = [x=>1., α => 1.] +@mtkcompile sys = System(eqs, t; defaults=def) -# @variables z(t) -# forward_subs = [t + α/(x+t) => z ] -# backward_subs = [ x => α/(z-t) - t] +@variables z(t) +forward_subs = [t + α/(x+t) => z ] +backward_subs = [ x => α/(z-t) - t] -# new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=true, t0=0.) -# # output should be equivalent to -# # t^2 + α - z^2 + 2 (but this simplification is not found automatically) +new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=true, t0=0.) +# output should be equivalent to +# t^2 + α - z^2 + 2 (but this simplification is not found automatically) -# tspan = (0., 1.) -# prob = ODEProblem(sys,[],tspan) -# new_prob = ODEProblem(new_sys,[],tspan) +tspan = (0., 1.) +prob = ODEProblem(sys,[],tspan) +new_prob = ODEProblem(new_sys,[],tspan) -# sol = solve(prob, Tsit5()) -# new_sol = solve(new_prob, Tsit5()) +sol = solve(prob, Tsit5()) +new_sol = solve(new_prob, Tsit5()) -# @test isapprox(sol[x][end], new_sol[x][end], rtol=1e-4) +@test isapprox(sol[x][end], new_sol[x][end], rtol=1e-4) -# # Linear transformation to diagonal system -# @independent_variables t -# @variables x(t)[1:3] -# x = reshape(x, 3, 1) -# D = Differential(t) -# A = [0. -1. 0.; -0.5 0.5 0.; 0. 0. -1.] -# right = A*x -# eqs = vec(D.(x) .~ right) +# Linear transformation to diagonal system +@independent_variables t +@variables x(t)[1:3] +x = reshape(x, 3, 1) +D = Differential(t) +A = [0. -1. 0.; -0.5 0.5 0.; 0. 0. -1.] +right = A*x +eqs = vec(D.(x) .~ right) -# tspan = (0., 10.) -# u0 = [x[1] => 1.0, x[2] => 2.0, x[3] => -1.0] +tspan = (0., 10.) +u0 = [x[1] => 1.0, x[2] => 2.0, x[3] => -1.0] -# @mtkcompile sys = System(eqs, t; defaults=u0) -# prob = ODEProblem(sys,[],tspan) -# sol = solve(prob, Tsit5()) +@mtkcompile sys = System(eqs, t; defaults=u0) +prob = ODEProblem(sys,[],tspan) +sol = solve(prob, Tsit5()) -# T = eigen(A).vectors -# T_inv = inv(T) +T = eigen(A).vectors +T_inv = inv(T) -# @variables z(t)[1:3] -# z = reshape(z, 3, 1) -# forward_subs = vec(T_inv*x .=> z) -# backward_subs = vec(x .=> T*z) +@variables z(t)[1:3] +z = reshape(z, 3, 1) +forward_subs = vec(T_inv*x .=> z) +backward_subs = vec(x .=> T*z) -# new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=true) +new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=true) -# new_prob = ODEProblem(new_sys, [], tspan) -# new_sol = solve(new_prob, Tsit5()) +new_prob = ODEProblem(new_sys, [], tspan) +new_sol = solve(new_prob, Tsit5()) -# # test RHS -# new_rhs = [eq.rhs for eq in equations(new_sys)] -# new_A = Symbolics.value.(Symbolics.jacobian(new_rhs, z)) -# A = diagm(eigen(A).values) -# A = sortslices(A, dims=1) -# new_A = sortslices(new_A, dims=1) -# @test isapprox(A, new_A, rtol = 1e-10) -# @test isapprox( new_sol[x[1],end], sol[x[1],end], rtol=1e-4) +# test RHS +new_rhs = [eq.rhs for eq in equations(new_sys)] +new_A = Symbolics.value.(Symbolics.jacobian(new_rhs, z)) +A = diagm(eigen(A).values) +A = sortslices(A, dims=1) +new_A = sortslices(new_A, dims=1) +@test isapprox(A, new_A, rtol = 1e-10) +@test isapprox( new_sol[x[1],end], sol[x[1],end], rtol=1e-4) # Change of variables for sde noise_eqs = ModelingToolkit.get_noise_eqs @@ -115,7 +115,7 @@ def = [x=>0., μ => 2., σ=>1.] @mtkcompile sys = System(eqs, t; defaults=def) forward_subs = [log(x) => y] backward_subs = [x => exp(y)] -new_sys = change_of_variable_SDE(sys, t, forward_subs, backward_subs) +new_sys = changeofvariables(sys, t, forward_subs, backward_subs) @test equations(new_sys)[1] == (D(y) ~ μ - 1/2*σ^2) @test noise_eqs(new_sys)[1] === value(σ) @@ -130,7 +130,7 @@ def = [x=>0., y=> 0., u=>0., μ => 2., σ=>1., α=>3.] @mtkcompile sys = System(eqs, t; defaults=def) forward_subs = [log(x) => z, y^2 => w, log(u) => v] backward_subs = [x => exp(z), y => w^.5, u => exp(v)] -new_sys = change_of_variable_SDE(sys, t, forward_subs, backward_subs) +new_sys = changeofvariables(sys, t, forward_subs, backward_subs) @test equations(new_sys)[1] == (D(z) ~ μ - 1/2*σ^2) @test equations(new_sys)[2] == (D(w) ~ α^2) @test equations(new_sys)[3] == (D(v) ~ μ - 1/2*(α^2 + σ^2)) From c0453f470f4f22af219d296ce5787a2c3bbe2c94 Mon Sep 17 00:00:00 2001 From: fchen121 Date: Mon, 9 Jun 2025 16:29:29 -0400 Subject: [PATCH 10/19] Change of variable for non-simplified SDE --- src/systems/diffeqs/basic_transformations.jl | 24 ++++++++++++++++---- test/changeofvariables.jl | 12 ++++++++-- 2 files changed, 29 insertions(+), 7 deletions(-) diff --git a/src/systems/diffeqs/basic_transformations.jl b/src/systems/diffeqs/basic_transformations.jl index e28149819e..c779d1b64b 100644 --- a/src/systems/diffeqs/basic_transformations.jl +++ b/src/systems/diffeqs/basic_transformations.jl @@ -98,9 +98,6 @@ function changeofvariables( sys::System, iv, forward_subs, backward_subs; simplify=true, t0=missing, isSDE=false ) - if !iscomplete(sys) - sys = mtkcompile(sys) - end t = iv old_vars = first.(backward_subs) @@ -110,7 +107,12 @@ function changeofvariables( # use: dY = (∂f/∂t + μ∂f/∂x + (1/2)*σ^2*∂2f/∂x2)dt + σ∂f/∂xdW old_eqs = equations(sys) neqs = get_noise_eqs(sys) - if neqs !== nothing + brownvars = brownians(sys) + + + if neqs === nothing && length(brownvars) === 0 + neqs = ones(1, length(old_eqs)) + elseif neqs !== nothing isSDE = true neqs = [neqs[i,:] for i in 1:size(neqs,1)] @@ -118,7 +120,19 @@ function changeofvariables( unwrap(only(@brownian $name)) end else - neqs = ones(1, length(old_eqs)) + isSDE = true + neqs = Vector{Any}[] + for (i, eq) in enumerate(old_eqs) + neq = Any[] + right = eq.rhs + for Bv in brownvars + lin_exp = linear_expansion(right, Bv) + right = lin_exp[2] + push!(neq, lin_exp[1]) + end + push!(neqs, neq) + old_eqs[i] = eq.lhs ~ right + end end # df/dt = ∂f/∂x dx/dt + ∂f/∂t diff --git a/test/changeofvariables.jl b/test/changeofvariables.jl index 1c0204e8fd..7f3ee5eb0b 100644 --- a/test/changeofvariables.jl +++ b/test/changeofvariables.jl @@ -127,9 +127,10 @@ new_sys = changeofvariables(sys, t, forward_subs, backward_subs) D = Differential(t) eqs = [D(x) ~ μ*x + σ*x*Bx, D(y) ~ α*By, D(u) ~ μ*u + σ*u*Bx + α*u*By] def = [x=>0., y=> 0., u=>0., μ => 2., σ=>1., α=>3.] -@mtkcompile sys = System(eqs, t; defaults=def) forward_subs = [log(x) => z, y^2 => w, log(u) => v] backward_subs = [x => exp(z), y => w^.5, u => exp(v)] + +@mtkcompile sys = System(eqs, t; defaults=def) new_sys = changeofvariables(sys, t, forward_subs, backward_subs) @test equations(new_sys)[1] == (D(z) ~ μ - 1/2*σ^2) @test equations(new_sys)[2] == (D(w) ~ α^2) @@ -139,4 +140,11 @@ new_sys = changeofvariables(sys, t, forward_subs, backward_subs) @test noise_eqs(new_sys)[2,1] === value(0) @test noise_eqs(new_sys)[2,2] === value(substitute(2*α*y, backward_subs[2])) @test noise_eqs(new_sys)[3,1] === value(σ) -@test noise_eqs(new_sys)[3,2] === value(α) \ No newline at end of file +@test noise_eqs(new_sys)[3,2] === value(α) + +# Test for Brownian instead of noise +@named sys = System(eqs, t; defaults=def) +new_sys = changeofvariables(sys, t, forward_subs, backward_subs; simplify=false) +@test simplify(equations(new_sys)[1]) == simplify((D(z) ~ μ - 1/2*σ^2 + σ*Bx)) +@test simplify(equations(new_sys)[2]) == simplify((D(w) ~ α^2 + 2*α*w^.5*By)) +@test simplify(equations(new_sys)[3]) == simplify((D(v) ~ μ - 1/2*(α^2 + σ^2) + σ*Bx + α*By)) \ No newline at end of file From 166b0f82c479cb856846c7653bf7ffcfdf843414 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Tue, 10 Jun 2025 09:43:05 +0000 Subject: [PATCH 11/19] Update test/changeofvariables.jl --- test/changeofvariables.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/test/changeofvariables.jl b/test/changeofvariables.jl index 7f3ee5eb0b..776d44f7e7 100644 --- a/test/changeofvariables.jl +++ b/test/changeofvariables.jl @@ -105,7 +105,7 @@ noise_eqs = ModelingToolkit.get_noise_eqs value = ModelingToolkit.value @independent_variables t -@brownian B +@brownians B @parameters μ σ @variables x(t) y(t) D = Differential(t) @@ -121,7 +121,7 @@ new_sys = changeofvariables(sys, t, forward_subs, backward_subs) #Multiple Brownian and equations @independent_variables t -@brownian Bx By +@brownians Bx By @parameters μ σ α @variables x(t) y(t) z(t) w(t) u(t) v(t) D = Differential(t) From 0884e59cc855691601be01da8d5d82a964843773 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 11 Jun 2025 13:12:31 +0000 Subject: [PATCH 12/19] Update src/systems/diffeqs/basic_transformations.jl --- src/systems/diffeqs/basic_transformations.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/systems/diffeqs/basic_transformations.jl b/src/systems/diffeqs/basic_transformations.jl index c779d1b64b..e2b49808be 100644 --- a/src/systems/diffeqs/basic_transformations.jl +++ b/src/systems/diffeqs/basic_transformations.jl @@ -117,7 +117,7 @@ function changeofvariables( neqs = [neqs[i,:] for i in 1:size(neqs,1)] brownvars = map([Symbol(:B, :_, i) for i in 1:length(neqs[1])]) do name - unwrap(only(@brownian $name)) + unwrap(only(@brownians $name)) end else isSDE = true From 70a5b3975cde4279f16dadb713531af58180032c Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 11 Jun 2025 22:40:53 +0000 Subject: [PATCH 13/19] Update Project.toml --- Project.toml | 2 -- 1 file changed, 2 deletions(-) diff --git a/Project.toml b/Project.toml index e9ea4cd50e..0f9789472b 100644 --- a/Project.toml +++ b/Project.toml @@ -164,8 +164,6 @@ StochasticDiffEq = "6.72.1" SymbolicIndexingInterface = "0.3.39" SymbolicUtils = "3.26.1" Symbolics = "6.40" -TermInterface = "2.0.0" -Test = "1.11.0" URIs = "1" UnPack = "0.1, 1.0" Unitful = "1.1" From 2d5c9e7128f904fe6d710177580be07ee1e7114b Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 11 Jun 2025 22:41:04 +0000 Subject: [PATCH 14/19] Update Project.toml --- Project.toml | 2 -- 1 file changed, 2 deletions(-) diff --git a/Project.toml b/Project.toml index 0f9789472b..67d54f1757 100644 --- a/Project.toml +++ b/Project.toml @@ -65,8 +65,6 @@ StochasticDiffEq = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0" SymbolicIndexingInterface = "2efcf032-c050-4f8e-a9bb-153293bab1f5" SymbolicUtils = "d1185830-fcd6-423d-90d6-eec64667417b" Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7" -TermInterface = "8ea1fca8-c5ef-4a55-8b96-4e9afe9c9a3c" -Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" URIs = "5c2747f8-b7ea-4ff2-ba2e-563bfd36b1d4" UnPack = "3a884ed6-31ef-47d7-9d2a-63182c4928ed" Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" From 2b32ad48e15e15abd408383f4d3b9c0afec6dbaf Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 11 Jun 2025 22:41:11 +0000 Subject: [PATCH 15/19] Update Project.toml --- Project.toml | 1 - 1 file changed, 1 deletion(-) diff --git a/Project.toml b/Project.toml index 67d54f1757..37bf660098 100644 --- a/Project.toml +++ b/Project.toml @@ -46,7 +46,6 @@ OffsetArrays = "6fe1bfb0-de20-5000-8ca7-80f57d26f881" OrderedCollections = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed" OrdinaryDiffEqCore = "bbf590c4-e513-4bbe-9b18-05decba2e5d8" -OrdinaryDiffEqDefault = "50262376-6c5a-4cf5-baba-aaf4f84d72d7" PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a" RecursiveArrayTools = "731186ca-8d62-57ce-b412-fbd966d074cd" Reexport = "189a3867-3050-52da-a836-e630ba90ab69" From 2b4864858aa9ba28059e9dabbb8bcbbd5a511537 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 11 Jun 2025 22:41:24 +0000 Subject: [PATCH 16/19] Update Project.toml --- Project.toml | 1 - 1 file changed, 1 deletion(-) diff --git a/Project.toml b/Project.toml index 37bf660098..201436fa19 100644 --- a/Project.toml +++ b/Project.toml @@ -60,7 +60,6 @@ SimpleNonlinearSolve = "727e6d20-b764-4bd8-a329-72de5adea6c7" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" -StochasticDiffEq = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0" SymbolicIndexingInterface = "2efcf032-c050-4f8e-a9bb-153293bab1f5" SymbolicUtils = "d1185830-fcd6-423d-90d6-eec64667417b" Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7" From 07affe6eec93afa8ff255280c8f3f782be4e08af Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 11 Jun 2025 22:41:32 +0000 Subject: [PATCH 17/19] Update Project.toml --- Project.toml | 1 - 1 file changed, 1 deletion(-) diff --git a/Project.toml b/Project.toml index 201436fa19..9645a412b0 100644 --- a/Project.toml +++ b/Project.toml @@ -44,7 +44,6 @@ NaNMath = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3" NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec" OffsetArrays = "6fe1bfb0-de20-5000-8ca7-80f57d26f881" OrderedCollections = "bac558e1-5e72-5ebc-8fee-abe8a469f55d" -OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed" OrdinaryDiffEqCore = "bbf590c4-e513-4bbe-9b18-05decba2e5d8" PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a" RecursiveArrayTools = "731186ca-8d62-57ce-b412-fbd966d074cd" From 7d5e9155d40f57dcd5bf4b72cfcca5efc21890ad Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Wed, 11 Jun 2025 22:42:17 +0000 Subject: [PATCH 18/19] Update src/ModelingToolkit.jl --- src/ModelingToolkit.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/ModelingToolkit.jl b/src/ModelingToolkit.jl index 294010e831..8aaecaa6a3 100644 --- a/src/ModelingToolkit.jl +++ b/src/ModelingToolkit.jl @@ -297,7 +297,7 @@ export isinput, isoutput, getbounds, hasbounds, getguess, hasguess, isdisturbanc hasunit, getunit, hasconnect, getconnect, hasmisc, getmisc, state_priority export liouville_transform, change_independent_variable, substitute_component, - add_accumulations, noise_to_brownians, Girsanov_transform, changeofvariables, change_of_variable_SDE + add_accumulations, noise_to_brownians, Girsanov_transform, changeofvariables export PDESystem export Differential, expand_derivatives, @derivatives export Equation, ConstrainedEquation From 96609b33c00da03cfc9f508aba03ed783fa5ef5e Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Thu, 12 Jun 2025 10:13:47 +0000 Subject: [PATCH 19/19] Update runtests.jl --- test/runtests.jl | 1 + 1 file changed, 1 insertion(+) diff --git a/test/runtests.jl b/test/runtests.jl index 31e4cc609f..75b23a8cda 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -47,6 +47,7 @@ end @safetestset "Error Handling" include("error_handling.jl") @safetestset "StructuralTransformations" include("structural_transformation/runtests.jl") @safetestset "Basic transformations" include("basic_transformations.jl") + @safetestset "Change of variables" include("changeofvariables.jl") @safetestset "State Selection Test" include("state_selection.jl") @safetestset "Symbolic Event Test" include("symbolic_events.jl") @safetestset "Stream Connect Test" include("stream_connectors.jl")