interpolated solution of ODE

Hi, all! How can I take an spline interpolated function that is evoluated solution of the differential equation? The first way is obtain solution as a vector of function values using ode45. After that take a spline using spline. I dislike this way because solving DE we already obtain function values and it derivatives hence it does not make sense to call spline.

Can I take 2 tables: function values and derivative values for the one call ode45 (or some else function...)?

I need to obtain the solution of a DE as a spline function, but I would not like to call spline function because it leads to big calculation error. I can explain my question on the following example (there is Scilab code):

//////////////////////////////////// 
// solve DE x' = -x^2 
function [dx] = rhs(t, x) 
dx = - x ^ 2; 
endfunction
//////////////////////////////////// 
// solving... 
t = [0:0.1:1]'; 
func_values = ode(1, 0, t, rhs)';
//////////////////////////////////// 
// take a spline 
deriv_values = rhs(t, func_values); 
deriv_values2 = splin(t, func_values);
function y = f1(q) 
y = (interp(q, t, func_values, deriv_values) - 1 ./ (1 + q))^2; 
endfunction
function y = f2(q) 
y = (interp(q, t, func_values, deriv_values2) - 1 ./ (1 + q))^2; 
endfunction
// calculate integral quadratic error 
disp(intg(0, 1, f1)); // err = 1.648D-12 
disp(intg(0, 1, f2)); // err = 9.577D-11

[Jan]

I do not understand, if you are looking for a spline method, or want to avoid it.

[cupuyc a]

I fixed my question. If I obtain spline coefficients myself, the calculation error is 2 orders smaller then I obtain spline coefficients using *spline* function. It is clear, I know coefficients of the spline because I know derivative values in the corresponding points (I'm able to call function rhs() for any point). But the function *spline* doesn't use this information...

Is any way exist to obtain the function values and coefficients of the corresponding spline?