How to find polynomial roots in Simulink?

Question

ozan eren 2020-5-8

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/523989-how-to-find-polynomial-roots-in-simulink

评论： ozan eren 2020-6-24

Hi,

I am trying to solve a 4th order polynomial equation in Simulink. I need to solve the equation by using Simulink blocks. The coefficients are calculated in Simulink blocks as well and I need to find the roots of this equations for each iteration.

To be more clear, I am using a for iterator block in Simulink and iterate 1000 times, in this for loop certain coefficients (C1,C2, C6 etc) are calculated from math expressions in each iteration. And I neet to find the roots of the polynomial whose coefficients are those calculated C values and I need find the roots for each iteration. I do not want to use script or built-in functions, but solve this problem with Simulink blocks if possible. Below I share a piece of m script which does exactly what I need to do in Simulink.

Any help will be appriciated, thanks in advance.

    r = roots([C5 (-C1+2*C4-C6) (3*C2-3*C5) (C1-2*C3+C6) -C2]);
    r = r(r==conj(r)); % discard complex-conjugate roots
    r = r(r>0); % discard negative roots

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

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

Answer 1

Walter Roberson 2020-6-22

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/523989-how-to-find-polynomial-roots-in-simulink#answer_455641

在 MATLAB Online 中打开

By using roots() on symbolic variables, you can get four closed form expressions for the roots. They occur in pairs, A+/-B and P+/-Q where B and Q are sqrt(), so by detecting whether the sqrt() involve imaginary quantities you can eliminate conjugate pairs as you wanted.

The closed form expressions are long, but if you are willing to run some code overnight, you can matlabFunction() writing to a 'File' with 'optimize' turned on. That produces a serquence of simple MATLAB code to calculate the expression.

You indicated that you do not want to use scripts, just Simulink blocks. Each of the statements from the optimized code is simple enough to make it practical to implement as a small number of Math blocks (or similar).

It would be a bit tedious to create these expressions manually, but it it only about 100 of them. And you would get exact solutions without iteration.

t2 = C2.*3.0;
t3 = C3.*2.0;
t4 = C4.*2.0;
t5 = C5.*3.0;
t9 = 1.0./C5;
t13 = sqrt(3.0);
t14 = sqrt(6.0);
t6 = -t3;
t7 = -t4;
t8 = -t5;
t10 = t9.^2;
t11 = t9.^3;
t15 = C2.*t9;
t12 = t10.^2;
t16 = C1+C6+t6;
t17 = C1+C6+t7;
t18 = t15.*1.2e+1;
t19 = t2+t8;
t20 = -t18;
t21 = t17.^2;
t22 = t17.^3;
t24 = t9.*t16;
t25 = t9.*t19;
t26 = (t9.*t17)./4.0;
t32 = t10.*t16.*t17.*3.0;
t35 = (t10.*t16.*t17)./4.0;
t37 = (t10.*t17.*t19)./2.0;
t23 = t21.^2;
t27 = t10.*t21.*(3.0./8.0);
t28 = (t11.*t22)./8.0;
t36 = -t35;
t38 = t11.*t19.*t21.*(3.0./4.0);
t39 = (t11.*t19.*t21)./1.6e+1;
t29 = -t27;
t30 = -t28;
t31 = t12.*t23.*(3.0./2.56e+2);
t33 = t12.*t23.*(9.0./6.4e+1);
t40 = -t39;
t34 = -t33;
t41 = t25+t29;
t47 = t24+t30+t37;
t53 = t15+t31+t36+t40;
t42 = t41.^2;
t43 = t41.^3;
t48 = t47.^2;
t54 = t53.^2;
t55 = t53.^3;
t58 = t41.*t53.*(4.0./3.0);
t59 = t41.*t53.*7.2e+1;
t44 = t42.^2;
t45 = t43.*2.0;
t46 = t43./2.7e+1;
t49 = t48.^2;
t50 = t48.*2.7e+1;
t52 = t48./2.0;
t56 = t55.*2.56e+2;
t57 = t43.*t48.*4.0;
t61 = t42.*t54.*1.28e+2;
t62 = t41.*t48.*t53.*1.44e+2;
t51 = t49.*2.7e+1;
t60 = t44.*t53.*1.6e+1;
t63 = t51+t56+t57+t60+t61+t62;
t64 = sqrt(t63);
t65 = t13.*t64.*3.0;
t66 = (t13.*t64)./1.8e+1;
t67 = t45+t50+t59+t65;
t69 = t46+t52+t58+t66;
t68 = sqrt(t67);
t70 = t69.^(1.0./3.0);
t72 = 1.0./t69.^(1.0./6.0);
t71 = t70.^2;
t74 = t41.*t70.*6.0;
t76 = t14.*t47.*t68.*3.0;
t73 = t71.*9.0;
t75 = -t74;
t77 = -t76;
t78 = t20+t32+t34+t38+t42+t73+t75;
t79 = sqrt(t78);
t80 = 1.0./t78.^(1.0./4.0);
t81 = t42.*t79;
t83 = t53.*t79.*1.2e+1;
t84 = t73.*t79;
t85 = t71.*t79.*-9.0;
t86 = (t72.*t79)./6.0;
t88 = t41.*t70.*t79.*1.2e+1;
t82 = -t81;
t87 = -t86;
t89 = -t88;
t90 = t76+t82+t83+t85+t89;
t91 = t77+t82+t83+t85+t89;
t92 = sqrt(t90);
t93 = sqrt(t91);
t94 = (t72.*t80.*t92)./6.0;
t95 = (t72.*t80.*t93)./6.0;
GG = [t26+t87-t94; t26+t87+t94; t26+t86-t95; t26+t86+t95];

