[singhharsh1708]

Adds benchmarks/PETScDiffEq with one work-precision page per problem: Lotka-Volterra (non-stiff), HIRES (stiff), the Robertson DAE in residual and mass-matrix form, and a 1000-state 1-D Brusselator with a sparse Jacobian. PETSc's TSRK, TSRosW, TSImplicit (BDF), TSDAE and TSARKIMEX run through PETScDiffEq 0.14.0 next to OrdinaryDiffEq and Sundials. The coefficient-matched pairs (TSRK 3bs/BS3 and 5dp/DP5, TSRosW ra34pw2/ROS34PW2 and rodas3/Rodas3, TSARKIMEX 4/KenCarp4) are drawn in one color, and on the Brusselator each PETSc method also runs with direct LU (-ksp_type preonly -pc_type lu) next to PETSc's default GMRES with ILU(0). The pages check the pairs with fixed steps, where the final states agree to between 2e-14 and 3e-11.

Each page checks its reference against a second solver; the two agree to between 7e-14 and 4e-10, which is 68 to 191 times below the smallest error plotted. Each page also measures the fixed cost of a PETScDiffEq solve with a one-step solve (0.34 to 1.4 ms here, against 0.75 µs to 0.46 ms for the OrdinaryDiffEq solvers measured) and states it next to the plot. Each page prints time-to-error tables, the matched pairs per step attempt, and every point with its step, Jacobian and Newton counts, so the plots can be checked in text. The shared helpers are in wpd_common.jl, included the way NonStiffODE includes enright_pryce.jl.

On an Apple M5 at the target errors, with PETScDiffEq 0.14.0, the fastest PETScDiffEq method takes 11.8 times as long as the fastest method on Lotka-Volterra (TSRK 5dp is 8.3 times DP5), 7.5 times on HIRES (TSRosW ra34pw2 is 4.8 times ROS34PW2), at least 17.5 times on Robertson, and 3.7 times on the Brusselator, where TSRosW with LU is 1.1 times ROS34PW2.

GKSwstype=100 julia --threads=auto --project=. benchmark.jl benchmarks/PETScDiffEq on Julia 1.12.7, the Manifest's version, weaves all four pages in 239 s (36 to 94 s each, numruns = 20), and the woven markdown has no errors or warnings. The Manifest keeps PETSc_jll 3.22.2, because 3.25 pulls in HDF5, MUMPS, SuperLU_DIST and the AWS libraries.

AI Disclosure: Claude assisted with the code. Testing and review by @singhharsh1708.

[singhharsh1708]

These pages show PETScDiffEq in a poor state for small and medium problems: it is 7 to 52 times slower than the fastest method on every problem, and 1.5 to 21 times slower than OrdinaryDiffEq even where both run the same tableau and take the same steps. The counts on the pages point at three causes in PETScDiffEq rather than in PETSc's methods: a fixed cost of 0.34 to 1.4 ms per solve, overhead per step and per right-hand-side call, and far more Jacobians than OrdinaryDiffEq forms on stiff problems (244 against FBDF's 8 on HIRES, 574 against KenCarp4's 21 for TSARKIMEX).

I have started profiling and fixing each of these in PETScDiffEq, and will rerun these pages and update the numbers here as the fixes land.

[singhharsh1708]

The pages now run on PETScDiffEq 0.14.0. The fixed cost per solve fell from about 0.4 ms to 0.04-0.1 ms, steps no longer allocate, small dense Newton systems are solved in Julia, PETSc no longer rebuilds work vectors or looks up viewer options in every implicit step, and banded sparse Jacobians use LU by default. Against the fastest method: Lotka-Volterra 34.9x to 11.8x, HIRES 19.4x to 7.5x, Robertson 51.9x to 17.5x, Brusselator 7.7x to 3.7x. Same-method pairs: TSRK 5dp vs DP5 21.5x to 8.3x, TSRosW vs ROS34PW2 12.6x to 4.8x on HIRES and 27.4x to 8.7x on Robertson. Results are unchanged apart from rounding.