Isoparametric Formulation

The isoparametric formulation is one of the most elegant ideas in FEM. The word comes from Greek: iso (same) + parametric (parameters). It means we use the same shape functions to interpolate both geometry and field variables.

The Core Idea

In an isoparametric element:

Geometry mapping:

$$x = \sum_{i=1}^{n} N_i(\xi, \eta) \, x_i$$

Why Isoparametric?

1. Unified Framework

One set of shape functions handles everything:

Element geometry (curved edges)
Displacement field
Strain computation
Stiffness integration

2. Handles Arbitrary Shapes

The parent element is always a simple shape (square, triangle). The physical element can be:

3. Standardized Integration

All integration happens in the parent element with limits $[-1, 1]$:

$$\int_{\Omega_e} f(x,y) \, dA = \int_{-1}^{1} \int_{-1}^{1} f(\xi,\eta) \, |J| \, d\xi \, d\eta$$

The Jacobian $|J|$ handles the coordinate transformation.

The Parent Element

The parent element (also called reference or master element) lives in natural coordinates $(\xi, \eta)$:

Element	Parent Shape	Coordinate Range
Quad	Unit square	$\xi, \eta \in [-1, 1]$
Triangle	Right triangle	$\xi, \eta \in [0, 1]$, $\xi + \eta \leq 1$
Hex	Unit cube	$\xi, \eta, \zeta \in [-1, 1]$
Tet	Unit tetrahedron	$\xi, \eta, \zeta \in [0, 1]$

Shape functions are defined on the parent element. They're the same regardless of the physical element's shape.

The Jacobian Matrix

The Jacobian relates derivatives in natural coordinates to derivatives in physical coordinates.

Definition

$$[J] = \begin{bmatrix} \frac{\partial x}{\partial \xi} & \frac{\partial y}{\partial \xi} \\[6pt] \frac{\partial x}{\partial \eta} & \frac{\partial y}{\partial \eta} \end{bmatrix}$$

Using the isoparametric mapping:

$$\frac{\partial x}{\partial \xi} = \sum_{i=1}^{n} \frac{\partial N_i}{\partial \xi} x_i$$

Inverse Transformation

To compute strains, we need derivatives with respect to $x$ and $y$:

$$\begin{Bmatrix} \frac{\partial N_i}{\partial x} \\[6pt] \frac{\partial N_i}{\partial y} \end{Bmatrix} = [J]^{-1} \begin{Bmatrix} \frac{\partial N_i}{\partial \xi} \\[6pt] \frac{\partial N_i}{\partial \eta} \end{Bmatrix}$$

The inverse Jacobian:

$$[J]^{-1} = \frac{1}{|J|} \begin{bmatrix} \frac{\partial y}{\partial \eta} & -\frac{\partial y}{\partial \xi} \\[6pt] -\frac{\partial x}{\partial \eta} & \frac{\partial x}{\partial \xi} \end{bmatrix}$$

Jacobian Determinant

The determinant $|J|$ represents the ratio of areas:

$$|J| = \frac{\partial x}{\partial \xi} \frac{\partial y}{\partial \eta} - \frac{\partial x}{\partial \eta} \frac{\partial y}{\partial \xi}$$

$|J| > 0$: Valid element (positive area)
$|J| = 0$: Degenerate element (collapsed)
$|J| < 0$: Inverted element (invalid!)

For a rectangular element aligned with axes:

$$|J| = \frac{a \cdot b}{4}$$

where $a$ and $b$ are the element dimensions.

Computing the B-Matrix

The strain-displacement matrix requires physical derivatives:

Step 1: Natural Derivatives

For Q4, the shape function derivatives are simple:

$$\frac{\partial N_1}{\partial \xi} = -\frac{1}{4}(1-\eta), \quad \frac{\partial N_1}{\partial \eta} = -\frac{1}{4}(1-\xi)$$

Step 2: Jacobian at Integration Point

Compute $[J]$ at each Gauss point $(\xi_g, \eta_g)$:

$$J_{11} = \sum_i \frac{\partial N_i}{\partial \xi} x_i, \quad J_{12} = \sum_i \frac{\partial N_i}{\partial \xi} y_i$$

Step 3: Physical Derivatives

$$\frac{\partial N_i}{\partial x} = J_{22}^{-1} \frac{\partial N_i}{\partial \xi} - J_{12}^{-1} \frac{\partial N_i}{\partial \eta}$$