4 个评论
显示 2更早的评论隐藏 2更早的评论

Walter Roberson 2020-6-23

在 MATLAB Online 中打开

Maple finds a more compact sequence:

t2 = C1-2*C4+C6;
t3 = 1/C5;
t5 = 1/4*t2*t3;
t6 = C2-C5;
t8 = t2^2;
t9 = C5^2; 
t10 = 1/t9;
t13 = 3*t6*t3-3/8*t8*t10;
t14 = t13^2;
t15 = C2*t3;
t17 = 3^(1/2);
t18 = t14*t13;
t20 = C1-2*C3+C6;
t24 = 1/t9/C5;
t30 = t20*t3-1/8*t8*t2*t24+3/2*t6*t2*t10;
t31 = t30^2;
t34 = t8^2;
t35 = t9^2;
t37 = t34/t35;
t40 = t20*t2*t10;
t43 = 3*t6*t8*t24;
t45 = t15+3/256*t37-1/4*t40-1/16*t43;
t46 = t45^2;
t49 = t31^2;
t53 = t14^2;
t60 = (144*t13*t31*t45+128*t14*t46+4*t18*t31+256*t45*t46+16*t45*t53+27*t49)^(1/2);
t61 = t17*t60;
t65 = t13*t45;
t67 = 1/18*t61+1/27*t18+1/2*t31+4/3*t65;
t68 = t67^(1/3);
t69 = t13*t68;
t71 = t68^2;
t76 = t14-12*t15-6*t69+9*t71-9/64*t37+3*t40+3/4*t43;
t77 = t76^(1/2);
t78 = t67^(1/6);
t79 = 1/t78;
t81 = 1/6*t77*t79;
t83 = 12*t45*t77;
t85 = 9*t71*t77;
t86 = t14*t77;
t88 = 12*t69*t77;
t89 = 6^(1/2);
t96 = (3*t61+2*t18+27*t31+72*t65)^(1/2);
t98 = 3*t89*t30*t96;
t100 = (t83-t85-t86-t88+t98)^(1/2);
t102 = t76^(1/4);
t103 = 1/t102;
t105 = 1/6*t100*t79*t103;
t109 = (t83-t85-t86-t88-t98)^(1/2);
t112 = 1/6*t109*t79*t103;
F(1) = t5-t81-t105;
F(2) = t5-t81+t105;
F(3) = t81+t5-t112;
F(4) = t81+t5+t112;

ozan eren 2020-6-24

Dear Walter,

I already tested and started to built the previous solution in Simulink, but this one is much more compact than earlier one as you mentioned and I will switch to this solution. I was stuck in this point almost a month, this is a huge help.

Once more thanks a lot for your effort and best regards.

请先登录，再进行评论。

Answer 2

Sai Teja Paidimarri 2020-6-18

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/523989-how-to-find-polynomial-roots-in-simulink#answer_453336

编辑：Walter Roberson 2020-6-22

在 MATLAB Online 中打开

Hi Ozan Eren,

You can find roots in Simulink and Matlab as well.

For a polynomial c5 x^4 + c4 x^3 +c3 x^2 +c2 x + c1, roots r can be found out using below steps

 x = c4 + r*c5;
 y = c3 + r*x;
 z = c2 + r*y; 

r value for which c1+r*z is zero is root of equation.

Above equations are from the concept of synthetic division method.

Counter, logical operator and constants can be used in case of Simulink implementation.

For or while loop can be used in case of Matlab Implementation.

3 个评论
显示 1更早的评论隐藏 1更早的评论

Sai Teja Paidimarri 2020-6-23

Hi ozan eren

First you need to find x,y,z because z depends on y , y depends upon x .Then condition associated with z.After that you can use r (using counter in simulink and any loop in matlab) to check the condition assocaited with z. if condition satisfies you can move r value in roots (declare roots as an array).

You can refer synthetic division method for understanding the above formulae.

ozan eren 2020-6-23

Dear Sai Teja,

Thanks a lot for your time and kind answers. But I think the above solution will be better for my application. In any case thank you for your time and best regards.

请先登录，再进行评论。