# Transformation auf das Referenzelement In der Finite-Elemente-Methode werden Basisfunktionen und Quadraturregeln nur einmal auf einem **Referenzelement** $T$ definiert. Jedes Element $T_i$ der Triangulierung wird über eine affine Abbildung $\sigma_i$ darauf zurückgeführt. ## Allgemeines Prinzip Gegeben sei ein nicht entartetes Simplex $T_i\subset\mathbb{R}^d$ mit Ecken $\vec{p}_{i,1},\dots,\vec{p}_{i,d+1}$ sowie das Referenzsimplex $$ T = \Bigl\{\vec{t}\in\mathbb{R}^d \;:\; t_k\ge 0,\ \textstyle\sum_{k=1}^d t_k \le 1\Bigr\} $$ mit den Ecken $\vec{0}$ und $\vec{e}_k$ ($k=1,\dots,d$). Die Koordinaten im Referenzelement heissen $\xi$ (1D), $(\xi,\eta)$ (2D) bzw. $(\xi,\eta,\zeta)$ (3D), zusammengefasst im Vektor $\vec{t}$. **Abbildung Referenz → Element.** Die affine Abbildung ```{math} :label: eq-ref-sigma \sigma_i : T \to T_i,\qquad \vec{x} = \sigma_i(\vec{t}\,) = \vec{p}_{i,1} + A_i\cdot\vec{t}, \qquad A_i = \bigl(\vec{p}_{i,2}-\vec{p}_{i,1},\ \vec{p}_{i,3}-\vec{p}_{i,1},\ \dots,\ \vec{p}_{i,d+1}-\vec{p}_{i,1}\bigr)\in\mathbb{R}^{d\times d} ``` bildet die Ecken des Referenzelements auf die Ecken $\vec{p}_{i,k}$ ab. Die Spalten der Jacobi-Matrix $A_i = D\sigma_i$ sind die von $\vec{p}_{i,1}$ ausgehenden Kantenvektoren. Die Jacobi-Determinante bezeichnen wir mit $$ J_i = \det A_i . $$ **Rücktransformation Element → Referenz.** Ist $J_i\neq 0$ (d.h. $T_i$ ist nicht entartet), so gilt ```{math} :label: eq-ref-sigma-inv \sigma_i^{-1} : T_i \to T,\qquad \vec{t} = \sigma_i^{-1}(\vec{x}) = A_i^{-1}\cdot\bigl(\vec{x}-\vec{p}_{i,1}\bigr), \qquad D\sigma_i^{-1}(\vec{x}) = A_i^{-1}. ``` **Wichtige Folgerungen.** - Volumen: $|T_i| = |J_i|\,|T|$ mit $|T| = 1/d!$. - Integrale: $\displaystyle\int_{T_i} f(\vec{x})\,d\vec{x} = \int_{T} f\bigl(\sigma_i(\vec{t}\,)\bigr)\,|J_i|\,d\vec{t}\ $. - Gradienten: Für $\varphi(\vec{x})=\tilde{\varphi}(\sigma_i^{-1}(\vec{x}))$ gilt $D\varphi(\vec{x}) = D\tilde{\varphi}(\sigma_i^{-1}(\vec{x}))\cdot A_i^{-1}$, also $\nabla\varphi(\vec{x}) = {A_i^{-1}}^T\cdot\nabla\tilde{\varphi}(\sigma_i^{-1}(\vec{x}))$. - Orientierung: $J_i>0$ genau dann, wenn $\vec{p}_{i,1},\dots,\vec{p}_{i,d+1}$ gleich orientiert sind wie die Referenzecken. - Die Referenzkoordinaten $\vec{t}$ sind eng mit den baryzentrischen Koordinaten von $\vec{x}$ verwandt, siehe [](#sec-baryzentrisch). (sec-baryzentrisch)= ## Baryzentrische Koordinaten **Definition.** Sei $T_i$ ein nicht entartetes Simplex mit Ecken $\vec{p}_{i,1},\dots,\vec{p}_{i,d+1}$. Jeder Punkt $\vec{x}\in\mathbb{R}^d$ lässt sich eindeutig als gewichtete Summe der Ecken schreiben, deren Gewichte sich zu eins summieren: ```{math} :label: eq-baryzentrisch \vec{x} = \sum_{k=1}^{d+1}\lambda_k\,\vec{p}_{i,k},\qquad \sum_{k=1}^{d+1}\lambda_k = 1 . ``` Die Zahlen $\lambda_1(\vec{x}),\dots,\lambda_{d+1}(\vec{x})$ heissen **baryzentrische Koordinaten** von $\vec{x}$ bezüglich $T_i$. **Zusammenhang mit den Referenzkoordinaten.** Setzt man $\lambda_1 = 1-\sum_{k=2}^{d+1}\lambda_k$ in {eq}`eq-baryzentrisch` ein, so erhält man $$ \vec{x}-\vec{p}_{i,1} = \sum_{k=2}^{d+1}\lambda_k\,\bigl(\vec{p}_{i,k}-\vec{p}_{i,1}\bigr) = A_i\cdot\begin{pmatrix}\lambda_2\\\vdots\\\lambda_{d+1}\end{pmatrix}. $$ Das ist genau die Gleichung $\vec{x}=\sigma_i(\vec{t}\,)$. Weil $A_i$ invertierbar ist, sind die $\lambda_k$ eindeutig bestimmt, und es gilt $$ \lambda_{k+1}(\vec{x}) = t_k \quad (k=1,\dots,d),\qquad \lambda_1(\vec{x}) = 1-\sum_{k=1}^{d} t_k, \qquad\text{mit}\quad \vec{t}=\sigma_i^{-1}(\vec{x}). $$ Im Dreieck also $\lambda_1 = 1-\xi-\eta$, $\lambda_2=\xi$, $\lambda_3=\eta$. **Eigenschaften.** - An den Ecken gilt $\lambda_k(\vec{p}_{i,j}) = \delta_{kj}$. Da die $\lambda_k$ affin in $\vec{x}$ sind, sind sie genau die linearen Basisfunktionen auf $T_i$: $\lambda_{k+1}(\vec{x}) = \tilde{\varphi}_k(\sigma_i^{-1}(\vec{x}))$ mit $\tilde{\varphi}_0 = 1-\xi-\eta$, $\tilde{\varphi}_1=\xi$, $\tilde{\varphi}_2=\eta$ (2D). - $\vec{x}\in T_i$ genau dann, wenn alle $\lambda_k(\vec{x})\ge 0$. Auf der Seite gegenüber von $\vec{p}_{i,k}$ ist $\lambda_k=0$. - Der Schwerpunkt von $T_i$ hat die baryzentrischen Koordinaten $\lambda_k=\frac{1}{d+1}$. **Geometrische Deutung.** Ersetzt man im Simplex $T_i$ die Ecke $\vec{p}_{i,k}$ durch den Punkt $\vec{x}$, so entsteht ein Teilsimplex $T_i^{(k)}(\vec{x})$. Es gilt $$ \lambda_k(\vec{x}) = \frac{|T_i^{(k)}(\vec{x})|}{|T_i|}, $$ wobei die Volumina orientiert zu nehmen sind (für $\vec{x}\notin T_i$ werden einzelne $\lambda_k$ negativ). Das folgt aus der Cramerschen Regel für $A_i\cdot\vec{t}=\vec{x}-\vec{p}_{i,1}$: Ersetzt man in $A_i$ eine Spalte durch $\vec{x}-\vec{p}_{i,1}$, so ist die Determinante dieser Matrix das $d!$-fache Volumen des entsprechenden Teilsimplex. Für $\lambda_1$ folgt die Aussage dann aus $\sum_k\lambda_k=1$, da sich die Teilvolumina zu $|T_i|$ summieren. Im Dreieck ist $\lambda_k$ also der Anteil des Teildreiecks, das der Ecke $\vec{p}_{i,k}$ gegenüberliegt, an der Fläche von $T_i$; im Tetraeder entsprechend der Volumenanteil eines Teiltetraeders. ```{figure} baryzentrisch_2d.svg :name: fig-baryzentrisch :width: 45% Baryzentrische Koordinaten im Dreieck: $\lambda_k(\vec{x})$ ist der Flächenanteil des Teildreiecks, das der Ecke $\vec{p}_{i,k}$ gegenüberliegt. Im Beispiel ist $\vec{x} = 0.2\,\vec{p}_{i,1} + 0.5\,\vec{p}_{i,2} + 0.3\,\vec{p}_{i,3}$. ``` ## 1D: Interval Mit $p_{i,1}=a$ und $p_{i,2}=b$ ist $A_i=J_i=b-a$ ein Skalar. Die Referenzkoordinate ist hier der Skalar $\vec{t}=\xi\in[0,1]$. $$ \begin{aligned} \text{Abbildung:}\quad & x = \sigma_i(\xi) = a + (b-a)\,\xi, \\ \text{Rücktransformation:}\quad & \xi = \sigma_i^{-1}(x) = \frac{x-a}{b-a}, \\ \text{Jacobi:}\quad & J_i=\frac{dx}{d\xi} = b-a,\qquad \frac{d\xi}{dx}=\frac{1}{b-a},\\ \text{Integral:}\quad & \int_a^b f(x)\,dx = (b-a)\int_0^1 f\bigl(a+(b-a)\xi\bigr)\,d\xi,\\ \text{Ableitung:}\quad & \varphi'(x) = \frac{1}{b-a}\,\tilde{\varphi}'\bigl(\sigma_i^{-1}(x)\bigr). \end{aligned} $$ **Beispiel.** $[a,b]=[1,4]$: $J_i=3$, $\xi=(x-1)/3$, $x = 1+3\xi$. Der Mittelpunkt $x=2.5$ wird auf $\xi = 0.5$ abgebildet. ```{figure} transformation_1d.svg :name: fig-trafo-1d :width: 80% Affine Abbildung zwischen $T_i=[a,b]$ und dem Einheitsintervall $T=[0,1]$. ``` ## 2D: Dreieck Seien $\vec{p}_{i,k}=(x_k,y_k)^T$ die Ecken von $T_i$. Das Referenzdreieck $T$ hat die Ecken $(0,0)$, $(1,0)$ und $(0,1)$. Ein Punkt in $T$ wird mit $\vec{t}=(\xi,\eta)^T$ bezeichnet. **Abbildung.** $$ \sigma_i(\vec{t}\,) = \sigma_i(\xi,\eta) = \vec{p}_{i,1} + A_i\cdot\begin{pmatrix}\xi\\\eta\end{pmatrix}, \qquad A_i = \begin{pmatrix} x_2-x_1 & x_3-x_1\\ y_2-y_1 & y_3-y_1\end{pmatrix}, $$ bzw. in baryzentrischer Form $\sigma_i(\xi,\eta) = (1-\xi-\eta)\,\vec{p}_{i,1} + \xi\,\vec{p}_{i,2} + \eta\,\vec{p}_{i,3}$. **Jacobi-Determinante.** $$ J_i = \det A_i = (x_2-x_1)(y_3-y_1)-(x_3-x_1)(y_2-y_1) = \pm 2\,|T_i| , \qquad dx\,dy = |J_i|\,d\xi\,d\eta . $$ $J_i$ ist genau dann positiv, wenn $\vec{p}_{i,1},\vec{p}_{i,2},\vec{p}_{i,3}$ im Gegenuhrzeigersinn angeordnet sind. **Inverse Jacobi-Matrix.** $$ A_i^{-1} = \frac{1}{J_i} \begin{pmatrix} y_3-y_1 & -(x_3-x_1)\\ -(y_2-y_1) & x_2-x_1\end{pmatrix}. $$ **Rücktransformation.** $$ \begin{pmatrix}\xi\\\eta\end{pmatrix} = \sigma_i^{-1}(x,y) = A_i^{-1}\cdot\left(\begin{pmatrix}x\\y\end{pmatrix}-\vec{p}_{i,1}\right), \qquad \begin{aligned} \xi &= \frac{(y_3-y_1)(x-x_1)-(x_3-x_1)(y-y_1)}{J_i},\\ \eta &= \frac{(x_2-x_1)(y-y_1)-(y_2-y_1)(x-x_1)}{J_i}. \end{aligned} $$ Nach [](#sec-baryzentrisch) ist $\xi=\lambda_2(\vec{x})$ und $\eta=\lambda_3(\vec{x})$: Der Zähler von $\xi$ ist das Doppelte der orientierten Fläche des Teildreiecks $(\vec{p}_{i,1},\vec{x},\vec{p}_{i,3})$, der Zähler von $\eta$ jenes von $(\vec{p}_{i,1},\vec{p}_{i,2},\vec{x})$. Geteilt durch $J_i = 2\,|T_i|$ ergeben sich die Flächenanteile. **Beispiel.** $\vec{p}_{i,1}=(1,1)^T$, $\vec{p}_{i,2}=(4,2)^T$, $\vec{p}_{i,3}=(2,4)^T$: $$ A_i=\begin{pmatrix}3&1\\1&3\end{pmatrix},\quad J_i = 8,\quad |T_i| = 4,\quad A_i^{-1}=\frac18\begin{pmatrix}3&-1\\-1&3\end{pmatrix}. $$ Probe: $\sigma_i^{-1}(\vec{p}_{i,2})=\tfrac18(3\cdot 3-1,\,-3+3)^T=(1,0)^T$ und $\sigma_i^{-1}(\vec{p}_{i,3})=(0,1)^T$. Der Schwerpunkt $(\tfrac73,\tfrac73)^T$ wird auf $(\tfrac13,\tfrac13)^T$ abgebildet. ```{figure} transformation_2d.svg :name: fig-trafo-2d :width: 95% Affine Abbildung zwischen einem beliebigen Dreieck $T_i$ und dem Einheitsdreieck $T$. Die grauen Linien $\xi=\text{const}$ und $\eta=\text{const}$ werden auf parallele Geraden abgebildet; die Spalten von $A_i$ sind die farbigen Kantenvektoren. ``` ## 3D: Tetraeder Das Referenztetraeder $T$ hat die Ecken $(0,0,0)$, $(1,0,0)$, $(0,1,0)$ und $(0,0,1)$. Ein Punkt in $T$ wird mit $\vec{t}=(\xi,\eta,\zeta)^T$ bezeichnet. Mit den Kantenvektoren $$ \vec{a} = \vec{p}_{i,2}-\vec{p}_{i,1},\qquad \vec{b} = \vec{p}_{i,3}-\vec{p}_{i,1},\qquad \vec{c} = \vec{p}_{i,4}-\vec{p}_{i,1} $$ gilt Folgendes. **Abbildung.** $$ \vec{x} = \sigma_i(\vec{t}\,) = \sigma_i(\xi,\eta,\zeta) = \vec{p}_{i,1} + A_i\cdot\begin{pmatrix}\xi\\\eta\\\zeta\end{pmatrix}, \qquad A_i = (\vec{a},\ \vec{b},\ \vec{c}), $$ bzw. $\sigma_i(\xi,\eta,\zeta) = (1-\xi-\eta-\zeta)\,\vec{p}_{i,1} + \xi\,\vec{p}_{i,2} + \eta\,\vec{p}_{i,3} + \zeta\,\vec{p}_{i,4}$. **Jacobi-Determinante (Spatprodukt).** $$ J_i = \det A_i = \vec{a}\cdot(\vec{b}\times\vec{c}) = \pm 6\,|T_i| , \qquad dx\,dy\,dz = |J_i|\,d\xi\,d\eta\,d\zeta . $$ $J_i$ ist genau dann positiv, wenn $(\vec{a},\vec{b},\vec{c})$ ein Rechtssystem bilden. **Inverse Jacobi-Matrix.** Die Zeilen von $A_i^{-1}$ sind (bis auf den Faktor $1/J_i$) Kreuzprodukte der Kantenvektoren: $$ A_i^{-1} = \frac{1}{J_i} \begin{pmatrix} (\vec{b}\times\vec{c})^T\\ (\vec{c}\times\vec{a})^T\\ (\vec{a}\times\vec{b})^T \end{pmatrix}. $$ Probe: $(\vec{b}\times\vec{c})\cdot\vec{a} = J_i$ und $(\vec{b}\times\vec{c})\cdot\vec{b} = (\vec{b}\times\vec{c})\cdot\vec{c}=0$. Geometrisch sind die Zeilen also Normalenvektoren der Seitenflächen durch $\vec{p}_{i,1}$. **Rücktransformation.** $$ \begin{pmatrix}\xi\\\eta\\\zeta\end{pmatrix} = \sigma_i^{-1}(\vec{x}) = A_i^{-1}\cdot\bigl(\vec{x}-\vec{p}_{i,1}\bigr), \qquad \begin{aligned} \xi &= \frac{(\vec{b}\times\vec{c})\cdot(\vec{x}-\vec{p}_{i,1})}{J_i},\\ \eta &= \frac{(\vec{c}\times\vec{a})\cdot(\vec{x}-\vec{p}_{i,1})}{J_i},\\ \zeta &= \frac{(\vec{a}\times\vec{b})\cdot(\vec{x}-\vec{p}_{i,1})}{J_i}. \end{aligned} $$ Nach [](#sec-baryzentrisch) ist $\xi=\lambda_2$, $\eta=\lambda_3$ und $\zeta=\lambda_4$. Zum Beispiel ist $(\vec{b}\times\vec{c})\cdot(\vec{x}-\vec{p}_{i,1}) = \det(\vec{x}-\vec{p}_{i,1},\ \vec{b},\ \vec{c})$ das Sechsfache des orientierten Volumens des Teiltetraeders $(\vec{p}_{i,1},\vec{x},\vec{p}_{i,3},\vec{p}_{i,4})$, in dem $\vec{p}_{i,2}$ durch $\vec{x}$ ersetzt ist. Geteilt durch $J_i = 6\,|T_i|$ ergibt sich der Volumenanteil. **Beispiel.** $\vec{p}_{i,1}=(1,1,0)^T$, $\vec{p}_{i,2}=(3,1,0)^T$, $\vec{p}_{i,3}=(1.5,3,0)^T$, $\vec{p}_{i,4}=(1.5,1.5,2)^T$: $$ A_i = \begin{pmatrix} 2 & 0.5 & 0.5\\ 0 & 2 & 0.5\\ 0 & 0 & 2\end{pmatrix},\qquad J_i = 8,\qquad |T_i| = \tfrac{8}{6} = \tfrac43, $$ $$ A_i^{-1} = \begin{pmatrix} \tfrac12 & -\tfrac18 & -\tfrac{3}{32}\\ 0 & \tfrac12 & -\tfrac18\\ 0 & 0 & \tfrac12\end{pmatrix} = \begin{pmatrix} 0.5 & -0.125 & -0.09375\\ 0 & 0.5 & -0.125\\ 0 & 0 & 0.5\end{pmatrix}. $$ Probe: $\sigma_i^{-1}(\vec{p}_{i,4}) = A_i^{-1}\cdot(0.5,\,0.5,\,2)^T = (0,0,1)^T$. Der Schwerpunkt $\tfrac14\sum_k\vec{p}_{i,k} = (1.75,\,1.625,\,0.5)^T$ wird auf $(\tfrac14,\tfrac14,\tfrac14)^T$ abgebildet. ```{figure} transformation_3d.svg :name: fig-trafo-3d :width: 95% Affine Abbildung zwischen einem beliebigen Tetraeder $T_i$ und dem Einheitstetraeder $T$ (verdeckte Kanten gestrichelt). Der Schwerpunkt wird auf den Schwerpunkt abgebildet. ``` ## Numerische Kontrolle Die Beispiele lassen sich für alle Dimensionen mit derselben Funktion nachrechnen: ```python import numpy as np from math import factorial def affine_map(pts): """Jacobi-Matrix A, sigma und sigma^{-1} für ein Simplex mit Ecken pts (Zeilen).""" pts = np.asarray(pts, dtype=float) p1 = pts[0] A = (pts[1:] - p1).T # Spalten = Kantenvektoren Ainv = np.linalg.inv(A) sigma = lambda t: p1 + A @ t # Referenz -> Element sigma_inv = lambda x: Ainv @ (x - p1) # Element -> Referenz return A, Ainv, sigma, sigma_inv for pts in ([[1], [4]], [[1, 1], [4, 2], [2, 4]], [[1, 1, 0], [3, 1, 0], [1.5, 3, 0], [1.5, 1.5, 2]]): A, Ainv, sigma, sigma_inv = affine_map(pts) d = len(pts) - 1 J = np.linalg.det(A) centroid = np.mean(pts, axis=0) print(f"d={d}: J = {J:.4g}, |T_i| = {abs(J)/factorial(d):.4g}, " f"sigma^-1(Schwerpunkt) = {sigma_inv(centroid)}") ``` ## Zusammenfassung | | 1D: Intervall | 2D: Dreieck | 3D: Tetraeder | |---|:---:|:---:|:---:| | $A_i$ | $b-a$ | $(\vec{p}_{i,2}-\vec{p}_{i,1},\ \vec{p}_{i,3}-\vec{p}_{i,1})$ | $(\vec{a},\ \vec{b},\ \vec{c})$ | | $J_i=\det A_i$ | $b-a$ | $\pm 2\vert T_i\vert$ | $\vec{a}\cdot(\vec{b}\times\vec{c})=\pm 6\vert T_i\vert$ | | $\sigma_i(\vec{t}\,)$ | $a+(b-a)\xi$ | $\vec{p}_{i,1}+A_i\cdot\vec{t}\,$ | $\vec{p}_{i,1}+A_i\cdot\vec{t}\,$ | | $\sigma_i^{-1}(\vec{x})$ | $\frac{x-a}{b-a}$ | $A_i^{-1}\cdot(\vec{x}-\vec{p}_{i,1})$ | $A_i^{-1}\cdot(\vec{x}-\vec{p}_{i,1})$ | | $\vert T\vert$ | $1$ | $\tfrac12$ | $\tfrac16$ | Für alle Fälle gilt $$ \begin{gathered} \int_{T_i} f(\vec{x})\,d\vec{x} = \int_{T} f\bigl(\sigma_i(\vec{t}\,)\bigr)\,|J_i|\,d\vec{t},\qquad \nabla\varphi(\vec{x}) = {A_i^{-1}}^T\cdot\nabla\tilde{\varphi}\bigl(\sigma_i^{-1}(\vec{x})\bigr),\\ \int_{T_i} \nabla\varphi_j(\vec{x})\cdot\nabla\varphi_k(\vec{x})\,d\vec{x} = \int_{T} \nabla\tilde{\varphi}_j(\vec{t}\,)^T\cdot C_i\cdot\nabla\tilde{\varphi}_k(\vec{t}\,)\;|J_i|\,d\vec{t}, \qquad C_i = A_i^{-1}\cdot{A_i^{-1}}^T . \end{gathered} $$ Da $\sigma_i$ affin ist, sind $A_i$, $J_i$ und $C_i$ auf ganz $T_i$ konstant. Deshalb genügt es, Basisfunktionen und Quadraturregeln nur einmal auf $T$ zu definieren.