1

I have a real algebraic variety given by the equation $x (z+v) - y (x+y)=0$ that I want to draw inside a probability simplex, so inside a tetrahedron given by the equations $x + y + z + v = 1$ and $0 < x, y, z, v < 1$.

My question is how to plot this using Mathematica.

I tried to simply substitute $v=1-x-y-z$ into the first equation and use counterplot or 3dplot or regionplot to achieve the result but this does not yield what I require: it does show some surface in 3D but not the simplex and strangely the z-coordinate is not at all where I thought it would be. So I guess the projection is somehow messed up by the naive substitution (I tried to replicate the problem as a projection from 3D to 2D, so onto a surface x+y+z=1 using Gram Schmidt and a proper orthogonal projection but I'm lacking intuition in order to find errors in my calculations.).

Below is a picture of what I imagine my surface should look like (roughly, just to get the idea).

enter image description here

I'm happy to do some more maths on paper in order to simplify my problem before throwing it at mathematica but at the moment I just don't know in what direction to go. Any help would be greatly appreciated!

Community
  • 1
howwouldIknow
  • 11
  • 1
  • 1
    Just a remark: in your image, the gray surface is a 2D object (even though it's supposed to represent a 3D object). In the end, to you want a 2D surface in 3D? Otherwise I don't think it's possible (except if you allow time as an additional parameter). – anderstood Nov 16 '16 at 20:36

2 Answers2

3

Using an approach similar to the 4D tetrahedron insphere here, here's a way to project the surface from 4D to 3D. I parametrized the surface over two of the variables (z, v).

proj0 = N@ Orthogonalize[IdentityMatrix[4]~Prepend~{1, 1, 1, 1}][[2 ;; 4]];
proj[pts_?VectorQ] := proj0.pts;
proj[pts_?MatrixQ] := Transpose[proj0.Transpose@pts];

dom = DiscretizeRegion[
   ImplicitRegion[v + z <= 1, {{z, 0, 1}, {v, 0, 1}}],
   MaxCellMeasure -> {"Length" -> 0.05}];

mr = MeshRegion[
   proj@Transpose[
     {x, y, z, v} /. 
       First@Solve[{x*(z + v) - y*(x + y) == 0, 
          x + y + z + v == 1}, {x, y}] /. 
      Thread[{z, v} -> Transpose@MeshCoordinates[dom]]],
   MeshCells[dom, 2]
   ];

Show[
 Graphics3D[
  {EdgeForm[Black], Opacity[0.], Tetrahedron[proj@IdentityMatrix[4]]}
  ],
 mr]

Mathematica graphics

Community
  • 1
Michael E2
  • 235,386
  • 17
  • 334
  • 747
  • You seem to have forgotten to include your definition dom in your answer? – Edmund Nov 17 '16 at 01:27
  • @Edmund Oops. Thanks! – Michael E2 Nov 17 '16 at 01:35
  • That's great, thank you! Is there a way to cancel the box around it and add some names to the vertices? – howwouldIknow Nov 17 '16 at 12:45
  • Add the option Boxed -> False to Graphics3D[] or to Show[]. You can add Text[..] as in Rahul's answer. The vertices are at the projections of the 4D vertices, e.g. Text[HoldForm[{1, 0, 0, 0}], proj[{1, 0, 0, 0}], offset], where offset controls how to offset the text from the point as described in the documentation for Text. In 3D graphics, it is sometimes better to add a vector to the point than to use offset. – Michael E2 Nov 17 '16 at 13:32
3

One way to do this using ContourPlot3D: Choose a mapping from 3D space $\{(p_x,p_y,p_z)\}$ to the 4D simplex $\{(x,y,z,v): x+y+z+v=1\}$, express the desired equation $x(z+v)-y(x+y)=0$ in terms of $(p_x,p_y,p_z)$, and draw its contour.

You can choose any four points to map the vertices of the simplex to, but I'll use the built-in regular tetrahedron from PolyhedronData.

vertices = PolyhedronData["Tetrahedron", "VertexCoordinates"]
xyzv[px_, py_, pz_] := 
 Evaluate@First@
   Solve[{Total[vertices {x, y, z, v}] == {px, py, pz}, 
     x + y + z + v == 1}, {x, y, z, v}]
min = Min /@ Transpose[vertices];
max = Max /@ Transpose[vertices];
Show[Graphics3D[{FaceForm[None], Tetrahedron[vertices], 
   Text["(1,0,0,0)", vertices[[1]], {-1, 0}], 
   Text["(0,1,0,0)", vertices[[2]], {0, 1}], 
   Text["(0,0,1,0)", vertices[[3]], {1, 0}], 
   Text["(0,0,0,1)", vertices[[4]], {0, -1}]}], 
 ContourPlot3D[
  Evaluate[(x (z + v) - y (x + y) == 0) /. xyzv[px, py, pz]], 
  {px, min[[1]], max[[1]]}, {py, min[[2]], max[[2]]}, {pz, min[[3]], max[[3]]}, 
  RegionFunction -> Function[{px, py, pz}, 
    And @@ Thread[({x, y, z, v} /. xyzv[px, py, pz]) >= 0]], 
  Mesh -> False], Boxed -> False, ViewPoint -> {-2, -2, 0}]

enter image description here

  • This looks exactly like I wanted! Great, thank so you much! Can you also tell me whether it's possible to generalise this solution to the case when I have more than one equation defining my surface inside the simplex? – howwouldIknow Nov 17 '16 at 12:46
  • If you have two equations, you would get a curve, not a surface, right? In that case, take a look at this question: http://mathematica.stackexchange.com/q/5968/484 ...Three equations would give you just a point (or a finite set of points), in which case you should just use Solve. –  Nov 17 '16 at 20:07