98 views (last 30 days)

I have the equation:

eqn = ((d^3*pi*((4251*d + 5951400)/(25*d))^(1/2))/(2*(d + 1400)))*(pi*(d^2)/4)-180 == 0

and want to find the real root/roots.

when i attempt to solve it with

x = vpasolve(eqn,d)

Matlab only return the imaginary part of the solution:

- 3.593452147062167996782934136112 - 1.3074862631155484137468498544529i

How do i find the real solution?

Roger Stafford
on 20 Dec 2016

This is an equation that can be manipulated so that d is one of the roots of a ninth degree polynomial equation of which only one of its nine roots is real.

The original equation is:

((d^3*pi*((4251*d+5951400)/(25*d))^(1/2))/(2*(d+1400)))*(pi*(d^2)/4)==180

Since 4251*d+5951400 = 4251*(d+1400), the d+1400 partially cancels with the same quantity in the denominator and we get the equation

pi^2/40*d^5*(4251/(d*(d+1400)))^(1/2)==180

or

pi^2/40*d^5*4251^(1/2) == 180*(d*(d+1400))^(1/2)

Squaring both sides and transposing gives

4251/1600*pi^4*d^10-32400*d^2-45360000*d == 0

One factor d can be factored out since d = 0 is clearly not a solution of the original equation and that finally leaves the polynomial equation

4251/1600*pi^4*d^9-32400*d-45360000 == 0

The ‘roots’ function can be used for this and it shows that there is only one real root, namely

d = 3.82617917662798

Massimo Zanetti
on 19 Dec 2016

Edited: Massimo Zanetti
on 19 Dec 2016

To force the vpasolve finding only real solutions you cannot use assume. If you know the real solution is only one a workaround is to set search range to [-Inf,Inf]:

x = vpasolve(eqn,d,[-Inf,Inf])

x =

3.8261791766279723703687735908663

Massimo Zanetti
on 20 Dec 2016

Walter Roberson
on 19 Dec 2016

Walter Roberson
on 20 Dec 2016

sols = solve(eqn, x);

sols(imag(sols)~=0) = []; %remove the ones that have a non-zero imaginary component.

However, in R2016b (and perhaps some other releases), this code can fail due to a bug in the Symbolic Toolbox (I notified them of the bug a few weeks ago.) The work-around for the moment is

sols(imag(vpa(sols))~=0) = [];

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

Start Hunting!
## 0 Comments

Sign in to comment.