previous up next SYLLABUS  Previous: 3.1 Mathematical background  Up: 3 FINITE ELEMENT METHOD  Next: 3.3 Numerical quadrature


3.2 An engineer's formulation


Slide : [ Weak form - By parts - time - space - Ax=b || VIDEO login]

After a section on what a physicist believes is a mathematician's view of the subject, it is time for an example. Using the format of a ``recipe'' that is applicable to a broad class of practical problems, we show here how the advection-diffusion equation (1.3.2#eq.2) is discretized using Galerkin linear finite elements (FEMs) and how it is implemented.

* Derive a weak variational form.
``Multiply'' your equation by $ \forall g\in \mathcal{C}^1(\Omega), \int_\Omega dx\; g^*(x)$ where the complex conjugate is necessary only if your equation(s) is (are) complex

$\displaystyle \forall g\in \mathcal{C}^1(\Omega),\hspace{1cm} \int_{x_L}^{x_R} ...
...+u\frac{\partial f}{\partial x} -D\frac{\partial^2 f}{\partial x^2} \right] = 0$ (1)

* Integrate by parts,
so as to avoid having to use quadratic basis functions to discretize the second order diffusion operator

$\displaystyle \int_{x_L}^{x_R} dx \;\;\; \left[ g\frac{\partial f}{\partial t} ...
...ft.Dg\frac{\partial f}{\partial x}\right\vert _{x_L}^{x_R} = 0 \;\;\; \forall g$ (2)

Assuming a periodic domain, the surface term has here been canceled, imposing natural boundary conditions.

* Discretize time
with a finite difference backward in time and using a partially implicit evaluation of the unknown $ f=(1-\theta)f^t + \theta f^{t+\Delta t}$ where $ \theta\in[1/2;1]$

$\displaystyle \int_{x_L}^{x_R} dx \left[ g\left(\frac{f^{t+\Delta t}-f^t}{\Delt...
...ial}{\partial x} \left[(1-\theta)f^t + \theta f^{t+\Delta t}\right] \right] = 0$    

$\displaystyle \int_{x_L}^{x_R}\!\!dx \left[ \frac{g}{\Delta t} +\theta ug\frac{...
...!\theta) D \frac{\partial g}{\partial x}\frac{\partial}{\partial x} \right] f^t$ (3)

$ \forall g$ . All the unknowns have been reassembled on the left. Re-scale by $ \Delta t$ , and

* Discretize space
using linear roof-tops and a Galerkin choice for the test functions

$\displaystyle f^t(x)=\sum_{j=1}^N f_j^t e_j(x), \hspace{2cm}\forall g\in\{e_i(x)\},\;\; i=1,N$ (4)


$\displaystyle \int_{x_L}^{x_R}\!\!$ $\displaystyle dx$ $\displaystyle \left[
e_i\sum_{j=1}^N f_j^{t+\Delta t} e_j
+\Delta t\theta e_i u...
...tial x}
\sum_{j=1}^N f_j^{t+\Delta t} \frac{\partial e_j}{\partial x}
\right] =$     
$\displaystyle \int_{x_L}^{x_R}\!\!$ $\displaystyle dx$ $\displaystyle \left[
e_i\sum_{j=1}^N f_j^{t } e_j
+\Delta t(\theta-1) e_i u \su...
... e_i}{\partial x}
\sum_{j=1}^N f_j^{t } \frac{\partial e_j}{\partial x}
\right]$  
     $\displaystyle \forall i=1,N$ (5)

Note how the condition $ \forall g\in \mathcal{C}^1([x_L; x_R])$ is here used to create as many independent equations $ i=1,N$ as there are unknowns $ \{f_j^{t+\Delta t}\}, j=1,N$ . All the essential boundary conditions are then imposed by allowing $ e_1(x)$ and $ e_N(x)$ to overlap in the periodic domain. Since only the basis and test functions $ e_i(x), e_j(x)$ and perhaps the problem coefficients $ u(x), D(x)$ remain space dependent, the discretized equations are all written in terms of inner products, involving overlap integrals of the form $ (e_i, u e_j^\prime)=\int_{x_L}^{x_R}\!\!dx \;\; u(x) e_i(x)e_j^\prime(x)$ . Reassembling them all in a matrix,

* Write a linear system
through which the unknown values from the next time step $ f_j^{t+\Delta t}$ can implicitly be obtained from the current values $ f_j^t$ by solving a linear system

$\displaystyle \sum_{j=1}^N\mathbf{A_{ij}} f_j^{t+\Delta t} = \sum_{j=1}^N\mathbf{B_{ij}} f_j^t$ (6)

To relate the Galerkin linear FEM scheme (3.2#eq.5) with the code that has been implemented in the JBONE applet, it is necessary now to evaluate the overlap integrals. This is usually performed numerically using a quadrature. In the case of a homogeneous mesh $ x_j=j\Delta x$ , a constant advection $ u(x)=u$ and diffusion coefficient $ D(x)=D$ , the coming section shows how explicit expressions can be obtained from the same rules to define directly the matrix elements.

SYLLABUS  Previous: 3.1 Mathematical background  Up: 3 FINITE ELEMENT METHOD  Next: 3.3 Numerical quadrature

      
back up next contents bibliography Copyright © Lifelong-learners at 09:23:40, November 18th, 2017