How to tell Matlab to give real cube roots instead of complex ones?

I am using symbolic math toolbox (Matlab2011b) to do some geometric computations to derive an expression for a geometric length. At some point of my calculations some cube roots are introduced and unfortunately Matlab seems to automatically simplify terms like (-1)^(1/3) to yield 0.5000 + 0.8660i.

I realize this is a correct solution, however i want my result to be real - is there a way to tell Matlab to always return the real cube root?

These imaginary ones are messing up my expressions...

(Also i wonder why matlab chooses to consider one of two complex roots instead of giving all three solutions or the real one...)

Answer 1

[Eric]

This is almost certainly occurring because MATLAB is doing the math by first converting to log space so that it can just multiply the exponent by the log of the argument and then convert it back. This is similar to how computers and calculators commonly do complicated arithmetic like multiplication (which becomes addition in log space) by first converting to log space, because adding logs is typically much easier, faster, and usually the only option available to the computer in the first place. The following is likely how MATLAB/the computer goes through the calculation:

 (-1)^(1/3) = exp(exp( log(log(-1)) - log(3) ))
 log(-1) = pi*i
 log(pi*i) = 1.1447 + 1.5708i
 log(3) = 1.0986
 log(log(-1)) - log(3) = (1.1447 + 1.5708i) - (1.0986) = 0.0461 + 1.5708i
 exp(0.0461 + 1.5708i) = exp(0.0461)*exp(1.5708i) = 1.0472i
 exp(1.0472i) = 0.5 + 0.866i % A complex cube root of -1, along with its conjugate

Since the log of a negative number (namely -1 in this case) is going to be complex, this is what (almost certainly) makes the output of the operation also complex and why there won't be a setting to change it to default to real. Having said that, it is fairly easy to get the real answer of a cube root by doing something like

 sign(x).*abs(x).^(1/3)

assuming your inputs are real. This takes advantage of what we already know to be true, namely that

 (-x)^(1/3) = -(x^(1/3)) % For x>=0

and should, in theory, work for Ruye's issue, though I am no expert on symbolic forms and am almost certainly using a later MATLAB release (R2016a).

At this point, it should be worth mentioning that MATLAB already has a function called

 nthroot(X,N)

which works for any root and will return a real result or bust. If you look at the code for nthroot(), you will see that it is essentially doing the above sign*abs method, but with fancy checking to see that the root will be real first.

(I know this is an old question, but I wanted to compile previous answers and give a plausible explanation for why it might give a complex number)

Answer 2

[Friedrich]

Hi,

seems more like MATLAB is converting your symbolic expression into a complex double value. So what are you doing with the symbolic expression? Or better, what is your actual code?

[Thomas]

Thanks for the response. Actually my example (-1)^(1/3) was only a minimal example, during my computations i get very long expressions with cube roots introduced as factors or summands at varius places. I managed to manually re-format one solution by eliminating a few (-1)^(1/3) factors so that Matlab does not replace these with complex numbers ruining the overall result. For example one expression yields:

    (1/(az*dx*dy*(2*ax*nx + dx*nx + 2*ay*ny + dy*ny + 2*az*nz)*sqrt(nx^2 + ny^2 + nz^2)))*(-2*alpha*ax^2*az*dx*dy*nx^2 - alpha*ax*az*dx^2*dy*nx^2 - 4*alpha*ax*ay*az*dx*dy*nx*ny - alpha*ay*az*dx^2*dy*nx*ny - alpha*ax*az*dx*dy^2*nx*ny - 2*alpha*ay^2*az*dx*dy*ny^2 - alpha*ay*az*dx*dy^2*ny^2 -    4*alpha*ax*az^2*dx*dy*nx*nz - alpha*az^2*dx^2*dy*nx*nz - 4*alpha*ay*az^2*dx*dy*ny*nz - alpha*az^2*dx*dy^2*ny*nz - 2*alpha*az^3*dx*dy*nz^2 - (-1)^(1/3)*az^(2/3)*dx^(2/3)*dy^(2/3)*(ax*nx + ay*ny + az*nz)^(1/3)*(2*ax*nx + dx*nx + 2*ay*ny + dy*ny + 2*az*nz)^(2/3)*(nx^2 + ny^2 + nz^2)^(1/3)*    ((2*dx*dy*ho*(3*az^2 - 3*az*ho + ho^2)*(ax*nx + dx*nx + ay*ny + az*nz)*(ax*nx + ay*ny + dy*ny + az*nz)*(ax*nx + dx*nx + ay*ny + dy*ny + az*nz))/az^2 + alpha^3*az*dx*dy*(ax*nx + ay*ny + az*nz)^2*(2*ax*nx + dx*nx + 2*ay*ny + dy*ny + 2*az*nz))^(1/3))

