ODE solvers solution accuracy at large values

I am comparing some simple ODE solvers, using sine or cosine, in order to get a feel for the relative accuracy at large values. I have previously used some FORTRAN code( Adams Method provided by LLNL) and have found typical losses of 1 order of accuracy per order of the solutions when comparing to the input values for Absolute tolerances and Relative tolerances.

I am obviously doing something incorrect as solutions though MATLAB are not anywhere close to this order of accuracy. For Example, solving for cos(t), I get something like {10^4} orders of accuracy loss each order of {10} value of t. If a larger range is chosen it makes relatively little difference in the order of accuracy of the output. I have tried 'refine' which makes no impact as well as adjusting the maximum steps size. Also I have used a few different ODE Solvers, ODE113 in particular. I also have created a loop shown below to solve at specific smaller times steps in order to compare to the other solvers.

Any help would be greatly appreciated, I also understand that the Relative Tolerance will be set to it's minimum value each run, but left as is, in case this is helping produce some error.

function main 

% a simple example to solve ODE's
% Uses ODE113 to solve 
% dy1 = -y2;
% dy2 = y1;
% set an error
options = odeset('RelTol',1e-15,'AbsTol', 1e-15);
% define the interval to solve over repeatedly
tstep = 0.01;
% final value of t to integrate to
tfinal = 10^4;
% number of steps to get to tfinal
N = tfinal/tstep;
% initial conditions for firsat run
Yo = [1;0];
for n=0:N

    StepN=n*tstep;
    %timespan
    tspan = [n*0.1,0.1+ n*0.1]; 
    %call the solver
    [T,Y] = ode113(@rigid113,tspan,Yo,options); 

return
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [dy] = rigid113(t,y) 
%a function which returns a rate of change vector
dy = zeros(2,1);    % a column vector
dy(1) = -y(2);
dy(2) = y(1);

Thank you for any help, I'd love to figure out what to do to improve the solution!

[Jan]

Absolute and relative tolerances in a magnitude of 1e-15 are extremely small and I would not expect, that an ODE solver can handle them with an sufficient accuracy.