(where $J_{ij}^{-1}$ are components of $[J]^{-1}$)

Step 4: Assemble B-Matrix

$$[B] = \begin{bmatrix} \frac{\partial N_1}{\partial x} & 0 & \frac{\partial N_2}{\partial x} & 0 & \cdots \\[4pt] 0 & \frac{\partial N_1}{\partial y} & 0 & \frac{\partial N_2}{\partial y} & \cdots \\[4pt] \frac{\partial N_1}{\partial y} & \frac{\partial N_1}{\partial x} & \frac{\partial N_2}{\partial y} & \frac{\partial N_2}{\partial x} & \cdots \end{bmatrix}$$

Stiffness Matrix Computation

Putting it all together:

$$[K_e] = \int_{-1}^{1} \int_{-1}^{1} [B]^T [D] [B] \, t \, |J| \, d\xi \, d\eta$$

Using Gauss quadrature:

$$[K_e] = \sum_{i=1}^{n_g} \sum_{j=1}^{n_g} w_i w_j \, [B(\xi_i, \eta_j)]^T [D] [B(\xi_i, \eta_j)] \, t \, |J(\xi_i, \eta_j)|$$

Algorithm Summary

for each Gauss point (ξᵢ, ηⱼ):
    1. Compute shape functions N and derivatives dN/dξ, dN/dη
    2. Compute Jacobian J and its determinant |J|
    3. Compute inverse Jacobian J⁻¹
    4. Transform: dN/dx = J⁻¹ × dN/dξ
    5. Build B-matrix from dN/dx, dN/dy
    6. K += wᵢwⱼ × Bᵀ × D × B × t × |J|

Subparametric and Superparametric

While isoparametric is most common, alternatives exist:

Subparametric

Lower-order shape functions for geometry
Higher-order for field variables
Example: Linear geometry, quadratic displacement

Superparametric

Higher-order shape functions for geometry
Lower-order for field variables
Rarely used (hard to justify)

Why Isoparametric Dominates

Simplest implementation (one set of functions)
Guaranteed patch test passage
Well-understood error behavior

Element Quality and Jacobian

The Jacobian varies within distorted elements. Quality metrics:

Jacobian Ratio

$$JR = \frac{\min(|J|)}{\max(|J|)}$$

$JR = 1$: Perfect (rectangular/equilateral)
$JR > 0.5$: Good
$JR < 0.3$: Poor (accuracy degraded)
$JR \leq 0$: Invalid element!

Aspect Ratio

$$AR = \frac{\text{longest edge}}{\text{shortest edge}}$$

Keep $AR < 5$ for good accuracy.

Skewness

For quads, interior angles should be near 90°:

Ideal: 90°
Acceptable: 45° to 135°
Poor: < 30° or > 150°

Common Issues

1. Negative Jacobian

Symptom: Solver fails or gives nonsense results Causes:

Inverted element (nodes in wrong order)
Severely distorted shape
Concave quadrilateral

Fix: Check mesh quality, reorder nodes, or remesh

2. Near-Zero Jacobian

Symptom: Ill-conditioned stiffness matrix Causes:

Nearly collapsed element
Very thin elements
Poor aspect ratio

Fix: Refine mesh in problem areas

3. Highly Variable Jacobian

Symptom: Poor stress accuracy Causes:

Severely distorted element
Strong curvature with low-order elements

Fix: Use higher-order elements or refine mesh

Practical Tips

Check Jacobian at all integration points before solving
Use quality metrics during mesh generation
Refine where Jacobian varies significantly
Higher-order elements handle curvature better
Align elements with expected stress gradients

Key Takeaways

Isoparametric: Same shape functions for geometry and field variables
Parent element: Standard shape in natural coordinates $(\xi, \eta)$
Jacobian $[J]$: Transforms derivatives between coordinate systems
Jacobian determinant $|J|$: Ratio of physical to natural areas; must be positive
B-matrix: Requires inverse Jacobian to get physical derivatives
Stiffness integral: $\int [B]^T[D][B] \, |J| \, d\xi \, d\eta$ over parent element
Element quality: Jacobian ratio, aspect ratio, skewness affect accuracy
Negative Jacobian = invalid element: Check mesh quality!

What's Next

With the mathematical framework complete, we turn to linear solvers — the algorithms that actually solve the system $[K]\{u\} = \{F\}$ for the unknown displacements.