Issue plotting function that is complex under certain domain

edit

I have been able to recreate my issue as the following simple problem.

Say I want to take a simple cubic expression like the following and solve it for x instead of f. This will no longer be a function since there are multiple answers for the same value of f (namely 3 answers). But I could simply take the 3 expressions and overlay their graphs on the same plot to get the graph I want.

var('x,y')
solution = solve( y==x\*(x-1)\*(x+1), x )
show( solution[0].rhs() )
plot( solution[0].rhs(), (y,0,1.5) )

I have been able to get the results I'd expect using Wolfram Alpha.

Now I believe that source of my issue is coming from the fact that over a certain domain the result is complex. I have been playing around with trying to learn about symbolic rings, but I haven't been able to make heads or tails out of how to use them properly, since every time I do I get a new error. Ideally I'd like to just be able to be able to plot only where the function is only real.