Skip to content

Commit

Permalink
update
Browse files Browse the repository at this point in the history
  • Loading branch information
ratnania committed Feb 27, 2024
1 parent 2cb3a70 commit 11dab4e
Show file tree
Hide file tree
Showing 7 changed files with 169 additions and 109 deletions.
1 change: 1 addition & 0 deletions _toc.yml
Original file line number Diff line number Diff line change
Expand Up @@ -33,6 +33,7 @@ parts:
- file: chapter1/space
- file: chapter1/rules
- file: chapter1/poisson
- file: chapter1/boundary-conditions
- caption: Linear Problems
chapters:
- file: chapter2/poisson
Expand Down
64 changes: 64 additions & 0 deletions chapter1/boundary-conditions.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,64 @@
# Boundary Conditions
*Author: Ahmed Ratnani*

SymPDE & Psydac allows you to use both strong and weak boundary conditions.
We start first by explaining how to identify a boundary in **NCube** domains, such as **Line**, **Square** and a **Cube**.

## Identifying a boundary in **NCube** domains

Since we're dealing with **NCube** domains, the best way to identify a boundary is to use the couple **(axis, extremity)**. The axis is usually defined as perpendicular to the boundary, while the **extremity** gives the lower or upper bound of the interval associated to the **axis**.

Let's see it through the following examples;

For a **Line**, the following figure gives the value to access a boundary.

![png](images/boundary-conditions/line-boundaries.png)

For a **Square**, we have,

![png](images/boundary-conditions/square-boundaries.png)

For a **Cube**, it is similar and would not be helpful to visualize it here.

## Essential boundary conditions

Essential boundary conditions are usually treated in two ways, either in the strong or weak form. The latter can be achieved through the use of Nitsche's method. You can check the associated examples in the sequel.

For the essential boundary conditions, usually, one needs one of the following expressions;

$$
\begin{align}
u &= 0, \quad \text{(homogeneous case)}, \\
u &= f, \quad \text{(inhomogeneous case)}, \\
\mathbf{v} &= 0, \quad \text{(homogeneous case)}, \\
\mathbf{v} &= \mathbf{g}, \quad \text{(homogeneous case)},
\end{align}
$$

where $u$ and $f$ are scalar functions while $\mathbf{v}$ and $\mathbf{g}$ denote vector functions.


Let `Gamma` denotes the boundary $\Gamma$.
The normal derivative is an available operator in SymPDE, it is provided as `Dn` operator. Alternatively, you can also define it manually,

```python
nn = NormalVector('nn')
dn = lambda a: dot(grad(a), nn)
```


we have the following equivalences,

