Re: Drawing polytopes
- To: mathgroup at smc.vnet.net
- Subject: [mg24270] Re: [mg24236] Drawing polytopes
- From: Andrzej Kozlowski <andrzej at tuins.ac.jp>
- Date: Wed, 5 Jul 2000 23:10:53 -0400 (EDT)
- Sender: owner-wri-mathgroup at wolfram.com
First a few corrections to some silly errors which have (as usual) crept in
into my earlier posting. What I carelessly referred to as "lines" meant of
course "line segments", while "2 dimensional linear subspaces" meant (3
dimensional) polytopes bounded by such.
The other point I wanted to make that the procedure below for obtaining the
"polytope" could quite easily be turned into a package (rather like the
existing 2-dimensional FilledPlot package.). The only somewhat tricky point
is to make sure that one carefully takes account of all the possible outputs
that LogicalExpand[InequalitySolve[]] (or alternatively,
LogicalExpand[Experimental`CylindricalAlgebraicDecomposition[]]) may return.
I have no time at present to try to do it tmyslef but it seems to me like a
worthwhile idea for a MathSource package.
on 7/5/00 4:29 PM, Andrzej Kozlowski at andrzej at tuins.ac.jp wrote:
> on 7/4/00 9:39 AM, Matthijs Sypkens Smit at matthijs at helena.tux.nu wrote:
>
>> I'm wondering if it's possible to draw (convex) polytopes defined by
>> systems of linear inequalities. For example the system:
>>
>> x + 2y + z <= 3
>> 2x - y >= 0
>> -x + 2y + z >= 0
>> y + z >= 0
>> 2x + 3y + 3z >= 1
>>
>> Can I get Mathematica to draw the polytope defined on R^3 by this system?
>>
>> I'm able to solve this system with
>> InequalitySolve[{x+2y+2z<=3,2x-y>=0,-x+2y+z>=0,y+z>=0,2x+3y+3z>=1},{x,y,z}]
>> but the data is hard to interpret this way.
>>
>> btw; the same problem, but a little more intuitive, would be to draw a
>> cube with the system defined below:
>>
>> {x >= -1, x <= 1, y >= -1, y <= 1, z >= -1, z <= 1}
>>
>> Any hints are appreciated,
>>
>
> There is no built in way to do this. I can think of two reasonable approaches,
> in one the code is simple but it is computation intensive and unlikely to
> produce an enlightening result, and the second in which the code is tedious to
> write out but which gives accurate results and does it fast. Let m estart with
> the second, (which is actually pretty obvious). Let me start with the second.
>
> First we load
> In[1]:=
> << Algebra`InequalitySolve`
>
> We define the inequalitis
>
> In[2]:=
> ineqs = {x + 2y + z <= 3,
> 2x - y >= 0, -x + 2y + z >= 0,
> y + z >= 0,
> 2x + 3y + 3z >= 1};
>
> and solve them
> In[3]:=
> ineqs2 = LogicalExpand[InequalitySolve[ineqs, {x, y, z}]]
> Out[3]=
> 3 3 3
> x == - && z == - - 2 y && y <= - ||
> 2 2 2
>
> 3
> y == 3 - x && z == 3 - x - 2 y && x < - && 1 <= x ||
> 2
>
> x <= -1 && x - 2 y <= z && y <= 2 x && z <= 3 - x - 2 y ||
>
> 1 1
> -1 < x && x < - && x - 2 y <= z && y <= - (-1 + 5 x) &&
> 2 3
>
> 1
> z <= 3 - x - 2 y || x < 1 && y < x && - <= x &&
> 2
>
> x - 2 y <= z && z <= 3 - x - 2 y ||
>
> 3
> x < - && y < x && 1 <= x && x - 2 y <= z &&
> 2
>
> 1
> z <= 3 - x - 2 y || -1 < x && x < - &&
> 2
>
> 1 1
> - (-1 + 5 x) < y && - (1 - 2 x - 3 y) <= z && y <= 2 x &&
> 3 3
>
> 1
> z <= 3 - x - 2 y || x < 1 && - <= x && x <= y &&
> 2
>
> -y <= z && y <= 2 x && z <= 3 - x - 2 y ||
>
> 3
> x < - && y < 3 - x && 1 <= x && x <= y && -y <= z &&
> 2
>
> z <= 3 - x - 2 y
>
> (Applying LogicalExpand makes the output clearer).
>
> Next the tedious part. First we shall define some global variables, that we
> shall use in the plot:
>
> In[4]:=
> Minx = -10; Miny = -10; Minz = -10; rx = 0.1; ry = 0.1;
>
> Next, we take each inequality and simply create the corresponding plot (using
> our parameters Minx,Miny and Minz whenever the lower bounds are not determined
> by the inequalities). We notice that the first two inequalities determine
> straight lines which can be ploted using ParametricPlot. The other
> inequalities determine 2 dimensional linear subspaces.
>
> In[7]:=
> plot[1] =
> ParametricPlot3D[{3/2, y, 3/2 - 2y}, {y, Miny, 3/2},
> DisplayFunction -> Identity];
>
> In[8]:=
> plot[2] =
> ParametricPlot3D[{x, 3 - x, 3 - x - 2}, {x, 1, 3/2},
> DisplayFunction -> Identity];
> In[9]:=
> plot[3] =
> Table[Line[{{x, y, x - 2y}, {x, y, 3 - x - 2y}}], {x, Minx, -1, rx}, {y,
> Miny, 2x, ry}];
> In[10]:=
> plot[4] =
> Table[Line[{{x, y, x - 2y}, {x, y, 3 - x - 2y}}], {x, -1, 0.5, rx}, {y,
> Miny, (5x - 1)/3, ry}];
> In[11]:=
> plot[5] =
> Table[Line[{{x, y, x - 2y}, {x, y, 3 - x - 2y}}], {x, 0.5, 1, rx}, {y,
> Miny, x, ry}];
> In[12]:=
> plot[6] =
> Table[Line[{{x, y, x - 2y}, {x, y, 3 - x - 2y}}], {x, 1, 3/2, rx}, {y,
> Miny, x, ry}];
> In[13]:=
> plot[7] =
> Table[Line[{{x, y, (1 - 2 x - 3 y)/3}, {x, y, 3 - x - 2 y}}], {x, -1, 0.5,
> rx}, {y, (5 x - 1)/3, 2 x, ry}];
> In[14]:=
> plot[8] =
> Table[Line[{{x, y, -y}, {x, y, 3 - x - 2 y}}], {x, 0.5, 1, rx}, {y, x, 2x,
> ry}];
>
> In[15]:=
> plot[9] =
> Table[Line[{{x, y, -y}, {x, y, 3 - x - 2 y}}], {x, 1, 3/2, rx}, {y, x,
> 3 - x, ry}];
>
> The polytope (which is really of course unbounded) can now be seen with:
>
> In[16]:=
> polytope = Show[Graphics3D[Table[{RGBColor[i/10, 0, 0], plot[i]}, {i, 3, 9}]]]
>
> (It is easier to see it properly if one loads in real-time 3d grpahics using
> <<RealTime3D`).
>
> The one dimensional pieces can be obtained in the same way:
>
> In[18]:=
> lines = Show[plot /@ {1, 2}, DisplayFunction -> $DisplayFunction]
>
> Of course one can display them together, though the result is not very
> impressive:
>
> In[19]:=
> Show[polytope, lines]
>
>
> Next, briefly, about the second approach. This method works very nicely with
> the cube example, but less well with the more complex one. I 'll describe it
> only for the cube:
>
> We only the package:
>
> In[1]:=
> << Graphics`ContourPlot3D`
>
> define the test function:
>
> In[2]:=
> test[x_, y_, z_] =
> x >= -1 && x <= 1 && y >= -1 && y <= 1 && z >= -1 && z <= 1;
>
> and the corresponding numerical function:
>
> In[3]:=
> f[x_, y_, z_] /; test[x, y, z] := 1;
> f[x_, y_, z_] := 0
>
> Unfortunately ContourPlot3D is useless here so we need to use
> ListContourPLot3D. We make a table of data
>
> In[5]:=
> data = Table[f[x, y, z],
> {z, -2, 2, .25},
> {y, -2, 2, .25},
> {x, -2, 2, .25}];
>
> now:
>
> In[6]:=
> ListContourPlot3D[data, Contours -> {1}, Lighting -> False, Axes -> True,
> ContourStyle -> {{RGBColor[0, 1, 0]}}]
>
> will produce a nice cube. Unfortunately when we apply the same method to th
> eother problem the computation takes far longer and the resulting polytope is
> very inaccurate.
>
--
Andrzej Kozlowski
Toyama International University
JAPAN
http://platon.c.u-tokyo.ac.jp/andrzej/
http://sigma.tuins.ac.jp/