# On integral of a function over a simplex

Help w/the following general calculation and references would be appreciated.

Let $ABC$ be a triangle in the plane. Then for any linear function of two variables $u$. $$\int_{\triangle}|\nabla u|^2=\gamma_{AB}(u(A)-u(B))^2+ \gamma_{AC}(u(A)-u(C))^2+\gamma_{BC}(u(B)-u(C))^2,$$ where $$\gamma_{AB}=\frac{1}{2}\cot(\angle C), \gamma_{AC}=\frac{1}{2}\cot(\angle B), \gamma_{BC}=\frac{1}{2}\cot(\angle A).$$

What is a good reference for the formula? Is it due to R. Duffin?

Is there generalization to linear functions of three variables? The number of parameters fits (in any $nD$), and it seems like one needs to calculate 6x6 Cayley-Menger like determinants, but it's difficult and the geometric interpretation is not clear. Thank you.

-
Yes, there is, using barycentric coordinates is the way to go. – Shuhao Cao May 18 '13 at 5:35
@Shuhao Cao:Is the system in my answer simple in the coordinates? – DVD May 19 '13 at 2:33
Yes, much simpler in Barycentric coordinates, it is basically a local coordinate system over simplex. I will derive the answer in 3D maybe later. – Shuhao Cao May 19 '13 at 2:36
I don't understand this —maybe you could explain the notation... If $u$ is a linear function, then $\nabla u$ is a constant vector, and the integral you are trying to compute is seems to be just the squared norm of that vector times the area of the triangle. – Mariano Suárez-Alvarez May 19 '13 at 5:51
@MarianoSuárez-Alvarez Yes, but I guess Daved wants to get a formula like in the post using trigonometric functions. – Shuhao Cao May 19 '13 at 5:57

It is always easier to work on an affine coordinate system, Barycentric coordinate system, when the integration is performed locally over a piecewise linear structure, than to work on the global Cartesian coordinates.

First let's replicate the formula you gave in a triangle (2-simplex). For the following triangle $T = \triangle ABC$:

denote the top vertex as $A$, lower left vertex as $B$, and lower right vertex as $C$. With slightly abuse of notation, $A$ also denotes the Cartesian coorinates: $A = (x_A,y_A)$, etc. The unit outward normal vector $\nu_{A}$ opposite to the vertex $A$, etc. The edge vector $e_A = C- B$, $e_B = A - C$, $e_C = B- A$, rotating counterclockwisely.

Then for any point $P$ inside or on the boundary can be written as the weight sum : $$P = A\lambda_A + B \lambda_B + C\lambda_C$$ where $\lambda_A(A)= 1$, and $\lambda_A(B)= \lambda_A(C) = 0$. This is, for any vertex $V$, $\lambda_V$ is a linear function valued 1 on the vertex $V$, 0 on two other vertices, also we have: $$\lambda_A(P) + \lambda_B(P) +\lambda_C(P) = 1.$$ Any linear function $u(x,y) = ax + by +c$ restricted on this triangle is a plane over this triangle:

Now We can rewrite $u(P)$ using the barycentric coordinate system (for you to check): $$u(P) =u\big(A\lambda_A(P) + B \lambda_B(P) + C\lambda_C(P) \big) = u(A)\lambda_A(P) + u(B) \lambda_B(P) + u(C)\lambda_C(P) .$$ in short we just write $$u = u(A)\lambda_A + u(B) \lambda_B + u(C)\lambda_C$$

