while solving a couple ordinary differential equation, I'm getting an error of Unable to perform assignment because the left and right sides have a different number of element

1 次查看（过去 30 天）
SAHIL SAHOO
SAHIL SAHOO 2022-7-6
回答： Johan 2022-7-6
output:-
Unable to perform assignment because the left and right sides have a different number of element.
Error in laser2 (line 52)
A1(i+1) = A1(i) + (h/6)*(k1+(2*(k2+k3))+k4);
input:-
format long
a = 0;
b = 30;
h = 0.1;
o = 3.5;
tc = 3.0;
tf = 2.3;
a1 = 0.1;
a2 = 0.1;
P1 = 0.2;
P2 = 0.2;
k= 0.9;
V = 1;
A1(1) = 10E-4;
A2(1) = 10E-5;
G1(1) = 10E-6;
G2(1) = 10E-6;
O(1) = 0;
n = (b-a)/h ;
t = ((a) + (0:n)*h).*10E-6;
func1 = @(G1,A1,A2,O) (1/tc).*(G1- a1).*(A1./V) +(k./tc).*(A2./V).*cos(O);
func3 = @(G2,A2,A1,O) (1/tc).*(G2- a2).*(A2./V) +(k./tc).*(A1./V).*cos(O);
func2 = @(G1) (1/tf).*(P1 - G1.*(A1.^2+1));
func4 = @(G2) (1/tf).*(P2 - G2.*(A2.^2+1));
func5 = @(A1,A2,O) o - (k./tc).*((A1./A2) + (A2./A1)).*sin(O);
for i = 1:n
k1 = feval(func1,G1(i),A1(i),A2(i),O(i));
l1 = feval(func2, G1(i));
m1 = feval(func3, G2(i),A1(i),A2(i),O(i));
n1 = feval(func4, G2(i));
q1 = feval(func5,A1(i),A2(i),O(i));
k2 = feval(func1, A2(i)+(m1/2)*h, A1(i)+(k1/2)*h, G1(i)+(l1/2)*h, O(i)+(q1/2)*h);
l2 = feval(func2, G1(i)+(l1*h/2));
m2 = feval(func3, A2(i)+(m1/2)*h, A1(i)+(k1/2)*h, G2(i)+(n1/2)*h, O(i)+(q1/2)*h);
n2 = feval(func4, G2(i)+(n1/2)*h);
q2 = feval(func5, A2(i)+(m1/2)*h, A1(i)+(k1/2)*h, O(i)+(q1/2)*h);
k3 = feval(func1, A2(i)+(m2/2)*h, A1(i)+(k2/2)*h, G1(i)+(l2/2)*h, O(i)+(q2/2)*h);
l3 = feval(func2, G1(i)+(l2*h/2));
m3 = feval(func3, A2(i)+(m2/2)*h, A1(i)+(k2/2)*h, G2(i)+(n2/2)*h, O(i)+(q2/2)*h);
n3 = feval(func4, G2(i)+(n2/2)*h);
q3 = feval(func5, A2(i)+(m2/2)*h, A1(i)+(k2/2)*h, O(i)+(q2/2)*h );
k4 = feval(func1, A2(i)+(m3*h), A1(i)+(k3*h), G1(i)+(l3*h), O(i)+(q3*h));
l4 = feval(func2, G1(i)+(l3*h));
m4 = feval(func3, A2(i)+(m3*h), A1(i)+(k3*h), G2(i)+(n3*h), O(i)+(q3*h));
n4 = feval(func4, G2(i)+(n3*h));
q4 = feval(func5, A2(i)+(m3*h), A1(i)+(k3*h), O(i)+(q3*h));
A1(i+1) = A1(i) + (h/6)*(k1+(2*(k2+k3))+k4);
G1(i+1) = G1(i) + (h/6)*(l1+(2*(l2+l3))+l4);
A2(i+1) = A2(i) + (h/6)*(m1+(2*(m2+m3))+m4);
G2(i+1) = G2(i) + (h/6)*(n1+(2*(n2+n3))+n4);
O(i+1) = O(i) + (h/6)*(q1+(2*(q2+q3))+q4);
end
subplot 211
plot(t,A1)
grid on
subplot 212
plot(t,O)
it maybe due to the func2 = @(G1) (1/tf).*(P1 - G1.*(A1.^2+1));
func4 = @(G2) (1/tf).*(P2 - G2.*(A2.^2+1));
this two code, but I couldn't understand what's the problem in it.

请先登录，再回答此问题。

采纳的回答

