Choral
Functions/Subroutines
diffusion_Neumann.f90 File Reference

Go to the source code of this file.

Functions/Subroutines

program diffusion_neumann
  SOLVES THE DIFFUSION PROBLEM with Neumann boundary conditions More...
 
real(rp) function u (x)
 exact solution 'u' More...
 
real(rp) function, dimension(3) grad_u (x)
 gradient of the exact solution More...
 
real(rp) function f (x)
 right hand side 'f' More...
 
real(rp) function rho (x)
 density '(x)' More...
 
real(rp) function bl_form_b (x, w1, w2)
 
real(rp) function theta (x)
 Angle of the principal direction with e_x. More...
 
real(rp) function, dimension(2, 2) b (x)
 Diffusion matrix B(x) at point x. More...
 
real(rp) function x_le_0 (x)
 To characterise {x=0}. More...
 
real(rp) function x_ge_1 (x)
 To characterise {x=1}. More...
 
real(rp) function y_le_0 (x)
 To characterise {x=0}. More...
 
real(rp) function y_ge_1 (x)
 To characterise {x=1}. More...
 
real(rp) function g_x_0 (x)
 g(x) on {x=0} More...
 
real(rp) function g_x_1 (x)
 g(x) on {x=1} More...
 
real(rp) function g_y_0 (x)
 g(x) on {y=0} More...
 
real(rp) function g_y_1 (x)
 g(x) on {y=1} More...
 

Function/Subroutine Documentation

◆ b()

real(rp) function, dimension(2,2) diffusion_neumann::b ( real(rp), dimension(3), intent(in)  x)

Diffusion matrix B(x) at point x.

Definition at line 432 of file diffusion_Neumann.f90.

Here is the call graph for this function:

◆ bl_form_b()

real(rp) function diffusion_neumann::bl_form_b ( real(rp), dimension(3), intent(in)  x,
real(rp), dimension(3), intent(in)  w1,
real(rp), dimension(3), intent(in)  w2 
)

Definition at line 404 of file diffusion_Neumann.f90.

Here is the call graph for this function:

◆ diffusion_neumann()

program diffusion_neumann ( )

SOLVES THE DIFFUSION PROBLEM with Neumann boundary conditions

\(~~~~~~~~~ -\dv(B(x) \nabla u) + \rho(x) u = f ~~~\) on \(~~~ \Omega= [0,1]^2 \)

\(~~~~~~~~~~ B(x)\nabla u.n = g ~~~\) on \(~~~\Gamma = \partial \Omega ~~~\) with

  • \( g \) = g_x_0 on \(\Gamma \cap \{x=0\} \)
  • \( g \) = g_x_1 on \(\Gamma \cap \{x=1\} \)
  • \( g \) = g_y_0 on \(\Gamma \cap \{y=0\} \)
  • \( g \) = g_y_1 on \(\Gamma \cap \{y=1\} \)



DIFFUSIVITY TENSOR \( B(x) = P(x) D ^TP(x) \) with \( D = \left(\begin{array}{cc} L_1 & 0 \\ 0 & L_2 \end{array}\right)\)
and \( P(x) \) is the rotation matrix with angle \( \theta = \theta(x)\).
See choral/maple/example_diffusion_Neumann.mw for the computation of the problem data.



VARIATIONAL FORMULATION

  • The finite element space is \( X_h\subset \Hu\).

Numerical probelm: find \(u\in X_h\) such that \( ~~~~\forall ~v \in X_h, ~~~~ \)

\[ \int_\Omega B(x) \nabla u\cdot\nabla v \,\dx ~+~ \int_\Omega \rho(x)\, u\, v \,\dx ~=~ \int_\Omega f \, v \,\dx ~+~ \int_\Gamma g \, v \,\dx \]



NUMERICAL RESOLUTION

