Consults - Githubissues

Hi,do have a detailed report?

Hi, I do not have a fully detailed report but I can summarize the main parts of the code and relate to the Finite Element theory.

Recall that we obtain the weak formulation for the Poisson (or Laplace if $f=0$) equation as follows:

\forall (x,y)\in\Omega\;,\;-\triangle u = -\nabla\cdot\nabla u = f

multiply along a test function $v$ such that $v=0$ on $\partial\Omega$ and integrate over $\Omega$
```
-\int_\Omega (\triangle u)v = \int_\Omega fv
```

use Green's first identity

\int_\Omega \nabla u \cdot \nabla v  - \int_{\partial\Omega} (v\nabla u)\cdot\boldsymbol{n} = \int_\Omega fv

Apply Dirichlet boundary condition ($u=u_D$ and $v=0$ on $\partial\Omega$)
```
\int_\Omega \nabla u \cdot \nabla v = \int_\Omega fv
```
Then we obtain a linear system to be solved using the finite elements discretization (introducing the shape functions $N$):
```
Ku = f
```
with
```
K_{ij} = \int_\Omega \nabla N_i \cdot \nabla N_j \text{ \quad and \quad } f_{i} = \int_\Omega N_i f
```

femFunctions.py: This file contains the definitions and methods for elements, mesh and solving the linear system

defineRefElement is used to set up the basic the reference element, which can be either a square (Type 0) or a triangle (Type 1). This reference element helps us to build the quadrature, the shape functions and their derivatives in an natural system of coordinates ($\xi,\eta$) and will then allow us to simply use isoparametric formulation to map to general element of the mesh.
Mesh creates the global structure of a cartesian mesh, providing the nodal and connectivity matrices.
elementMatrices compute the elemental matrix $K_i$ and vector $f_i$ associated to a given element using the isoparametric mapping (compute $\partial_x N$ and $\partialy N$ from $\partial\xi N$ and $\partial_\eta N$ using the Jacobian of the transformation $(x,y)=f(\xi,\eta )$ ).
findSolution performs system reduction by removing the columns/rows corresponding to $u$ being on $\partial\Omega$ (hence we already know its value) and then solve the linear system.

computeError computes the error of Finite Elements approximation $u^h$ in $\mathcal{L}^2$ and $\mathcal{H}^1$ norms:

E_{\mathcal{L}^2} = \left.\sqrt{\int_\Omega(u-u^h)^2} \middle/ \sqrt{\int_\Omega u^2}\right.

E_{\mathcal{H}^1} = \left.\sqrt{\int_\Omega(\nabla u-\nabla u^h)^2} \middle/ \sqrt{\int_\Omega (\nabla u)^2}\right.

to check the order of convergence depending on the element type.

plot_mesh and plot_errors --> post-processing tools.

FEM_2D.py: Main script to run the example, note that the assembly of elemental matrices $K_i$ and vectors $f_i$ into global stiffness matrix $K$ and force vector $f$ is performed here.

font{
    line-height: 1.6;
}
ul,ol{
    padding-left: 20px;
    list-style-position: inside;
}

Thank you!

                            2074375758

                                ***@***.***

---- Replied Message ----

     From 

        Hugo ***@***.***>

     Date 

    06/13/2024 17:24

     To 

        ***@***.***>

     Cc 

        ***@***.***>
        ,

        ***@***.***>

     Subject 

          Re: [hvalayer/fem_laplace_2D] Consults (Issue #1)

Hi, I do not have a fully detailed report but I can summarize the main parts of the code and relate to the Finite Element theory. Recall that we obtain the weak formulation for the Poisson (or Laplace if $f=0$) equation as follows: $$\forall (x,y)\in\Omega\;,\;\triangle u = \nabla\cdot\nabla u = f$$

multiply along a test function $v$ such that $v=0$ on $\partial\Omega$ and integrate over $\Omega$

$$\int\Omega v \nabla\cdot\nabla u = \int\Omega vf$$

use integration by parts

$$\int\Omega \nabla u \cdot \nabla v = \int\Omega vf - v \int_{\partial\Omega} \nabla u$$

Apply Dirichlet boundary condition ($u=u_D$ and $v=0$ on $\partial\Omega$)

$$\int\Omega \nabla u \cdot \nabla v = \int\Omega vf$$ Then we obtain a linear system to be solved using the finite elements discretization (introducing the shape functions $N$): $$Ku = f$$ with $$K{ij} = \int\Omega \nabla N_i \cdot \nabla Nj \text{ \quad and \quad } f{i} = \int_\Omega N_i f$$ femFunctions.py: This file contains the definitions and methods for elements, mesh and solving the linear system

defineRefElement is used to set up the basic the reference element, which can be either a square (Type 0) or a triangle (Type 1). This reference element helps us to build the quadrature, the shape functions and their derivatives in an natural system of coordinates ($\xi,\eta$) and will then allow us to simply use isoparametric formulation to map to general element of the mesh.

Mesh creates the global structure of a cartesian mesh, providing the nodal and connectivity matrices.

elementMatrices compute the elemental matrix $K_i$ and vector $f_i$ associated to a given element using the isoparametric mapping (compute $\partial_x N$ and $\partialy N$ from $\partial\xi N$ and $\partial_\eta N$ using the Jacobian of the transformation $(x,y)=f(\xi,\eta )$ ).

findSolution performs system reduction by removing the columns/rows corresponding to $u$ being on $\Omega_D$ (hence we already know its value) and then solve the linear system.

computeError computes the error of Finite Elements approximation $u^h$ in $\mathcal{L}^2$ and $\mathcal{H}^1$ norms:

$$E{\mathcal{L}^2} = \left.\sqrt{\int\Omega(u-u^h)^2} \middle/ \sqrt{\int\Omega u^2}\right.$$ $$E{\mathcal{H}^1} = \left.\sqrt{\int\Omega(\nabla u-\nabla u^h)^2} \middle/ \sqrt{\int\Omega (\nabla u)^2}\right.$$ to check the order of convergence depending on the element type.

plot_mesh and plot_errors --> post-processing tools.

FEM_2D.py: Main script to run the example, note that the assembly of elemental matrices $K_i$ and vectors $f_i$ into global stiffness matrix $K$ and force vector $f$ is performed here.

—Reply to this email directly, view it on GitHub, or unsubscribe.You are receiving this because you authored the thread.Message ID: @.***>

hvalayer / fem_laplace_2D

Consults #1