The gradient of $u$ over this triangle is: $$\nabla u = u(A)\nabla\lambda_A + u(B) \nabla\lambda_B + u(C)\nabla\lambda_C .$$ For vertex $A$, the gradient lies in the direction of descent of a plane 1 at $A$, 0 at $B$ and $C$, and is inversely proportional to the height $h_A$ on $e_A$: $$\nabla \lambda_A = -\frac{1}{h_A} \nu_A = -\frac{|e_A|}{2|T|}\nu_A. \tag{1}$$ Then \begin{aligned} &|\nabla u |^2 = \nabla u \cdot \nabla u \\ =& (u(A)\nabla\lambda_A + u(B) \nabla\lambda_B + u(C)\nabla\lambda_C)\cdot (u(A)\nabla\lambda_A + u(B) \nabla\lambda_B + u(C)\nabla\lambda_C) \\ =& u(A)^2 |\nabla\lambda_A|^2 + u(B)^2 |\nabla\lambda_B|^2 + u(C)^2 |\nabla\lambda_C|^2 \\ &+ 2u(A) u(B)\nabla\lambda_A\cdot \nabla\lambda_B + 2u(B) u(C)\nabla\lambda_B\cdot \nabla\lambda_C + 2u(C) u(A)\nabla\lambda_C\cdot \nabla\lambda_A \end{aligned} By (1): $$\nabla\lambda_A\cdot \nabla\lambda_B = \frac{|e_A|}{2|T|}\nu_A\cdot \frac{|e_B|} {2|T|}\nu_B = -\frac{|e_A||e_B|}{4|T|^2} \cos(\angle C)$$ we have: \begin{aligned} |\nabla u |^2 &= u(A)^2\frac{|e_A|^2}{4|T|^2} + u(B)^2\frac{|e_B|^2}{4|T|^2} + u(C)^2\frac{|e_C|^2}{4|T|^2} \\ &\quad - 2u(A) u(B) \frac{|e_A||e_B|}{4|T|^2} \cos\angle C - 2u(B) u(C) \frac{|e_B||e_C|}{4|T|^2} \cos\angle A \\ &\quad - 2u(C) u(A) \frac{|e_C||e_A|}{4|T|^2} \cos\angle B.\tag{2} \end{aligned} To get the formula, denote the height $h_A$ touches the base $e_A$ at $A'$, similar for height $h_B$ and $h_C$. Notice $u(A)^2$'s coefficient times area $|T|$ is: $$\frac{|e_A|^2}{4|T|} = \frac{|e_A|}{2 h_A} = \frac{|BA'|+|A'C|}{2 h_A} = \frac{1}{2} (\cot\angle B + \cot\angle C) .$$ similar for other two. The $u(A)u(B)$'s coefficient times area $|T|$ is: $$\frac{|e_A||e_B|}{4|T|} \cos\angle C = \frac{|e_A||e_B|}{2 |e_A| h_A} \frac{|A'C|}{|e_B|} = \frac{|A'C|}{2h_A} =\frac{1}{2} \cot \angle C$$ Therefore, (2) becomes: \begin{aligned} |\nabla u |^2 &= u(A)^2\frac{1}{2|T|} (\cot\angle B + \cot\angle C)\\ &\quad + u(B)^2\frac{1}{2|T|} (\cot\angle C + \cot\angle A) \\ &\quad + u(C)^2\frac{1}{2|T|} (\cot\angle A + \cot\angle B) \\ &\quad - 2u(A) u(B) \frac{1}{2|T|} \cot \angle C \\ &\quad - 2u(B) u(C) \frac{1}{2|T|} \cot \angle A \\ &\quad - 2u(C) u(A) \frac{1}{2|T|} \cot \angle B.\tag{3} \end{aligned} Finally for we are integrating a constant number on this triangle $T$: $$\int_{T} |\nabla u|^2 = |\nabla u|^2 |T|,$$ and (3) multiplying with the area $|T|$ is the formula you gave.

For a $u(x,y,z) = ax + by +cz +d$, the integration of $\nabla u$ on a 3-simplex $T = V_1 V_2 V_3 V_4$ (a tetrahedron), denote the barycentric coordinate of vertex $V_i$ as $\lambda_i$, $\lambda_i(V_j) = \delta_{ij}$. Also let $u(V_i) = u_i$, then: $$\nabla u = \sum_{i=1}^4 u_i \nabla \lambda_i.$$ And $\nabla \lambda_i$ is pointed in the opposite direction of the unit outer normal vector $\nu_i$ to the face $F_i$(a triangle) opposite to vertex $V_i$, also inversely proportional to the height: $$\nabla \lambda_i = -\frac{1}{h_i}\nu_i = -\frac{|F_i|}{3|T|} \nu_i.$$ In the cross terms of $|\nabla u|^2$: $$\nu_i\cdot \nu_j = \cos \angle F_i F_j$$ where $\angle F_i F_j$ is the dihedral angle (angle between two planes, notice the angle is defined using normal vectors' inner product, for higher dimensions, I am afraid we have to use normals instead of trigonometric functions), working the same routine as in 2 dimensional case, we will get a similar formula: $$\int_T |\nabla u|^2 = \sum^4_{i=1}u_i^2\frac{|F_i|^2}{9|T|} - 2\sum_{i,j} u_i u_j\frac{|F_i||F_j|}{9|T|}\cos\angle F_i F_j, \tag{\dagger}$$ Consider $|F_1|^2/(9|T|)$: $$\frac{|F_1|^2}{9|T|} = \frac{|F_1|}{3 h_1}.$$ We draw the height of base $F_1$ from vertex $V_1$, denote $P_1$ as the point where this height $h_1$ touches the base $F_1$. $F_1$ has three sides: $e_{12}$, $e_{13}$, and $e_{14}$, where $e_{ij}$ is the common edge of $F_i$ and $F_j$.

From $P_1$ draw three perpendicular line segment $m_{12}$, $m_{13}$, and $m_{14}$ to $e_{12}$, $e_{13}$, and $e_{14}$ (this happens on the triangle $F_1$), then we can decompose the area $|F_1|$ into three parts: \begin{aligned} \frac{|F_1|}{3 h_1} &= \frac{|m_{12}| |e_{12}| + |m_{13}| |e_{13}| + |m_{14}| |e_{14}|}{6 h_1} \\ &= \frac{1}{6}( \cot\angle F_1 F_2 |e_{12}| + \cot\angle F_1 F_3 |e_{13}| + \cot\angle F_1 F_4 |e_{14}|) \end{aligned}\tag{5} we can decompose the term for other three vertices also in this fashion.

Now for coefficient $|F_1||F_2|\cos\angle F_1 F_2 /(9|T|)$ (we have six this types of coefficient corresponding to six edges of a tetrahedron):

$$\frac{|F_1||F_2|}{9|T|}\cos\angle F_1 F_2 = \frac{|F_1||F_2|}{3|F_2| h_2}\cos\angle F_1 F_2 = \frac{|m_{12}||e_{12}|}{6 h_1} = \frac{1}{6} \cot\angle F_1 F_2 |e_{12}| \tag{6}$$ doing this for other 5 edges also. Now plugging (5) and (6) into $(\dagger)$ we have: $$\int_T |\nabla u|^2 = \sum_{i<j} \big(u_i - u_j \big)^2 \frac{1}{6} |e_{ij}|\cot\angle F_i F_j, \tag{\ddagger}$$ which is a formula in three dimension. Just a reminder for notation: $u_i$ is the value of $u$ at vertex $V_i$, $F_i$ is the face opposite to $F_i$, $\angle F_i F_j$ is the dihedral angle between face $F_i$ and $F_j$, $e_{ij}$ is the common edge of $F_i$ and $F_j$, not the edge from vertex $i$ to $j$.

The result can be generalized to the integration the gradient square of a linear $u:\mathbb{R}^n\to \mathbb{R}$ on an $n$-simplex, $\triangle^n$, which is the convex hull of $(n+1)$ points.

Things are similar to the case in 3-simplex, $F_i$ is the $(n-1)$-face (codimension is 1) "opposite" to $V_i$, $\nu_i$ is unit outward normal to this $(n-1)$-face. We should be able to get a formula like $(\dagger)$, but I am afraid I could not see a trigonometrical interpretation like $(\ddagger)$.

-
So, it's $1/L$, where $L$ is the length of the interval, in 1D. It's $\cot(\angle)/2$ in 2D. Having trouble seeing the geometric pattern for nD... Please, help. – DVD May 19 '13 at 4:49
@Daved I got a formula in 3D, however in higher dimension I don't think I could get a formula with such degree of trigonometrical meaning. – Shuhao Cao May 19 '13 at 5:39
It is a nice solution but it is too long. – xpaul May 22 '13 at 19:49

The solution by @Shuhao Cao is nice. However it is too long to follow. There is a short way to do this. It is easy to check that $\int_{\Delta}|\nabla u|^2$ is invariant under translation and rotation $$x\to x+h,y\to y+k, x\to x\cos\theta-y\sin\theta,y\to x\sin\theta+y\cos\theta.$$ Thus we can take the coordinates of $A,B,C$ as $$A(0,0), B(x_B,0), C(x_C,y_C)$$ respectively, from which we obtain $$\cot A=\frac{x_C}{y_C}, \cot B=\frac{x_B-x_C}{y_C}, |\Delta|=\frac{1}{2}x_By_C.$$ From the identity $$\cot A\cot B+\cot B\cot C+\cot C\cot A=1$$ we have $$\cot C= \frac{-x_Bx_C+x_C^2+y_C^2}{x_By_C}$$ Let $u=\alpha x+\beta y+\gamma$ and hence $|\nabla u|^2=\alpha^2+\beta^2$. So we have to show the following identity $$(\alpha^2+\beta^2)|\Delta|=\frac{1}{2}\alpha^2 x_B^2\cot C+\frac{1}{2}(\alpha x_C-\beta y_C)^2\cot B+\frac{1}{2}[\alpha (x_B-x_C)-\beta y_C]^2\cot A.$$ Note \begin{eqnarray*} RHS&=&\frac{\alpha^2}{2}[x_B^2\cot C+x_C^2\cot B+(x_B-x_C)^2\cot A]\\ &&+\frac{\alpha^2}{2}y_C^2(\cot B+\cot A)+\alpha\beta[-x_Cy_C\cot B+(x_B-x_C)\cot A]. \end{eqnarray*} It is easy to check $$y_C^2(\cot B+\cot C)=x_By_C,-x_Cy_C\cot B+(x_B-x_C)\cot A=0$$ and $$x_B^2\cot C+x_C^2\cot B+(x_B-x_C)^2\cot A=x_By_C.$$ Thus $$RHS=(\alpha^2+\beta^2)|\Delta|=LHS.$$

-
But what about 3D or nD?.. – DVD May 22 '13 at 22:51
@Daved: for 3D or nD, the proof is the same since the integral is invariant under translation and O(n). So what you need to do is just to work in 2D plane. – xpaul May 23 '13 at 1:51
There are $n^2$ order parameters/unknowns in the problem. The translation and rotation invariance decreases their number only by order $n^2/2$. So, the problem is still complicated in nD... – – DVD May 28 '13 at 3:34

Let $(x_k, y_k,z_k)$ be four points in 3D. Then for a linear function

$$u(x,y,z)=Ax+By+Cz+D,$$ and $$x_{kl}=x_k-x_l, \dots$$ matching coefficients in the identity in question, one gets the following system of equations:

$$\begin{pmatrix} x_{12}^2 & x_{13}^2 & \dots & x_{34}^2 \\ y_{12}^2 & y_{13}^2 & \dots & y_{34}^2 \\ z_{12}^2 & \dots & \dots & \dots \\ x_{12}y_{12} & x_{13}y_{13} & \dots & \dots \\ x_{12}z_{12} & \dots & \dots & \dots \\ y_{12}z_{12} & y_{13}z_{13} & \dots & y_{34}z_{34} \\ \end{pmatrix} \begin{pmatrix} \gamma_{12} \\ \gamma_{13} \\ \dots \\ \dots \\ \dots \\ \gamma_{34} \end{pmatrix} = Vol_{\triangle} \begin{pmatrix} 1 \\ 1 \\ 1 \\ 0 \\ 0 \\ 0 \end{pmatrix},$$ where $$x_{kl}+x_{lm}+x_{mk}=0, \dots$$ Can these expression be simplified? Is there geometric interpretation of the solution?

-
If you are intending to ask a question, doing so in an answer box is not a good idea. Please ask the question as a new question! – Mariano Suárez-Alvarez May 19 '13 at 5:49
Okay, thank you, got the answer from @Shuhao Cao. – DVD May 19 '13 at 6:06