| Mathematical Expression | SymPDE Expression | Example |
| --------------------------------------------------------- | ------------------------------ | ---------------------------------- |
| $u = 0$ on $\Gamma$ | `EssentialBC(u, 0, Gamma)` | |
| $u = f$ on $\Gamma$ | `EssentialBC(u, f, Gamma)` | |
| $\partial_n u = 0$ on $\Gamma$ | `EssentialBC(dn(u), 0, Gamma)` | |
| $\partial_n u = f$ on $\Gamma$ | `EssentialBC(dn(u), f, Gamma)` | |
| $\mathbf{v} = 0$ on $\Gamma$ | `EssentialBC(v, 0, Gamma)` | |
| $\mathbf{v} = \mathbf{g}$ on $\Gamma$ | `EssentialBC(v, g, Gamma)` | |
| $\mathbf{v} \cdot \mathbf{n} = 0$ on $\Gamma$ | ? | |
| $\mathbf{v} \times \mathbf{n} = 0$ on $\Gamma$ | ? | |
| $\mathbf{v} \cdot \mathbf{n} = f,\mathbf{g}$ on $\Gamma$ | ? | |
| $\mathbf{v} \times \mathbf{n} = f,\mathbf{g}$ on $\Gamma$ | ? | |

Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
93 changes: 56 additions & 37 deletions chapter2/poisson-nitsche.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -18,8 +18,8 @@
"\n",
"$$\n",
"\\begin{align}\n",
" - \\nabla^2 u = f \\quad \\text{in~$\\Omega$}, \\quad \\quad \n",
" u = g \\quad \\text{on~$\\partial \\Omega$}.\n",
" - \\nabla^2 u = f \\quad \\text{in $\\Omega$}, \\quad \\quad \n",
" u = g \\quad \\text{on $\\partial \\Omega$}.\n",
"\\end{align}\n",
"$$"
]
Expand Down Expand Up @@ -54,12 +54,12 @@
},
{
"cell_type": "code",
"execution_count": 2,
"execution_count": 18,
"id": "ecd1330a",
"metadata": {},
"outputs": [],
"source": [
"from sympde.expr import BilinearForm, LinearForm, integral\n",
"from sympde.expr import BilinearForm, LinearForm, integral, Norm\n",
"from sympde.expr import find, EssentialBC\n",
"from sympde.topology import ScalarFunctionSpace, Square, element_of\n",
"from sympde.calculus import grad, dot, Dn\n",
Expand All @@ -80,16 +80,19 @@
"u,v = [element_of(V, name=i) for i in ['u', 'v']]\n",
"\n",
"# bilinear form\n",
"a = BilinearForm((u,v), integral(domain, dot(grad(v), grad(u)) + kappa*u*v - u*Dn(v) - v*Dn(u)))\n",
"a = BilinearForm((u,v), integral(domain, dot(grad(v), grad(u))))\n",
"an = BilinearForm((u,v), integral(domain.boundary, kappa*u*v - u*Dn(v) - v*Dn(u)))\n",
"\n",
"# linear form\n",
"g = sin(pi*x)*sin(pi*y)\n",
"f = 2*pi**2*sin(pi*x)*sin(pi*y)\n",
"ue = sin(pi*x)*sin(pi*y)\n",
"g = ue\n",
"f = 2*pi**2*sin(pi*x)*sin(pi*y)\n",
"\n",
"l = LinearForm(v, integral(domain, f*v + kappa*g*v - g*Dn(v)))\n",
"l = LinearForm(v, integral(domain, f*v))\n",
"ln = LinearForm(v, integral(domain.boundary, kappa*g*v - g*Dn(v)))\n",
"\n",
"# Variational problem\n",
"equation = find(u, forall=v, lhs=a(u, v), rhs=l(v))"
"equation = find(u, forall=v, lhs=a(u,v) + an(u,v), rhs=l(v) + ln(v))"
]
},
{
Expand All @@ -102,7 +105,7 @@
},
{
"cell_type": "code",
"execution_count": 3,
"execution_count": 19,
"id": "3d86674b",
"metadata": {},
"outputs": [],
Expand All @@ -113,7 +116,7 @@
},
{
"cell_type": "code",
"execution_count": 6,
"execution_count": 20,
"id": "94c83267",
"metadata": {},
"outputs": [],
Expand All @@ -125,7 +128,7 @@
"Vh = discretize(V, domain_h, degree=degree)\n",
"\n",
"# Discretize equation\n",
"#equation_h = discretize(equation, domain_h, [Vh, Vh])"
"equation_h = discretize(equation, domain_h, [Vh, Vh])"
]
},
{
Expand All @@ -138,53 +141,69 @@
},
{
"cell_type": "code",
"execution_count": 7,
"execution_count": 21,
"id": "8f1b2d68",
"metadata": {},
"outputs": [],
"source": [
"#uh = equation_h.solve()"
"equation_h.set_solver('gmres', info=False, tol=1e-8)"
]
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 22,
"id": "e96f321c",
"metadata": {},
"outputs": [],
"source": []
"source": [
"uh = equation_h.solve(kappa=1e3)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "c7ae0b9e",
"cell_type": "markdown",
"id": "d87fd79f",
"metadata": {},
"outputs": [],
"source": []
"source": [
"## Computing the error norm"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "9ad386a6",
"cell_type": "markdown",
"id": "08ee5cac",
"metadata": {},
"outputs": [],
"source": []
"source": [
"### Computing the $L^2$ norm"
]
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 23,
"id": "068bdb95",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
"id": "aa88a918",
"metadata": {},
"outputs": [],
"source": []
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"0.00021402394834116976\n"
]
}
],
"source": [
"u = element_of(V, name='u')\n",
"\n",
"# create the formal Norm object\n",
"l2norm = Norm(u - ue, domain, kind='l2')\n",
"\n",
"# discretize the norm\n",
"l2norm_h = discretize(l2norm, domain_h, Vh)\n",
"\n",
"# assemble the norm\n",
"l2_error = l2norm_h.assemble(u=uh)\n",
"\n",
"# print the result\n",
"print(l2_error)"
]
}
],
"metadata": {
Expand Down
Loading

0 comments on commit 11dab4e

Please sign in to comment.