updated nipdes

Mathematics/4th/Numerical_integration_of_partial_differential_equations/Numerical_integration_of_partial_differential_equations.tex

@@ -553,20 +553,18 @@
     $$
   \end{proof}
   \subsection{Introduction to finite element methods}
-  \begin{definition}
-    The \emph{finite element method} is one of the most popular, general,
-    powerful and elegant approaches for approximating the solutions of
-    PDEs. Unlike finite difference methods, it naturally handles complicated domains (useful for engines and aeroplanes) and minimally
-    regular data (such as discontinuous forcing terms).
+  The \emph{finite element method} is one of the most popular, general,
+  powerful and elegant approaches for approximating the solutions of
+  PDEs. Unlike finite difference methods, it naturally handles complicated domains (useful for engines and aeroplanes) and minimally
+  regular data (such as discontinuous forcing terms).
 
-    There are four basic ingredients in the finite element method:
-    \begin{enumerate}
-      \item A variational formulation of the problem in an infinite-dimensional space $V$.
-      \item A variational formulation of the problem in a finite-dimensional space $V_h\subset V$.
-      \item The construction of a basis for $V_h$.
-      \item The assembly and solution of the resulting linear system of equations.
-    \end{enumerate}
-  \end{definition}
+  There are four basic ingredients in the finite element method:
+  \begin{enumerate}
+    \item A variational formulation of the problem in an infinite-dimensional space $V$.
+    \item A variational formulation of the problem in a finite-dimensional space $V_h\subset V$.
+    \item The construction of a basis for $V_h$.
+    \item The assembly and solution of the resulting linear system of equations.
+  \end{enumerate}
   \subsubsection{The variational formulation}
   \begin{definition}
     Let $\Omega\subseteq \RR^n$ be an open bounded connected set such that $\Fr{U}$ is of class $\mathcal{C}^1$, $f\in\mathcal{C}(\Omega)$ and $g\in\mathcal{C}(\Fr{\Omega})$. Consider the following Dirichlet problem of finding $u\in\mathcal{C}^2(\Omega)\cap \mathcal{C}(\overline{\Omega})$ such that:
@@ -578,24 +576,31 @@
     \end{equation}
     Let
     $$V:=\{v:\Omega\rightarrow\RR:\norm{v}_{L^2(\Omega)}+\norm{\grad{v}}_{L^2(\Omega)}<\infty, v|_{\Fr{\Omega}}=0\}$$
-    The \emph{variational formulation} of the problem is to find $u\in V$ such that:
+    The \emph{variational formulation} (or \emph{weak formulation}) of the problem is to find $u\in V$ such that:
     \begin{equation}\label{NIPDE:varDirichlet}
-      \int_\Omega \grad{u}\cdot\grad{v}\dd{x}=\int_\Omega fv\dd{x}\quad\forall v\in V
+      \int_\Omega \grad{u}\cdot\grad{v}\dd{\vf{x}}=\int_\Omega fv\dd{\vf{x}}\quad\forall v\in V
     \end{equation}
   \end{definition}
   \begin{remark}
-    The variational formulation can be obtained by multiplying \cref{NIPDE:Dirichlet} by $v$ and the use the \mnameref{PDE:greenidentities}.
+    The variational formulation can be obtained by multiplying \cref{NIPDE:Dirichlet} by $v$ and using the \mnameref{PDE:greenidentities}.
   \end{remark}
-  \begin{lemma}
-    If $u\in V$ is a solution to \cref{NIPDE:varDirichlet}, then $u$ is a solution to \cref{NIPDE:Dirichlet}.
-  \end{lemma}
   \begin{theorem}
     If $f\in \mathcal{C}(\Omega)$, then the solutions to \cref{NIPDE:varDirichlet} are $\mathcal{C}^2(\Omega)$.
   \end{theorem}
