Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I have a polynomial in four variables x,a,b and c. The number of roots of the polynomial in x depends of the choice of a, b and c. I would like to have a 3D-Plot with a, b and c on the axes, while the number of roots >0 at a point (a,b,c) is symbolized by different colours.

The most important to me is to see where exactly the transitions are.

I think I need CountRoots[Polynom,{x,0,Infinity}], but as I am new to Mathematica I can't figure out how to do this. Thanks a lot for your help!

The function is

-a b (1 + 
    c) (a^8 (1 + c)^4 (6 - 56 c + 68 c^2 - 56 c^3 + 6 c^4 + 
       3 b (-1 + c)^4 (1 + c)) + 
    a^7 b (1 + c)^4 (9 b (-1 + c)^4 (1 + c) - 256 c (1 - c + c^2)) + 
    3 a^6 b^3 (-1 + c)^4 (1 + c)^5 - 
    8 a^6 b^2 (1 + c)^4 (15 + 52 c - 22 c^2 + 52 c^3 + 15 c^4) - 
    15 a^5 b^4 (-1 + c)^4 (1 + c)^5 - 
    128 a^5 b^3 (1 + c)^4 (3 + 2 c + 4 c^2 + 2 c^3 + 3 c^4) + 
    3 a^3 b^6 (1 + c) (-1 + c^2)^4 - 
    128 a^3 b^5 (1 + c)^4 (3 + 2 c + 4 c^2 + 2 c^3 + 3 c^4) - 
    15 a^4 b^5 (-1 + c)^4 (1 + c)^5 - 
    20 a^4 b^4 (1 + c)^4 (27 + 4 c + 50 c^2 + 4 c^3 + 27 c^4) + 
    3 a b^8 (-1 + c)^4 (1 + c)^5 - 
    256 a b^7 c (1 + c)^4 (1 - c + c^2) + 
    9 a^2 b^7 (-1 + c)^4 (1 + c)^5 - 
    8 a^2 b^6 (1 + c)^4 (15 + 52 c - 22 c^2 + 52 c^3 + 15 c^4) + 
    2 b^8 (1 + c)^4 (3 - 28 c + 34 c^2 - 28 c^3 + 3 c^4)) - 
 a b (1 + c) (21 a^6 b (-1 + c)^4 (1 + c)^3 - 
    30 a^6 (1 + c)^2 (5 - 12 c + 30 c^2 - 12 c^3 + 5 c^4) + 
    21 a^5 b^2 (-1 + c)^4 (1 + c)^3 - 
    60 a^5 b (1 + c)^2 (7 - 4 c + 42 c^2 - 4 c^3 + 7 c^4) - 
    42 a^3 b^4 (-1 + c)^4 (1 + c)^3 - 
    120 a^3 b^3 (1 + c) (1 + 37 c + 42 c^2 + 42 c^3 + 37 c^4 + c^5) - 
    42 a^4 b^3 (-1 + c)^4 (1 + c)^3 - 
    30 a^4 b^2 (1 + c)^2 (11 + 76 c + 66 c^2 + 76 c^3 + 11 c^4) + 
    21 a b^6 (-1 + c)^4 (1 + c)^3 - 
    60 a b^5 (1 + c)^2 (7 - 4 c + 42 c^2 - 4 c^3 + 7 c^4) + 
    21 a^2 b^5 (-1 + c)^4 (1 + c)^3 - 
    30 a^2 b^4 (1 + c)^2 (11 + 76 c + 66 c^2 + 76 c^3 + 11 c^4) - 
    30 b^6 (1 + c)^2 (5 - 12 c + 30 c^2 - 12 c^3 + 5 c^4)) x - 
 a b (1 + c) (54 a^3 b^2 (-1 + c)^4 (1 + c) - 
    168 a^3 b (1 + c) (1 + 11 c + 11 c^2 + c^3) - 
    54 a^4 b (-1 + c)^4 (1 + c) - 
    6 a^4 (29 - 4 c + 286 c^2 - 4 c^3 + 29 c^4) - 
    54 a b^4 (-1 + c)^4 (1 + c) - 
    168 a b^3 (1 + c)^2 (1 + 10 c + c^2) + 
    54 a^2 b^3 (-1 + c)^4 (1 + c) + 
    12 a^2 b^2 (1 - 340 c - 330 c^2 - 340 c^3 + c^4) - 
    6 b^4 (29 - 4 c + 286 c^2 - 4 c^3 + 29 c^4)) x^2 - 
 a b (1 + c) (-20 a b (9 + 50 c + 9 c^2) - 
    10 a^2 (9 + 50 c + 9 c^2) - 10 b^2 (9 + 50 c + 9 c^2)) x^3 + 
 72 a b (1 + c) x^4
share|improve this question
    
Can you give us the polynomial explicitly..? –  PlatoManiac Sep 16 '12 at 23:50
    
I did not at first because it is really ugly. –  Apatura Sep 16 '12 at 23:55
    
"Number of roots" as in "number of real roots"? A cubic always has three roots, for instance... –  J. M. Sep 17 '12 at 0:02
    
Number of real roots >0. –  Apatura Sep 17 '12 at 0:04
add comment

