3 views (last 30 days)

Show older comments

syms x

eq4=x^3*(464250790297358007291860.3258421-10669886405219997721570.253037505i) + x^2*(4435793796620487051084389.1847765 + 338982820194409710227262107.98938i) - x*(82245532791487483177653457826.664 - 3925722719589374275368064048.4952i) - (434627852636126382139310330371.17 + 6697207434326522126273421132204.8i)

d=solve(eq4,x);

b=double(d)

Solutions are obtained but when I'm substituting the values back into the equation. It doesn't satisfies

David Goodmanson
on 7 Mar 2021

Hello Varun,

Your code gives for solutions

b = 30.857 - 175.74i

-94.032 - 211.56i

70.398 - 342.71i

Rather than substitute those values directly into eq4, the code below converts the polynomial coefficients to doubles and then uses polyval. One reason, but not the only one, is that who can tell the sizes of those coefficients just by looking at the opaque format of eq4? Certainly not me.

format short g

c = double([(464250790297358007291860.3258421-10669886405219997721570.253037505i);

+ (4435793796620487051084389.1847765 + 338982820194409710227262107.98938i);

- (82245532791487483177653457826.664 - 3925722719589374275368064048.4952i);

- (434627852636126382139310330371.17 + 6697207434326522126273421132204.8i)])

R = polyval(c,b)

c =

4.6425e+23 - 1.067e+22i

4.4358e+24 + 3.3898e+26i

-8.2246e+28 + 3.9257e+27i

-4.3463e+29 - 6.6972e+30i

R = % values at roots should be small

-4.2221e+14 + 0i

1.8296e+15 + 1.1259e+15i

1.2666e+15 - 3.3777e+15i

Denoting the value of the polynomial by P and the value at the roots by R, certainly the values of R are huge, nowhere near 0. But the question is, huge compared to what? The c coefficents are so large that they affect what is considered large and small.

If the coefficients of the polynomial are all divided by c(1), the polynomial takes a standard form with the coefficient of the highest power equal to 1 (same roots).

cnew = c/c(1)

Rnew = polyval(cnew,b)

cnew =

1 + 0i

-7.223 + 730.01i

-1.7726e+05 + 4382.1i

-6.0432e+05 - 1.444e+07i

Rnew = % values at root should be small

-2.3283e-10 - 3.7253e-09i

3.0268e-09 + 0i

-3.9581e-09 + 5.5879e-09i

This looks a lot better. However, although rescaling has a good psychological effect, it has not fundamentally changed anything. There still needs to be some kind of comparison to get a better idea of large vs small.

The following contour plot of abs(P) indicates the positions of the roots and the size of P in that region.

[xx yy] = meshgrid(linspace(-120,100,200),linspace(-400,-100,200));

P = polyval(cnew,xx+i*yy);

contour(xx,yy,abs(P),32)

colorbar

Using data tips in the area in between the three roots shows a maximum abs(P) of about P0 = 1e6. With this as the benchmark, you can see that the values of Rnew are actually very small,

abs(Rnew)/P0 ~~ 1e-15.

Walter Roberson
on 7 Mar 2021

Your final coefficient is only accurate to +/- 0.05 . However if you convert your coefficients as if they were exact fractions instead of being partly unknown, and you use the exact formula for solutions to cubic equations, and then take the given solution + 0.05 (the margin of error in your coefficients), the result will be on the order of 10^26. The range of your coefficients is high and small errors make a big difference in the calculations. Indeed if you take the exact solutions under the assumption that all of the provided coefficients are exact, then the error in the solution x cannot exceed 1e-29 in order for the error in back substitution to be less than 1.

It is a category error to use solve() with in-exact coefficients. solve() is for indefinitely precise calculations, but it is meaningless to get an indefinitely precise result when your coefficients are only a few decimal places accurate. Garbage In, Garbage Out. You should be using vpasolve() instead of solve(), and acknowledging that imprecision in the coefficients leads to imprecision in (finite-length) solutions.

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

Start Hunting!