+  \begin{lemma}
+    If $u\in V$ is a solution to \cref{NIPDE:varDirichlet}, then $u$ is a solution to \cref{NIPDE:Dirichlet}.
+  \end{lemma}
+  \begin{proof}
+    Note that $\grad u\cdot\grad v=\div(v\grad u) - v\laplacian u$. Thus, using the \mnameref{DG:divergenceRn} we have:
+    $$
+      \int_\Omega \grad{u}\cdot\grad{v}-fv\dd{\vf{x}}=\int_\Omega v(-\laplacian u-f)\dd{\vf{x}}+\int_{\Fr{\Omega}}v\grad{u}\cdot\vf{n}\dd{\vf{s}}=0
+    $$
+    because $v=0$ on $\Fr{\Omega}$. Now using the \mnameref{PDE:fundamentallemma}, we conclude that we must have $-\laplacian u=f$ in $\Omega$.
+  \end{proof}
   \begin{definition}[Galerkin approximation]
     Let $V_h\subset V$ be a finite-dimensional subspace of $V$. The \emph{Galerkin approximation} of \cref{NIPDE:varDirichlet} is to find $u_h\in V_h$ such that:
     \begin{equation}\label{NIPDE:Galerkin}
-      \int_\Omega \grad{u_h}\cdot\grad{v_h}\dd{x}=\int_\Omega fv_h\dd{x}\quad\forall v_h\in V_h
+      \int_\Omega \grad{u_h}\cdot\grad{v_h}\dd{\vf{x}}=\int_\Omega fv_h\dd{\vf{x}}\quad\forall v_h\in V_h
     \end{equation}
   \end{definition}
   \subsubsection{Construction of function spaces}
@@ -605,20 +610,31 @@
       \item $\Int(K_i) \cap \Int(K_j) = \varnothing$ for all $i\neq j$.
       \item $\bigcup_{i=1}^N K_i = \overline{\Omega}$.
     \end{enumerate}
-    The cells are usually chosen to be $n$-simplices or $n$-parallelepipeds.
+    The cells are usually chosen to be $n$-simplexes or $n$-parallelepipeds.
   \end{definition}
   \begin{definition}
-    The \emph{finite element method} is a particular choice of Galerkin approximation, where the discrete function space
+    The \emph{finite element method} (\emph{FEM}) is a particular choice of Galerkin approximation, where the discrete function space
     $V_h$ is:
     \begin{multline*}
       V_h:=\{v\in\mathcal{C}(\Omega):v\text{ is piecewise linear when restricted}\\\text{to a cell}\}
     \end{multline*}
-    Note that the functions in $V_h$ are uniquely determined by its values at the vertices of the cell. The vertices of the cells are called \emph{nodes}.
+    Note that the functions in $V_h$ are uniquely determined by its values at the vertices of the cell because of the unicity of the interpolating polynomial. The vertices of the cells are called \emph{nodes}.
   \end{definition}
   \begin{definition}
-    Given the location $x_i$ of $M$ \emph{nodes} in $\Int\Omega$, we define the \emph{nodal basis} $(\phi_1,\ldots,\phi_M)$ as the functions $\phi_i$ such that:
-    $$\phi_i(x_j)=\delta_{ij}$$
+    Given the locations $\vf{x}_i$ of the $M$ \emph{nodes} in $\Int\Omega$, we define the \emph{nodal basis} $(\phi_1,\ldots,\phi_M)$ as the functions $\phi_i$ such that:
+    $$\phi_i(\vf{x}_j)=\delta_{ij}$$
   \end{definition}
+  \begin{lemma}
+    The nodal basis is indeed a basis of $V_h$.
+  \end{lemma}
+  \begin{proof}
+    Let $v\in V_h$. Then, $v$ can be written as:
+    $$
+      v=\sum_{i=1}^M v(\vf{x}_i)\phi_i
+    $$
+    Since it is uniquely determined by its values at the nodes, the equality holds.
+    So, $\langle \phi_1, \ldots, \phi_M\rangle=V_h$. Furthermore, if we have $\sum_{i=1}^M c_i\phi_i=0$, then evaluating at $\vf{x}_j$ we have $c_j=0$ $\forall j=1,\ldots,M$.
+  \end{proof}
   \subsubsection{Linear algebraic formulation}
   \begin{proposition}
     Given a mesh of $\Omega$, consider the space $V_h\subset V$ and its associate nodal basis. Suppose:
