# How to plot a 2D region when there are 3 free parameters?

How could I obtain the region of payments, with the parameters a, b, c, d between 0 and 1 and the sum of all equal to 1?

After the comment/recommendation of belisarius, I edited the original question, and ask you a question in the last line of code.

P = {{a, b}, {c, d}}

A = {{2, -1}, {0, 1}}
B = {{-1, 1}, {2, -1}}

A[[1]].P[[1]] + A[[2]].P[[2]]

f1[a_, b_, c_, d_] = A[[1]].P[[1]] + A[[2]].P[[2]]
f2[a_, b_, c_, d_] = B[[1]].P[[1]] + B[[2]].P[[2]]

ParametricPlot[{f1[a, b, c, d], f2[a, b, c, d]}, {  a, 0, 1}, {b, 0,
1}, {c, 0, 1}, {d, 0, 1} ]


In the last line, how can I code to consider only a, b, c, d in $[0,1]$ with $a+b+c+d=1$?

It´s similar to obatining all the convex combinations of the points $(a(i,j),b(i,j))$

• At least for me, it isn't clear what you're trying to do Commented Dec 20, 2013 at 16:11
• @belisarius Yes It´s not very clear! :-) I want to draw the pairs (f1[a,b,c,d],f2[a,b,c,d]), with a,b,c,d between 0 and 1, and verifyng a+b+c+d=1 Commented Dec 20, 2013 at 16:17
• It´s possible that it´s difficult to use Manipulate in this case, because there are 4 parameters (a,b,c,d) interconnected and depending each other. Commented Dec 20, 2013 at 16:18

