Logo
Finite Element Embedded Library and Language in C++
Feel++ Feel++ on Github Feel++ on Travis-CI Feel++ on Twitter Feel++ on YouTube Feel++ community
 All Classes Namespaces Files Functions Variables Typedefs Pages
Quickstart Laplacian

As an introduction to the aim and the way to do with Feel++, we provide a sort of Hello World program to evaluate the library. We solve for the laplacian with Dirichlet boundary conditions.

Problem setup

We want to solve the following problem:

\[ - \Delta u = 1,\quad u_{|\partial \Omega} = g \]

where \(\Omega \in \mathbb{R}^n, n\in{1,2,3}\).

The variational form reads: looks for \(u \in H^1_g\left( \Omega \right) = \{ v \in H^1(\Omega) | v = g \mbox{ on } \partial \Omega \}\) such that

\[a\left( u,v \right) = l\left( v \right)\quad \forall v \in H^1\left( \Omega \right).\]

with:

\[a\left( u,v \right) = \int_{\Omega} \nabla u \cdot \nabla v ,\quad l\left( v \right) = \int_{\Omega} v \]

The aim of Feel++ is to provide the simplest way to write the \(a\) and \(l\) forms.

From a discrete point of view, we introduce \(V_h\subset H^1\left( \Omega \right)\) such that:

\( \begin{aligned} V_h = \left\{ v \in C^0\left( \Omega \right), \forall K\in \mathcal{T}_h, \right.v\left|_K \in P_1\left( K \right) \right\}, \end{aligned} \)


where \(\mathcal{T}_h\) is the set of element \(K\) forming the mesh of \(\Omega\).
We now look for \(u_h \in V_h\) such that:

\( \begin{aligned} \forall v_h\in V_h, a\left( u_h,v_h \right)=l\left( v_h \right). \end{aligned} \)


About the code

This section is here to declare that we want to use the namespace Feel++, to passe the command line options to the created environnement and add some informations (basics Feel++ options, application name).

using namespace Feel;
Environment env( _argc=argc, _argv=argv,
_desc=feel_options(),
_directory=".",
_about=about(_name="qs_laplacian",
_author="Feel++ Consortium",
_email="feelpp-devel@feelpp.org"));

We have to define the mesh, the approximation space and our test and trial functions.

auto mesh = unitSquare();
auto Vh = Pch<1>( mesh );
auto u = Vh->element();
auto v = Vh->element();

We create now our bilinear and linear forms, we add the homogeneous Dirichlet conditions and solve the discretized (linear) system.

auto l = form1( _test=Vh );
l = integrate(_range=elements(mesh),
_expr=id(v));
auto a = form2( _trial=Vh, _test=Vh );
a = integrate(_range=elements(mesh),
_expr=gradt(u)*trans(grad(v)) );
a+=on(_range=boundaryfaces(mesh), _rhs=l, _element=u,
_expr=constant(0.) );
a.solve(_rhs=l,_solution=u);

Feel++ provides the possibility to save the results:

auto e = exporter( _mesh=mesh );
e->add( "u", u );
e->save();
return 0;

top


First execution & visualization

To test that part of code, please go to:

  cd $FeelppBuildDir/quickstart

and execute the code, by:

  ./feelpp_qs_laplacian

This will produce several files in directory $HOME/feel/qs_laplacian/np_1/:

  feel.geo
  feel.msh
  qs_laplacian-1_0.case
  qs_laplacian-1_0.geo
  qs_laplacian-1.sos
  qs_laplacian-paraview-1.sos
  qs_laplacian.timeset
  qs_laplacian.u-1_0.001

You can visualize the results using any Ensight file reader, such as Paraview, opening qs_laplacian-paraview-1.sos.

  paraview qs_laplacian-paraview-1.sos

You may have a look to the options provided by the application

  ./feelpp_qs_laplacian --help

as well as the options provided by the Feel++ library

  ./feelpp_qs_laplacian --help-lib