change the time range of ode45 makes the integral of postprocess different

i want to solve the ode and get y(t),y'(t) and p1=c2*(y'(t))^2/r and pa=int(p1,40,40+t0),it means the int time range is the same but when i chnage the time range of ode45,for example i change it from 100 to 50,or other value ,why the pa changed and how to correct it?

And also when I change ode45 to ode23 or ode23s it is also different

And when I change tspan = [0 50] to tspan = [0:0.01:50] it is also different

And when I change the options = odeset('RelTol',1e-4,'AbsTol',[1e-4 1e-5]) ii is also different

So why all of them influence and how should I change and which answer is the most accurate?

    %
    function dy = fun1(t,y)
    dy = zeros(2,1);
    dy(1) =cos(t);
    dy(2) = -sin(t);
    %

and

    %
    function [p1] = my_p1( x1t )
    c2=32000;
    r=0.2;

    p1=c2*x1t;

    end
    %

and

    %main program
    %to solve the ode
    tspan = [0 50];
    options = odeset('RelTol',1e-4,'AbsTol',[1e-4 1e-5]);
    y0 = [0;1];
     [T,Y] = ode45('fun1',tspan,y0,options);

    figure(1)
    plot(T,Y(:,1),'-')
    %to postprocess the ode45's results and p1=c2*(y'(t))^2
    %and pa=int(p1,40,40+T0)/T0 
    dx=Y(:,2);

    p1=zeros(length(dx),1);   
    for i=1:1:length(dx);
      p1(i)=my_p1(dx(i));
    end

    figure(3)      
    plot(T,p1,'-r');

    T0=2*pi;

    t=linspace(40,40+T0,length(p1));
    t=t';

    pa=1/T0*trapz(t,p1)
%when tspan = [0 50],it is 964.6038 while when tspan = [0 100],it is  -138.9761

Answer 1

[Jan]

 when I change ode45 to ode23 or ode23s it is also different

Of course it is different. These are 3 different integrators. Do not interchange solvers for stiff and non-stiff problems.

 when I change tspan = [0 50] to tspan = [0:0.01:50] it is also different

This concerns the trapz command at the end? The question is not clear to me. On the top of the question your mentioned "when i chnage the time range of ode45,for example i change it from 100 to 50". So where does which time range matter? And is it really surprising, that the output changes, if you integrate over a different range?

 when I change the options = odeset('RelTol',1e-4,'AbsTol',[1e-4 1e-5]) it is also different

Of course, otherwise the tolerances would be useless. They are used to determine the _optimal_step size, which is a balance between the local discretization error and the accumulated rounding error.

 when I change tspan = [0 50] to tspan = [0:0.01:50] it is also different

Please post the code to reproduce, what you mean here. Differences caused by rounding are expected.

 which answer is the most accurate

It depends on how you define "accurate". To run a numerical intergration, you specify a tolerance to control the step size. Then the trajectory is kept inside a certain range around the estimated correct solution. Reducing the step size will decrease the local discretization error (roughly: The error caused by approximating the trajectory by a polynomial locally), but increase the accumulated rounding errors. In consequence the computed trajectory depends on the tolerances and it is not meaningful to assess the "accuracy" without considering the tolerance.

Of course you can compare the result computed by a numerical integration with the symbolic solution (if there is one), but this remains a comparison between apples and fish. This gets immediately clear, if you process an ODE, which does not have an analytical solution. Then asking "what is the accurate solution" becomes meaningless. Even with a deterministic analytical solution, a numerical integration does have a limited power concerning the "accuracy": Simulate the throwing of a dice by integrating the simple equations of motion. You will never get a reliable "accurate" solution, which predicts the trajectory in a general case.

[dcydhb dcydhb]

but don't you think it is queer that i just change the tspan of ode45 and i have not change the int range,so why the result is dufferent?

[Torsten]

Please post your code to reproduce the differences. Are they visible to the naked eye ?

[dcydhb dcydhb]

i mean when tspan is from 0 to 50,the pa is diffetent from that tspan is from 0 to 100,it is queer because the int range to get pa is never changed

[dcydhb dcydhb]

is it not clear that when you change the tspan then you can see that the pa changed?

[Torsten]

The range of integration for pa does not change, but p1 used in the call to "trapz" changes.

[dcydhb dcydhb]

yes it is,can you say it more specificly and how to correct it to make the result changeless

[Jan]

@dcydhb dcydhb: I still do not understand, what you change where. Please post some code, which reproduces your expectations. I'm convinced that this is not a bug of Matlab, such that there is either a bug in your code or a misunderstanding of the mathematical meaning. Maybe there are rounding effects in addition.

fun1  = @(t,y) [cos(t); -sin(t)];
my_p1 = @(x) 32000 * x;
tspan = [0 50];
options = odeset('RelTol',1e-4,'AbsTol',[1e-4 1e-5]);
[T, Y] = ode45(fun1, tspan, [0;1], options);

dx = Y(:,2);
p1 = my_p1(dx);

T0 = 2*pi;
t  = linspace(40, 40+T0, length(p1)).';
pa = 1 / T0 * trapz(t, p1)

Now you change what and observe what? If you change tspan, your trapz is applied to different data. Of course you get a different result. You do not crop out the part of the trajectory Y for the times [40, 40+T0], but you have completely different dx and in consequence p1 and let trapz only assume, that the data are available for the same time.

t   = linspace(40, 40+T0, length(p1)).';
pa1 = 1 / T0 * trapz(t, p1)
t   = linspace(21789, 21789+T0, length(p1)).';
pa2 = 1 / T0 * trapz(t, p1)

The same result is expected.

[dcydhb dcydhb]

For example

    fun1  = @(t,y) [cos(t); -sin(t)];
    my_p1 = @(x) 32000 * x;
    tspan = [0 50];% here is the difference
    options = odeset('RelTol',1e-4,'AbsTol',[1e-4 

    1e-5]);
    [T, Y] = ode45(fun1, tspan, [0;1], options);
    dx = Y(:,2);
    p1 = my_p1(dx);
    T0 = 2*pi;
    t  = linspace(40, 40+T0, length(p1)).';
    pa = 1 / T0 * trapz(t, p1)

it runs out that pa= 964.6038

but when

    fun1  = @(t,y) [cos(t); -sin(t)];
    my_p1 = @(x) 32000 * x;
    tspan = [0 100];% here is the difference

    options = odeset('RelTol',1e-4,'AbsTol',[1e-4 

    1e-5]);
    [T, Y] = ode45(fun1, tspan, [0;1], options);
    dx = Y(:,2);
    p1 = my_p1(dx);
    T0 = 2*pi;
    t  = linspace(40, 40+T0, length(p1)).';
    pa = 1 / T0 * trapz(t, p1)

it runs out pa=-138.9761

I mean this and why?

[dcydhb dcydhb]

and why

 If you change tspan, your trapz is applied to different data.

[dcydhb dcydhb]

p1 is of the original data and t is of the original time range and why

   your trapz is applied to different data.

[dcydhb dcydhb]

for example i use ode45 to solve an ode and get x(t) and i want to int x(t) in some time range,for example the fixed time range 0-10s,so if i change the tspan of ode45 will change this value of int?

it is queer