Variation of initial values and slope calculation

Hi everyone,

I'd like to plot a function with different initial values, and calculate the slope(v0) between the 2. and 3. point, but it doesn't work. Could you please help me.

Thank you!!

Hier is my code:

% Main program
n=15;
s00=logspace(1,3,n);
for k=1:n;
    s0=s00(k);  
    v0(k)=(x(4)(3)-x(4)(2))/(t(3)-t(2));
end
e0=2;
c0=0;
p0=0;
t=linspace(0,1000,100);
x0=[s0 c0 e0 p0];

[t,x]=ode45(@(t,x) num_lsg1(t,x), t, x0);  % function call
%figure
plot(t,x(:,4),'-*')
hold all;
xlabel('Time [s]');
ylabel('Product');
%%Definition function
function dx=num_lsg1(t,x) 
% Definition parameters
k1=0.01;
km1=10^(-6);
k2=1;
dx=zeros(4,1);
dx(1) = km1*x(2)-k1*x(1)*x(3);
dx(2) = k1*x(1)*x(3)-(km1+k2)*x(2);
dx(3) = (km1+k2)*x(2)-k1*x(1)*x(3);
dx(4) = k2*x(2);
end

[Stephan]

What do you want to achieve with this (invalid) syntax:

(x(4)(3)-x(4)(2))

in line 6?

[Tatjana Henning]

I want to calculate the slope of the graph from x(4) (dx(4) = k2*x(2)) between the second and third point.

Answer 1

[Stephan]

Hi,

it appears that there are at least 2 problems in your code. Due to this (and the fact that i do have to understand fully how exactly you want to calculate the slope), i suggest to solve this question stepwise.

1. Calculate with different initial values

%% Main program
% Constant values
e0=2;
c0=0;
p0=0;
timesteps = 100;
t=linspace(0,1000,timesteps);
n=15;
s0=logspace(1,3,n);
result = zeros(timesteps,5,n);%Preallocate a 3D array to store ALL results

% Calculate the results for all s0
hold on
for k=1:n
    x0=[s0(k) c0 e0 p0];
    [t,x]=ode45(@(t,x) num_lsg1(t,x), t, x0);  % function call
    result(:,1,k) = t; % All t values in the first column in every "plane"
    result(:,2:5,k) = x(:,1:4); %All results in the columns 2-5 in every "plane"
    plot(result(:,1,k),result(:,4,k),'-*')
end

%figure details
xlabel('Time [s]');
ylabel('Product');
hold off

%Part 2 - Calculate slope
%v0(k)=(x(4)(3)-x(4)(2))/(t(3)-t(2));

%%Definition function
function dx=num_lsg1(~,x) 
% Definition parameters
k1=0.01;
km1=10^(-6);
k2=1;
dx=zeros(4,1);
dx(1) = km1*x(2)-k1*x(1)*x(3);
dx(2) = k1*x(1)*x(3)-(km1+k2)*x(2);
dx(3) = (km1+k2)*x(2)-k1*x(1)*x(3);
dx(4) = k2*x(2);
end

The code so far gives this result:

2. Calculate the slope

Now - for part 2 i have to understand how exactly you want to calculate slope - can you try to explain in detail? The only problem in your code is, that you use an invalid syntax and this should be possible to fix easily, once we have understood what you want to do. As far as i understand i would say you want the slope of all the values of x(4) at the times t=2 and t=3 for all 15 curves. Is that correct?

If we assume that this is correct, you can calculate the slope this way (no for loop is needed for this):

%Part 2 - Calculate slope
v0 = squeeze(result(3,5,:)-result(2,5,:))./squeeze(result(3,1,:)-result(2,1,:));

In words:

Take the 3. entry in result column 5 (is equal to x(4) at timstep=3) for all 15 values minus the corresponding 2. entry in result column 5 (is equal to x(4) at timstep=2) for all 15 values... - This is how indexing woorks in Matlab.

The only "special" thing to know is squeeze, to get a 15x1 result for v0.

The result then is:

v0 =

    0.1430
    0.1945
    0.2625
    0.3504
    0.4610
    0.5953
    0.7510
    0.9218
    1.0978
    1.2677
    1.4216
    1.5537
    1.6619
    1.7474
    1.8133

Best regards

Stephan

[madhan ravi]

+1

[Tatjana Henning]

Hi Stephan,

Thank you very much for the solution of my problem!! It works. I only changed the line 18 to plot(result(:,1,k),result(:,5,k),'-*'). And yes, the 2. question is correct. :)

[Stephan]

nice, yes i see - i plotted x(3) over t.