Why does DifferentialEquations.jl gives me an error when running a code from the official documentation?

Question

I'm trying to run an example from the DifferentialEquations.jl documentation and it is not working.

The example is the Simple Harmonic Oscillator (https://docs.sciml.ai/DiffEqDocs/stable/examples/classical_physics/#Simple-Harmonic-Oscillator), please see the code below:

# Simple Harmonic Oscillator Problem
using OrdinaryDiffEq, Plots

#Parameters
ω = 1

#Initial Conditions
x₀ = [0.0]
dx₀ = [π / 2]
tspan = (0.0, 2π)

ϕ = atan((dx₀[1] / ω) / x₀[1])
A = √(x₀[1]^2 + dx₀[1]^2)

#Define the problem
function harmonicoscillator(ddu, du, u, ω, t)
    ddu .= -ω^2 * u
end

#Pass to solvers
prob = SecondOrderODEProblem(harmonicoscillator, dx₀, x₀, tspan, ω)
sol = solve(prob, DPRKN6())

#Plot
plot(sol, vars = [2, 1], linewidth = 2, title = "Simple Harmonic Oscillator",
    xaxis = "Time", yaxis = "Elongation", label = ["x" "dx"])
plot!(t -> A * cos(ω * t - ϕ), lw = 3, ls = :dash, label = "Analytical Solution x")
plot!(t -> -A * ω * sin(ω * t - ϕ), lw = 3, ls = :dash, label = "Analytical Solution dx")

The error seems to be related to the plotting of the results. The error message is:

ERROR: BoundsError: attempt to access 1-element Vector{Float64} at index [2]
Stacktrace:
  [1] getindex
    @ .\essentials.jl:13 [inlined]
  [2] _ode_interpolant(Θ::Float64, dt::Float64, y₀::RecursiveArrayTools.ArrayPartition{…}, y₁::RecursiveArrayTools.ArrayPartition{…}, k::Vector{…}, cache::OrdinaryDiffEq.DPRKN6Cache{…}, idxs::Int64, T::Type{…}, differential_vars::Nothing)
    @ OrdinaryDiffEq c:\Users\jdoe\.julia\packages\OrdinaryDiffEq\uQcyc\src\dense\interpolants.jl:2936
  [3] ode_interpolant(Θ::Float64, dt::Float64, y₀::RecursiveArrayTools.ArrayPartition{…}, y₁::RecursiveArrayTools.ArrayPartition{…}, k::Vector{…}, cache::OrdinaryDiffEq.DPRKN6Cache{…}, idxs::Int64, T::Type{…}, differential_vars::Nothing)
    @ OrdinaryDiffEq c:\Users\jdoe\.julia\packages\OrdinaryDiffEq\uQcyc\src\dense\generic_dense.jl:594
  [4] _evaluate_interpolant(f::Function, Θ::Float64, dt::Float64, timeseries::Vector{…}, i₋::Int64, i₊::Int64, cache::OrdinaryDiffEq.DPRKN6Cache{…}, idxs::Int64, deriv::Type, ks::Vector{…}, ts::Vector{…}, p::Int64, differential_vars::Nothing)
    @ OrdinaryDiffEq c:\Users\jdoe\.julia\packages\OrdinaryDiffEq\uQcyc\src\dense\generic_dense.jl:300
  [5] evaluate_interpolant(f::Function, Θ::Float64, dt::Float64, timeseries::Vector{…}, i₋::Int64, i₊::Int64, cache::OrdinaryDiffEq.DPRKN6Cache{…}, idxs::Int64, deriv::Type, ks::Vector{…}, ts::Vector{…}, id::OrdinaryDiffEq.InterpolationData{…}, p::Int64, differential_vars::Nothing)
    @ OrdinaryDiffEq c:\Users\jdoe\.julia\packages\OrdinaryDiffEq\uQcyc\src\dense\generic_dense.jl:334
  [6] (::OrdinaryDiffEq.var"#742#743"{…})(j::Int64)
    @ OrdinaryDiffEq c:\Users\jdoe\.julia\packages\OrdinaryDiffEq\uQcyc\src\dense\generic_dense.jl:371
  [7] iterate
    @ .\generator.jl:47 [inlined]
  [8] _collect(c::Vector{…}, itr::Base.Generator{…}, ::Base.EltypeUnknown, isz::Base.HasShape{…})
    @ Base .\array.jl:854
  [9] collect_similar
    @ .\array.jl:763 [inlined]
 [10] map
    @ .\abstractarray.jl:3285 [inlined]
 [11] ode_interpolation(tvals::Vector{…}, id::OrdinaryDiffEq.InterpolationData{…}, idxs::Int64, deriv::Type{…}, p::Int64, continuity::Symbol)
    @ OrdinaryDiffEq c:\Users\jdoe\.julia\packages\OrdinaryDiffEq\uQcyc\src\dense\generic_dense.jl:353
 [12] InterpolationData
    @ c:\Users\jdoe\.julia\packages\OrdinaryDiffEq\uQcyc\src\interp_func.jl:169 [inlined]
 [13] (::ODESolution{…})(t::Vector{…}, ::Type{…}, idxs::Int64, continuity::Symbol)
    @ SciMLBase c:\Users\jdoe\.julia\packages\SciMLBase\m3AcC\src\solutions\ode_solutions.jl:175
 [14] #_#440
    @ c:\Users\jdoe\.julia\packages\SciMLBase\m3AcC\src\solutions\ode_solutions.jl:147 [inlined]
 [15] AbstractODESolution
    @ c:\Users\jdoe\.julia\packages\SciMLBase\m3AcC\src\solutions\ode_solutions.jl:145 [inlined]
 [16] solplot_vecs_and_labels(dims::Int64, vars::Vector{…}, plott::Vector{…}, sol::ODESolution{…}, plot_analytic::Bool, plot_analytic_timeseries::Nothing)
    @ SciMLBase c:\Users\jdoe\.julia\packages\SciMLBase\m3AcC\src\solutions\solution_interface.jl:443
 [17] diffeq_to_arrays(sol::ODESolution{…}, plot_analytic::Bool, denseplot::Bool, plotdensity::Int64, tspan::Nothing, vars::Vector{…}, tscale::Symbol, plotat::Nothing)
    @ SciMLBase c:\Users\jdoe\.julia\packages\SciMLBase\m3AcC\src\solutions\solution_interface.jl:339
 [18] macro expansion
    @ c:\Users\jdoe\.julia\packages\SciMLBase\m3AcC\src\solutions\solution_interface.jl:185 [inlined]
 [19] apply_recipe(plotattributes::AbstractDict{Symbol, Any}, sol::SciMLBase.AbstractTimeseriesSolution)
    @ SciMLBase c:\Users\jdoe\.julia\packages\RecipesBase\BRe07\src\RecipesBase.jl:300
 [20] _process_userrecipes!(plt::Any, plotattributes::Any, args::Any)
    @ RecipesPipeline c:\Users\jdoe\.julia\packages\RecipesPipeline\BGM3l\src\user_recipe.jl:38
 [21] recipe_pipeline!(plt::Any, plotattributes::Any, args::Any)
    @ RecipesPipeline c:\Users\jdoe\.julia\packages\RecipesPipeline\BGM3l\src\RecipesPipeline.jl:72
 [22] _plot!(plt::Plots.Plot, plotattributes::Any, args::Any)
    @ Plots c:\Users\jdoe\.julia\packages\Plots\HyyIK\src\plot.jl:223
 [23] plot(args::Any; kw...)
    @ Plots c:\Users\jdoe\.julia\packages\Plots\HyyIK\src\plot.jl:102
 [24] top-level scope
    @ c:\data\JuliaLanguage\Teste_DE.jl:25

I tried to adjust the plot parameters, without success. Also I can run the other examples of the DifferentialEquations.jl without problems.