Johan
Johan 2022-7-6
Ran in:
func2 and func4 uses A1 and A2 respectively which means they are using the whole arrays, either call A1 and A2 in the func definition (see example below) or set it to A1(1) and A2(1) if you only want to use the first term.
On a side note it's considered bad practice to add 1 2 3 .... to your variable names. It makes your code hard to understand and debug there is a post that explain why here : https://fr.mathworks.com/support/search.html/answers/304528-tutorial-why-variables-should-not-be-named-dynamically-eval.html
%change some variable names so that they make sense
time_start = 0;
time_end = 30;
time_step = 0.1;
o = 3.5;
tc = 3.0;
tf = 2.3;
a1 = 0.1;
a2 = 0.1;
P1 = 0.2;
P2 = 0.2;
k= 0.9;
V = 1;
n_step = (time_end-time_start)/time_step ;
time = ((time_start) + (0:n_step)*time_step).*10E-6;
func1 = @(G1,A1,A2,O) (1/tc).*(G1- a1).*(A1./V) +(k./tc).*(A2./V).*cos(O);
func3 = @(G2,A2,A1,O) (1/tc).*(G2- a2).*(A2./V) +(k./tc).*(A1./V).*cos(O);
func2 = @(G1,A1) (1/tf).*(P1 - G1.*(A1.^2+1));
func4 = @(G2,A2) (1/tf).*(P2 - G2.*(A2.^2+1));
func5 = @(A1,A2,O) o - (k./tc).*((A1./A2) + (A2./A1)).*sin(O);
%Initialize arrays
A1 = zeros(n_step,1);
A2 = A1;
G1 = A1;
G2 = A1;
O = A1;
A1(1) = 10E-4;
A2(1) = 10E-5;
G1(1) = 10E-6;
G2(1) = 10E-6;
O(1) = 0;
for i_loop = 1:n_step
k1 = func1(G1(i_loop),A1(i_loop),A2(i_loop),O(i_loop));
l1 = func2( G1(i_loop),A1(i_loop));
m1 = func3( G2(i_loop),A1(i_loop),A2(i_loop),O(i_loop));
n1 = func4( G2(i_loop),A2(i_loop));
q1 = func5(A1(i_loop),A2(i_loop),O(i_loop));
k2 = func1( A2(i_loop)+(m1/2)*time_step, A1(i_loop)+(k1/2)*time_step, G1(i_loop)+(l1/2)*time_step, O(i_loop)+(q1/2)*time_step);
l2 = func2( G1(i_loop)+(l1*time_step/2), A1(i_loop)+(k1/2)*time_step);
m2 = func3( A2(i_loop)+(m1/2)*time_step, A1(i_loop)+(k1/2)*time_step, G2(i_loop)+(n1/2)*time_step, O(i_loop)+(q1/2)*time_step);
n2 = func4( G2(i_loop)+(n1/2)*time_step, A2(i_loop)+(m1/2)*time_step);
q2 = func5( A2(i_loop)+(m1/2)*time_step, A1(i_loop)+(k1/2)*time_step, O(i_loop)+(q1/2)*time_step);
k3 = func1( A2(i_loop)+(m2/2)*time_step, A1(i_loop)+(k2/2)*time_step, G1(i_loop)+(l2/2)*time_step, O(i_loop)+(q2/2)*time_step);
l3 = func2( G1(i_loop)+(l2*time_step/2), A1(i_loop)+(k2/2)*time_step);
m3 = func3( A2(i_loop)+(m2/2)*time_step, A1(i_loop)+(k2/2)*time_step, G2(i_loop)+(n2/2)*time_step, O(i_loop)+(q2/2)*time_step);
n3 = func4( G2(i_loop)+(n2/2)*time_step, A2(i_loop)+(m2/2)*time_step);
q3 = func5( A2(i_loop)+(m2/2)*time_step, A1(i_loop)+(k2/2)*time_step, O(i_loop)+(q2/2)*time_step );
k4 = func1( A2(i_loop)+(m3*time_step), A1(i_loop)+(k3*time_step), G1(i_loop)+(l3*time_step), O(i_loop)+(q3*time_step));
l4 = func2( G1(i_loop)+(l3*time_step), A1(i_loop)+(k3*time_step));
m4 = func3( A2(i_loop)+(m3*time_step), A1(i_loop)+(k3*time_step), G2(i_loop)+(n3*time_step), O(i_loop)+(q3*time_step));
n4 = func4( G2(i_loop)+(n3*time_step),A2(i_loop)+(m3*time_step));
q4 = func5( A2(i_loop)+(m3*time_step), A1(i_loop)+(k3*time_step), O(i_loop)+(q3*time_step));
A1(i_loop+1) = A1(i_loop) + (time_step/6)*(k1+(2*(k2+k3))+k4);
G1(i_loop+1) = G1(i_loop) + (time_step/6)*(l1+(2*(l2+l3))+l4);
A2(i_loop+1) = A2(i_loop) + (time_step/6)*(m1+(2*(m2+m3))+m4);
G2(i_loop+1) = G2(i_loop) + (time_step/6)*(n1+(2*(n2+n3))+n4);
O(i_loop+1) = O(i_loop) + (time_step/6)*(q1+(2*(q2+q3))+q4);
end
subplot 211
plot(time,A1)
grid on
subplot 212
plot(time,O)

