Tutorial/Laplace equation with Neumann boundary condition

We propose here to solve a first very simple academic example with GMSH and GetDP. The problem is the following: \begin{equation} \begin{cases} \Delta u + u = f & \text{in } \Omega\\ \displaystyle{\frac{\partial u}{\partial \mathbf{n}} = 0} & \text{on }\partial\Omega, \end{cases} \end{equation} where $\Omega$ is the unit square, $\partial\Omega$ its boundary and $f$ the two-variables function defined by $$ f(x,y) = (1+2\pi^2)\cos(\pi x)\cos(\pi y) $$