Example 2: The DofHandler concept

This example is meant to be a practical application of the theoretical preliminaries given in the section called “Defining finite element spaces in BETL”. This tutorial example is directly build upon the previous example stated in the section called “Example 1: Processing input files and exporting vtu-files” and its complete source code can be found at $BETL_ROOT/tutorial/example_2/main.cpp.

In detail, this section deals with

For the creation of the finite element space the two following header files need to be invoked:

#include "febasis/febasis.hpp"
#include "dof_handler/dof_handler.hpp"

We will create two bases for two different finite element spaces. One space will be made up the lowest order Lagrangian based functions. The other space will be based upon the lowest order edge functions. While the first space will be approximated discontinuously the latter one will represent a continuous approximation. The declarations of the bases are

// define a constant, discontinuous lagrange febasis
typedef FEBasis< element_t, 
                 CONSTANT,
                 Discontinuous,
                 LagrangeTraits > feb_e3_const_disc_lagr_t;
// define a liner, discontinuous div-febasis
typedef FEBasis< element_t, 
                 LINEAR,
                 Continuous,
                 EdgeDivTraits > feb_e3_lin_cont_div_t;

where element_t corresponds to the flat triangular element type betl::Element<3>. With the declaration of the FE basis the declaration and instantiation of the dof handler is straightforward

// define two dofhandler objects
typedef DoFHandler< feb_e3_const_disc_lagr_t > dofhandler_lagr_t;
typedef DoFHandler< feb_e3_lin_cont_div_t    > dofhandler_div_t;
// instantiate dofhandler objects
dofhandler_lagr_t dofhandler_lagr;
dofhandler_div_t  dofhandler_div;

For the complete setup of the finite element space it remains to distribute the degrees of freedom on a particular range of elements. Here, the degrees of freedom are distributed on the complete mesh

// get element iterators
typedef mesh_t::const_element_iterator const_iterator;
const_iterator begin = mesh.e_begin( );
const_iterator end   = mesh.e_end( );
// distribute dofs
dofhandler_lagr.distributeDoFs( begin, end );
dofhandler_div.distributeDoFs( begin, end );

That's it. The method distributeDoFs finalises the construction of the finite element space.

From time to time it might be helpful to display some debug data. This debug data might either be displayed directly on screen or via a file-stream into a textfile. In BETL, this output is usually performed by calling the overloaded stream operator (operator<<). Its use is illustrated in the following example:

// debug output (on screen)
std::cout << dofhandler_lagr << std::endl;
// debug output (to file)
const std::string filename = basename + "_DH_DIV_DEBUG.dat";
std::ofstream dh_div_out( filename.c_str( ) );
dh_div_out << dofhandler_div << std::endl;
dh_div_out.close( );

As for the former example this example can be build via make betl_ex2. When calling betl_ex2 tetrahedron_4.msh the additional on-screen output should be:

<<<-------------------------------------------------------------------
                              DoFHandler                              
----------------------------------------------------------------------
  distribute dofs: 4 dofs created.
------------------------------------------------------------------->>>

<<<-------------------------------------------------------------------
                              DoFHandler                              
----------------------------------------------------------------------
  distribute dofs: 6 dofs created.
------------------------------------------------------------------->>>
==============================================================
       LIST OF DEGREES OF FREEDOM AND THEIR SUPPORT           
             (represented by element indices)                 
==============================================================
Index = 0, Support = 0  
Index = 1, Support = 1  
Index = 2, Support = 2  
Index = 3, Support = 3  
==============================================================
  LIST OF ELEMENTS WITH THEIR ASSOCIATED DEGREES OF FREEDOM   
             (represented by dof indices)                     
==============================================================
Ele-Idx = 0, Dofs = 0 
Ele-Idx = 1, Dofs = 1 
Ele-Idx = 2, Dofs = 2 
Ele-Idx = 3, Dofs = 3 

The output between the <<<- and the >>>-signs features the standard messages generated by the distributeDoFs-calls. For the given mesh, four degrees of freedom have been created by the Lagrangian dof handler type while the edge-based finite element space consists of six degrees of freedom. The remainder of the on-screen output is nothing but the Lagrangian dof handler's debug output. The debug output for the edge dof handler can be found in the generated textfile tetrahedron_4_DH_DIV_DEBUG.dat. Its content is given below.

==============================================================
       LIST OF DEGREES OF FREEDOM AND THEIR SUPPORT           
             (represented by element indices)                 
==============================================================
Index = 0, Support = 0  3  
Index = 1, Support = 0  2  
Index = 2, Support = 0  1  
Index = 3, Support = 1  2  
Index = 4, Support = 1  3  
Index = 5, Support = 2  3  
==============================================================
  LIST OF ELEMENTS WITH THEIR ASSOCIATED DEGREES OF FREEDOM   
             (represented by dof indices)                     
==============================================================
Ele-Idx = 0, Dofs = 0 1 2 
Ele-Idx = 1, Dofs = 3 4 2 
Ele-Idx = 2, Dofs = 3 1 5 
Ele-Idx = 3, Dofs = 4 5 0 

The output above concludes this section. With the dof handler at hand the general setup of the discretisation scheme is finished. Now we are ready to directly dive into the more boundary element specific topics. This will be done within the following examples.