更多回答（2 个）

KSSV
KSSV 2022-7-6
This function:
func2 = @(G1) (1/tf).*(P1 - G1.*(A1.^2+1));
should be:
func2 = @(G1,A1) (1/tf).*(P1 - G1.*(A1.^2+1));
Walter Roberson
Walter Roberson 2022-7-6
Ran in:
The posted code executes for me
format long
a = 0;
b = 30;
h = 0.1;
o = 3.5;
tc = 3.0;
tf = 2.3;
a1 = 0.1;
a2 = 0.1;
P1 = 0.2;
P2 = 0.2;
k= 0.9;
V = 1;
A1(1) = 10E-4;
A2(1) = 10E-5;
G1(1) = 10E-6;
G2(1) = 10E-6;
O(1) = 0;
n = (b-a)/h ;
t = ((a) + (0:n)*h).*10E-6;
func1 = @(G1,A1,A2,O) (1/tc).*(G1- a1).*(A1./V) +(k./tc).*(A2./V).*cos(O);
func3 = @(G2,A2,A1,O) (1/tc).*(G2- a2).*(A2./V) +(k./tc).*(A1./V).*cos(O);
func2 = @(G1) (1/tf).*(P1 - G1.*(A1.^2+1));
func4 = @(G2) (1/tf).*(P2 - G2.*(A2.^2+1));
func5 = @(A1,A2,O) o - (k./tc).*((A1./A2) + (A2./A1)).*sin(O);
for i = 1:n
k1 = feval(func1,G1(i),A1(i),A2(i),O(i));
l1 = feval(func2, G1(i));
m1 = feval(func3, G2(i),A1(i),A2(i),O(i));
n1 = feval(func4, G2(i));
q1 = feval(func5,A1(i),A2(i),O(i));
k2 = feval(func1, A2(i)+(m1/2)*h, A1(i)+(k1/2)*h, G1(i)+(l1/2)*h, O(i)+(q1/2)*h);
l2 = feval(func2, G1(i)+(l1*h/2));
m2 = feval(func3, A2(i)+(m1/2)*h, A1(i)+(k1/2)*h, G2(i)+(n1/2)*h, O(i)+(q1/2)*h);
n2 = feval(func4, G2(i)+(n1/2)*h);
q2 = feval(func5, A2(i)+(m1/2)*h, A1(i)+(k1/2)*h, O(i)+(q1/2)*h);
k3 = feval(func1, A2(i)+(m2/2)*h, A1(i)+(k2/2)*h, G1(i)+(l2/2)*h, O(i)+(q2/2)*h);
l3 = feval(func2, G1(i)+(l2*h/2));
m3 = feval(func3, A2(i)+(m2/2)*h, A1(i)+(k2/2)*h, G2(i)+(n2/2)*h, O(i)+(q2/2)*h);
n3 = feval(func4, G2(i)+(n2/2)*h);
q3 = feval(func5, A2(i)+(m2/2)*h, A1(i)+(k2/2)*h, O(i)+(q2/2)*h );
k4 = feval(func1, A2(i)+(m3*h), A1(i)+(k3*h), G1(i)+(l3*h), O(i)+(q3*h));
l4 = feval(func2, G1(i)+(l3*h));
m4 = feval(func3, A2(i)+(m3*h), A1(i)+(k3*h), G2(i)+(n3*h), O(i)+(q3*h));
n4 = feval(func4, G2(i)+(n3*h));
q4 = feval(func5, A2(i)+(m3*h), A1(i)+(k3*h), O(i)+(q3*h));
A1(i+1) = A1(i) + (h/6)*(k1+(2*(k2+k3))+k4);
G1(i+1) = G1(i) + (h/6)*(l1+(2*(l2+l3))+l4);
A2(i+1) = A2(i) + (h/6)*(m1+(2*(m2+m3))+m4);
G2(i+1) = G2(i) + (h/6)*(n1+(2*(n2+n3))+n4);
O(i+1) = O(i) + (h/6)*(q1+(2*(q2+q3))+q4);
end
subplot 211
plot(t,A1)
grid on
subplot 212
plot(t,O)

类别

在 Help Center 和 File Exchange 中查找有关 Particle & Nuclear Physics 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by