@@ -629,20 +645,60 @@
     $$\vf{Au}=\vf{b}$$
     where $\vf{A}=(a_{ij})$ and $\vf{b}=(b_i)$ are defined as:
     \begin{equation*}
-      a_{ij} =\int_\Omega \grad{\phi_i}\cdot\grad{\phi_j}\dd{x} \qquad b_i    =\int_\Omega f\phi_i\dd{x}
+      a_{ij} =\int_\Omega \grad{\phi_i}\cdot\grad{\phi_j}\dd{\vf{x}} \qquad b_i    =\int_\Omega f\phi_i\dd{\vf{x}}
     \end{equation*}
+    The matrix $\vf{A}$ is usually called the \emph{stiffness matrix} and $\vf{b}$ the \emph{load vector}.
   \end{proposition}
   \begin{proof}
     Since, $u_h\in V_h$, and using the linearity of the integral we have:
     \begin{align*}
-      \int_\Omega \grad{u_h}\cdot\grad{v_h}\dd{x}                     & =\int_\Omega f v_h\dd{x}                    \\
-      \sum_{i=1}^M v_i \int_\Omega \grad{u_h}\cdot\grad{\phi_i}\dd{x} & =\sum_{i=1}^M v_i \int_\Omega f\phi_i\dd{x}
+      \int_\Omega \grad{u_h}\cdot\grad{v_h}\dd{\vf{x}}                     & =\int_\Omega f v_h\dd{\vf{x}}                    \\
+      \sum_{i=1}^M v_i \int_\Omega \grad{u_h}\cdot\grad{\phi_i}\dd{\vf{x}} & =\sum_{i=1}^M v_i \int_\Omega f\phi_i\dd{\vf{x}}
     \end{align*}
     As this holds for all $v_h\in V_h$, we have that this is equivalent to
-    $$\int_\Omega \grad{u_h}\cdot\grad{\phi_i}\dd{x} =\int_\Omega f\phi_i\dd{x}$$
-    which implies:
-    $$\sum_{j=1}^M u_j \int_\Omega \grad{\phi_j}\cdot\grad{\phi_i}\dd{x} =\int_\Omega f\phi_i\dd{x}$$
+    $$\int_\Omega \grad{u_h}\cdot\grad{\phi_i}\dd{\vf{x}} =\int_\Omega f\phi_i\dd{\vf{x}}$$
+    for $i=1,\ldots,M$, which implies:
+    $$\sum_{j=1}^M u_j \int_\Omega \grad{\phi_j}\cdot\grad{\phi_i}\dd{\vf{x}} =\int_\Omega f\phi_i\dd{\vf{x}}$$
   \end{proof}