The basis functions are \( (v_i)_{1\le i\le N} \) (basis of \(X_h\)).

  • Computation of the stiffness matrix

    \[ S =[s_{i,\,j}]_{1\le i,\,j\le N}, \quad \quad s_{i,\,j} = \int_\Omega B(x) \nabla v_i\cdot\nabla v_j \,\dx \]

  • Computation of the mass matrix

    \[ M =[s_{i,\,j}]_{1\le i,\,j\le N}, \quad \quad m_{i,\,j} = \int_\Omega \rho(x)\, v_i\, v_j \,\dx \]

  • Computation of the right hand side for the volumoc source term \( f \)

    \[ F = (f_i)_{1\le i\le N}, \quad \quad f_i = \int_\Omega f \, v_i \,\dx \]

  • Computation of the right hand side for the surfacic source term \( g \)

    \[ G = (g_i)_{1\le i\le N}, \quad \quad g_i = \int_\Gamma g \, v_i \,\dx \]

  • Resolution of the (symmetric positive definite) system

    \[ (M+S) U_h = F + G \]



POST TREATMENT

  • Computation of the numerical error between the exact solution \(u\) and the numerical solution \(u_h\),

    \[ \int_\Omega |u-u_h|^2 \dx~,\quad\quad \int_\Omega |\nabla u-\nabla u_h|^2 \dx \]

  • Graphical display of the numerical solution with gmsh

Charles PIERRE, October 2019

Definition at line 84 of file diffusion_Neumann.f90.

Here is the call graph for this function:

◆ f()

real(rp) function diffusion_neumann::f ( real(rp), dimension(3), intent(in)  x)

right hand side 'f'

Definition at line 356 of file diffusion_Neumann.f90.

Here is the call graph for this function:

◆ g_x_0()

real(rp) function diffusion_neumann::g_x_0 ( real(rp), dimension(3), intent(in)  x)

g(x) on {x=0}

Definition at line 499 of file diffusion_Neumann.f90.

Here is the call graph for this function:

◆ g_x_1()

real(rp) function diffusion_neumann::g_x_1 ( real(rp), dimension(3), intent(in)  x)

g(x) on {x=1}

Definition at line 516 of file diffusion_Neumann.f90.

Here is the call graph for this function:

◆ g_y_0()

real(rp) function diffusion_neumann::g_y_0 ( real(rp), dimension(3), intent(in)  x)

g(x) on {y=0}

Definition at line 527 of file diffusion_Neumann.f90.

Here is the call graph for this function:

◆ g_y_1()

real(rp) function diffusion_neumann::g_y_1 ( real(rp), dimension(3), intent(in)  x)

g(x) on {y=1}

Definition at line 544 of file diffusion_Neumann.f90.

Here is the call graph for this function:

◆ grad_u()

real(rp) function, dimension(3) diffusion_neumann::grad_u ( real(rp), dimension(3), intent(in)  x)

gradient of the exact solution

Definition at line 343 of file diffusion_Neumann.f90.

◆ rho()

real(rp) function diffusion_neumann::rho ( real(rp), dimension(3), intent(in)  x)

density '(x)'

Definition at line 392 of file diffusion_Neumann.f90.

◆ theta()

real(rp) function diffusion_neumann::theta ( real(rp), dimension(3), intent(in)  x)

Angle of the principal direction with e_x.

Definition at line 422 of file diffusion_Neumann.f90.

◆ u()

real(rp) function diffusion_neumann::u ( real(rp), dimension(3), intent(in)  x)

exact solution 'u'

Definition at line 332 of file diffusion_Neumann.f90.

◆ x_ge_1()

real(rp) function diffusion_neumann::x_ge_1 ( real(rp), dimension(3), intent(in)  x)

To characterise {x=1}.

Definition at line 467 of file diffusion_Neumann.f90.

◆ x_le_0()

real(rp) function diffusion_neumann::x_le_0 ( real(rp), dimension(3), intent(in)  x)

To characterise {x=0}.

Definition at line 456 of file diffusion_Neumann.f90.

◆ y_ge_1()

real(rp) function diffusion_neumann::y_ge_1 ( real(rp), dimension(3), intent(in)  x)

To characterise {x=1}.

Definition at line 488 of file diffusion_Neumann.f90.

◆ y_le_0()

real(rp) function diffusion_neumann::y_le_0 ( real(rp), dimension(3), intent(in)  x)

To characterise {x=0}.

Definition at line 477 of file diffusion_Neumann.f90.