On a separate note: using a finite difference to calculate the reaction rate is not ideal — you are throwing away a lot of the accuracy of the solution.
You can get the “exact” time derivative (or at least, as accurate as your solution) just by plugging the solution into your right-hand-side function crnn!. But even more simply you can do e.g.
dmass_r1_dt = getindex.(sol1.(time, Val{1}), 1) # Val{1} = first derivative
or replace the getindex call with
dmass_r1_dt = sol1.(time, Val{1}, idxs=1)
or
dmass_r1_dt = sol1(time, Val{1}, idxs=1).u
See also DifferentialEquations and derivatives of solutions? and “Interpolations and Calculating Derivatives” in the manual.