Being this the first answer at the site using the Region context, I think it deserves an entry on its own:
   AppendTo[$ContextPath, "Region"] z = 0; RegionPlot[ RegionProperty[ TransformedRegion[ ParametricRegion[{{p1, p2, p4}, 0 < p1 < 1 && 0 < p2 < 1 - p1 && 0 < p4 < 1 - p1 - p2}], {2 #1 - #2 + #3, 2 - 3 #1 - #2 - 3 #3} &], {x, y, z}, ImplicitDescription], {x, -3, 3}, {y, -3, 3}]  • Beautiful, many thanks! Btw the ParametricRegion can be defined in a more straightforward manner using just the problem statement : ParametricRegion[{{p1,p2,p3,p4},0<p1<1 && 0<p2<1 && 0<p3<1&& 0<p4<1 && p1+p2+p3+p4==1}],{2 #1 - #2 + #4, - #1 + #2 +2 #3 - #4}(...). – A.G. Commented Dec 22, 2013 at 1:09 To start with, I will work in 4 dimensions to accommodate a,b,c,d. P = {a, b, c, d}; M = {{2, -1, 0, 1}, {-1, 1, 2, -1}};  The points of your 2D set are the linear transforms, of the form M.P (* {2 a - b + d, -a + b + 2 c - d} *)  where P lies in some 4D set which is defined by$a+b+c+d=1$and all positive, that is some polytope in 4D. Since you get your 2D points by a linear transform on the 4D points (M.P), and since the linear transform of a polytope is the convex hull of the transform of the extreme points of the polytope, your set will be the convex hull of the points M.P for all P being an extreme point of the 4D polytope. Not, this 4D polytope is a simple one, you can get its extreme points by all combinations {a,b,c,d} of 0's and 1's that have a sum <= 1. You can generate all the combinations and then Select the right ones : Flatten[Table[{a, b, c, d}, {a, 0, 1}, {b, 0, 1}, {c, 0, 1}, {d, 0, 1}], 3] extremePoints4D = Select[%, Total[#] <= 1 &] (* {{0, 0, 0, 0}, {0, 0, 0, 1}, {0, 0, 1, 0}, {0, 1, 0, 0}, {1, 0, 0, 0}} *)  (yes, you could have guessed that just fine!) Now, transform these to their 2D transforms : Points2D = M.# & /@ extremePoints4D (* {{0, 0}, {1, -1}, {0, 2}, {-1, 1}, {2, -1}} *)  (not all of these will end up extreme points of the 2D region but that's OK, we will take the convex hull and interior points do not matter) Now, plot the convex hull : << ComputationalGeometry Graphics[Polygon[Points2D[[ConvexHull[Points2D]]]], Frame -> True]  Note: This approach needs to generate all extreme points of a polytope; in higher dimensions (4D -> nD) and with constraints not so simple as$a+b+c+d=1$, getting those extreme points gets trickier. About the bounding box @Silvia got me thinking about the bounding box of the 2D region. With the above method it is quite obvious (min and max x and y over extreme points) but I thought I would suggest another way to do this that has the advantage of working in higher dimensions. The idea is to get the min and max x and y over all possible (x,y). Linear programming will do just that : (* Min x (or -x, y, -y) subject to linear constraints a + b + c + d = 1 2a - b - d -x = 0 -a + b + 2c - d - y = 0 0 <= a, b, c, d <= 1 -infinity < x, y < infinity with a 6D vector {a,b,c,d,x,y}. *) obj = {{0, 0, 0, 0, 1, 0}, {0, 0, 0, 0, -1, 0}, {0, 0, 0, 0, 0, 1}, {0, 0, 0, 0, 0, -1}}; (* obj will be used to min/max y/y *) M = {{1, 1, 1, 1, 0, 0}, {2, -1, 0, 1, -1, 0}, {-1, 1, 2, -1, 0, -1}}; (* M.{a,b,c,d,x,y}=0 generates all 3 constraints *) bounds = {{0, 1}, {0, 1}, {0, 1}, {0, 1}, {-Infinity, Infinity}, {-Infinity, Infinity}} (* Now get the 4 solutions to min/max x/y: *) sol = LinearProgramming[#, M, {{1, 0}, {0, 0}, {0, 0}}, bounds] & /@ obj (* and keep the {x,y}'s sol[[1 ;; 4, 5 ;; 6]] (* all downhill from here *) xmin = Min[points[[1 ;; 4, 1]]] xmax = Max[points[[1 ;; 4, 1]]] ymin = Min[points[[1 ;; 4, 2]]] ymax = Max[points[[1 ;; 4, 2]]] ListPlot[points, AxesOrigin -> {xmin - 1, ymin - 1}, Epilog -> {Orange, Line[{{xmin, ymin}, {xmax, ymin}, {xmax, ymax}, {xmin, ymax}, {xmin, ymin}}]}]  • Hi A.G., I guess you missed an image by carelessness. I added it back, hope you don't mind. Commented Dec 22, 2013 at 3:10 You have 3 free parameters, since the imposed condition over their sum takes one away. So, let's plot in 3D the region for those three free parms: RegionPlot3D[(0 < p1 < 1 && 0 < p2 < 1 && 0 < p4 < 1 && p1 + p2 + p4 < 1), {p1, 0, 1}, {p2, 0, 1}, {p4, 0, 1}, PlotPoints -> 30]  Which is equivalent to: Graphics3D[gc = GraphicsComplex[{{0, 0, 0}, {1, 0, 0}, {0, 1, 0}, {0, 0, 1}}, Polygon@{{1, 2, 3}, {2, 3, 4}, {1, 2, 4}, {1, 3, 4}}], Axes -> True]  Now your functions f1 and f2 can bee written as: f[p_, x_] := Total[p x, 2] p = {{p1, p2}, {p3, p4}}; a = {{2, -1}, {0, 1}}; b = {{-1, 1}, {2, -1}}; {f[p, a], f[p, b]} (* {2 p1 - p2 + p4, -p1 + p2 + 2 p3 - p4} *)  Which is your target transform. Now we find a transformation to map any {p1, p2, p3, p4} (minus one parm set by the sum restriction) into your target: r = ToRules@ Reduce[ForAll[{p1, p2, p4}, TransformationFunction[{{a1, a2, a3, a4}, {b1, b2, b3, b4}, {c1, c2, c3, c4}, {d1, d2, d3, d4}}][{p1, p2, p4}] == {f[p, a], f[p, b], 0} /. p3 -> 1 - p1 - p2 - p4], {a1, a2, a3, a4, b1, b2, b3, b4, c1, c2, c3, c4, d1, d2, d3, d4}]; t = TransformationFunction[{{a1, a2, a3, a4}, {b1, b2, b3, b4}, {c1, c2, c3, c4}, {d1, d2, d3, d4}}] /. r /. a1 -> 1  And now we apply that transformation to our graphics complex: Graphics3D[{FaceForm[Black], EdgeForm[None], GeometricTransformation[gc, t]}, Axes -> True]  And we can show it in 2D as: Graphics[Normal@ GeometricTransformation[gc, t] /. {x_, y_, 0} :> {x, y} /. (VertexNormals -> x_) :> Sequence[], Frame -> True]  • Your first RegionPlot3D could be simpler: RegionPlot3D[a + b + c <= 1, {a, 0, 1}, {b, 0, 1}, {c, 0, 1}] – A.G. Commented Dec 21, 2013 at 3:09 Reduce can be used to turn the system {x == 2 a - b + d, y == -a + b + 2 c - d, 0 <= a <= 1, 0 <= b <= 1, 0 <= c <= 1, 0 <= d <= 1, a + b + c + d == 1}  into a system of inequalities of the form (via a CylindricalDecomposition) (x1 < x <= x2 && f1[x] < y < f2[x] && ..) || ..  where the omitted inequalities in the first line constrain the parameters a, b, c, and d in terms of x and y. Only the first pair of each && expression is needed to plot the {x, y} region. From there, the functions f1 and f2 form the boundary of the region, so it can be plotted. The command that will convert the system is basically Reduce[system, {x, y, a, b, c, d}, Reals, Backsubstitution -> True]  A key is putting x and y first. Then we can extract the x and y parts with Reduce[system, {x, y, a, b, c, d}, Reals, Backsubstitution -> True][[All, {1, 2}]]  Some of the entries will be strict equalities. These represent boundary points or lines that bound the inequalities and may be removed. Finally, some of the cases may be combined. Reduce initially divides the region according to the x coordinates of points of intersection and then it subdivides these. If we gather the ones with the same x inequality (with Split below), Reduce can then be used to combine them. Here is this first step: Clear[x, y, a, b, c, d]; P = {{a, b}, {c, d}}; A = {{2, -1}, {0, 1}}; B = {{-1, 1}, {2, -1}}; redSystem = Reduce[{ {x, y} == Inner[Dot, {Hold @@@ A, Hold @@@ B}, Hold @@@ P, Plus] /. Hold -> List, 0 <= a <= 1, 0 <= b <= 1, 0 <= c <= 1, 0 <= d <= 1, a + b + c + d == 1}, {x, y, a, b, c, d}, Reals, Backsubstitution -> True]; removeBoundaries[ineqs_] := Select[ineqs, ! MemberQ[#, _Equal] &]; pieces = removeBoundaries @ redSystem[[All, {1, 2}]]; combinedPieces = Reduce[#, {x, y}, Reals] & /@ Or @@@ SplitBy[List @@ pieces, #[[1, 1]] &] (* {-1 < x <= 0 && -x < y < 2 + x, 0 < x <= 5/7 && -x < y < 1/2 (4 - 3 x), 5/7 < x <= 1 && -x < y < 1/2 (4 - 3 x), 1 < x < 2 && -1 < y < 1/2 (4 - 3 x)} *)  There are two ways to process these inequalities. If the image is a polygon, as in this case, then we can get its convex hull. If the region is bounded by curves, then we can use Plot. Convex Hull method vertices[ineqx_ && ineqy_] := With[{ineq1 = List @@ ineqx, ineq2 = List @@ ineqy}, Thread[{x, ineq2[[{1, -1}]]}] /. ({x -> #} & /@ #) &@ ineq1[[{1, -1}]]]; points = DeleteDuplicates @ Flatten[vertices /@ combinedPieces, 2] (* {{-1, 1}, {0, 0}, {0, 2}, {5/7, -(5/7)}, {5/7, 13/14}, {1, -1}, {1, 1/2}, {2, -1}} *) Needs["ComputationalGeometry"]; Graphics[{ EdgeForm[ColorData[1][1]], RGBColor[0.368417, 0.506779, 0.709798], Opacity[0.2], Polygon[points[[ComputationalGeometryConvexHull @ points]]]}, Frame -> True]  Plot method Clear[rplot]; rplot[_[x1_, ___, x2_] && _[y1_, ___, y2_], opts : OptionsPattern[Plot]] := Plot[{y1, y2}, {x, x1, x2}, Filling -> {1 -> {2}}, opts] Show[ rplot[#, PlotStyle -> ColorData[1][1], FillingStyle -> Directive[Opacity[1], Lighter[RGBColor[0.368417, 0.506779, 0.709798], 0.8]]] & /@ combinedPieces, PlotRange -> All, Frame -> True, Axes -> False, AspectRatio -> Automatic ]  Inspired by belisarius and Simon's answers, I came up with a similar method, but without burdensome RegionPlot. GraphicsRegionRegionInit[]; parareg = TransformedRegion[ ParametricRegion[{ {p1, p2, p4}, 0 < p1 < 1 && 0 < p2 < 1 - p1 && 0 < p4 < 1 - p1 - p2 }], {2 #1 - #2 + #3, 2 - 3 #1 - #2 - 3 #3} &] meshreg = MeshRegion@RegionToMeshObject[Evaluate@{ {x, y}, Prepend[ (* Thanks to belisarius, we have a good way to obtain the boundingbox: *) 1.1 RegionProperty[parareg, "BoundingBox"], RegionProperty[parareg, {x, y}, "ImplicitDescription"] ] }] GeometryPlot[meshreg, PlotRange -> 1.1 RegionProperty[parareg, "BoundingBox"]]  • Why doesn't just GeometryPlot@RegionConvert@parareg show the result? Commented Dec 22, 2013 at 1:03 • @belisarius It looks like RegionConvert will truncate the object at the boundary of a$2\times 2\$ square. Commented Dec 22, 2013 at 1:16
• Yep, GeometryPlot@ RegionConvert@ TransformedRegion[parareg, ScalingTransform[1/10 {1, 1}]] ... but there must be a way to tell RegionConvert to expand the domain Commented Dec 22, 2013 at 1:23
• @belisarius I Trace-ed the RegionConvert and eventually arrived RegionToMeshObject, where the {{-1,1},{-1,1}} boundary is just inserted "suddenly". I think more internal trace will be needed for revealing the key point. Commented Dec 22, 2013 at 1:32
• @Silvia Try RegionProperty[parareg, "BoundingBox"] Commented Jun 10, 2014 at 13:03