fix: handle observeds in get_all_timeseries_indexes for split = false systems

src/systems/abstractsystem.jl

@@ -770,11 +770,19 @@ for traitT in [
                    is_variable(sys, operation(s)(get_iv(sys)))
                 # DDEs case, to detect x(t - k)
                 push!(ts_idxs, ContinuousTimeseries())
-            elseif has_index_cache(sys) && (ic = get_index_cache(sys)) !== nothing
-                if (ts = get(ic.observed_syms_to_timeseries, s, nothing)) !== nothing
-                    union!(ts_idxs, ts)
-                elseif (ts = get(ic.dependent_pars_to_timeseries, s, nothing)) !== nothing
-                    union!(ts_idxs, ts)
+            else
+                if has_index_cache(sys) && (ic = get_index_cache(sys)) !== nothing
+                    if (ts = get(ic.observed_syms_to_timeseries, s, nothing)) !== nothing
+                        union!(ts_idxs, ts)
+                    elseif (ts = get(ic.dependent_pars_to_timeseries, s, nothing)) !==
+                           nothing
+                        union!(ts_idxs, ts)
+                    end
+                else
+                    # for split=false systems
+                    if has_observed_with_lhs(sys, sym)
+                        push!(ts_idxs, ContinuousTimeseries())
+                    end
                 end
             end
         end

test/odesystem.jl

@@ -1521,6 +1521,6 @@ end
     @parameters p[1:2] = [1.0, 2.0]
     @mtkbuild sys = ODESystem([D(x) ~ x, y^2 ~ x + sum(p)], t)
     prob = DAEProblem(sys, [D(x) => x, D(y) => D(x) / 2y], [], (0.0, 1.0))
-    sol = solve(prob, DFBDF(), abstol=1e-8, reltol=1e-8)
+    sol = solve(prob, DFBDF(), abstol = 1e-8, reltol = 1e-8)
     @test sol[x]≈sol[y^2 - sum(p)] atol=1e-5
 end

test/symbolic_indexing_interface.jl

@@ -213,3 +213,14 @@ end
     prob = ODEProblem(sys, [], (0.0, 1.0))
     @test prob.ps[b] == prob.ps[:b]
 end
+
+@testset "`get_all_timeseries_indexes` with non-split systems" begin
+    @variables x(t) y(t) z(t)
+    @parameters a
+    @named sys = ODESystem([D(x) ~ a * x, y ~ 2x, z ~ 0.0], t)
+    sys = structural_simplify(sys, split = false)
+    for sym in [x, y, z, x + y, x + a, y / x]
+        @test only(get_all_timeseries_indexes(sys, sym)) == ContinuousTimeseries()
+    end
+    @test isempty(get_all_timeseries_indexes(sys, a))
+end