Lab 23: Finding Areas of Arbitrary Polygons

Polygons

You may remember from elementary or middle school the formulas for special polygons such as triangles, trapezoids, and regular hexagons. In this lab we will use the theorems of Green and Stokes to find the area of arbitrary polygons in the plane and in 3-space. In theory we could use Monte Carlo integration to find the area, but for an arbitrary polygon it can be hard to determine whether a point is inside or outside the polygon. Recall Green’s Theorem which says:

Green’s Theorem: Let \(C\) be a positively oriented, piecewise-smooth, simple closed curve in the plane and let \(D\) be the region bounded by \(C\). If \(P\) and \(Q\) have continuous partial derivatives on an open region that contains \(D\), then

\[\int_C P \, dx + Q \, dy = \iint_D \left(\frac{\partial Q}{\partial x} - \frac{\partial P}{\partial y}\right) \, dA.\]

Since we know that \(\iint_D 1 \, dA\) gives the area of a region \(D\) if we can find a \(P\) and \(Q\) such that \(\frac{\partial Q}{\partial x} - \frac{\partial P}{\partial y} = 1\) then we could use Green’s theorem to find the area. There are several ways to do this but one way would be to set \(P(x,y) = -\frac 12y\) and \(Q(x,y) = \frac 12 x\), then \(A = \frac 12 \oint_C x \, dy - y \, dx\). A few moments of thought (and recalling problem 25 in Section 16.4) will show that if \(C\) is the line segment connecting the point \((x_1, y_1)\) to the point \((x_2,y_2)\) then

\[\int_C x \, dy - y \, dx = x_1 y_2 - x_2 y_1.\]

Task 1

Write a function, area(x), that takes an array of points in \(\mathbb R^2\) in counter-clockwise order (positively oriented), x, and computes the area of the polygon they define by using Green’s Theorem.

Polyhedra

Suppose now you wish to find the area of the face of a polyhedron in \(\mathbb R^3\). One way you could accomplish this would be to rotate and translate the polyhedron so that the face lies in the \(xy\) plane. Alternatively you could use Stokes’ theorem, which we state below.

Stokes’ Theorem: Let \(S\) be an oriented piecewise-smooth surface that is bounded by a simple, closed piecewise-smooth boundary curve \(C\) with positive orientation. Let \(\mathbf F\) be a vector field whose components have continuous partial derivatives on an open region in \(\mathbb R^3\) that contains \(S\). Then

\[\int_C \mathbf F \cdot d \mathbf R = \iint_S \mathrm{curl} \mathbf F \cdot d\mathbf S = \iint_S \mathrm{curl} \mathbf F \cdot\mathbf n dS.\]

Recalling that the surface area of \(S\) is equal to \(\iint_S 1 \, dS\), we need only find appropriate values for \(\mathbf F\) and \(\mathbf n\) so that \(\mathrm{curl} \mathbf F \cdot \mathbf n = 1\). Let’s begin with finding an appropriate value for \(\mathbf n\). Since we are interested in finding the surface area of a face of a polyhedron, we know that the face is a plane. The normal to the plane (and hence the surface) can be determined by finding the equation of the plane.

Task 2

Write a function, plane(x), that takes in three points (as an array) in \(\mathbb R^3\), x, and calculates the the coefficients defining the plane that contains all three points i.e. return a, b, c, d such that a * x + b * y + c * z == d for x[0, :], x[1, :], x[2, :].

The Anti-curl

Now that we have found \(\mathbf n = \langle a,b,c \rangle\) it remains to find an appropriate value of \(\mathbf F\) so that \(\mathrm{curl} \mathbf F \cdot \mathbf n = 1.\) In general, you cannot always find the anti-curl. That is to say, given a vector field \(\mathbf G\) it is not always possible to find a vector field \(\mathbf F\) such that \(\mathrm{curl} \mathbf F = \mathbf G.\) However, for this specific problem it can be done.

Suppose that

\[\mathbf G = \left \langle \frac{1}{3a}, \frac{1}{3b}, \frac{1}{3c} \right\rangle.\]

A straightforward calculation shows that \(\mathbf G \cdot \mathbf n = 1.\) Since \(\mathbf G\) is divergence-free (i.e. the divergence of \(\mathbf G\) equals zero), we can find an \(\mathbf F\) such that \(\mathrm{curl} \mathbf F = \mathbf G.\) Unsurprisingly, this \(\mathbf F\) is not unique and one such example of an \(\mathbf F\) is

\[\mathbf F = \frac 16 \left\langle \frac zb - \frac yc, \frac xc - \frac za, \frac ya - \frac xb \right\rangle.\]

(Check this! It’s good practice for your final exam.)

If \(C\) is the line segment connection the point \((x_1, y_1, z_1)\) to the point \((x_2, y_2, z_2)\) then

\[\int_C \mathbf F \cdot d\mathbf r = \frac 16 \left( \frac{x_1y_2-x_2y_1}{c} + \frac{x_2z_1-x_1z_2}{b} + \frac{y_1z_2-y_2z_1}{a} \right).\]

Task 3

Use your code from the previous exercises to write a function, oriented_area(x), that takes an array of points in \(\mathbb R^3\) in counter-clockwise order (positively oriented), x, and computes the area of the polygon they define by using Stokes’ Theorem. Make sure to first verify that all points lie in the same plane (by finding the plane and using np.isclose or np.allclose to verify each point satisfies the plane equation) and raise ValueError("all points must lie in the same plane") otherwise.

Equation of a plane

Suppose we want to find the equation of the plane that passes through the points \(A = (3, 0, -1)\), \(B = (-2, -2, 3)\) and \(C = (7,1,-4)\). The first thing is to create two vectors from the points. We determine the vectors \(\mathbf{AB} = \langle -5, -2, 4\rangle\) and \(\mathbf{AC} = \langle 4, 2, -3\rangle\). These are two vectors that lie on the plane. The next step is to find a vector orthogonal to these two vectors, which can be done by finding the cross product. Hence \(\langle a,b,c \rangle = \mathbf n = \mathbf{AB} \times \mathbf{AC} = \langle -2, 1, -2\rangle.\) Once we have found this normal vector we can use the equation of the plane

\[a(x-x_0) + b(y-y_0) + c(z-z_0) = 0,\]

where \((x_0, y_0, z_0)\) is a point on the plane. Hence the equation of our example plane is

\[-2(x-3) + 1(y-0) - 2(z-(-1)) = 0, \quad \text{ i.e. } \quad 2x - y + 2z = 4.\]