+  \begin{remark}
+    Solving this system of linear equations we obtain the approximation by finite elements of the Dirichlet problem for the Poisson equation (\mcref{NIPDE:Dirichlet}). Note that the approximate solution is a piecewise linear function which may not be differentiable at the vertices of the cells. Even so, the approximate solution converges to the exact solution as the mesh is refined.
+  \end{remark}
+  \begin{remark}
+    On the computation of the coefficients $a_{ij}$ we should proceed as follows:
+    $$
+      a_{ij}=\sum_{m=1}^N\int_{K_m} \grad{\phi_i}\cdot\grad{\phi_j}\dd{\vf{x}}
+    $$
+    Note, however, that many of these integrals will be zero as if $P_i\notin K_m$, then $\varphi_i=0$ on the nodes of $K_m$, and therefore $\varphi_i=0$ and $\grad{\varphi_i}=0$ on $K_m$. Thus, we only need to compute the integrals for $K_m$ such that $P_i, P_j\in K_m$. For these (a priori) non-zero integrals, we use a reference $n$-simplex to compute them. In the following proposition we expose the case $n=2$.
+  \end{remark}
+  \begin{proposition}
+    Let $S$ be an $n$-simplex with vertices at $Q_0=\vf{0}$, $Q_i=\vf{e}_i$ (thought as a point), $i=1,\ldots,n$, where $\vf{e}_i$ is the $i$-th vector of the canonical basis of $\RR^n$. Consider the FEM method for the \mcref{NIPDE:Dirichlet}. Then:
+    \begin{align*}
+      \int_{K_m}\grad\varphi_{K_m,\ell}\cdot \grad\varphi_{K_m,k}\dd{\vf{x}} & =\frac{d_m}{n!}{\grad\psi_\ell}{\left(\transpose{\vf{D\sigma}_m}\vf{D\sigma}_m\right)}^{-1}\transpose{\grad\psi_k}
+    \end{align*}
+    where $\sigma_m$ is the affine transformation that carries the reference simplex $S$ onto $K_m$, $d_m=\abs{\det\vf{D\sigma}_m}$, $\phi_{K_m,\ell}$ denote that basis function such that evaluates to 1 at the $\ell$-th vertex of $K_m$ (with an ordering fixed), $\ell =0,\ldots,n$, and:
+    $$
+      \psi_k(\vf{x})=\begin{cases}
+        1-\sum_{i=1}^n x_i & k=0          \\
+        x_k                & k=1,\ldots,n
+      \end{cases}
+    $$
+  \end{proposition}
+  \begin{proof}
+    Note $\psi_k(Q_k)=\delta_{ij}$ and so by the unicity of the interpolation we have $\varphi_{K_m,\ell}\circ \sigma_m=\psi_\ell$, $\ell=0,\ldots,n$. Thus, by the chain rule, $\grad\psi_\ell=\transpose{\vf{D\sigma}_m}\grad\varphi_{K_m,\ell}$, and so:
+    \begin{align*}
+      \int_{K_m} & \grad\varphi_{K_m,\ell}\cdot \grad\varphi_{K_m,k}\dd{\vf{y}}  =\int_S\grad\varphi_{K_m,\ell}\cdot \transpose{\grad\varphi_{K_m,k}}d_m\dd{\vf{x}} \\
+                 & =\int_S\grad\psi_\ell{\left(\vf{D\sigma}_m\right)}^{-1}\transpose{{\left(\vf{D\sigma}_m\right)}^{-1}}\transpose{\grad\psi_k}d\dd{\vf{x}}         \\
+                 & =\frac{d_m}{n!}\grad\psi_\ell{\left(\transpose{\vf{D\sigma}_m}\vf{D\sigma}_m\right)}^{-1}\transpose{\grad\psi_k}
+    \end{align*}
+    where we used that the volume of the $n$-simplex $S$ is $1/n!$ and all inside the integral is constant.
+  \end{proof}
+  \begin{remark}
+    With the same idea, the integrals $b_i$ can be computed as:
+    $$
+      \int_{K_m}f\varphi_{K_m,\ell}=d_m\int_Sf\circ\sigma_m\psi_\ell
+    $$
+    and we use a quadrature formula to approximate over a triangle.
+  \end{remark}
   %   Let $\Omega\subset\RR^2n$ be a bounded domain with a polygonal boundary. Thus, $\Omega$ can be exactly covered by a finite number of triangles. It will be assumed that any pair of triangles in a triangulation of $\Omega$ intersect along a complete edge, at a vertex, or not at all, as shown in Fig. 2.3. We will denote by $h_K$ the diameter (longest side) of the triangle $K$, and we define $h = \max_K h_K$. We define the associate basis functions 
   % \end{definition}
 \end{multicols}