I couldn't understand the problem in this code.
2 views (last 30 days)
Show older comments
SAHIL SAHOO
on 30 Dec 2022
ti = 0;
tf = 70E-8;
tspan=[ti tf];
k = (1E-3);
KC = k;
here I want to make the k and yita_mn time dependent.
y0= [(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
((-3.14).*rand(20,1) + (3.14).*rand(20,1));
0];
yita_mn = [
0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1;
1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 1 0 1 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 0 1 0 1 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 0 0 1 0 1 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 0 0 0 1 0 1 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 0 0 0 0 1 0 1 0 0 0 0 0 0 0;
0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 0 0 0 0 0;
0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 0 0 0 0;
0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 0 0 0;
0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 0 0;
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 0;
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0;
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1;
1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0
].*(KC);
N = 20;
tp = 1E-12;
o = sort(10e4*rand(1,20),'ascend');
[T,Y]= ode45(@(t,y) rate_eq(t,y,yita_mn,N,o,k),tspan./tp,y0);
r = ((1/20).*( exp(i.*Y(:,3)) + exp(i.*Y(:,6)) + exp(i.*Y(:,9)) + exp(i.*Y(:,12)) + exp(i.*Y(:,15)) ...
+exp(i.*Y(:,18)) +exp(i.*Y(:,21)) +exp(i.*Y(:,24)) + exp(i.*Y(:,27)) + exp(i.*Y(:,30)) + exp(i.*Y(:,33)) ...
+ exp(i.*Y(:,36)) + exp(i.*Y(:,39)) +exp(i.*Y(:,42)) + exp(i.*Y(:,45)) + exp(i.*Y(:,48)) + exp(i.*Y(:,51)) + exp(i.*Y(:,54))+ exp(i.*Y(:,57)) + exp(i.*Y(:,60))));
figure(1)
plot(T,(Y(:,61)),'linewidth',0.8);
hold on
plot(T,(Y(:,62)),'linewidth',2);
plot(T,(Y(:,63)),'linewidth',2);
plot(T,(Y(:,64)),'linewidth',2);
plot(T,(Y(:,65)),'linewidth',2);
plot(T,(Y(:,66)),'linewidth',2);
plot(T,(Y(:,67)),'linewidth',2);
plot(T,(Y(:,68)),'linewidth',2);
plot(T,(Y(:,69)),'linewidth',2);
plot(T,(Y(:,70)),'linewidth',2);
plot(T,(Y(:,71)),'linewidth',2);
plot(T,(Y(:,72)),'linewidth',2);
plot(T,(Y(:,73)),'linewidth',2);
plot(T,(Y(:,74)),'linewidth',2);
plot(T,(Y(:,75)),'linewidth',2);
plot(T,(Y(:,76)),'linewidth',2);
plot(T,(Y(:,77)),'linewidth',2);
plot(T,(Y(:,78)),'linewidth',2);
plot(T,(Y(:,79)),'linewidth',2);
plot(T,(Y(:,80)),'linewidth',2);
hold off
grid on
xlabel("time")
ylabel("phase difference")
set(gca,'fontname','times New Roman','fontsize',18,'linewidth',1.8);
figure(2)
plot(T,(Y(:,2)),'linewidth',2);
hold on
plot(T,(Y(:,5)),'linewidth',2);
plot(T,(Y(:,8)),'linewidth',2);
plot(T,(Y(:,11)),'linewidth',2);
plot(T,(Y(:,14)),'linewidth',2);
plot(T,(Y(:,17)),'linewidth',2);
plot(T,(Y(:,20)),'linewidth',2);
plot(T,(Y(:,23)),'linewidth',2);
plot(T,(Y(:,26)),'linewidth',2);
plot(T,(Y(:,29)),'linewidth',2);
plot(T,(Y(:,32)),'linewidth',2);
plot(T,(Y(:,35)),'linewidth',2);
plot(T,(Y(:,38)),'linewidth',2);
plot(T,(Y(:,41)),'linewidth',2);
plot(T,(Y(:,44)),'linewidth',2);
plot(T,(Y(:,47)),'linewidth',2);
plot(T,(Y(:,50)),'linewidth',2);
plot(T,(Y(:,53)),'linewidth',2);
plot(T,(Y(:,56)),'linewidth',2);
plot(T,(Y(:,59)),'linewidth',2);
hold off
grid on
xlabel("time")
ylabel("amplitude")
set(gca,'fontname','times New Roman','fontsize',18,'linewidth',1.8);
figure(3)
plot(T,abs(r))
xlabel("time")
ylabel("orderparameter")
grid on
set(gca,'fontname','times New Roman','fontsize',18,'linewidth',1.8);
function dy = rate_eq(t,y,yita_mn,N,o,k)
dy = zeros(4*N + 2,1);
dGdt = zeros(N,1);
dAdt = zeros(N,1);
dOdt = zeros(N,1);
P = 0.5;
a = 0;
T = 2000;
tp = 1E-12;
Gt = y(1:3:3*N-2);
At = y(2:3:3*N-1);
Ot = y(3:3:3*N-0);
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) + o(1,i).*tp;
for j = 1:N
dAdt(i) = dAdt(i) + yita_mn(i,j)*(At(j))*cos(Ot(j)-Ot(i));
dOdt(i) = dOdt(i) + yita_mn(i,j)*((At(j)/At(i)))*sin(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:20)';
n2 = circshift(n1,-1);
n61 = n1 +60;
n62 = circshift(n61,-1);
n80 = circshift(n61,1);
j2 = 3*(1:20)-1;
j5 = circshift(j2,-1);
j8 = circshift(j2,-2);
j59 = circshift(j2,1);
dy(n61) = (o(1,n2).' - o(1,n1).').*tp + a.*(Gt(n2) - Gt(n1)) - (k).*(y(j2)./y(j5)).*sin(y(n61)) - (k).*(y( j5)./y(j2)).*sin(y(n61)) + (k).*(y(j8)./y(j5)).*sin(y(n62)) + (k).*(y(j59)./y(j2)).*sin(y(n80));
dy(81) = k;
dy(82) = yita_mn;
end
0 Comments
Accepted Answer
Sulaymon Eshkabilov
on 30 Dec 2022
Here is the corrected code. There were a number of unreadable white spaces (due to incorrect copy and paste), i (mean to be imaginary number and thus, it was substituted with 1i) and one missing IC for y0, etc.;
ti = 0;
tf = 70E-8;
tspan=[ti tf];
k = (1E-3);
KC = k;
y0= [(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
(1e-2)*rand(2,1); ((-3.14).*rand(1,1) + (3.14).*rand(1,1));
((-3.14).*rand(20,1) + (3.14).*rand(20,1)); 0; 0];
yita_mn = [ 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1;
1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 1 0 1 0 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 1 0 1 0 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 0 1 0 1 0 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 0 0 1 0 1 0 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 0 0 0 1 0 1 0 0 0 0 0 0 0 0;
0 0 0 0 0 0 0 0 0 0 1 0 1 0 0 0 0 0 0 0;
0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 0 0 0 0 0;
0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 0 0 0 0;
0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 0 0 0;
0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 0 0;
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 0;
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0;
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1;
1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 ].*(KC);
N = 20;
tp = 1E-12;
o = sort(10e4*rand(1,20),'ascend');
[T,Y]= ode45(@(t,y) rate_eq(t,y,yita_mn,N,o,k),tspan./tp,y0);
r = ((1/20).*( exp(1i*Y(:,3)) + exp(1i*Y(:,6)) + exp(1i.*Y(:,9)) + exp(1i.*Y(:,12)) + exp(1i.*Y(:,15)) ...
+exp(1i.*Y(:,18)) +exp(1i.*Y(:,21)) +exp(1i.*Y(:,24)) + exp(1i.*Y(:,27)) + exp(1i.*Y(:,30)) + exp(1i.*Y(:,33)) ...
+ exp(1i.*Y(:,36)) + exp(1i.*Y(:,39)) +exp(1i.*Y(:,42)) + exp(1i.*Y(:,45)) + exp(1i.*Y(:,48)) + exp(1i.*Y(:,51)) + exp(1i.*Y(:,54))+ exp(1i.*Y(:,57)) + exp(1i.*Y(:,60))));
figure(1)
plot(T,(Y(:,61)),'linewidth',0.8);
hold on
plot(T,(Y(:,62)),'linewidth',2);
plot(T,(Y(:,63)),'linewidth',2);
plot(T,(Y(:,64)),'linewidth',2);
plot(T,(Y(:,65)),'linewidth',2);
plot(T,(Y(:,66)),'linewidth',2);
plot(T,(Y(:,67)),'linewidth',2);
plot(T,(Y(:,68)),'linewidth',2);
plot(T,(Y(:,69)),'linewidth',2);
plot(T,(Y(:,70)),'linewidth',2);
plot(T,(Y(:,71)),'linewidth',2);
plot(T,(Y(:,72)),'linewidth',2);
plot(T,(Y(:,73)),'linewidth',2);
plot(T,(Y(:,74)),'linewidth',2);
plot(T,(Y(:,75)),'linewidth',2);
plot(T,(Y(:,76)),'linewidth',2);
plot(T,(Y(:,77)),'linewidth',2);
plot(T,(Y(:,78)),'linewidth',2);
plot(T,(Y(:,79)),'linewidth',2);
plot(T,(Y(:,80)),'linewidth',2);
hold off
grid on
xlabel("time")
ylabel("phase difference")
set(gca,'fontname','times New Roman','fontsize',18,'linewidth',1.8);
figure(2)
plot(T,(Y(:,2)),'linewidth',2);
hold on
plot(T,(Y(:,5)),'linewidth',2);
plot(T,(Y(:,8)),'linewidth',2);
plot(T,(Y(:,11)),'linewidth',2);
plot(T,(Y(:,14)),'linewidth',2);
plot(T,(Y(:,17)),'linewidth',2);
plot(T,(Y(:,20)),'linewidth',2);
plot(T,(Y(:,23)),'linewidth',2);
plot(T,(Y(:,26)),'linewidth',2);
plot(T,(Y(:,29)),'linewidth',2);
plot(T,(Y(:,32)),'linewidth',2);
plot(T,(Y(:,35)),'linewidth',2);
plot(T,(Y(:,38)),'linewidth',2);
plot(T,(Y(:,41)),'linewidth',2);
plot(T,(Y(:,44)),'linewidth',2);
plot(T,(Y(:,47)),'linewidth',2);
plot(T,(Y(:,50)),'linewidth',2);
plot(T,(Y(:,53)),'linewidth',2);
plot(T,(Y(:,56)),'linewidth',2);
plot(T,(Y(:,59)),'linewidth',2);
hold off
grid on
xlabel("time")
ylabel("amplitude")
set(gca,'fontname','times New Roman','fontsize',18,'linewidth',1.8);
figure(3)
plot(T,abs(r))
xlabel("time")
ylabel("orderparameter")
grid on
set(gca,'fontname','times New Roman','fontsize',18,'linewidth',1.8);
function dy = rate_eq(t,y,yita_mn,N,o,k)
dy = zeros(4*N + 2,1);
dGdt = zeros(N,1);
dAdt = zeros(N,1);
dOdt = zeros(N,1);
P = 0.5;
a = 0;
T = 2000;
tp = 1E-12;
Gt = y(1:3:3*N-2);
At = y(2:3:3*N-1);
Ot = y(3:3:3*N-0);
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) + o(1,i).*tp;
for j = 1:N
dAdt(i) = dAdt(i) + yita_mn(i,j)*(At(j))*cos(Ot(j)-Ot(i));
dOdt(i) = dOdt(i) + yita_mn(i,j)*((At(j)/At(i)))*sin(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:20)';
n2 = circshift(n1,-1);
n61 = n1 +60;
n62 = circshift(n61,-1);
n80 = circshift(n61,1);
j2 = 3*(1:20)-1;
j5 = circshift(j2,-1);
j8 = circshift(j2,-2);
j59 = circshift(j2,1);
dy(n61) = (o(1,n2).' - o(1,n1).').*tp + a.*(Gt(n2) - Gt(n1)) - (k).*(y(j2)./y(j5)).*sin(y(n61)) - (k).*(y( j5)./y(j2)).*sin(y(n61)) + (k).*(y(j8)./y(j5)).*sin(y(n62)) + (k).*(y(j59)./y(j2)).*sin(y(n80));
dy(81) = k;
dy(82) = 0;
end
0 Comments
More Answers (1)
Stephen23
on 30 Dec 2022
Somehow you have managed to include NO BREAK SPACEs in the your code:
rather than the normal, required ASCII SPACE character:
most likely this is due to either:
- copy pasting from some unreliable source (website, formatted document, etc)
- third-party editor adding these characters as a "feature".
Because they all appear to be leading characters, you can probably use the MATLAB editor to remove them by selecting all text and then pressing ctrl+i. In future, use a more reliable editor or a better code source.
I used Notepad++ to highlight in magenta the NO BREAK SPACE characters in your code:
0 Comments
See Also
Categories
Find more on Entering Commands 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!