## 1. Solve a partial differential equation

With all the previously notions we approach, the definition of a partial differential equation and boundary conditions are our next step. More details on these aspects can be retrieve at this page.

#### 1.1. Variational formulation

This example refers to a laplacian problem, define by

Strong formulation
$-\Delta u=1 \text{ in } \Omega=[0,1]^2, \quad u=1 \text{ on } \partial \Omega$

After turning the Strong formulation into its weak form, we have

$\int_\Omega \nabla u \cdot \nabla v = \int_\Omega v,\quad u=1 \text{ on } \partial \Omega$

where $u$ is the unknown and $v$ a test function. The left side is known as the bilinear form $a$ and the right side is the linear form $l$.

$a(u,v) = l(v), \quad u=1 \text{ on } \partial \Omega$

#### 1.2. Implementation

The steps to implement this problem are

• Loading a 2D mesh, creating the function space $V_h$, composed of piecewise polynomial functions of order 2, and its associated elements

`````` auto mesh = loadMesh(_mesh=new Mesh<Simplex<2>>);
auto Vh = Pch<2>( mesh );
auto u = Vh->element();
auto v = Vh->element();``````
• Define the linear form $l$ with test function space $V_h$

``````auto l = form1( _test=Vh );
l = integrate(_range=elements(mesh),
_expr=id(v));``````
• Define the bilinear form $a$ with $V_h$ as test and trial function spaces

``````auto a = form2( _trial=Vh, _test=Vh);
a = integrate(_range=elements(mesh),
_expr=gradt(u)*trans(grad(v)) );``````

`form1` and `form2` are used to define respectively the left and right side of our partial differential equation.

• Add Dirichlet boundary condition on u

``````a+=on(_range=boundaryfaces(mesh),
_rhs=l, _element=u, _expr=cst(0.) );``````

We impose, in this case, $u=0$ on $\partial\Omega$, with the keyword `on`.

• Solving the problem

``a.solve(_rhs=l,_solution=u);``
• Exporting the solution

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

The complete code reads as follows :

``````#include <feel/feel.hpp>

int main(int argc, char**argv )
{
using namespace Feel;

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

auto mesh = loadMesh(_mesh=new Mesh<Simplex<2>>);
auto Vh = Pch<2>( mesh );
auto u = Vh->element();
auto v = Vh->element();

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=cst(0.) );
a.solve(_rhs=l,_solution=u);

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

and the corresponding config file

``````[gmsh]
hsize=1e-1``````