Hi there, Before i state what my actual question is, let me give some background on some work i have done on this. I have been trying to solve the 0-D PRKE transient problem for a step and a ramp input numerically using couple of high precision methods. 1) Using MATLAB's inbuilt ode23s solver. 2) Using a custom written adaptive solver My problem is that i want to find the order of accuracy or maybe even just the exact error in my adaptive method in comparison to the actual true solution. Without the exact analytical solution, i would have to assume that the solution from MATLAB is exact and compute the error. But i am reluctant to do that because i want to test few embedded numerical methods with higher orders of accuracy than matlab's 2-3 method. So now, i did derive the analytical solution for PRKE with a step input (constant reactivity) and tested my method with the solution i computed analytically. It perfectly matches ! But when the reactivity is a function of time (Ramp or say parabolic change), i end up with an integro-differential form which does require approximating the integral numerically, which in turn will introduce truncation errors in the calculation. Now is there any way to get around this ? OR have i completely missed out some other alternative ? I would appreciate any and all help that someone can provide. Thanks in advance !