Unable to perform assignment because the indices on the left side are not compatible with the size of the right side.

n1=1.5; Vn=5.56; beta1=(table2array(beta))

for i=1:length(beta1)

	a=1; b(i)=2*Vn*cos(beta1(i)); c=((Vn)^2); d(i)=-(sin(beta1(i)))^(2*n1/(n1-2));

	p1(i)=[a b(i) c 0 0 0 0 0 d(i)];

	r1(i)=roots(p1(i));
	
	r1(i)=r1(r1(i)==real(r1(i))&r1(i)>0);
end

Dear all:

I am going to find the roots of polynomial. Since the one parameter (beta1) of the polynomial has different values (6843 value in total, array form), so I consider to use for loop to solve the polynomial. For expample, if beta 1=36 degree, I wish to find the roots of p1; then if the beta 1=50 degree, I wish to find the roots of p1 in this circumstance, the loop will be continued as so on. The beta1 is an array and the data is attached as the mat file.

I am not sure why it always shows the indices of LHS do not match with the indices at RHS. Could our friend help me fix the error? I spend all night but I still can't solve it.

By the way, I also wish to save only the positive real number answer.

Regards

[Steven Lord]

Why do you assume that there is one and only one positive real root? A polynomial with purely real coefficients can still have only purely imaginary roots.

syms x
s = expand((x-1i)*(x+1i)*(x-2i)*(x+2i))

s = 

By inspection the roots of s are . Each of those four values makes one of the terms inside the expand call equal to 0. Let's check numerically.

p = sym2poly(s)

p = 1x5

1 0 5 0 4

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

r = roots(p)

r = 
   0.0000 + 2.0000i
   0.0000 - 2.0000i
   0.0000 + 1.0000i
   0.0000 - 1.0000i

The roots function confirms the result we received by inspection. [Okay, to be fair there's some very small (on the order of eps) non-zero real part of some of the elements of r. That's noise.]

So what would you expect to be stored if the polynomial whose roots you're computing turned out to be something like p in the code above? None of the elements in r satisfy even the first of the conditions you're using.

r(r == real(r))

ans =

  0x1 empty double column vector

[Kaiyuan]

Thanks for the detailed explaination with patience, Steven.

What you have explained is also corrected as well for general cases. I am quiet lucky in my case that I have got the results which meet my expectation.

Answer 1

[Torsten]

Only a single number can be saved in an array element.

But you try to save a vector ([a b(i) c 0 0 0 0 0 d(i)]) in a single element (p1(i)) in this line of your code:

p1(i)=[a b(i) c 0 0 0 0 0 d(i)];

Use

n1=1.5; Vn=5.56; 
u = load("beta.mat");
beta1 = table2array(u.RouteWindData);
size(beta1)

ans = 1x2

6843 1

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

beta1 = beta1(~isnan(beta1));
size(beta1)

ans = 1x2

5817 1

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

for i=1:length(beta1)

	a=1; b=2*Vn*cos(beta1(i)); c=((Vn)^2); d=-(sin(beta1(i)))^(2*n1/(n1-2));

	r=roots([a b c 0 0 0 0 0 d]);
	
	r1{i}=r(abs(imag(r))<1e-6 & real(r)>0);
end
% Fortunately, each cell contains only one root such that r1 
% can be converted to a usual array
r1 = cell2mat(r1)  

r1 = 1x5817

30.5607 9.9333 2.1891 1.1534 0.6205 0.5661 0.9622 0.5850 1.1534 9.9333 30.5607 30.5607 30.5607 30.5607 30.5607 30.5607 30.5607 30.5607 30.5607 30.5607 5.9370 0.8463 0.7423 0.5668 0.8566 8.4736 2.5872 9.0579 0.6520 1.9065

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

[Kaiyuan]

Thanks for the very professional and detailed explaination, Torsten.

I have learnt from it and the results are exactly what I expect!