# Computational Biomechanics

Dr. Kewei Li

## Numerical Integration

### One-dimensional Integration

Numerical integration is a method to approximate a definite integral up to a given degree of accuracy.

$I = \int_{-1}^{1} f(x) \mathrm{d} x \approx \sum_{i=1}^{p} f(x_i) w_i$

where $$w_i$$ are integration weights and $$x_i$$ are integration points. It is also called Gauss point. In the following table, you can find Gauss-Legendre quadrature points and weights.

$$p$$ $$x_i$$ $$w_i$$
1 $$x_1 = 0.0$$ $$w_1=2.0$$
2 $$x_1 = -1/\sqrt{3}$$
$$x_2 = 1/\sqrt{3}$$
$$w_1=1.0$$
$$w_2=1.0$$
3 $$x_1=-\sqrt{0.6}$$
$$x_2 = 0.0$$
$$x_3 = \sqrt{0.6}$$
$$w_1=5.0/9.0$$
$$w_2 = 8.0/9.0$$
$$w_3=5.0/9.0$$

Thus, if the two-point integration rule is used, then

$I \approx f(x_1 = -\frac{1}{\sqrt{3}} ) + f(x_2 =\frac{1}{\sqrt{3}})$

For example, if $$f(x) = x^2$$, then

$\int_{-1}^{1} x^2 \mathrm{d} x \approx f(x_1 = -\frac{1}{\sqrt{3}} ) + f(x_2 =\frac{1}{\sqrt{3}}) = \frac{2}{3}$

The exact solution is

$\int_{-1}^{1} x^2 \mathrm{d} x = \left[ \frac{x^3}{3} \right]^{1}_{-1} = \frac{2}{3}$

For higher order functions, you may need to use more integration points, see the Gaussian Quadrature page on Wikipedia for more information.

### Two-dimensional Integration

For two-dimensional domain $$\{(x, y) \mid -1 \le x \le 1, -1 \le y \le 1\}$$, we usually use the $$2 \times 2$$ product rule with $$p=4$$ integration points, thus

$I = \int_{-1}^{1}\int_{-1}^{1} f(x,y) \mathrm{d} x \mathrm{d} y \approx \sum_{i=1}^{p} f(x_i, y_i) w_i$

The integration points and weights are given in the following table.

$$i$$ $$x_i$$ $$y_i$$ $$w_i$$
1 $$-1/\sqrt{3}$$ $$-1/\sqrt{3}$$ $$1.0$$
2 $$1/\sqrt{3}$$ $$-1/\sqrt{3}$$ $$1.0$$
3 $$1/\sqrt{3}$$ $$1/\sqrt{3}$$ $$1.0$$
4 $$-1/\sqrt{3}$$ $$1/\sqrt{3}$$ $$1.0$$

This integration scheme is often used in the calculation of stiffness matrix for four-node quadrilateral elements.

### Three-dimensional Integration

For a three-dimensional domain such as,

$\{(x, y, z) \mid -1 \le x \le 1, -1 \le y \le 1, -1 \le z \le 1\}$

we usually use the $$2 \times 2 \times 2$$ product rule with $$p=8$$ integration points, thus

$I = \int_{-1}^{1}\int_{-1}^{1}\int_{-1}^{1} f(x,y, z) \mathrm{d} x \mathrm{d} y \mathrm{d} z \approx \sum_{i=1}^{p} f(x_i, y_i, z_i) w_i$

The integration points and weights are given in the following table.

$$i$$ $$x_i$$ $$y_i$$ $$z_i$$ $$w_i$$
1 $$-1/\sqrt{3}$$ $$-1/\sqrt{3}$$ $$-1/\sqrt{3}$$ $$1.0$$
2 $$1/\sqrt{3}$$ $$-1/\sqrt{3}$$ $$-1/\sqrt{3}$$ $$1.0$$
3 $$1/\sqrt{3}$$ $$1/\sqrt{3}$$ $$-1/\sqrt{3}$$ $$1.0$$
4 $$-1/\sqrt{3}$$ $$1/\sqrt{3}$$ $$-1/\sqrt{3}$$ $$1.0$$
5 $$-1/\sqrt{3}$$ $$-1/\sqrt{3}$$ $$1/\sqrt{3}$$ $$1.0$$
6 $$1/\sqrt{3}$$ $$-1/\sqrt{3}$$ $$1/\sqrt{3}$$ $$1.0$$
7 $$1/\sqrt{3}$$ $$1/\sqrt{3}$$ $$1/\sqrt{3}$$ $$1.0$$
8 $$-1/\sqrt{3}$$ $$1/\sqrt{3}$$ $$1/\sqrt{3}$$ $$1.0$$

This integration scheme is often used in the calculation of stiffness matrix for eight-node hexahedral elements.

Go Back