Restrict ModelingToolkitTearing to v1.19 to fix crash

[hersle]

With v1.20 this is unstable:

using SymBoltz
lmax = 5
M = ΛCDM(K = nothing; lmax) # flat
pars = parameters_Planck18(M)
prob = CosmologyProblem(M, pars)
sol = solve(prob, 20.0)

Cosmology solution for model ΛCDM
Stages:
  Background: return code Terminated; solved with Rodas5P; 795 points
  Perturbations: return codes Unstable; solved with Rodas5P; 149-149 points; x1 k ∈ [20.0, 20.0] H₀/c (interpolation in-between)