how to rid of this warning and get correct solution?
    4 views (last 30 days)
  
       Show older comments
    
    SAHIL SAHOO
 on 14 Oct 2022
  
    
    
    
    
    Answered: Gokul Nath S J
    
 on 18 Oct 2022
            ti = 0;
tf = 1E-7;
tspan=[ti tf];
k = 5E-3;
h = 10E-2;
y0= [(h)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
     (h)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
     (h)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1)); 
     (h)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1)); 
     (h)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
      ((-3.14).*rand(5,1) + (3.14).*rand(5,1))];
yita_mn = [
    0 1 0 0 1;
    1 0 1 0 0; 
    0 1 0 1 0;
    0 0 1 0 1;
    1 0 0 1 0;
      ]*(k);
N = 5;
tp = 1E-12;
[T,Y]= ode45(@(t,y) rate_eq(t,y,yita_mn,N),tspan./tp,y0);
figure(1)
plot(T,(Y(:,16)),'linewidth',0.8);
hold on 
for m = 16:20
    plot(T,(Y(:,m)),'linewidth',0.8);
end 
hold off
grid on
xlabel("time")
ylabel("phase difference")
set(gca,'fontname','times New Roman','fontsize',18,'linewidth',1.8);
function dy = rate_eq(t,y,yita_mn,N,o)
dy = zeros(4*N,1);
dGdt = zeros(N,1);
dAdt = zeros(N,1);
dOdt = zeros(N,1);
P = 1.25;
a = 5;
T = 2000;
Gt = y(1:3:3*N-2);
At = y(2:3:3*N-1);
Ot = y(3:3:3*N-0);
k = 5E-5;
for i = 1:N
    dGdt(i) = (P - Gt(i) - (1 + 2.*Gt(i)).*(At(i))^2)./T ; 
    dAdt(i) = (Gt(i).*(At(i)));
    dOdt(i) = -a.*(Gt(i));  
    for j = 1:N
        dAdt(i) = dAdt(i)+yita_mn(i,j).*(At(j))*sin(Ot(j)-Ot(i));
        dOdt(i) = dOdt(i)+yita_mn(i,j).*((At(j)/At(i)))*cos(Ot(j)-Ot(i));
  end
end
dy(1:3:3*N-2) = dGdt;
dy(2:3:3*N-1) = dAdt;
dy(3:3:3*N-0) = dOdt;
n1 = (1:5)';
n2 = circshift(n1,-1);
n16 = n1 + 15;
n17 = circshift(n16,-1);
n20 = circshift(n16,1);
j2 =  3*(1:5)-1;
j5 = circshift(j2,-1);
j8 = circshift(j2,-2);
j19 = circshift(j2,1);
dy(n16) =  -a.*(Gt(n2)-Gt(n1)) + (k).*(y(j2)./y(j5)).*cos(y(n16)) - (k).*(y( j5)./y(j2)).*cos(y(n16)) + (k).*(y(j8)./y(j5)).*cos(y(n17)) - (k).*(y(j19)./y(j2)).*cos(y(n20));  
end
0 Comments
Accepted Answer
  Gokul Nath S J
    
 on 18 Oct 2022
        Dear Sahil Sahoo, 
I have tried running your code on my machine. By setting the relative and absolute error values to 1e-2, the warning messages can be avoided.  
options = odeset('RelTol',1e-2,'AbsTol',1e-2); 
The code should be included before calling ode45 (line 22) 
However, at the same time, it is essential to check the validity of the eventual answer.  
Also, I have found similar cases that might be helpful. Refer to the link below for the same. 
0 Comments
More Answers (0)
See Also
Categories
				Find more on Numerical Integration and Differential Equations in Help Center and File Exchange
			
	Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!