2 Answers

up vote 10 down vote accepted

Defining your polynomial as

p[x_, a_, b_, c_] := the formula

I'd do this :

Manipulate[ CountRoots[ p[x, a, b, c], {x, 0, Infinity}],
            {a, -100, 100}, {b, -100, 100}, {c, -100, 100}]

enter image description here

If a plot is needed one can proceed this way :

ListContourPlot3D[ Table[ CountRoots[ p[x, a, b, c], {x, 0, Infinity}],
                          {a, -10, 10}, {b, -10, 10}, {c, -10, 10} ], Contours -> 3,
                   Mesh -> None, ContourStyle -> {Red, Yellow, Lighter @ Blue},              
                   DataRange -> {{-10, 10}, {-10, 10}, {-10, 10}}]

enter image description here

Another way (the most expensive) is to make use of RegionPlot3D. Here are 3d regions where a given polynomial has respectively : at least 1, 2 and 3 roots for x > 0

GraphicsRow[ 
    Table[ RegionPlot3D[ CountRoots[ p[x, a, b, c], {x, 0, Infinity}] >= k,
                         {a, -10, 10}, {b, -10, 10}, {c, -10, 10}, 
                         PlotStyle -> Directive[ Orange, Opacity[0.5], Specularity[White, 30]]],
           {k,  3}]]

enter image description here

One can see that for for higher k we need a better resolution, nevertheless it appears very expensive to compute regions using higher PlotPoints and MaxRecursion options.

 RegionPlot3D[ CountRoots[ p[x, a, b, c], {x, 0, Infinity}] >= 2, 
               {a, -10, 10}, {b, -10, 10}, {c, -10, 10}, 
               PlotStyle -> Directive[ Orange, Opacity[0.5], Specularity[White, 30]],
               PlotPoints -> 40, MaxRecursion -> 2]

enter image description here

share|improve this answer
    
Thanks a lot, that is great for me to get to know the function! But I really need a plot I can print on paper. –  Apatura Sep 17 '12 at 0:29
    
@Apatura I updated my answer with a ListContourPlot3D. –  Artes Sep 17 '12 at 0:45
    
Cute plot work! –  belisarius Sep 17 '12 at 4:53
    
@Artes: RegionPlot does exactly what I want, thank you very much! Unfortunately, when I try a different function, it gives back a strange error. What is wrong with RegionPlot3D[ CountRoots[\[Sigma]q33[\[Omega], V1, V2, V3], {\[Omega], 0, Infinity}] == 1, {V1, 0, 10}, {V2, 0, 10}, {V3, 0, 10}, PlotStyle -> Directive[Gray, Opacity[0.5]], Mesh -> None] if [Sigma]q33 is a polynomial in omega, Exp[-V1],Exp[-V2] and Exp[-V3]? –  Apatura Sep 17 '12 at 14:52
    
@Apatura You are welcome. Note that when I set PlotPoints -> 60, MaxRecursion -> 3 it took more than 5 minutes to evaluate. If you have an error most likely it is because of a wrong definition of your function \[Sigma]q33. Look carefully at the documentation of RegionPlot3D. –  Artes Sep 17 '12 at 14:57
show 1 more comment

Here's one approach. Generate data (for a simple polynomial):

data = Table[CountRoots[a*x^2 + b*x + c, x],
  {a, -5, 5, 0.5}, {b, -5, 5, 0.5}, {c, -5, 5, 0.5}];

Then we'll display a collection of Cuboids with Opacity based on this data.

Graphics3D[{Opacity[0.2], EdgeForm[],
  MapIndexed[Which[
    # == 0, {Red, Cuboid[#2, #2 + {1, 1, 1}]},
    # == 1, {White, Cuboid[#2, #2 + {1, 1, 1}]},
    # == 2, {Blue, Cuboid[#2, #2 + {1, 1, 1}]}] &,
  data, {3}]}]

enter image description here

Now it appears to turn out, that your more complicated functions all have zero, one, or two roots. Changing the With statement accordingly, we get the following.

enter image description here

You might also look into CUDAVolumetricRender, if you have a good NVidia graphics card.

share|improve this answer
    
Might also be a good idea to generate a separate plot for the different root counts... –  J. M. Sep 17 '12 at 0:15
    
@J.M. I'm not sure I follow - that image is based on root counts. I turns out that all polynomials in that family have on zero, two, or four roots - hence, the three colors. –  Mark McClure Sep 17 '12 at 0:19
    
Hmm, I thought you were doing something like ContourPlot3D[(* root-counting function *), {a, ...}, {b, ...}, {c, ...}, Contours -> {0, 2, 4}]... –  J. M. Sep 17 '12 at 0:22
    
I'd have used Switch[] instead of Which[] myself... something like {Switch[#, 0, Red, 1, White, 2, Blue], Cuboid[#2, #2 + {1, 1, 1}]} &. –  J. M. Sep 17 '12 at 0:39
    
@MarkMcClure will you please explain a bit how to use CUDAVolumetricRender in this context? –  PlatoManiac Sep 17 '12 at 8:38
show 2 more comments

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.