Question 2 of 7: Gauss quadrature over a rectangular domain

Question 2: Gauss quadrature over a rectangular domain (20 marks)

Question text not reproduced: the examination questions are © Engineers and Geoscientists BC. Open the official past paper (linked at the top of this page) to read the question, then follow the worked solution below.

Given. The integrand and the rectangular integration domain listed below, together with the four bilinear shape functions that map the rectangle onto the parent square.

Find. (a) the Gauss-quadrature estimate of g, using a rule of sufficient order; (b) an explanation of how that estimate compares with the exact value.

Approach. Map the rectangle onto the parent square with the bilinear shape functions, evaluate the constant Jacobian determinant, choose the Gauss order from the polynomial degree of the mapped integrand, and sum the weighted sample values.

1. Build the isoparametric map. With corner coordinates $x_{i} = (4,6,6,4)$ and $y_{i} = (2,2,8,8)$ for nodes 1–4, the geometric interpolation $x = \sum N_{i}x_{i}$, $y = \sum N_{i}y_{i}$ collapses to the familiar mid-point-plus-half-span form for a rectangle:$$x = \frac{4+6}{2} + \frac{6-4}{2}\,\xi = 5 + \xi, \qquad y = \frac{2+8}{2} + \frac{8-2}{2}\,\eta = 5 + 3\eta$$
2. Evaluate the Jacobian. Because the map is affine, the Jacobian matrix is constant over the element:$$[J] = \begin{bmatrix} \partial x/\partial \xi & \partial y/\partial \xi \\\partial x/\partial \eta & \partial y/\partial \eta \end{bmatrix}= \begin{bmatrix} 1 & 0 \\ 0 & 3 \end{bmatrix}, \qquad |J| = 3$$ so $dx\,dy = |J|\,d\xi\,d\eta = 3\,d\xi\,d\eta$.
3. Transform the integral. Substituting the map and the Jacobian gives an integral over the parent square:$$g = \int_{-1}^{1}\!\!\int_{-1}^{1} (5+\xi)^{2}(5+3\eta)^{3}\,(3)\,d\xi\,d\eta$$
4. Choose the quadrature order. An n-point Gauss–Legendre rule integrates a polynomial of degree $2n-1$ exactly. The mapped integrand is degree 2 in $\xi$ and degree 3 in $\eta$, so$$2n-1 \ge 3 \;\Longrightarrow\; n = 2 \ \text{in each direction}$$ A $2 \times 2$ rule with stations $\xi,\eta = \pm 1/\sqrt{3} = \pm 0.57735$ and unit weights is therefore sufficient (and, as shown below, exact).
5. Sum the weighted sample values. Writing $\Phi(\xi,\eta) = 3(5+\xi)^{2}(5+3\eta)^{3}$, the rule is$$g \approx \sum_{i=1}^{2}\sum_{j=1}^{2} W_{i}W_{j}\,\Phi(\xi_{i},\eta_{j}),\qquad W_{i} = W_{j} = 1$$ The four contributions are $\Phi(-0.57735,-0.57735) = 2\,047.92$, $\Phi(+0.57735,-0.57735) = 3\,256.89$, $\Phi(-0.57735,+0.57735) = 17\,903.11$ and $\Phi(+0.57735,+0.57735) = 28\,472.08$.
6. Add the four terms. Summing with unit weights:$$\boxed{\;g \approx 2\,047.92 + 3\,256.89 + 17\,903.11 + 28\,472.08 = 51\,680\;}$$

(b) Comparison with the exact value. The exact integral separates because the integrand is a product of a function of x and a function of y:$$g_{\text{exact}} = \int_{4}^{6} x^{2}dx \int_{2}^{8} y^{3}dy= \frac{6^{3}-4^{3}}{3}\cdot\frac{8^{4}-2^{4}}{4} = \frac{152}{3}\cdot 1020 = 51\,680$$ The quadrature result is not merely close to this value, it is identical. There is no discretization error at all.

The agreement is exact for two reasons acting together. The geometric map is affine, so the Jacobian is a constant that can be carried outside the sampling and introduces no rational function of the natural coordinates; and the mapped integrand is a polynomial of degree 2 and 3 in the two natural directions, which is precisely the degree a two-point rule reproduces exactly ($2n-1 = 3$). Gauss quadrature is not an approximation here, it is an algebraic identity.

The contrast with a lower-order rule makes the point sharper. A single-point rule samples only the centroid, giving $g \approx 4\,\Phi(0,0) = 37\,500$, an error of 27.4 per cent, because a one-point rule reproduces only linear variation and the integrand is cubic in $\eta$. A $3 \times 3$ rule returns 51 680 again: once the rule is exact, adding points buys nothing but arithmetic. In practice this is why element stiffness integrations are sized from the polynomial degree of $\mathbf{B}^{T}\mathbf{D}\mathbf{B}\,|J|$ rather than being over-integrated.