Choral
Data Types | Functions/Subroutines
fespace_mod Module Reference

DERIVED TYPE feSpace: finite element spaces More...

Data Types

interface  assemble
 
interface  clear
 
interface  fespace
 Derived type feSpace: finite element space on a mesh. More...
 
interface  interp_vect_func
 
interface  print
 print a short description More...
 
interface  set
 
interface  valid
 
interface  write
 write a finite element space to a file More...
 

Functions/Subroutines

subroutine fespace_clear (X_h)
 Destructor for feSpace type. More...
 
type(fespace) function fespace_create (m)
  Constructor for the type feSpace More...
 
subroutine fespace_set (X_h, ft)
 Apply the finite element method 'ft' to all cells with comptible geometry. More...
 
subroutine fespace_print (X_h)
 Print a short description. More...
 
logical function fespace_valid (X_h)
 Check if the structure content is correct. More...
 
subroutine fespace_write (X_h, fileName, fileFormat)
 Write feSpace to file 'fileName' with format 'fileFormat'. More...
 
subroutine gmsh_write (X_h, fileName)
 Write feSpace m to file 'fileName' gmsh ASCII format. More...
 
subroutine, public gmsh_addview (X_h, v, fileName, vname, time, idx)
 add a view to X_h output file gmesh file format More...
 
subroutine, public interp_scal_fe (u2, X_h2, u1, X_h1)
 Interpolate a scalar finite element function u1 on a second finite element space X_h2. More...
 
subroutine, public eval_scal_fe_at_points (val, P, u, X_h)
 u is a finite element function on the finite element space X_h More...
 
subroutine locate_nodes (T, P, X_h, ndToNd)
 Locate each point P(:,i) inside a cell of X_hmesh associated to a non void finite element. More...
 
integer function closest_cell (x, X_h)
 finds the cell in X_h with non void feType the closest from node x More...
 
subroutine fespace_assemble (X_h)
 sets cell dof local indexing –> global dof indexing connectivity More...
 
subroutine x_h_dofcoord (X_h)
 
subroutine initialise_doftab (dofTab, cl_row, X_h)
 
subroutine scan_dof (dofTab, X_h, wdt)
 Checks every local dof to say wether it is equal to dof 'dof2' of cell 'cl2' for cl2 < cell. More...
 
subroutine build_cltodof (dofTab, cl_row, X_h)
 
subroutine get_perm_3 (p, t1, t2)
 t2 is a permutation of t1 More...
 
subroutine get_perm_2 (p, t1, t2)
 t2 is a permutation of t1 More...
 
subroutine, public interp_scal_func (uh, u, X_h)
 Interpolation of a scalar function to a scalar finite element function. More...
 
subroutine fespace_interp_vect_func (phi_h, phi, X_h)
 Interp : vector function case. More...
 
subroutine, public flag_x_h_dof (dof_flag, X_h, dim, f)
 OUTPUT: flag(:) array of logical with size X_hnbDof. More...
 

Detailed Description

DERIVED TYPE feSpace: finite element spaces

Given a mesh \( \mathcal{T}\) of a domain \( \Omega \):

Let 'X_h' of type feSpace:

Tutorial examples:

Author
Charles Pierre

Function/Subroutine Documentation

◆ build_cltodof()

subroutine fespace_mod::build_cltodof ( integer, dimension(:,:), intent(inout)  dofTab,
integer, dimension(:), intent(in)  cl_row,
type(fespace), intent(inout)  X_h 
)
private

Definition at line 1507 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ closest_cell()

integer function fespace_mod::closest_cell ( real(rp), dimension(3), intent(in)  x,
type(fespace), intent(in)  X_h 
)
private

finds the cell in X_h with non void feType the closest from node x

distance = sum_i^nbVtx |x-Ai|_1/nbVtx, A_i vertex of cell cl

Definition at line 902 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ eval_scal_fe_at_points()

subroutine, public fespace_mod::eval_scal_fe_at_points ( real(rp), dimension(:), allocatable  val,
real(rp), dimension(:,:), intent(in)  P,
real(rp), dimension(:), intent(in)  u,
type(fespace), intent(in)  X_h 
)

u is a finite element function on the finite element space X_h

P is a set of points

OUTPUT : val = values of the function u at the points P: val(i) = u( P(:,i))

INPUT : u = finite element function X_h = finite element space for u1 P = set of points

It is assumed the points in P are in the mesh domain !

Definition at line 605 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ fespace_assemble()

subroutine fespace_mod::fespace_assemble ( type(fespace), intent(inout)  X_h)
private

sets cell dof local indexing –> global dof indexing connectivity

Definition at line 952 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ fespace_clear()

subroutine fespace_mod::fespace_clear ( type(fespace), intent(inout)  X_h)
private

Destructor for feSpace type.

Definition at line 180 of file feSpace_mod.F90.

◆ fespace_create()

type(fespace) function fespace_mod::fespace_create ( type(mesh), target  m)
private

Constructor for the type feSpace

  • OUTPUT:
    • X_h = feSpace, finite element space

Definition at line 206 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ fespace_interp_vect_func()

subroutine fespace_mod::fespace_interp_vect_func ( real(rp), dimension(:), allocatable  phi_h,
procedure(r3tor3)  phi,
type(fespace), intent(in)  X_h 
)
private

Interp : vector function case.

To be reconsidered

Definition at line 1675 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ fespace_print()

subroutine fespace_mod::fespace_print ( type(fespace), intent(in)  X_h)
private

Print a short description.

