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}\label{eq:problemU} \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 $$ \forall x,y\in [0,1]^2,\qquad f(x,y) = (1+2\pi^2)\cos(\pi x)\cos(\pi y) $$ Then, the unique solution of problem \eqref{eq:problemU} is $$ \forall x,y\in[0,1]^2, \qquad u(x,y) = \cos(\pi x)\cos(\pi y). $$