\definecolor{dgreen}{rgb}{0,.8,.4}
\definecolor{dblue}{rgb}{0,.4,.4}
-\def\todo#1{\textcolor{dred}{#1}}
-\def\why#1{\textcolor{blue}{#1}}
-\def\Matlab{{\sc Matlab}}
+
+
\def\q{\Q}
\let\mod\relax
\DeclareMathOperator{\sgn}{sgn}
\DeclareMathOperator{\arsinh}{arsinh}
+\DeclareMathOperator{\diam}{diam}
+\DeclareMathOperator{\dist}{dist}
+\DeclareMathOperator{\len}{len}
+
\def\Ta{$T_a$}
\def\Tb{$T_b$}
\newcommand{\ul}[1]{\underline{#1}}
\newcommand{\dif}[2]{\frac{\partial#1}{\partial#2}}
-
-% \newcommand{\enorm}[1]{\left| \! \left| \! \left|#1\right| \! \right| \! \right|}
\newcommand{\enorm}[1]{| \! | \! |#1| \! | \! |}
+\newcommand{\Enorm}[1]{\left| \! \left| \! \left|#1\right| \! \right| \! \right|}
\newcommand{\norm}[1]{\lVert#1\rVert}
+\newcommand{\Norm}[1]{\left\lVert#1\right\rVert}
\newcommand{\abs}[1]{\lvert#1\rvert}
+\newcommand{\Abs}[1]{ \left\lvert#1\right\rvert}
+
+\newcommand{\todo}[1]{\textcolor{dred}{#1}}
+\newcommand{\why}[1]{\textcolor{dblue}{#1}}
+\newcommand{\Matlab}{{\sc Matlab}}
\newcommand{\showMesh}[2][]{\begin{figure}[ht]
\caption{#1}
\newtheorem{sat}[defi]{Satz}
\newtheorem{bew}[defi]{Beweis}
\newtheorem{alg}[defi]{Algorithmus}
+\newcommand{\beweis}{{\it Beweis. }}
% \numberwithin{lem}{section}
\numberwithin{defi}{section}
\section{Einleitung}
\todo{
\begin{itemize}
- \item kurz + praegnant, worum es geht (vgl. Mayr-Bakk, wo das exakt
-eine Seite ist)
+\item kurz + praegnant, worum es geht (vgl. Mayr-Bakk, wo das exakt eine Seite ist)
\end{itemize}
}
\noindent
In dieser Arbeit beschäftigen wir uns mit der Randelementmethode für die homogene Laplace-Gleichung mit Dirichlet-Randbedingungen
\begin{align}
- - \varDelta u &= 0 &&\text{ in } \Omega \subset \R^3,&&\\
- u &= g &&\text{ auf }\Gamma,&&
+- \varDelta u &= 0 &&\text{ in } \Omega \subset \R^3,&&\\
+u &= g &&\text{ auf }\Gamma,&&
\end{align}
wobei $\varDelta u := \partial_x^2u+\partial_y^2u$ den Laplace-Operator bezeichnet und $\Omega \subset \R^3$ eine beschränkte Teilmenge von $\R^3$ mit Lipschitz-Rand $\Gamma := \partial \Omega$ ist.\\
\section{Randelementmethode}
\todo{
\begin{itemize}
- \item kuerzer als in Ferraz-DA \cite{fer:errbem}
- \item insbesondere fastregulaere Triangulierung + K-Gitter (vgl. Ferraz-DA \cite{fer:errbem})
+\item kuerzer als in Ferraz-DA \cite{fer:errbem}
+\item insbesondere fastregulaere Triangulierung + K-Gitter (vgl. Ferraz-DA \cite{fer:errbem})
\end{itemize}
}
\subsection{Netze}
Für die Diskretisierung des Problems wollen wir nun einige Begriffe definieren.
-\begin{defi}[Rechteck] Sei $\bs x \in \R^3$, $\bs a,\bs b \in \{(1,0,0)^T,(0,1,0)^T,(0,0,1)^T\}$ mit $\bs a\neq \bs b$ und $a,b \in \R$ mit $a,b > 0$. Dann heißt
-\begin{align}
- T := \{\bs x + \lambda_1 a {\bs a} + \lambda_2 b {\bs b} ~|~ \lambda_1,\lambda_2 \in[0,1]\}
-\end{align}
+\begin{defi}\label{math:def:T} Sei $\bs v \in \R^3$, $\bs a,\bs b \in \{(1,0,0)^T,(0,1,0)^T,(0,0,1)^T\}$ mit $\bs a\neq \bs b$ und $a,b \in \R$ mit $a,b > 0$. Dann heißt
+\begin{align*}
+T := \{\bs v + \lambda_1 a {\bs a} + \lambda_2 b {\bs b} ~|~ \lambda_1,\lambda_2 \in[0,1]\}
+\end{align*}
achsenorientiertes Rechteck.
\end{defi}
\begin{defi}
- Sei $T$ ein Rechteck, Dann heißt
- \begin{align}
- \gamma := [0,1]^2 \to T: \lambda_1,\lambda_2 \mapsto \bs x + \lambda_1 a {\bs a} + \lambda_2 b {\bs b}
- \end{align}
+Sei $T$ ein Rechteck, Dann heißt
+\begin{align*}
+ \gamma := [0,1]^2 \to T: \lambda_1,\lambda_2 \mapsto \bs v + \lambda_1 a {\bs a} + \lambda_2 b {\bs b}
+\end{align*}
die zu $T$ zugehörige Parametrisierung.
\end{defi}
-
-
-\begin{defi}[hängender Knoten]
- Ist $c \in T_j \cap T_k$ ein Eckpunkt von $T_j$ aber nicht von $T_k$, so nennen wir $c$ einen hängenden Knoten.
+\begin{defi}Sei $a,b\in\R$ für $T$ definiert wie in Def. \ref{math:def:T}, dann heißt
+\begin{align*}
+ \diam (T) &= (a^2+b^2)^{1/2}
+\end{align*}
+Durchmesser von T. Weiterhin nennen wir
+\begin{align*}
+ \diam_{|\bs a} (T) = a
+\end{align*}
+die Seitenlänge des Rechtecks $T$ in Richtung $\bs a$. Die Seitenlänge in Richtung $\bs b$ definieren wir analog.
\end{defi}
-
-
-\begin{defi}[Partition]
- Wir nennen $\T_{\ell} = \{T_1,T_2,\dots,T_N\}$ eine Partition von $\Gamma$ falls:
+\begin{defi}Der Abstand zweier Rechtecke $T_j$ und $T_k$ sei definiert durch
+\begin{align*}
+ \dist(T_j,T_k) &= \inf \{\abs{\bs x-\bs y} ~|~ \bs x \in T_j,\bs y \in T_k\}.
+\end{align*}
+\end{defi}
+\todo{Unterscheiden zwischen parallel und orthogonal und dafür $\dist$ definieren}
+\begin{defi}
+Ist $c \in T_j \cap T_k$ ein Eckpunkt von $T_j$ aber nicht von $T_k$, so nennen wir $c$ einen hängenden Knoten.
+\end{defi}
+\begin{defi}
+Wir nennen $\T_{\ell} = \{T_1,T_2,\dots,T_N\}$ eine Partition von $\Gamma$ falls:
\begin{itemize}
- \item $\overline{\Gamma} = \bigcup_{j=1}^NT_j$
- \item alle Elemente aus $T_{\ell}$ sind abgeschlossene achsenorientierte Rechtecke
- \item $\abs{T_j \cap T_k}=0$ für den Schnitt zweier Elemente $T_j,T_k\in\T_{\ell}$ mit $T_j\neq T_k$
+\item $\overline{\Gamma} = \bigcup_{j=1}^NT_j$
+\item alle Elemente aus $T_{\ell}$ sind abgeschlossene achsenorientierte Rechtecke
+\item $\abs{T_j \cap T_k}=0$ für den Schnitt zweier Elemente $T_j,T_k\in\T_{\ell}$ mit $T_j\neq T_k$
\end{itemize}
Hier bezeichne $\abs{\cdot}$ das 2-Dimensionale Oberflächenmaß für Rechtecke. Weiterhin sei der Schnitt $T_j \cap T_k$ zweier Elemente $T_j,T_k\in\T_{\ell}$ mit $T_j\neq T_k$
- \begin{itemize}
+\begin{itemize}
\item leer,
\item Knoten von $T_j$ und $T_k$,
\item Kante von $T_j$ und $T_k$,
\item Knoten von $T_j$ und $T_k$ und o.B.d.A. Kante von $T_j$, also nur einen hängenden Knoten pro Kante.
- \end{itemize}
+\end{itemize}
\end{defi}
\begin{figure}[ht]
\centering
\label{fig:net}
\end{figure}
\begin{defi}
- Außerdem wollen wir zu jedem Netz $\T$ die Menge der Seiten $\S$ definieren, wobei zu jeder Seite $e_T \in \S$ eine bezüglich des Elements $T$ gegenüberliegende Seite $\tilde e_T$ existiert.
+Außerdem wollen wir zu jedem Netz $\T$ die Menge der Seiten $\S$ definieren, wobei zu jeder Seite $e_T \in \S$ eine bezüglich des Elements $T$ gegenüberliegende Seite $\tilde e_T$ existiert.
\end{defi}
-
-\todo{
-\subsection{Abstand und Größe der Elemente}
-\begin{defi}[dist]
- \end{defi}
-\begin{defi}[diam]
- \end{defi}
-\begin{defi}[hmin]
- \end{defi}
-}
-
\subsection{Verfeinern}
\begin{defi}[Lokale Verfeinerung]
Ein Element $T \in \T$ wird isotrop in vier Elemente $T_1,\ldots,T_4$ geteilt, wenn $T = \bigcup_{n=1}^4 T_n$ gilt und die Seitenlängen der Elemente $T_1,\ldots,T_4$ gleich groß sind. Weiterhin wird ein Element $T \in \T$ anisotrop in zwei Elemente $T_1,T_2$ geteilt, wenn ebenfalls $T = T_1 \cup T_2$ gilt und $T_1,T_2$ gleich große Seitenlängen haben. Hierbei kann $T$ entweder horizontal oder vertikal geteilt werden, wie in Abb.~\ref{fig:refType} gezeigt ist.
\end{figure}
\begin{defi}
- Für eine Partition $\T_{\ell}$ bezeichnen wir mit $\widehat \T_{\ell}$ jene Partition die entsteht, wenn alle Elemente isotrop verfeinert werden.
+Für eine Partition $\T_{\ell}$ bezeichnen wir mit $\widehat \T_{\ell}$ jene Partition die entsteht, wenn alle Elemente isotrop verfeinert werden.
\end{defi}
\begin{alg}[Verfeinern] \label{alg:refine} Sei $\T_{\ell}$ eine Partition und $\S_{\ell}^{(0)}$ eine Menge markierter Kanten, wobei zu jedem $e_T \in \S_{\ell}^{(0)}$ auch $\tilde e_T \in \S_{\ell}^{(0)}$ ist. Nun sei $i=0$ und gehe so vor:
- \begin{enumerate}
+\begin{enumerate}
\renewcommand{\theenumi}{(\roman{enumi})}
\item \label{alg:refine:first}
Definiere $\S_{\ell}^{(i+1)} := \{ e_T, \tilde e_T ~|~ T \in \T_{\ell}$ mit $\exists$ hängender-Knoten auf $T$ und $\exists T' \in\T_{\ell}$ sodass $T \cap T' = e_{T'},\in \S_{\ell}^{(i)} \}$
\item
Falls $\S_{\ell}^{(i)} \subsetneq \S_{\ell}^{(i+1)}$, $i = i+1$ und gehe zu Schritt \ref{alg:refine:first}
\item Teile alle Elemente aus $\T_{\ell}$ bezüglich der markierten Kanten $\S_{\ell}^{(0,1,\ldots,i)}$, wie in der lokalen Verfeinerung vorgegeben
- \end{enumerate}
+\end{enumerate}
\end{alg}
\clearpage
\section{Analytische und Semi-analytische Berechnung}
\todo{
\begin{itemize}
- \item Modellproblem $A_{kj} = \int_{Tj}\int_{Tk} \kappa(x,y) ds_y ds_x$
-mit $\kappa(.,.)$ asymptotisch glatt und $Tj, Tk$ achsenorientierte
+\item Modellproblem $A_{kj} = \int_{Tj}\int_{Tk} \kappa(x,y) ds_y ds_x$
+mit $\kappa(.,.)$ asymptotisch glatt und $Tj, Tk$ achsenorientiertes
Rechtecke.
\item Das ist das Kernstueck der Arbeit (Uebertragung der Mayr-Bakk auf 3D).
\item Warum erwarten man Instabilitaet bei analytischer Berechnung?
\begin{align}\label{math:gal:frac}
A_{jk} &=\frac{1}{4\pi} \int_{T_j} \int_{T_k} \frac{1}{|\bs x- \bs y|} ds_{\bs y} ds_{\bs x}.
\end{align}
-unter bestimmten Vorraussetzungen an die affinen Randstücke $T_j,T_k$ und den asymptotisch glatten Integranden $\kappa : \R^3 \times \R^3 \to \R$.
+unter bestimmten Voraussetzungen an die affinen Randstücke $T_j,T_k$ und den asymptotisch glatten Integranden $\kappa : \R^3 \times \R^3 \to \R$.
\subsection{Interpolation}
Zunächst wollen wir die Chebyshev-Interpolation auf dem Intervall [0,1] definieren. Hierzu sei $\P_p$ die Menge aller Polynome vom Grad $p$.
\begin{defi}
- Seien $x_0,\ldots,x_p \in [0,1]$ paarweise verschiedene Knoten und $y_0,\ldots,y_p$ die zugehörigen Funktionswerte. Die Lagrange'sche Interpolationsaufgabe lautet dann folgendermaßen: Finde ein Polynom $q$ vom Grad $p-1$, sodass
- \begin{align*}
+Seien $x_0,\ldots,x_p \in [0,1]$ paarweise verschiedene Knoten und $y_0,\ldots,y_p$ die zugehörigen Funktionswerte. Die Lagrange'sche Interpolationsaufgabe lautet dann folgendermaßen: Finde ein Polynom $q$ vom Grad $p-1$, sodass
+\begin{align*}
q(x_i) &=y_i \quad \text{für alle} i=0,\ldots,p.
- \end{align*}
+\end{align*}
\end{defi}
Ferner definieren wir die zu den Knoten $x_0,\ldots,x_p$ zugehörigen Lagrange-Polynome durch
\begin{align}
- L_j(x) &= \prod_{i=1,i\neq j}^p \frac{x-x_i}{x_j-xi}.
+L_j(x) &= \prod_{i=1,i\neq j}^p \frac{x-x_i}{x_j-xi}.
\end{align}
\todo{cite} zeigt, dass die Lagrange'sche Interpolationsaufgabe eine eindeutige Lösung besitzt, die gegeben ist durch
\begin{align}
- q&=\sum_{j=0}^py_jL_j.
+q&=\sum_{j=0}^py_jL_j.
\end{align}
Wir definieren den Interpolationsoperator $\I_pu$ durch
\begin{align}
- \I_pu &:= \sum_{j=0}^p u(x_j)L_j.
+\I_pu &:= \sum_{j=0}^p u(x_j)L_j.
\end{align}
Außerdem definieren wir die Lebesgue-Konstante durch
\begin{align}
- \Lambda_p &:= \max_{x\in[0,1]} \sum_{j=0}^p\abs{L_j(x)}.
+\Lambda_p &:= \max_{x\in[0,1]} \sum_{j=0}^p\abs{L_j(x)}.
\end{align}
\todo{cite} stellt einen Bezug zwischen der Bestapproximation und der Polynominterpolation her und ist später von Bedeutung.
-\begin{sat}
- Sei $\Lambda_p$ die Lebesgue-Konstante zu $p$ paarweise verschiedener Knoten $x_1,\ldots,x_p$, die alle im Intervall $[0,1]$ liegen. Dann gilt für jedes $u \in \C[0,1]$
- \begin{align}
+\begin{sat}\label{math:ipol}
+Sei $\Lambda_p$ die Lebesgue-Konstante zu $p$ paarweise verschiedener Knoten $x_1,\ldots,x_p$, die alle im Intervall $[0,1]$ liegen. Dann gilt für jedes $u \in \C[0,1]$
+\begin{align}
\norm{\I_pu-u}_{\infty} &\leq (1 + \Lambda_p) \min_{v\in \P_p} \norm{u-v}_{\infty},
- \end{align}
+\end{align}
wobei $\P_p$ die Menge aller Polynome vom Grad $p$ bezeichnet.
\hfill $\square$
\end{sat}
\Q_n(f) &:= \sum_{k=0}^n w_kp(x_k)
\end{align}
gilt, wobei $p$ das Interpolationspolynom von $f$ vom Grad $n$ zu den Knoten $x_0,\ldots,x_n$ ist.
+\todo{\subsection{Bedingung}
+}
\begin{align}
- \Q(f) := \sum_{k=0}^n w_kf(x_k) \approx \int_0^s f(x) dx
-\end{align}
-
-\subsection{Bedingung}
-\begin{align}
- dist (T, \tilde T)&\geq \zeta \min\{ diam(T) , diam(\tilde T)\}\\
- dist (I_1, J_1) &\geq \zeta \max\{ len(I_2) , len(J_2)\}
+ \dist (T_j, T_k)&\geq \zeta \min\{ \diam(T_j) , \diam(T_k)\}\\
+ \dist (T_j,T_k) &\geq \zeta_1 \max\{ \diam_{|\bs a}(T_j) , \diam_{|\bs a}(T_k)\}
\end{align}
\subsection{Semianalytisch}
+% Das Integral über $T_j$ kann durch die analytische Doppelintegral Funktion ersetzt werden. Die Quadratur über $T_k$ wird nun mittels einer doppelten Gauss-Quadratur berechnet.
\begin{align}
- \int_{T} \int_{\tilde T} \approx \q_{T} \int_{\tilde T}
-\end{align}
-Das Integral über \Tb~ kann durch die analytische Doppelintegral Funktion ersetzt werden. Die Quadratur über \Ta~ wird nun mittels einer doppelten Gauss-Quadratur berechnet.
-\begin{align}
- \int_{T} \int_{\tilde T} \approx \q_{I_1} \q_{J_1} \int_{I_2}\int_{J_2}
-\end{align}
-\begin{align}
- \int_{T} \int_{\tilde T} \approx \q_{I_1} \int_{I_2} \int_{J_1} \int_{J_2}
+ \int_{T_j} \int_{T_k} &\approx \q_{T_j} \int_{T_k}\\
+ \int_{T_j} \int_{T_k} &\approx \q_{I_1} \q_{J_1} \int_{I_2}\int_{J_2}\\
+ \int_{T_j} \int_{T_k} &\approx \q_{I_1} \int_{I_2} \int_{J_1} \int_{J_2}\\
+ \int_{T_j} \int_{T_k} &\approx \q_{T_j} \q_{J_1} \int_{J_2}
\end{align}
-
-
\begin{defi}
- Die Kern-Funktion $\kappa(\bs x,\bs y)$ heißt asymptotisch glatt, falls sie glatt ist für $\bs x \neq \bs y$ und Konstanten $c_1,c_2>0$ und eine Ordnung der Singularität $s \in \R$ existieren, sodass
- \begin{align}
- \abs{\partial_{\bs x}^{\alpha}\partial_{\bs y}^{\beta}\kappa(\bs x, \bs y)} &\leq c_1(c_2(\abs{\bs x - \bs y})^{-(\abs{\alpha}-\abs{\beta}+s)}(\alpha + \beta)!
- \end{align}
-für alle Multiindizes $\alpha,\beta \in \N_0^d$ mit $\abs{\alpha} + \abs{\beta}\geq 1$ gilt.
+ Die Kern-Funktion $\kappa(\bs x,\bs y)$ heißt asymptotisch glatt, falls sie glatt ist für $\bs x \neq \bs y$ und Konstanten $c_1,c_2>0$ und eine Ordnung der Singularität $s \in \R$ existieren, sodass
+% \begin{align}
+% \abs{\partial_{\bs x}^{\alpha}\partial_{\bs y}^{\beta}\kappa(\bs x, \bs y)}
+% &\leq c_1(c_2(\abs{\bs x - \bs y})^{-(\abs{\alpha}-\abs{\beta}+s)}(\alpha + \beta)!
+% \end{align}
+ \begin{align}\label{math:sem:glatt}
+ \abs{\partial_{x_a}^{\alpha_a}\partial_{x_b}^{\alpha_b}\partial_{y_a}^{\beta_a}\partial_{y_b}^{\beta_b}\kappa(\bs x, \bs y)}
+ &\leq c_1(c_2(\abs{\bs x - \bs y})^{-(\abs{\alpha_a}+\abs{\alpha_b}+\abs{\beta_a}+\abs{\beta_b}+s)}(\alpha_a + \alpha_b + \beta_a + \beta_b)!
+ \end{align}
+ für alle Multiindizes $\alpha_a, \alpha_b, \beta_a, \beta_b \in \N_0^d$ mit $\abs{\alpha_a}+\abs{\alpha_b}+\abs{\beta_a}+\abs{\beta_b}\geq 1$ gilt.
\end{defi}
Wie wir im Folgenden sehen werden, lassen sich asymptotisch glatte Kernfunktionen besonders gut durch Polynome interpolieren. Hierzu benötigen wir das Lemma \todo{cite}:
\begin{lem}
- Sei $J \subseteq \R$ ein abgeschlossenes und beschränktes Intervall und angenommen die Funktion $u$ erfüllt für Konstanten $C_u,\gamma_u > 0$
- \begin{align}
- \norm{u^{(n)}}_{\infty,J} &\leq C_u \gamma_u^nn!\text{~~~für alle} n \in \N_0.
- \end{align}
+Sei $J \subseteq \R$ ein abgeschlossenes und beschränktes Intervall und angenommen die Funktion $u$ erfüllt für Konstanten $C_u,\rho_u > 0$
+\begin{align}\label{math:sem:ipol}
+ \norm{u^{(n)}}_{\infty,J} &\leq C_u \rho_u^nn!\text{~~~für alle} n \in \N_0.
+\end{align}
Dann gilt für alle $k\in \N_0$
\begin{align}
- \min_{v \in \P_p} \norm{u-v}_{\infty,J} &\leq C_u 4e(1+\gamma_u\abs{J})(p+1)\left(1+\frac{2}{\gamma_u\abs{J}}\right)^{-(p+1)}.
+\min_{v \in \P_p} \norm{u-v}_{\infty,J} &\leq C_u 4e(1+\rho_u\abs{J})(p+1)\left(1+\frac{2}{\rho_u\abs{J}}\right)^{-(p+1)}.
\end{align}
\end{lem}
\hfill$\square$
+\begin{lem}
+ Seien $T_j,T_k$ zwei $\zeta_1$-zulässige, affine Randstücke mit $\diam(T_j) \leq \diam(T_k)$. Weiterhin sei $\kappa$ eine asymptotisch glatte Kernfunktion mit den Konstanten $c_1,c_2>0$ und Singularitätsordung $s\geq0$.
+ Dann gilt mit
+ \begin{align}
+ C_{\zeta,j,\kappa} &:= 4 e\frac{c_1}{(c_2 \todo{\dist(T_j,T_k)})^s}\left(1+\frac{\zeta}{c_2}\right)
+ \end{align}
+ die Abschätzung
+ \begin{align*}
+ \sup_{\bs y \in T_k,x_b} \norm{\kappa(\cdot,\bs y)-\I_p\kappa(\cdot,\bs y)}_{\infty,T_{j| x_{\bs b}}}
+ &\leq C_{\zeta,j,\kappa}(1+\Lambda_p)(p+1)\left( 1+ \frac{2c_1}{\zeta} \right)^{-(p+1)}.
+\end{align*}
+\end{lem}
+\beweis Wir definieren die Konstanten
+\begin{align}
+ C_{\kappa} &:= \frac{c_1}{(c_2 \todo{\dist(T_j,T_k)})^s}\geq 0 \quad \text{und} \quad
+ \rho_{\kappa} := \frac{\zeta}{c_2}\geq 0
+\end{align}
+und können die Konstante $C_{\zeta,j,\kappa}$ dann schreiben als:
+\begin{align}
+ C_{\zeta,j,\kappa} &= C_{\kappa} 4 e(1+\rho_{\kappa}).
+\end{align}
+Sei $\bs y \in T_k$ fest gewählt. Ist $\gamma_{j | x_{\bs b}}$ die Parametrisierung von $T_{j| x_{\bs b}}$, so gilt aufgrund $\abs{\gamma_{j|x_{\bs b}}'}=\abs{a_j\bs a_j}=\diam_{|\bs a}(T_j)$ und der Kettenregel
+\begin{align*}
+\abs{\partial_{x_{\bs a}}^n\kappa}
+&= \abs{\partial_{x_{\bs a}}^n(\kappa(\gamma_{j | x_{\bs b}}(x_{\bs a}),\bs y))}
+= (\diam_{|\bs a}(T_j))^n\abs{(\partial_{x_{\bs a}}^n\kappa)(\gamma_{j | x_{\bs b}}(x_{\bs a}),\bs y)}.
+\end{align*}
+Mithilfe von \eqref{math:sem:glatt} und den Konstanten $C_{\kappa},\rho_k$ gilt die Abschätzung\\
+\todo{Ist $\abs{\gamma_{j | x_{\bs b}}(x_{\bs a})-\bs y}$ wirklich $\abs{\bs x -\bs y} \Rightarrow \dist(T_j,T_k)$}
+\begin{align*}
+ (\diam_{|\bs a}(T_j))^n&\abs{(\partial_{x_{\bs a}}^n\kappa)(\gamma_{j | x_{\bs b}}(x_{\bs a}),\bs y)}\\
+ &\leq \diam_{|\bs a}(T_j))^n c_1(c_2 \abs{\gamma_{j | x_{\bs b}}(x_{\bs a})-\bs y})^{-n-s}n!\\
+ &= \diam_{|\bs a}(T_j))^n c_1(c_2 \todo{\abs{\bs x-\bs y}})^{-n-s}n!\\
+ &= \frac{c_1}{(c_2 \abs{\bs x-\bs y})^s} \left(\frac{\diam_{|\bs a}(T_j))}{c_2 \abs{\bs x-\bs y}} \right)^nn!\\
+ &\leq C_{\kappa}\rho_{\kappa}^nn!
+\end{align*}
+Da $T_j$ und $T_k$ zulässig sind ist $\kappa$ auf $T_j \times T_k$ glatt. Daher sind die Voraussetzungen für \eqref{math:sem:ipol} erfüllt und es gilt
+\begin{align*}
+ \min_{v \in \P^p} \norm{\kappa(\gamma_{j | x_{\bs b}}(\cdot),\bs y)-v}_{\infty,[0,1]}
+ &\leq C_{\kappa} 4 e(1+\rho_{\kappa})(p+1)\left( 1+ \frac{2}{\rho_{\kappa}} \right)^{-(p+1)}\\
+ &= C_{\zeta,j,\kappa}(p+1)\left( 1+ \frac{2}{\rho_{\kappa}} \right)^{-(p+1)}
+\end{align*}
+für alle $p \in \N_0$, wobei $\P^p$ die Menge aller Polynome vom Grad $p$ ist. Daraus folgt also
+\begin{align*}
+ \min_{v \in \P^p} \norm{\kappa(\gamma_{j | x_{\bs b}}(\cdot),\bs y)-v}_{\infty,[0,1]}
+ &\leq C_{\zeta,j,\kappa}(p+1)\left( 1+ \frac{2c_1}{\zeta} \right)^{-(p+1)}
+\end{align*}
+und mithilfe von Satz \ref{math:ipol} gilt
+\begin{align*}
+ \norm{\kappa(\cdot,\bs y)-\I_p\kappa(\cdot,\bs y)}_{\infty,T_{j| x_{\bs b}}}
+ &\leq C_{\zeta,j,\kappa}(1+\Lambda_p)(p+1)\left( 1+ \frac{2c_1}{\zeta} \right)^{-(p+1)}.
+\end{align*}
+Da $\bs y$ und $x_{\bs b}$ beliebig war und die rechte Seite unabhängig von $\bs y,x_{\bs b}$ ist, folgt die Behauptung für das Supremum.\hfill$\square$
+\begin{sat}
+Seien $Tj,Tk$ zwei $\zeta_1$-zulässige affine Randstücke mit zugehörigen Parametrisierungen $\gamma_j,\gamma_k$ und $\diam(T_j) \leq \diam(T_k)$. Sei weiterhin $\kappa : \R^3\times\R^3\to \R : (\bs x, \bs y) \mapsto \kappa(\bs x,\bs y)$ eine asymptotisch glatte Kernfunktion mit Konstanten $c_1,c_2,s$. Sei
+\begin{align}
+ \tilde C_{\zeta,j,\kappa}.
+\end{align}
+Dann gilt für das Integral
+\begin{align}
+ A_{jk}
+ &= \int_{T_j} \int_{T_k} \kappa(\bs x,\bs y) ds_{\bs y} ds_{\bs x}\\
+ &=\diam(T_j) \int_0^1 \int_{T_{j | x_{\bs a}}} \int_{T_k} \kappa(\gamma_{j | x_{\bs b}}(x_{\bs a}),\bs y) ds_{\bs y} ds_{x_{\bs b}}
+\end{align}
+und für den durch die Gauss-Quadratur von $A_{jk}$ zum Grad $p$ entstehenden Term
+\begin{align}
+ (A_p)_{jk}
+ &=\diam(T_j) \sum_{i=0}^p w_i \int_{T_{j | x_{\bs a}}} \int_{T_k} \kappa(\gamma_{j | x_{\bs b}}(x_{\bs a,i}),\bs y) ds_{\bs y} ds_{x_{\bs b}}
+\end{align}
+die Abschätzung
+\begin{align}
+.
+\end{align}
+\end{sat}
+\beweis Da die Gauss-Quadratur ist interpolatorisch mit dem Exaktheitsgrad $2p+1$ ist gilt:
+\begin{align*}
+ (A_p)_{jk}
+ &=\diam(T_j) \sum_{i=0}^p w_i \int_{T_{j | x_{\bs a}}} \int_{T_k} \kappa(\gamma_{j | x_{\bs b}}(x_{\bs a,i}),\bs y) ds_{\bs y} ds_{x_{\bs b}}\\
+ &=\diam(T_j) \int_0^1 \int_{T_{j | x_{\bs a}}} \int_{T_k} \I_{2p+1} \kappa(\gamma_{j | x_{\bs b}}(x_{\bs a}),\bs y),
+\end{align*}
+wobei $\I_{2p+1}$ den Chebyshev-Interpolationsoperator vom Grad $2p+1$ in $x_{\bs a}$ bezeichnet. Wegen der Additivität des Integrals und durch Hineinziehen des Betrags erhält man
+\begin{align*}
+ \abs{A_{jk} &- (A_p)_{jk}}\\
+ &=\diam(T_j) \Abs{ \int_0^1 \int_{T_{j | x_{\bs a}}} \int_{T_k}
+ \kappa(\gamma_{j | x_{\bs b}}(x_{\bs a}),\bs y)- \I_{2p+1} \kappa(\gamma_{j | x_{\bs b}}(x_{\bs a}),\bs y)ds_{\bs y} ds_{x_{\bs b}}}\\
+ &\leq \diam(T_j) \int_0^1 \int_{T_{j | x_{\bs a}}} \int_{T_k}
+ \abs{\kappa(\gamma_{j | x_{\bs b}}(x_{\bs a}),\bs y)- \I_{2p+1} \kappa(\gamma_{j | x_{\bs b}}(x_{\bs a}),\bs y)} ds_{\bs y} ds_{x_{\bs b}}\\
+ &=\why{\sup_{\bs y \in T_k,x_b} \norm{\kappa(\gamma_{j | x_{\bs b}}(\cdot),\bs y)-\I_{2p+1}\kappa(\gamma_{j | x_{\bs b}}(\cdot),\bs y)}_{\infty,[0,1]}}
+\end{align*}
+
+\todo{\hfill$\square$\\}
+
+\begin{lem}
+ Seien $Tj,T_k$ zwei Randstücke, $\kappa$ eine asymptotisch glatte Kernfunktion und
+ \begin{align}
+ A_{jk}
+ &= \int_{T_j} \int_{T_k} \kappa(\bs x,\bs y) ds_{\bs y} ds_{\bs x}.
+ \end{align}
+ Dann gilt
+ \begin{align}
+ A_{jk} = A_{k_j}.
+ \end{align}
+\todo{\beweis
+\hfill$\square$}
+
+\end{lem}
-\section{Analytische Berechnung vom Integral im Fall des
-Einfachschichtpotentials}
+
+\clearpage
+
+\section{Analytische Berechnung}
\todo{
\begin{itemize}
- \item Zusammenfassung des Maischak-Papers \cite{mai:3dbem}
- \item Ergebnisse ohne Beweise
- \item Betonen, wo Sie im Maischak-Paper einen Fehler gefunden haben
+\item Zusammenfassung des Maischak-Papers \cite{mai:3dbem}
+\item Ergebnisse ohne Beweise
+\item Betonen, wo Sie im Maischak-Paper einen Fehler gefunden haben
\end{itemize}
}
\noindent
\begin{lem}
Für das Integral
\begin{align*}
- g(p;y;x;\lambda) &:= \int \{(x-y)^2 + \lambda^2 \}^p dy
+g(p;y;x;\lambda) &:= \int \{(x-y)^2 + \lambda^2 \}^p dy
\end{align*}
mit $\lambda = 0$ gilt die Formel mit beliebigen Skalaren $x,y,p \in \R$
\begin{align*}
- (2p+1)g(p;y;x;0) = \begin{cases} \sgn(y-x) \log \abs{y-x}^{2p} & p=-1/2 \\(y-x) \abs{y-x}^{2p} & sonst \end{cases}
+(2p+1)g(p;y;x;0) = \begin{cases} \sgn(y-x) \log \abs{y-x}^{2p} & p=-1/2 \\(y-x) \abs{y-x}^{2p} & sonst \end{cases}
\end{align*}
falls für $p\leq-/2$ $x$ außerhalb des betrachteten Integrationsbereichs liegt.
Für beliebige $x,p,y, \lambda \in \R$ mit $\lambda \neq 0$ gilt die Rekursionsformel:
\begin{align*}
- (2p+1)g(p;y;x;\lambda) &= (y-x) \{(y-x)^2+\lambda^2\}^p + 2p\lambda g(p-1;y;x;\lambda)
+(2p+1)g(p;y;x;\lambda) &= (y-x) \{(y-x)^2+\lambda^2\}^p + 2p\lambda g(p-1;y;x;\lambda)
\end{align*}
Im weiteren werden wir noch die Formeln für $p \in \{1/2 , 0 , -1/2 , -1, -3/2\}$
\begin{align*}
- g(1/2;y;x;\lambda) &= \frac{y-x}{2}\{(y-x)^2+\lambda^2\}^{1/2} + \frac{\lambda^2}{2} \arsinh \frac{y-x}{|\lambda|}\\
- g(0;y;x;\lambda) &= y-x\\
- g(-1/2;y;x;\lambda) &= \arsinh \frac{y-x}{|\lambda|}\\
- g(-1;y;x;\lambda) &= \frac 1 {\abs{\lambda}} \arctan \frac{y-x}{|\lambda|}\\
- g(-3/2;y;x;\lambda) &= \frac{y-x}{\lambda^2}\{(y-x)^2+\lambda^2\}^{-1/2}
+g(1/2;y;x;\lambda) &= \frac{y-x}{2}\{(y-x)^2+\lambda^2\}^{1/2} + \frac{\lambda^2}{2} \arsinh \frac{y-x}{|\lambda|}\\
+g(0;y;x;\lambda) &= y-x\\
+g(-1/2;y;x;\lambda) &= \arsinh \frac{y-x}{|\lambda|}\\
+g(-1;y;x;\lambda) &= \frac 1 {\abs{\lambda}} \arctan \frac{y-x}{|\lambda|}\\
+g(-3/2;y;x;\lambda) &= \frac{y-x}{\lambda^2}\{(y-x)^2+\lambda^2\}^{-1/2}
\end{align*}
benötigen, welche alle in \cite[Seite 2-3]{mai:3dbem} bewiesen wurden.
\end{lem}
Hierbei wollen wir darauf hinweisen, dass wir für $p=-1$ nicht die in \cite[Seite 3]{mai:3dbem} vorgeschlagene Formel verwenden, sondern die selbst hergeleitete:
\begin{align*}
- \int \frac 1 {(x-y)^2 +\lambda^2} dy &= \int \frac 1 {\lambda^2((\frac{y-x}{\abs{\lambda}})^2 +1)} dy
- = \left|
- \begin{array}{c}
- z = \frac{y-x}{\abs{\lambda}}\\
- dz = \frac 1 {\abs{\lambda}} dy
- \end{array}
- \right|
- = \frac 1 {\abs{\lambda}} \int \frac 1 {z^2+1} dz\\
- &= \frac 1 {\abs{\lambda}} \arctan (z) +c = \frac 1 {\abs{\lambda}} \arctan (\frac{y-x}{\abs{\lambda}}) +c
+\int \frac 1 {(x-y)^2 +\lambda^2} dy &= \int \frac 1 {\lambda^2((\frac{y-x}{\abs{\lambda}})^2 +1)} dy
+= \left|
+\begin{array}{c}
+z = \frac{y-x}{\abs{\lambda}}\\
+dz = \frac 1 {\abs{\lambda}} dy
+\end{array}
+\right|
+= \frac 1 {\abs{\lambda}} \int \frac 1 {z^2+1} dz\\
+&= \frac 1 {\abs{\lambda}} \arctan (z) +c = \frac 1 {\abs{\lambda}} \arctan (\frac{y-x}{\abs{\lambda}}) +c
\end{align*}
\begin{lem}
Des Weiteren werden wir das Integral
\begin{align*}
- G(p;y_1;y_2;x_1;x_2;\lambda) &:= \int \int \{(x_1-y_1)^2 + (x_2-y_2)^2 + \lambda^2\}^p dy_2 dy_1
+G(p;y_1;y_2;x_1;x_2;\lambda) &:= \int \int \{(x_1-y_1)^2 + (x_2-y_2)^2 + \lambda^2\}^p dy_2 dy_1
\end{align*}
benötigen, welches sich für den den Fall $p = -3/2$ und $\lambda = 0 \wedge \lambda \neq 0$ schreiben lässt als:
\begin{align*}
- G(-3/2;y_1,y_2;x_1,x_2,0) &= - \frac{\sqrt{(y_1-x_1)^2+(y_2-x_2)^2}}{(y_1-x_1)(y_2-x_2)}\\
- G(-3/2;y_1,y_2;x_1,x_2,\lambda) &=- \frac{\sgn\{(y_1-x_1)(y_2-x_2)\}}{2|\lambda|}\\
- &\cdot \arccos\left( \frac{-2(y_1-x_1)^2(y_2-x_2)^2}{\{(y_1-x_1)^2+\lambda^2\}\{(y_2-x_2)^2+\lambda^2\}} +1 \right).
+G(-3/2;y_1,y_2;x_1,x_2,0) &= - \frac{\sqrt{(y_1-x_1)^2+(y_2-x_2)^2}}{(y_1-x_1)(y_2-x_2)}\\
+G(-3/2;y_1,y_2;x_1,x_2,\lambda) &=- \frac{\sgn\{(y_1-x_1)(y_2-x_2)\}}{2|\lambda|}\\
+&\cdot \arccos\left( \frac{-2(y_1-x_1)^2(y_2-x_2)^2}{\{(y_1-x_1)^2+\lambda^2\}\{(y_2-x_2)^2+\lambda^2\}} +1 \right).
\end{align*}
Für alle weiteren relevante Fälle $p = -3/2 +k : k \in \N\backslash0$ können wir folgende Rekursionsformel aufstellen:
\begin{align*}
- (2p+2)G(p;y_1,y_2;x_1,x_2,\lambda) &= 2p\lambda^2 G(p-1;y_1,y_2;x_1,x_2,\lambda)\\
- &+ (y_1-x_1)g(p;y_2;x_2;\sqrt{(y_1-x_1)^2+\lambda^2})\\
- &+ (y_2-x_2)g(p;y_1;x_1;\sqrt{(y_2-x_2)^2+\lambda^2}),
+(2p+2)G(p;y_1,y_2;x_1,x_2,\lambda) &= 2p\lambda^2 G(p-1;y_1,y_2;x_1,x_2,\lambda)\\
+&+ (y_1-x_1)g(p;y_2;x_2;\sqrt{(y_1-x_1)^2+\lambda^2})\\
+&+ (y_2-x_2)g(p;y_1;x_1;\sqrt{(y_2-x_2)^2+\lambda^2}),
\end{align*} welche ebenfalls in \cite[Seite 5-7]{mai:3dbem} bewiesen wurden.
\end{lem}
\subsubsection{Parallele Elemente}
Liegen die beiden Elemente parallel zueinander, können wir sie wie in \cite[Seite 13]{mai:3dbem} gezeigt, darstellen als:
\begin{align*}
- T_j &= \v + [(0,s_1) \times (0,s_2) \times \{0\}]\\
- T_k &= \tilde \v + [(0,\tilde s_1) \times (0,\tilde s_2) \times \{0\}], \text{ mit } \v,\tilde \v \in \R^3.
+T_j &= \v + [(0,s_1) \times (0,s_2) \times \{0\}]\\
+T_k &= \tilde \v + [(0,\tilde s_1) \times (0,\tilde s_2) \times \{0\}], \text{ mit } \v,\tilde \v \in \R^3.
\end{align*}
Wir setzen $\bs{\delta} = \tilde \v - \v \in \R^3$, dann können wir das gesuchte Integral
\begin{align*}
- \int_{T_j} \int_{T_k} \frac{1}{|\bs x- \bs y|} ds_{\bs y} ds_{\bs x}\\
-\end{align*} umformen zu:
+\int_{T_j} \int_{T_k} \frac{1}{|\bs x- \bs y|} ds_{\bs y} ds_{\bs x}
+\end{align*}
+umformen zu:
\begin{align*}
&= \int_{T_j} \int_{T_k}
\left ( (x_1-y_1)^2+(x_2-y_2)^2+(x_3-y_3)^2 \right)^{-1/2} ds_{\bs y} ds_{\bs x}\\
\end{align*}
Dies können wir als Stammfunktion
\begin{align*}
-F_{par}(x_1,x_2,y_1,y_2,\bs \delta) &:= \int \int \int \int \left\{(x_1-y_1-\delta_1)^2+(x_2-y_2-\delta_2)^2+\delta_3^2\right\}^{-1/2}dy_2 dy_1 dx_2 dx_1
+F_{par}&(x_1,x_2,y_1,y_2,\bs \delta) :=\\ &\int \int \int \int \left\{(x_1-y_1-\delta_1)^2+(x_2-y_2-\delta_2)^2+\delta_3^2\right\}^{-1/2}dy_2dy_1 dx_2 dx_1
\end{align*}
-schreiben, welche das Integral über zwei parallele Elemente löst.
-
-\begin{align*}
+schreiben, welche das Integral über zwei parallele Elemente löst:
+\begin{align}
F_{par}(x_1&,x_2,y_1,y_2,\bs \delta)\\
- =& (x_1-y_1-\delta_1)(x_2-y_2-\delta_2)G(-1/2;x_1,x_2;y_1+\delta_1,y_2+\delta_2,\delta_3)\\
- &-(x_1-y_1-\delta_1) g(1/2;x_1;y_1+\delta_1;\{(x_2-y_2-\delta_2)^2 + \delta_3^2\}^{1/2})\\
- &-(x_2-y_2-\delta_2) g(1/2;x_2;y_2+\delta_2;\{(x_1-y_1-\delta_1)^2 + \delta_3^2\}^{1/2})\\
- &+\frac 1 3 \{(x_1-y_1-\delta_1)^2+(x_2-y_2-\delta_2)^2+\delta_3^2\}^{3/2}
-\end{align*}
+=& (x_1-y_1-\delta_1)(x_2-y_2-\delta_2)G(-1/2;x_1,x_2;y_1+\delta_1,y_2+\delta_2,\delta_3)\\
+&-(x_1-y_1-\delta_1) g(1/2;x_1;y_1+\delta_1;\{(x_2-y_2-\delta_2)^2 + \delta_3^2\}^{1/2})\\
+&-(x_2-y_2-\delta_2) g(1/2;x_2;y_2+\delta_2;\{(x_1-y_1-\delta_1)^2 + \delta_3^2\}^{1/2})\\
+&+\frac 1 3 \{(x_1-y_1-\delta_1)^2+(x_2-y_2-\delta_2)^2+\delta_3^2\}^{3/2}
+\end{align}
\subsubsection{Orthogonale Elemente}
Analog können wir orthogonal liegende Elemente wie in \cite[Seite 14]{mai:3dbem} gezeigt, schreiben als:
\begin{align*}
- T_j &= \v + [(0,s_1) \times (0,s_2) \times \{0\}]\\
- T_k &= \tilde \v + [ \{0\} \times (0,\tilde s_2) \times (0,\tilde s_3)], \text{ mit } \v,\tilde \v \in \R^3
+T_j &= \v + [(0,s_1) \times (0,s_2) \times \{0\}]\\
+T_k &= \tilde \v + [ \{0\} \times (0,\tilde s_2) \times (0,\tilde s_3)], \text{ mit } \v,\tilde \v \in \R^3
\end{align*} und setzen wir $\bs{\delta} = \tilde \v - \v \in \R^3$.
Dann können wir
\begin{align*}
- \int_{T_j} \int_{T_k} \frac{1}{|\bs x- \bs y|} ds_{\bs y} ds_{\bs x}\\
-\end{align*} umformen zu:
+\int_{T_j} \int_{T_k} \frac{1}{|\bs x- \bs y|} ds_{\bs y} ds_{\bs x}
+\end{align*}
+umformen zu:
\begin{align*}
&= \int_{T_j} \int_{T_k}
\left ( (x_1-y_1)^2+(x_2-y_2)^2+(x_3-y_3)^2 \right)^{-1/2} ds_{\bs y} ds_{\bs x}\\
\end{align*}
Dies können wir ebenfalls als Stammfunktion
\begin{align*}
-F_{ort}(x_1,x_2,y_2,y_3,\bs \delta) &:= \int \int \int \int \left\{(x_1-\delta_1)^2+(x_2-y_2-\delta_2)^2+(-y_3-\delta_3)^2\right\}^{-1/2}
+F_{ort}&(x_1,x_2,y_2,y_3,\bs \delta) :=\\ &\int \int \int \int \left\{(x_1-\delta_1)^2+(x_2-y_2-\delta_2)^2+(-y_3-\delta_3)^2\right\}^{-1/2}
dy_3 dy_2 dx_2 dx_1
\end{align*}
-für orthogonal liegende Elemente schreiben.
-
+für orthogonal liegende Elemente schreiben:
\begin{align*}
2F_{ort}(x_1&,x_2,y_2,y_3,\bs \delta)\\
=&-G(1/2;y_3,x_1;-\delta_3,\delta_1,x_2-y_2-\delta_2)\\
\clearpage
+
\section{fasträguläre Partionierung in \Matlab}
\todo{
\begin{itemize}
- \item Datenstruktur
- \item Netzverfeinerung
- \item Code + kurze Dokumentation des Codes (was passiert wie/wo)
+\item Datenstruktur
+\item Netzverfeinerung
+\item Code + kurze Dokumentation des Codes (was passiert wie/wo)
\end{itemize}
}
Für die Implementierung in \Matlab~und C++ wollen wollen wir eine einheitliche Datenstruktur einführen.
Die für die Partition $\T_{\ell} = \{T_1\ldots T_N\}$ benötigten Knoten $\K_{\ell} = \{C_1\ldots C_M\}$ stellen wir in einer $M \times 3$ Matrix $COO$ dar. Dabei enthält die $j$-te Zeile die Koordinaten des Knoten $C_j$, d.h. :
\begin{align}
- COO[j,1:3] = C_j := (x_j,y_j,z_j)^{T} \text{ wobei } x_j,y_j,z_j \in \R.
+COO[j,1:3] = C_j := (x_j,y_j,z_j)^{T} \text{ wobei } x_j,y_j,z_j \in \R.
\end{align}
Die Elemente $\T_{\ell}$ werden wir ebenfalls Zeilenweise in einer $N \times 4$ Matrix $ELE$ abspeichern. Dabei soll die $i$-te Zeile den Indizes der Knoten $\{C_j,C_k,C_{\ell},C_m\}$ des Elements $T_i$ entsprechen, also:
\begin{align}
- ELE[i,1:4] = T_i := (j,k,l,m).
+ELE[i,1:4] = T_i := (j,k,l,m).
\end{align}
Die Knoten wollen wir gegen den Uhrzeigersinn anordnen \todo{und der Knoten $C_j$ soll der kleinste bezüglich der Koordinaten sein}.
Für die bessere Handhabung der Elemente beim Verfeinern der Partition, wollen wir auch die Nachbarschaftsrelationen geeignet abspeichern. Dazu überlegen wir uns, dass wir Aufgrund der Netzstabilität maximal zwei Nachbarn pro Kante eines Elements zulassen wollen.
Wir legen also eine $M \times 8$ Matrix für die Indizes der Nachbarelemente an, wobei die $i$-te Zeile die Nachbarelemente $\{T_{n_1},\ldots,T_{n_8}\}$ zum Element $T_i$ enthält.
\begin{align}
- NEI[i,1:8] = N_i := (n_1,\ldots,n_8)
+NEI[i,1:8] = N_i := (n_1,\ldots,n_8)
\end{align}
Offensichtlich ist $i \notin N_i$. Wir wollen uns aber noch genauer eine geeignete Anordnung für die Nachbarelemente überlegen. Hierbei bezeichnen wir die Seite $[j,k]$ eines Elements als Seite 1 und gegen den Uhrzeigersinn alle weiteren mit $\{2,3,4\}$. Für den einfacheren Zugriff auf die Elemente einer Seite $s$, wollen wir die Nachbarelemente zur Seite $s$ unter den Indizes $n_s$ und $n_{s+4}$ abspeichern. Die Nachbarelemente $\{T_{n_s}, T_{n_{s+4}}\}$ liegen also an der Seite $s$ des Elements. Für Seiten die nur einen Nachbarn $T_{n_s}$ besitzen setzen wir $n_{s+4}=0$ und für Seiten mit keinem Nachbarn setzen wir $n_s = n_{s+4} = 0$. Daraus folgt unmittelbar, dass für $n_s =0$ auch $n_{s+4} = 0$ gilt, die Seite $s$ also keine Nachbarelemente besitzt und umgekehrt folgt aus $n_{s+4} \neq 0$ $n_s \neq 0$, womit die Seite $s$ genau zwei Nachbarelemente hat.
(Siehe Abb.:\ref{exmpl13:nei:part})
\caption{Nachbarschaftsrelationen Element 4 aus Abb.\ref{exmpl13}}
\label{exmpl13:nei:part}
\centering
- \subfloat[Lage]{\includegraphics{fig/Net_Neigh}}
- \subfloat[Nachbarn]{\input{fig/exmpl3_nei_part}}
+\subfloat[Lage]{\includegraphics{fig/Net_Neigh}}
+\subfloat[Nachbarn]{\input{fig/exmpl3_nei_part}}
\end{figure}
\noindent
Diese Funktion sei die Implementierung von Algorithmus \ref{alg:refine}.
Da wir im weiteren Verlauf sowohl adaptive also auch anisotrope Netzverfeinerung zulassen wollen, ist es sinnvoll eine Verfeinerungsfunktion zu implementieren, die alle Möglichen Teilungsprozesse auf einem Element unterstützt. Dabei sind vier nur wirklich relevant:
\begin{enumerate}
- \item keine Teilung
- \item volle Teilung in vier gleich große Elemente
- \item halbe Teilung in zwei gleichgroße horizontal liegende Elemente
- \item halbe Teilung in zwei gleichgroße vertikal liegende Elemente
+\item keine Teilung
+\item volle Teilung in vier gleich große Elemente
+\item halbe Teilung in zwei gleichgroße horizontal liegende Elemente
+\item halbe Teilung in zwei gleichgroße vertikal liegende Elemente
\end{enumerate}
Zusätzlich wurde auch Typ 5. belegt, welcher als Ergebnis eine volle Teilung 2. ausführt, diese aber schrittweise durch eine 3. Teilung und zwei 4. Teilungen. Aus Sicherheitsgründen wird auch jede vierte volle 2. Teilung durch eine 5. Teilung ausgeführt, da sonst kurzzeitig Seiten mit mehr als zwei Elementen auftreten könnten.
Zur anisotropen Verfeinerung wird weiterhin berechnet:
\begin{align*}
\begin{pmatrix}
- C_j^{(1)}\\ C_j^{(2)}\\C_j^{(3)}\\ C_j^{(4)}
+C_j^{(1)}\\ C_j^{(2)}\\C_j^{(3)}\\ C_j^{(4)}
\end{pmatrix}
&= \frac 1 4
\begin{pmatrix}
1 & -1 & -1 & 1\\
\end{pmatrix}
\begin{pmatrix}
- x_j^{(1)}\\ x_j^{(2)}\\x_j^{(3)}\\ x_j^{(4)}
+x_j^{(1)}\\ x_j^{(2)}\\x_j^{(3)}\\ x_j^{(4)}
\end{pmatrix}
\end{align*}
$\theta \in (0,1),i =0$
\begin{enumerate}
- \renewcommand{\theenumi}{(\roman{enumi})}
- \item Verfeinere $T_{\ell}^{(i)}$ um $\hat T_{\ell}$ zu erhalten
- \item Berechne die Galerkinlösung $\hat \phi_{\ell} \in P^0(\\hat T_{\ell})$
- \item Berechne Fehlerschätzer $\tilde \mu_{i} := \norm{\varrho^{\ell/2}(\hat \phi_{\ell} - \Pi_{\ell} \hat \phi_{\ell} )}$
- \item Wähle $M_{\ell} \subseteq T_{\ell}^{(i)}$ mit minimaler Kardinalität, so dass
+\renewcommand{\theenumi}{(\roman{enumi})}
+\item Verfeinere $T_{\ell}^{(i)}$ um $\hat T_{\ell}$ zu erhalten
+\item Berechne die Galerkinlösung $\hat \phi_{\ell} \in P^0(\\hat T_{\ell})$
+\item Berechne Fehlerschätzer $\tilde \mu_{i} := \norm{\varrho^{\ell/2}(\hat \phi_{\ell} - \Pi_{\ell} \hat \phi_{\ell} )}$
+\item Wähle $M_{\ell} \subseteq T_{\ell}^{(i)}$ mit minimaler Kardinalität, so dass
\begin{align}
- \theta \sum_{T\in \T^{(i)}_{\ell}} \tilde\mu_{i}^2 & \leq \sum_{T\in M_{\ell}} \tilde\mu_{i}^2
+\theta \sum_{T\in \T^{(i)}_{\ell}} \tilde\mu_{i}^2 & \leq \sum_{T\in M_{\ell}} \tilde\mu_{i}^2
\end{align}
- \item Verfeinere die Markierten Elemente $M_{\ell}$ um $\T_{\ell}^{(i+1)}$ zu erhalten
- \item $i \mapsto i+1$, gehe zu $(i)$
+\item Verfeinere die Markierten Elemente $M_{\ell}$ um $\T_{\ell}^{(i+1)}$ zu erhalten
+\item $i \mapsto i+1$, gehe zu $(i)$
\end{enumerate}
-Zum Plotten (Abb.\ref{exmplAA_2DQuad})werden noch folgende Schritte ausgeführt
+Zum Plotten (Abb.\ref{fig:exmplAA_2DQuad})werden noch folgende Schritte ausgeführt
\begin{itemize}
- \item Berechne Galerkinlösung $\phi_{l} \in P^0(\T_{\ell}^{(i)})$
+\item Berechne Galerkinlösung $\phi_{l} \in P^0(\T_{\ell}^{(i)})$
% \item $\enorm{\hat \phi_{\ell}}$
% \item $\enorm{\phi_{\ell}}$
- \item $error_{i} = \sqrt{\enorm{\phi}^2 - \enorm{\phi_{l}}^2}$
- \item $\mu_{i} = \norm{\varrho^{1/2}(\hat \phi_{\ell} - \phi_{l} )}$
- \item $\eta_{i} = \enorm{\hat \phi_{\ell} - \phi_{l}}$
- \item $\kappa_{i} = \enorm{\hat \phi_{\ell}^{(i)}}-\enorm{\hat \phi_{\ell}^{(i-1)}}$
- \item $\kappa2_{i} = \enorm{\phi_{l}^{(i)}}-\enorm{\phi_{l}^{(i-1)}}$
- \item $\kappa3_{i} = \enorm{\hat \phi_{\ell}^{(i)}-\hat \phi_{\ell}^{(i-1)}}$
+\item $error_{i} = \sqrt{\enorm{\phi}^2 - \enorm{\phi_{l}}^2}$
+\item $\mu_{i} = \norm{\varrho^{1/2}(\hat \phi_{\ell} - \phi_{l} )}$
+\item $\eta_{i} = \enorm{\hat \phi_{\ell} - \phi_{l}}$
+\item $\kappa_{i} = \enorm{\hat \phi_{\ell}^{(i)}}-\enorm{\hat \phi_{\ell}^{(i-1)}}$
+\item $\kappa2_{i} = \enorm{\phi_{l}^{(i)}}-\enorm{\phi_{l}^{(i-1)}}$
+\item $\kappa3_{i} = \enorm{\hat \phi_{\ell}^{(i)}-\hat \phi_{\ell}^{(i-1)}}$
\end{itemize}
\clearpage
+
\section{Numerische Experimente}
\todo{
\begin{itemize}
- \item Zusammenfassung der h-h/2 Strategie (Ergebnisse aus Ferraz-DA und
-Paper), interessant sind nur eta und tilde-mu
-\item adaptiver Algorithmus als Algorithmus und MATLAB-Code
-\item Warum braucht man Adaptivitaet? Warum reicht isotrope Verfeinerung
-nicht aus? ggf. anhand numerischem Beispiel
-\item Beispiele, bei denen anisotrop gerechnet wird, bis Instabilitaet
-auftritt, verschiedene Strategien (voll-analytisch + ...)
+ \item Zusammenfassung der h-h/2 Strategie (Ergebnisse aus Ferraz-DA und Paper), interessant sind nur eta und tilde-mu
+ \item adaptiver Algorithmus als Algorithmus und MATLAB-Code
+\item Warum braucht man Adaptivitaet? Warum reicht isotrope Verfeinerung nicht aus? ggf. anhand numerischem Beispiel
+\item Beispiele, bei denen anisotrop gerechnet wird, bis Instabilitaet auftritt, verschiedene Strategien (voll-analytisch + ...)
\item in Plots: Fehler/eta/tilde-mu
-\item genau beschreiben: Was ist das Beispiel? Was sieht man im Plot?
-Wie interpretiert man das?
+\item genau beschreiben: Was ist das Beispiel? Was sieht man im Plot? Wie interpretiert man das?
\end{itemize}
}
\begin{defi}Es bezeichne $\phi$ die Lösung von Formel , $\phi_{\ell}$ die Galerkinlösung auf dem Gitter $\T_{\ell}$ und $\hat \phi_{\ell}$ die Galerkinlösung auf dem uniform verfeinerten Gitter $\hat \T_{\ell}$. Dann gilt, der Schätzer
\begin{align}
- \eta_{\ell} &:= \enorm{\hat \phi_{\ell} - \phi_{\ell}}
+\eta_{\ell} &:= \enorm{\hat \phi_{\ell} - \phi_{\ell}}
\end{align}
ist effizient
\begin{align}
- \eta_{ell} &\leq \enorm{\phi - \phi_{\ell}}
+\eta_{ell} &\leq \enorm{\phi - \phi_{\ell}}
\end{align}
und unter der Saturationsannahme
\begin{align}
- \norm{\phi -\hat \phi_{\ell}} &\leq C_{sat} \cdot \norm{\phi - \phi_{\ell}} & 0 < C_{sat} < 1
+\norm{\phi -\hat \phi_{\ell}} &\leq C_{sat} \cdot \norm{\phi - \phi_{\ell}} & 0 < C_{sat} < 1
\end{align}
zuverlässig
\begin{align}
- \enorm{\phi - \phi_{\ell}}&\leq \frac 1{\sqrt{1-C_{sat}}}\eta_{\ell}.
+\enorm{\phi - \phi_{\ell}}&\leq \frac 1{\sqrt{1-C_{sat}}}\eta_{\ell}.
\end{align}
\end{defi}
\begin{sat}[A-posteriori Fehlerschätzer] Seien also:
\begin{align}
- \eta_{\ell} &= \enorm{\hat \phi_{\ell} - \phi_{\ell}}\\
- \tilde\eta_{\ell} &= \enorm{\hat \phi_{\ell} - \Pi_{\ell}\hat \phi_{\ell}}\\
- \mu_{\ell} &= \norm{\varrho^{1/2}(\hat \phi_{\ell} - \phi_{\ell})}_{L^2(\Gamma)}\\
- \tilde\mu_{\ell} &= \norm{\varrho^{1/2}(\hat \phi_{\ell} - \Pi_{\ell}\hat \phi_{\ell})}_{L^2(\Gamma)}
+\eta_{\ell} &= \enorm{\hat \phi_{\ell} - \phi_{\ell}}\\
+\tilde\eta_{\ell} &= \enorm{\hat \phi_{\ell} - \Pi_{\ell}\hat \phi_{\ell}}\\
+\mu_{\ell} &= \norm{\varrho^{1/2}(\hat \phi_{\ell} - \phi_{\ell})}_{L^2(\Gamma)}\\
+\tilde\mu_{\ell} &= \norm{\varrho^{1/2}(\hat \phi_{\ell} - \Pi_{\ell}\hat \phi_{\ell})}_{L^2(\Gamma)}
\end{align}
wobei $\Pi_{\ell}$ die $L_2$ Projektion auf $P^0(\T_{\ell})$ ist.
\end{sat}
$\mu_{\ell}$ ist da $\mu_{\ell} \approx \eta_{\ell}$ gilt noch immer zuverlässig und effizient.\\
Dann gilt auf isotropen Netzen:
\begin{itemize}
- \item Schätzer sind equivalent\\ $\tilde \mu_{\ell} \leq \mu_{\ell} \leq \sqrt{2}C_3 \eta_{\ell}$
- \item sie sind effizient
- \item sie sind unter Saturationsannahme auch zuverlässig
+\item Schätzer sind equivalent\\ $\tilde \mu_{\ell} \leq \mu_{\ell} \leq \sqrt{2}C_3 \eta_{\ell}$
+\item sie sind effizient
+\item sie sind unter Saturationsannahme auch zuverlässig
\end{itemize}
- Siehe S.F. Paper $\mapsto$ THM 3.2 \& 3.4
+Siehe S.F. Paper $\mapsto$ THM 3.2 \& 3.4
\end{align*}
-
-
-
-
\todo{Viele Bunte Bilder und Testergebnisse}
\begin{itemize}
- \item gerechnet wurde $V\phi = 1$
- \item Beispiele für verschiedene Figuren
- \begin{itemize}
- \item 2D LShape
- \item 2D Rechteck
- \item 3D Würfel
- \item Fischer Würfel
- \item L Figur
- \end{itemize}
- \item Beispiele für verschiedene Berechnungsarten
- \begin{itemize}
- \item VollAnalytisch
- \item SemiAnalytisch über Element
- \item SemiAnalytisch über Achse
- \item SemiAnalytisch über Seite
- \end{itemize}
- \item Fehlerschätzer untersuchen und prüfen ob sie sich wie erwartet Verhalten
+\item gerechnet wurde $V\phi = 1$
+\item Beispiele für verschiedene Figuren
+\begin{itemize}
+\item 2D LShape
+\item 2D Rechteck
+\item 3D Würfel
+\item Fischer Würfel
+\item L Figur
+\end{itemize}
+\item Beispiele für verschiedene Berechnungsarten
+\begin{itemize}
+\item VollAnalytisch
+\item SemiAnalytisch über Element
+\item SemiAnalytisch über Achse
+\item SemiAnalytisch über Seite
+\end{itemize}
+\item Fehlerschätzer untersuchen und prüfen ob sie sich wie erwartet Verhalten
\end{itemize}
\begin{figure}[ht]
\caption{Objekt Beispiele}
\centering
-\label{objects}
+\label{fig:objects}
\subfloat[2D L Shape]{\includegraphics[width=0.5\textwidth]{fig/exmpl_2DLShape_ref}}
\subfloat[2D Quad]{\includegraphics[width=0.5\textwidth]{fig/exmpl_2DQuad_ref}}\\
\subfloat[3D Cube]{\includegraphics[width=0.5\textwidth]{fig/exmpl_3DCube_ref}}
\begin{figure}[ht]
\caption{2D Quad adaptiv anisotrop vollanalytisch $V\phi = 1$}
\centering
-\label{exmplAA_2DQuad}
+\label{fig:exmplAA_2DQuad}
\subfloat[Fehler]{\includegraphics[width=0.5\textwidth]{fig/exmplAA_2DQuad_error}}
\subfloat[Energie Norm]{\includegraphics[width=0.5\textwidth]{fig/exmplAA_2DQuad_norm}}
\end{figure}
\begin{figure}[ht]
\caption{12t05n05 2DQuad 31}
\centering
-\label{12t05n05_2DQuad_31}
+\label{fig:12t05n05_2DQuad_31}
\subfloat[Fehler]{\includegraphics[width=0.5\textwidth]{fig/12t05n05_2DQuad_31_error}}
\subfloat[Energie Norm]{\includegraphics[width=0.5\textwidth]{fig/12t05n05_2DQuad_31_norm}}\\
\subfloat[Kondition]{\includegraphics[width=0.5\textwidth]{fig/12t05n05_2DQuad_31_cond}}
\begin{figure}[ht]
\caption{13t05n05 2DQuad 30}
\centering
-\label{13t05n05_2DQuad_30}
+\label{fig:13t05n05_2DQuad_30}
\subfloat[Fehler]{\includegraphics[width=0.5\textwidth]{fig/13t05n05_2DQuad_30_error}}
\subfloat[Energie Norm]{\includegraphics[width=0.5\textwidth]{fig/13t05n05_2DQuad_30_norm}}\\
\subfloat[Kondition]{\includegraphics[width=0.5\textwidth]{fig/13t05n05_2DQuad_30_cond}}
\todo{
Die wichtigsten Funktionen
\begin{enumerate}
- \item compute
- \item refine
- \item mark
- \item build\_V
- \item plot
+\item compute
+\item refine
+\item mark
+\item build\_V
+\item plot
\end{enumerate}
}
% \lstinputlisting[language=oC++]{../src/slpRectangle.cpp}
projects / bacc.git / commitdiff