typos

Author: oscardssmith

lib/OrdinaryDiffEqRosenbrock/src/rosenbrock_perform_step.jl

@@ -1230,7 +1230,7 @@ end
 
     # Initialize ks
     num_stages = size(A,1)
-    du = f(u, p, t)
+    du = f(uprev, p, t)
     linsolve_tmp = @.. du + dtd[1] * dT
     k1 = _reshape(W \ -_vec(linsolve_tmp), axes(uprev))
     # constant number for type stability make sure this is greater than num_stages
@@ -1246,27 +1246,27 @@ end
         OrdinaryDiffEqCore.increment_nf!(integrator.stats, 1)
 
         # Compute linsolve_tmp for current stage
-        linsolve_tmp1 = zero(du)
+        linsolve_tmp = zero(du)
         if mass_matrix === I
             for i in 1:stage-1
-                linsolve_tmp1 = @.. linsolve_tmp1 + dtC[stage, i] * ks[i]
+                linsolve_tmp = @.. linsolve_tmp + dtC[stage, i] * ks[i]
             end
         else
             for i in 1:stage-1
-                linsolve_tmp1 = @.. linsolve_tmp1 + dtC[stage, i] * ks[i]
+                linsolve_tmp = @.. linsolve_tmp + dtC[stage, i] * ks[i]
             end
-            linsolve_tmp1 = mass_matrix * linsolve_tmp1
+            linsolve_tmp = mass_matrix * linsolve_tmp
         end
-        linsolve_tmp = @.. du + dtd[stage] * dT + linsolve_tmp1
+        linsolve_tmp = @.. du + dtd[stage] * dT + linsolve_tmp
 
         ks = Base.setindex(ks, _reshape(W \ -_vec(linsolve_tmp), axes(uprev)), stage)
         integrator.stats.nsolve += 1
     end
     #@show ks
-    u = u .+ ks[end]
+    u = u .+ ks[num_stages]
 
     if integrator.opts.adaptive
-        atmp = calculate_residuals(ks[end], uprev, u, integrator.opts.abstol,
+        atmp = calculate_residuals(ks[num_stages], uprev, u, integrator.opts.abstol,
             integrator.opts.reltol, integrator.opts.internalnorm, t)
         integrator.EEst = integrator.opts.internalnorm(atmp, t)
     end

lib/OrdinaryDiffEqRosenbrock/src/stiff_addsteps.jl

@@ -322,7 +322,7 @@ function _ode_addsteps!(k, t, uprev, u, dt, f, p, cache::Rodas4ConstantCache,
         num_stages = size(A,1)
         du = f(u, p, t)
         linsolve_tmp = @.. du + dtd[1] * dT
-        k1 = _reshape(W \ -_vec(linsolve_tmp), axes(uprev))
+        k1 = _reshape(W \ _vec(linsolve_tmp), axes(uprev))
         # constant number for type stability make sure this is greater than num_stages
         ks = ntuple(Returns(k1), 10)
         # Last stage doesn't affect ks
@@ -335,18 +335,18 @@ function _ode_addsteps!(k, t, uprev, u, dt, f, p, cache::Rodas4ConstantCache,
             du = f(u, p, t + c[stage] * dt)
 
             # Compute linsolve_tmp for current stage
-            linsolve_tmp = @.. du + dtd[stage] * dT
+            linsolve_tmp = zero(du)
             if mass_matrix === I
                 for i in 1:stage-1
                     linsolve_tmp = @.. linsolve_tmp + dtC[stage, i] * ks[i]
                 end
             else
                 for i in 1:stage-1
-                    linsolve_tmp1 = mass_matrix * @..dtC[stage, i] * ks[i]
-                    linsolve_tmp = @.. linsolve_tmp + linsolve_tmp1
+                    linsolve_tmp = @.. linsolve_tmp + dtC[stage, i] * ks[i]
                 end
+                linsolve_tmp = mass_matrix * linsolve_tmp
             end
-
+            linsolve_tmp = @.. du + dtd[stage] * dT + linsolve_tmp
             ks = Base.setindex(ks, _reshape(W \ _vec(linsolve_tmp), axes(uprev)), stage)
         end