Parameter calculation by using loops in ODE45

For any initial values(Ofcourse with limits) I used it either excutes on the 4th iteration or it would be beyond the MAX iteration.

I am afraid something is wrong and I need a help.

zspan=[0,400];
v0mat = [1 0.05 1];
tol = 10^-4; % Treshold
MAX = 1000;

v2 = 0.4;
interval = [0.01 0.2]; % alpha interval 
a = interval(1);
b = interval(2);
alpha = (a+b)/2;
zsol = {};
v1sol = {};
v2sol = {};
v3sol = {};

for k=1:size(v0mat,1)
  v0=v0mat(k,:);
  [z,v]=ode45(@(z,v)rhs(z,v,alpha),zspan,v0);
  zsol{k}=z;
  v1sol{k}=v(:,1);
  v2sol{k}=v(:,2);
  v3sol{k}=v(:,3);
end
iter = 1;
while(abs((v(:, 2)) - v2) > tol)
 alpha= (a + b)/2;
[z,v]=ode45(@(z,v)rhs(z,v,alpha),interval,v0);
 if(abs(v(:,2))-v2 > 0)
b = alpha;
 else
 a = alpha;
 end
 
 iter = iter + 1;
 if(iter > MAX)
 return;
 end
end

for k=1:size(v0mat,1)
    figure(1)
     plot(v2sol{k},zsol{k},'g')
      xlabel('Velocity,w')
       ylabel('Height, z')
      grid on       
end

function parameters=rhs(z,v,alpha)
            Nsqr = 0.001;
            db= 2*alpha-(v(1).*v(3))./(2*v(2).^2);
            dw= (v(3)./v(2))-(2*alpha*v(2)./v(1));
            dgmark= -Nsqr-(2*alpha*v(3)./v(1)); 
   parameters=[db;dw;dgmark];
end

[darova]

It is unclear what are you trying to do here. Are you trying to check the last value?

if(abs(v(:,2))-v2 > 0)

[Dereje]

My main goal here is to get the value of the parameter alpha in the above ODE 45 equation. I just randomly pick v2 as a reference so that I would do the iterations. Is there any better way that I can do this?

[darova]

What is condition that alpha depends on?

[Dereje]

It doesn’t depend on anything. Also not a constant.

[darova]

How do you understand this line? What does it mean to you

if(abs(v(:,2))-v2 > 0)

[Dereje]

It subtracts v2 from the second value in the row (in this case the first intial value is 0.05). I am

expecting the loop to exciute this until I get a solution i.e. when v(2) is zero.

[darova]

But v(:,2) has many values. Which one should be greater than v2

[Dereje]

Since there is only one row, we have only a signle value not many values I think.

[darova]

I believe it's a column

[Walter Roberson]

There will be one row of v for every timepoint.

[Dereje]

Yes its a 2nd column.

[Star Strider]

@Dereje — Parameter Estimation for a System of Differential Equations might be informative. It is not possible to determine that from the posted code.

[Walter Roberson]

It is a column with one row for each timepoint. That makes (abs(v(:,2))-v2 > 0) a vector result that might be true for some of the time entries but false for other times. Which one time do you need it to be true for? If you want to be testing for the last time then

if (abs(v(end,2))-v2 > 0)

[Dereje]

@Walter: I really don’t understand what you meant by time point? You mean the time of the iterations when you said the last time?

[Walter Roberson]

In your context, timepoint is z value. Your ode45 call will process over z in zspan, producing one row for each of a number of z values. (abs(v(:,2))-v2 > 0) is a vector result that might be true for some of the z entries but false for other times. You need to decide which z the condition must hold for.

[Dereje]

I think It doesn't matter if the condition holds for , for example the last time. But the problem is here still this line of code

if (abs(v(end,2))-v2 > 0)

doesn't hold(True).