Definition at line 256 of file feSpace_mod.F90.

◆ fespace_set()

subroutine fespace_mod::fespace_set ( type(fespace), intent(inout)  X_h,
integer, intent(in)  ft 
)
private

Apply the finite element method 'ft' to all cells with comptible geometry.

Definition at line 225 of file feSpace_mod.F90.

◆ fespace_valid()

logical function fespace_mod::fespace_valid ( type (fespace), intent(in)  X_h)
private

Check if the structure content is correct.

Definition at line 286 of file feSpace_mod.F90.

◆ fespace_write()

subroutine fespace_mod::fespace_write ( type(fespace), intent(in)  X_h,
character(len=*), intent(in)  fileName,
character(len=*), intent(in)  fileFormat 
)
private

Write feSpace to file 'fileName' with format 'fileFormat'.

Definition at line 308 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ flag_x_h_dof()

subroutine, public fespace_mod::flag_x_h_dof ( logical, dimension(:), allocatable  dof_flag,
type(fespace), intent(in)  X_h,
integer, intent(in)  dim,
procedure(r3tor), optional  f 
)

OUTPUT: flag(:) array of logical with size X_hnbDof.

for the dof 'i' in the finire element space 'X_h':
flag(i) = .TRUE. iif:

  • the dof 'i' belongs to a cell 'cl' of dimension 'dim' in the mesh 'X_hmesh'
  • f(x) >=0 at all node 'x' of the cell (if the argument 'f' is provided)

Definition at line 1840 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ get_perm_2()

subroutine fespace_mod::get_perm_2 ( integer, intent(out)  p,
integer, dimension(2), intent(in)  t1,
integer, dimension(2), intent(in)  t2 
)
private

t2 is a permutation of t1

returns p such that :

t1(1) = t2( p )

Definition at line 1624 of file feSpace_mod.F90.

◆ get_perm_3()

subroutine fespace_mod::get_perm_3 ( integer, dimension(2), intent(out)  p,
integer, dimension(3), intent(in)  t1,
integer, dimension(3), intent(in)  t2 
)
private

t2 is a permutation of t1

returns p(:) such that :

t1(i) = t2( p(i) ) for i=1 .. 2

Definition at line 1586 of file feSpace_mod.F90.

◆ gmsh_addview()

subroutine, public fespace_mod::gmsh_addview ( type(fespace), intent(in)  X_h,
real(rp), dimension(:), intent(in)  v,
character(len=*), intent(in)  fileName,
character(len=*), intent(in)  vname,
real(rp), intent(in)  time,
integer, intent(in)  idx 
)

add a view to X_h output file gmesh file format

Definition at line 412 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ gmsh_write()

subroutine fespace_mod::gmsh_write ( type(fespace), intent(in)  X_h,
character(len=*), intent(in)  fileName 
)
private

Write feSpace m to file 'fileName' gmsh ASCII format.

Definition at line 328 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ initialise_doftab()

subroutine fespace_mod::initialise_doftab ( integer, dimension(:,:), intent(inout), allocatable  dofTab,
integer, dimension(:), intent(out)  cl_row,
type(fespace), intent(inout)  X_h 
)
private

Definition at line 1044 of file feSpace_mod.F90.

◆ interp_scal_fe()

subroutine, public fespace_mod::interp_scal_fe ( real(rp), dimension(:), intent(out)  u2,
type(fespace), intent(in)  X_h2,
real(rp), dimension(:), intent(in)  u1,
type(fespace), intent(in)  X_h1 
)

Interpolate a scalar finite element function u1 on a second finite element space X_h2.

OUTPUT : u2 = finite element function on the finite element space X_h2

INPUT : X_h1 = finite element space u1 = finite element function on X_h1 X_h = second finite element space

OUTPUT: u2 = finite element function on X_h2 interpolant of u1 at the dof of X_h2

It is assumed that X_h1 and X_h2 are Lagrangian spaces here. Consequently, u2 is given by the values of u1 at the nodes of X_h2.

Definition at line 468 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ interp_scal_func()

subroutine, public fespace_mod::interp_scal_func ( real(rp), dimension(:), allocatable  uh,
procedure(r3tor)  u,
type(fespace), intent(in)  X_h 
)

Interpolation of a scalar function to a scalar finite element function.

Definition at line 1641 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ locate_nodes()

subroutine fespace_mod::locate_nodes ( integer, dimension(:), intent(inout)  T,
real(rp), dimension(:,:), intent(in)  P,
type(fespace), intent(in)  X_h,
type(graph), intent(in)  ndToNd 
)
private

Locate each point P(:,i) inside a cell of X_hmesh associated to a non void finite element.

initially, T(ii) is the closest node of X_hmesh obtained with :

call closest_node(T, dist, P, X_hmesh, ndToNd)

Input = ndToNd = node to node connectivity of the mesh X_hmesh

Definition at line 732 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ scan_dof()

subroutine fespace_mod::scan_dof ( integer, dimension(:,:), intent(inout)  dofTab,
type(fespace), intent(inout)  X_h,
integer, intent(in)  wdt 
)
private

Checks every local dof to say wether it is equal to dof 'dof2' of cell 'cl2' for cl2 < cell.

In that case dofTab(1:2, dof) = (/cl2, dof2/)

Definition at line 1091 of file feSpace_mod.F90.

Here is the call graph for this function:

◆ x_h_dofcoord()

subroutine fespace_mod::x_h_dofcoord ( type(fespace), intent(inout)  X_h)
private

Definition at line 990 of file feSpace_mod.F90.