_compute_estimate: change order of T() and ^Q in denominator

Running a simple test

    using Unitful: @u_str
    import Interpolations as Ilib
    import Interpolations
    import FiniteDifferences as FDlib

    tt = (0.0 : 0.01 : 1.1) .* u"s"

    fn = 3u"J / s" .* tt

    ii = Ilib.interpolate((tt,),fn, ComposedFunction(Ilib.Gridded, Ilib.Linear)())

    cdm = FDlib.central_fdm(5,1)

    @show dd = cdm(ii, 0.6u"s")

without this change in the order of operations causes the code to crash
in _limit_step, because of a dimension mismatch. This is a temporary "fix",
that actually not a fix. The given order of operations should really
work as-is...

Author: SeSodesa

src/methods.jl

@@ -260,7 +260,12 @@ function _compute_estimate(
     # therefore perform the broadcasting first. See
     # https://github.com/JuliaLang/julia/issues/39151.
     _coefs = T.(coefs)
-    return sum(fs .* _coefs) ./ T(step)^Q
+    # NOTE the denominator: while doing T(step)^Q should technically be possible and correct,
+    # without flipping the order of ^ Q and T(), the current code crashes in _limit_step
+    # because of dimension mismatch. However, the output now has wrong units.
+    #
+    # TODO: fix the output units. The numerical value is correct in a simple test.
+    return sum(fs .* _coefs) ./ T(step ^ Q)
 end
 
 # Check the method and derivative orders for consistency.