Project file on GitLab
## Motivation of this test case

## Problem description

## Weak form

## Analytical solution

## Results

The aim of this test is:

- to provide a simple introductory example of how Python BCs can be used
- to show that Python BCs can be used to easily prescribe BCs based on analytical formulas
- to check whether both essential and natural BCs are implemented correctly in OpenGeoSys’ Python BC.

We solve Laplace’s Equation in 2D on a $1 \times 1$ square domain. The boundary conditions will be specified below.

Laplace’s equation is $$ \begin{equation} - \mathop{\mathrm{div}} (a \mathop{\mathrm{grad}} u) = 0 \end{equation} $$ The weak form is derived as usual by multiplying with a test function $v$ and integrating over the domain $\Omega$: $$ \begin{equation} - \int_{\Omega} v \mathop{\mathrm{div}} (a \mathop{\mathrm{grad}} u) , \mathrm{d}\Omega = 0 ,, \end{equation} $$ which can be transformed further to $$ \begin{equation} \int_{\Omega} a \mathop{\mathrm{grad}} v \cdot \mathop{\mathrm{grad}} u , \mathrm{d}\Omega = \int_{\Omega} \mathop{\mathrm{div}} (v a \mathop{\mathrm{grad}} u) , \mathrm{d}\Omega = \int_{\Gamma_{\mathrm{N}}} v a \mathop{\mathrm{grad}} u \cdot n , \mathrm{d}\Gamma ,, \end{equation} $$

where in the second equality Gauss’s theorem has been applied.
As usual, the domain boundary $\partial\Omega = \Gamma_{\mathrm{D}} \cup \Gamma_{\mathrm{N}}$ is subdivided
into the dirichlet and the Neumann boundary and $v$ vanishes on
$\Gamma_{\mathrm{D}}$.
The r.h.s. of the above equation is the total flux associated with $u$ flowing
**into** the domain $\Omega$ through $\Gamma_{\mathrm{N}}$:
$-a \mathop{\mathrm{grad}} u$ is the flux density and $-n$ is the inwards directed surface
normal.

The weak form just derived is implemented (after FEM discretization) in the groundwater flow process in OpenGeoSys. Note that for the application of Neumann boundary conditions, it is necessary to know whether the flux has to be applied with a positive or a negative sign!

The coefficient $a$ of Laplace’s equation is taken to be unity. By differentiation it can be easily checked that $$ \begin{equation} u(x, y) = \sin(bx) \sinh(by) \end{equation} $$ solves Laplace’s equation inside $\Omega$ for any $b$. In this example we set $b = \tfrac 23 \pi$.

As boundary conditions we apply Dirichlet BCs at the top, left and bottom of the domain with values from $u(x,y)|_{\Gamma_{\mathrm{D}}}$. On the right boundary of the domain a Neumann BC is applied. There $n = (1, 0)$, which implies that $a \mathop{\mathrm{grad}} u \cdot n = a , \partial u / \partial x$.

The numerical result obtained from OpenGeoSys is:

The absolute difference between the analytical and numerical solutions is smaller than $4 \cdot 10^{-4}$:

This article was written by Christoph Lehmann. If you are missing something or you find an error please let us know.
Generated with Hugo 0.96.0. Last revision: February 16, 2022 *
Commit: explanation stress concentration and running model 5c94789
| Edit this page on
*