diff --git a/Project.toml b/Project.toml
index da06bd7..ee7a7f1 100644
--- a/Project.toml
+++ b/Project.toml
@@ -24,6 +24,7 @@ Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7"
 RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"
 RotationalComponents = "824fa3ed-d7dc-4365-8bdf-ac77a1ac2011"
 RuntimeGeneratedFunctions = "7e49a35a-f44a-4d26-94aa-eba1b4ca6b47"
+SymbolicIndexingInterface = "2efcf032-c050-4f8e-a9bb-153293bab1f5"
 SynchCompiler = "5d9dccf6-a926-4748-b7e2-6521ccc431d1"
 SynchJulia = "a1b2c3d4-5e6f-7a8b-9c0d-e1f2a3b4c5d6"
 SynchToolkit = "f500ffa8-9682-42ac-8d34-0ca7926c2e94"
@@ -33,6 +34,7 @@ TranslationalComponents = "42acc370-819b-4877-bdd6-eb3294461f5d"
 [compat]
 ControlSystemsBase = "1"
 ControlSystemsMTK = "2.6"
+DiscreteComponents = "0.3"
 DyadInterface = "7.1.0"
 ElectricalComponents = "2.2.1"
 HydraulicComponents = "2.0.3"
@@ -43,9 +45,10 @@ OrdinaryDiffEqDefault = "2"
 Printf = "1"
 RecipesBase = "1.3"
 RotationalComponents = "2.0.0"
-SynchCompiler = "0.4.0"
-SynchJulia = "0.4.0"
-SynchToolkit = "0.4.8"
+SymbolicIndexingInterface = "0.3.51"
+SynchCompiler = "0.4.3"
+SynchJulia = "0.4.4"
+SynchToolkit = "0.5"
 ThermalComponents = "2.0.5"
 TranslationalComponents = "2.4.2"
 
diff --git a/src/program.jl b/src/program.jl
index b8b9b16..1e947a0 100644
--- a/src/program.jl
+++ b/src/program.jl
@@ -21,6 +21,9 @@ using ModelingToolkit
 using ModelingToolkit: t_nounits as t
 using DiscreteComponents: PeriodicClock
 using OrderedCollections: OrderedDict
+# `default_values` (a system's bindings, initial conditions and observed equations merged
+# into one map) is owned by SymbolicIndexingInterface; ModelingToolkit no longer re-exports it.
+using SymbolicIndexingInterface: default_values
 
 export ProgramLog, ProgramRuntime, run_program!, read_log
 
@@ -51,8 +54,7 @@ analysis forward its own value into it.
 function resolve_tunables(sys, syms::AbstractDict)
     SymT = ModelingToolkit.SymbolicT
     keysyms = [ModelingToolkit.unwrap(k) for k in keys(syms)]
-    dv = Dict{SymT, SymT}(
-        ModelingToolkit.default_values(ModelingToolkit.expand_connections(sys)))
+    dv = Dict{SymT, SymT}(default_values(ModelingToolkit.expand_connections(sys)))
     ModelingToolkit.evaluate_varmap!(dv, keysyms)
     out = OrderedDict{Symbol, Any}()
     for (k, field) in syms
diff --git a/test/Project.toml b/test/Project.toml
index 94fbf19..31f78e7 100644
--- a/test/Project.toml
+++ b/test/Project.toml
@@ -23,5 +23,6 @@ RobustAndOptimalControl = "21fd56a4-db03-40ee-82ee-a87907bee541"
 SeeToDee = "1c904df7-48cd-41e7-921b-d889ed9a470c"
 StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
 Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
+SymbolicIndexingInterface = "2efcf032-c050-4f8e-a9bb-153293bab1f5"
 SynchToolkit = "f500ffa8-9682-42ac-8d34-0ca7926c2e94"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
diff --git a/test/runtests.jl b/test/runtests.jl
index 6ac3ed2..f39efa0 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -4,6 +4,7 @@ using Test
 # backend and C export used below require Clang_unified_jll to be loaded.
 using Clang_unified_jll
 using ModelingToolkit
+using SymbolicIndexingInterface: default_values
 using MultibodyComponents
 # using DiscreteComponents
 using SynchToolkit
@@ -510,7 +511,7 @@ import DyadCompilerPasses
         # the whole friction law as one expression in `w`; compile that and evaluate it.
         expr = S.substitute(S.substitute(rhs("tau_f(t)"),
                                          Dict(sym("sw(t)") => rhs("sw(t)"))),
-                            Dict(ModelingToolkit.default_values(f)))
+                            Dict(default_values(f)))
         tau = S.build_function(expr, sym("w(t)"); expression = Val(false))
 
         # Odd in w: friction opposes the motion whichever way the axis turns.
@@ -535,7 +536,7 @@ import DyadCompilerPasses
             rhs2(n) = eqs2[findfirst(e -> string(e.lhs) == n, eqs2)].rhs
             e2 = S.substitute(S.substitute(rhs2("tau_f(t)"),
                                           Dict(sym2("sw(t)") => rhs2("sw(t)"))),
-                              Dict(ModelingToolkit.default_values(fi)))
+                              Dict(default_values(fi)))
             tq = S.build_function(e2, sym2("w(t)"); expression = Val(false))
             ws = range(1.0, 40.0, length = 200)          # the measured range
             @test q.kv > 0