All variables are defined as real syms, e.g.

syms ax ay az nx ny nz alpha dx dy ho real;

The result of this expression should be a length, so no not real number. When copy-pasting the above expression matlab instantly reformats it and adds imaginary terms... (except when i remove the (-1)^(1/3), however i'd like matlab to autmoatically select real simplifications instead of complex ones).

[Friedrich]

which MATLAB version are you using? For me in 13a MATLAB doesnt add imaginary terms:

>> syms ax ay az nx ny nz alpha dx dy ho real;
>>  (1/(az*dx*dy*(2*ax*nx + dx*nx + 2*ay*ny + dy*ny + 2*az*nz)*sqrt(nx^2 + ny^2 + nz^2)))*(-2*alpha*ax^2*az*dx*dy*nx^2 - alpha*ax*az*dx^2*dy*nx^2 - 4*alpha*ax*ay*az*dx*dy*nx*ny - alpha*ay*az*dx^2*dy*nx*ny - alpha*ax*az*dx*dy^2*nx*ny - 2*alpha*ay^2*az*dx*dy*ny^2 - alpha*ay*az*dx*dy^2*ny^2 -    4*alpha*ax*az^2*dx*dy*nx*nz - alpha*az^2*dx^2*dy*nx*nz - 4*alpha*ay*az^2*dx*dy*ny*nz - alpha*az^2*dx*dy^2*ny*nz - 2*alpha*az^3*dx*dy*nz^2 - (-1)*az^(2/3)*dx^(2/3)*dy^(2/3)*(ax*nx + ay*ny + az*nz)^(1/3)*(2*ax*nx + dx*nx + 2*ay*ny + dy*ny + 2*az*nz)^(2/3)*(nx^2 + ny^2 + nz^2)^(1/3)*    ((2*dx*dy*ho*(3*az^2 - 3*az*ho + ho^2)*(ax*nx + dx*nx + ay*ny + az*nz)*(ax*nx + ay*ny + dy*ny + az*nz)*(ax*nx + dx*nx + ay*ny + dy*ny + az*nz))/az^2 + alpha^3*az*dx*dy*(ax*nx + ay*ny + az*nz)^2*(2*ax*nx + dx*nx + 2*ay*ny + dy*ny + 2*az*nz))^(1/3))

ans =

-(2*alpha*az^3*dx*dy*nz^2 + alpha*ax*az*dx^2*dy*nx^2 + 2*alpha*ax^2*az*dx*dy*nx^2 + alpha*ay*az*dx*dy^2*ny^2 + 2*alpha*ay^2*az*dx*dy*ny^2 + alpha*az^2*dx^2*dy*nx*nz + alpha*az^2*dx*dy^2*ny*nz + alpha*ax*az*dx*dy^2*nx*ny + 4*alpha*ax*az^2*dx*dy*nx*nz + alpha*ay*az*dx^2*dy*nx*ny + 4*alpha*ay*az^2*dx*dy*ny*nz - az^(2/3)*dx^(2/3)*dy^(2/3)*(alpha^3*az*dx*dy*(ax*nx + ay*ny + az*nz)^2*(2*ax*nx + 2*ay*ny + 2*az*nz + dx*nx + dy*ny) + (2*dx*dy*ho*(3*az^2 - 3*az*ho + ho^2)*(ax*nx + ay*ny + az*nz + dx*nx)*(ax*nx + ay*ny + az*nz + dy*ny)*(ax*nx + ay*ny + az*nz + dx*nx + dy*ny))/az^2)^(1/3)*(nx^2 + ny^2 + nz^2)^(1/3)*(ax*nx + ay*ny + az*nz)^(1/3)*(2*ax*nx + 2*ay*ny + 2*az*nz + dx*nx + dy*ny)^(2/3) + 4*alpha*ax*ay*az*dx*dy*nx*ny)/(az*dx*dy*(nx^2 + ny^2 + nz^2)^(1/2)*(2*ax*nx + 2*ay*ny + 2*az*nz + dx*nx + dy*ny))

[Thomas]

Oh, you're right. I accidentally pasted my manually corrected version. I've edited my post above and pasted the uncorrected version which basically has a ^(1/3) after the (-1)... would be nice to know if Matlab 2013 deals better with this, if you could try once again.

[Friedrich]

When using the corrected code I get an "i" part in the result:

-(2*alpha*az^3*dx*dy*nz^2 + alpha*ax*az*dx^2*dy*nx^2 + 2*alpha*ax^2*az*dx*dy*nx^2 + alpha*ay*az*dx*dy^2*ny^2 + 2*alpha*ay^2*az*dx*dy*ny^2 + alpha*az^2*dx^2*dy*nx*nz + alpha*az^2*dx*dy^2*ny*nz + alpha*ax*az*dx*dy^2*nx*ny + 4*alpha*ax*az^2*dx*dy*nx*nz + alpha*ay*az*dx^2*dy*nx*ny + 4*alpha*ay*az^2*dx*dy*ny*nz + 4*alpha*ax*ay*az*dx*dy*nx*ny + az^(2/3)*dx^(2/3)*dy^(2/3)*((3^(1/2)*i)/2 + 1/2)*(alpha^3*az*dx*dy*(ax*nx + ay*ny + az*nz)^2*(2*ax*nx + 2*ay*ny + 2*az*nz + dx*nx + dy*ny) + (2*dx*dy*ho*(3*az^2 - 3*az*ho + ho^2)*(ax*nx + ay*ny + az*nz + dx*nx)*(ax*nx + ay*ny + az*nz + dy*ny)*(ax*nx + ay*ny + az*nz + dx*nx + dy*ny))/az^2)^(1/3)*(nx^2 + ny^2 + nz^2)^(1/3)*(ax*nx + ay*ny + az*nz)^(1/3)*(2*ax*nx + 2*ay*ny + 2*az*nz + dx*nx + dy*ny)^(2/3))/(az*dx*dy*(nx^2 + ny^2 + nz^2)^(1/2)*(2*ax*nx + 2*ay*ny + 2*az*nz + dx*nx + dy*ny))

[Thomas]

Too bad. Thanks for testing anyway.

[Sriman]

How to find in cubic root

Answer 3

[Carlos]

 Try this sign(x).*abs(x.^(1/3))

>> sign(-1).*abs(-1.^(1/3))

ans = -1

[Thomas]

Yeah, i thought about this too, but i cannot do it manually for all terms in my computations (at least i don't want to to save time in the long term)... but thanks for the idea anyway ;)

Answer 4

[Ruye Wang]

I have the same question. I need to evaluate a function in symbolic form, essentially like the following:

syms x;

f=x^(1/3);

subs(f,-1)

The output is one of the three roots: 0.5000 + 0.8660i

But what I need is the real cube root of x=-1 (which is -1).

I know I could get the real cube root if I explicitly hard code it like this:

nthroot(-1,3)

However, how do I do this in the form of symbolic function, which is needed in the program?

[Thomas]

I found no direct answer to the problem. I switched to use Mathematica, which is better for symbolic computations anyway (faster, easier input and better simplification routines imho). There i can bracket expressions with "Assuming[...]" to tell Mathematica that variables or subexpressions are real.

It think as og 2012 Matlab has an "assume" as well, however i am not sure if "assume(imag(x)=0)" or the like would work, as i still have Matlab 2011... maybe you could look at this.

Answer 5

[Steven Lord]

See this Answer for an explanation of why something like (-8)^(1/3) returns a complex result and an alternative function you can use.

[Jason Duvall]

This link is no lonber available, but the issue persists. Does anyone have an update?

[Walter Roberson]

MATLAB defines the .^ operator on numeric values in terms of log and exp, and that will give you complex results for negative base numbers except for integer powers.