Create polygon using edge lengths and area

Question

Four-sided convex land plots are usually denoted with edge lengths and area.

How do I create a Polygon Object with these parameters in Mathematica ?

Answer 1

I'll illustrate one way with an example. The idea is to fix the first two points at origin and on positive x axis, set up equations for the lengths and the area, solve them, then get rid of some symmetric solutions by insisting the third point have positive y value.

I cribbed the area formula from some random Internet site so I'm certain it's correct.

lens = {3, 1, 4, 2};
area = 3;

p1 = {0, 0};
p2 = {First[lens], 0};
p3 = {x3, y3};
p4 = {x4, y4};
{a, b, c, d} = lens;
psq = (p3 - p1).(p3 - p1);
qsq = (p4 - p2).(p4 - p2);
polys = {(p3 - p2).(p3 - p2) - lens[[2]]^2,
   (p4 - p3).(p4 - p3) - lens[[3]]^2,
   p4.p4 - lens[[4]]^2, 
   area^2 - 1/16 (4*psq*qsq - (a^2 + c^2 - b^2 - d^2)^2)};

sol = Solve[polys == 0];
Select[sol, (y3 /. #) >= 0 &] // N

(* Out[198]= {{y4 -> -1.88729833462, x4 -> 0.661895003862, 
  y3 -> 0.536554646919, x3 -> 3.8438655763}, {y4 -> 1.11270166538, 
  x4 -> -1.66189500386, y3 -> 0.733275958395, x3 -> 2.32006884993}} *)

--- edit --- Motivated by the response from @Michael E2, I packaged this as a Module and used NSolve (which, under the hood, also computes a Groebner basis over InexactNumbers).

quad[lens_, area_, x_, y_] := Module[
  {p1, p2, p3, p4, a, b, c, d, psq, qsq, polys, sol},
  p1 = {0, 0};
  p2 = {First[lens], 0};
  p3 = {x[3], y[3]};
  p4 = {x[4], y[4]};
  {a, b, c, d} = lens;
  psq = (p3 - p1).(p3 - p1);
  qsq = (p4 - p2).(p4 - p2);
  polys = {(p3 - p2).(p3 - p2) - lens[[2]]^2,
    (p4 - p3).(p4 - p3) - lens[[3]]^2, p4.p4 - lens[[4]]^2,
    area^2 - 1/16 (4*psq*qsq - (a^2 + c^2 - b^2 - d^2)^2)};
  sol = NSolve[polys];
  Select[sol, With[{y3 = y[3] /. #}, Im[y3] == 0 && y3 >= 0] &]
  ]

I cannot say much about the relative speeds since I have not succeeded in getting both methods to work. My choice of an area formula might lead to suboptimal performance. I show the example from the @Michael E2 response; one of my answers is the same as the one indicated therein.

qsols = quad[{2, 3, 1, 2}, Pi, x, y]

Out[8]= {{y[4] -> 1.72915590151, x[4] -> -1.0049974469, 
  y[3] -> 2.14622307243, 
  x[3] -> -0.0961217816172}, {y[4] -> 1.95571758855, 
  x[4] -> 0.418531616287, y[3] -> 2.80744298296, 
  x[3] -> 0.942520025049}}

--- end edit ---

Answer 2

I had roughly the same idea as Daniel Lichtblau, but with a playful twist that turned out to improve the speed. (I was surprised.) So on that basis I feel it's worth sharing.

We'll set up some polynomials that will define the desired quadrilateral. I made the origin an explicit argument, even though later the quadrilateral function will use it only as {0, 0}. The Total[..] line creates equations that specify the distance of three sides to be s2, s3,s4. (The first side from the origin topwill be preassigned top = {s1, 0}`. One could add the equation here, but we have the freedom to determine the position of one side of the quadrilateral.) To these are added the equation for the area.

area[orig_, p_, q_, r_] := Det[{p - orig, r - orig}]/2 + Det[{r - q, p - q}]/2;

quadIdeal[s1_, s2_, s3_, s4_, A_, orig_, p_, q_, r_] := 
  Join[
   Total[Subtract @@@ Partition[{p, q, r, orig}, 2, 1]^2, {2}] - {s2, s3, s4}^2,
   {area[orig, p, q, r] - A}];

We use GroebnerBasis to reduce these polynomials. It puts the polynomials in an order such that their roots can be solved sequentially -- a nice setup for Fold. We apply Reduce and convert the equations it outputs to Set, to assign the value to the coordinate variable in the equation. In GroebnerBasis we set CoefficientDomain -> InexactNumbers, so that any numeric arguments may be processed. One caveat: if the area is too great, the figure is impossible and complex roots will be returned.

quadrilateral[s1_?NumericQ, s2_?NumericQ, s3_?NumericQ, s4_?NumericQ, A_?NumericQ] := 
  Polygon @
   Block[{x2, y2, x3, y3},
    Fold[
      Reduce[#2 == 0, Variables[#2]] /. Equal -> Set &,
      0,
      GroebnerBasis[#, Reverse@Variables[#], CoefficientDomain -> InexactNumbers] &@
        quadIdeal[s1, s2, s3, s4, A, {0, 0}, {s1, 0}, {x2, y2}, {x3, y3}]];
    {{0, 0}, {s1, 0}, {x2, y2}, {x3, y3}}
    ];

Example -- "circle" (π) quadrature :-)

N @ quadrilateral[2, 3, 1, 2, Pi]

Polygon[{{0., 0.}, {2., 0.}, {0.94252, 2.80744}, {0.418532, 1.95572}}]

Graphics[quadrilateral[2, 3, 1, 2, Pi], Frame -> True]

Surprisingly using Fold etc. this way is faster than Solve or NSolve. Edit - updated timings. On a fresh kernel everything ran a little faster today, so I'll report them with the other updates. Also GroebnerBasis with InexactNumbers uses very high precision arithmetic apparently (perhaps a hundred digits). One can save a little if at least one machine precision number is passed to quadrilateral:

Do[quadrilateral[2, 3, 1, 3, N @ Pi], {100}] // Timing // First

Do[quadrilateral[2, 3, 1, 3, Pi], {100}] // Timing // First

Do[Solve[quadIdeal[2, 3, 1, 3, Pi, {0, 0}, {2, 0}, {x2, y2}, {x3, y3}] == 0,
         {x2, y2, x3, y3}], {100}] // Timing // First

Do[NSolve[quadIdeal[2, 3, 1, 3, Pi, {0, 0}, {2, 0}, {x2, y2}, {x3, y3}] == 0,
          {x2, y2, x3, y3}], {100}] // Timing // First

0.202163
0.220306
1.774986
1.049178

Using @DanielLichtblau's quad (also a little faster with machine precision input):

Do[quad[{2, 3, 1, 2}, N @ Pi, x, y], {100}] // Timing // First
Do[quad[{2, 3, 1, 2}, Pi, x, y], {100}] // Timing // First

1.548769
1.581755

---

To see how it works, let's examine the steps. Here are the initial polynomials.

quadIdeal[1, 4, 2, 5, 6, {0, 0}, {x1, 0}, {x2, y2}, {x3, y3}]

{-16 + (1 - x2)^2 + y2^2,  -4 + (x2 - x3)^2 + (y2 - y3)^2,
 -25 + x3^2 + y3^2,  -6 + y3/2 + 1/2 (y2 - x3 y2 - y3 + x2 y3)}

GroebnerBasis reduces the set to a simpler one that has the same roots. Because of the nature of this problem and the way GroebnerBasis works, it turns out to be easy to solve for the roots.

GroebnerBasis[%, Reverse@Variables[%]]

{21 - 26 x2 + 5 x2^2,  33 - 20 x2 + 17 x3, -9 + x2 + 2 y2,  -89 + 5 x2 + 17 y3}

The root of the first polynomial yields x2; here Reduce returns two equations, and when they are converted to Set, the second assignment overrides the first. The second polynomial, with x2 being found, yields x3. Finally y2 and y3 are set, again x2 having been set earlier.