diff --git a/src/model/westergaard.cpp b/src/model/westergaard.cpp
index a96f248..7e340d9 100644
--- a/src/model/westergaard.cpp
+++ b/src/model/westergaard.cpp
@@ -1,221 +1,225 @@
 /**
  * @file
  *
  * @author Lucas Frérot <lucas.frerot@epfl.ch>
  *
  * @section LICENSE
  *
  * Copyright (©)  2017 EPFL  (Ecole Polytechnique  Fédérale de
  * Lausanne)  Laboratory (LSMS  -  Laboratoire de  Simulation  en Mécanique  des
  * Solides)
  *
  * Tamaas is free  software: you can redistribute it and/or  modify it under the
  * terms  of the  GNU Lesser  General Public  License as  published by  the Free
  * Software Foundation, either version 3 of the License, or (at your option) any
  * later version.
  *
  * Tamaas is  distributed in the  hope that it  will be useful, but  WITHOUT ANY
  * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
  * A  PARTICULAR PURPOSE. See  the GNU  Lesser General  Public License  for more
  * details.
  *
  * You should  have received  a copy  of the GNU  Lesser General  Public License
  * along with Tamaas. If not, see <http://www.gnu.org/licenses/>.
  *
  */
 /* -------------------------------------------------------------------------- */
 #include "westergaard.hh"
 #include "fft_plan_manager.hh"
 #include "influence.hh"
 #include "loop.hh"
 #include "static_types.hh"
 #include "types.hh"
 /* -------------------------------------------------------------------------- */
 __BEGIN_TAMAAS__
 /* -------------------------------------------------------------------------- */
 
 template <model_type type>
 Westergaard<type>::Westergaard(Model* model) : BEEngine(model) {
   // Copy sizes
   std::array<UInt, trait::dimension> sizes;
   std::copy(model->getDiscretization().begin(),
-            model->getDiscretization().end(), sizes.begin());
+	    model->getDiscretization().end(), sizes.begin());
   auto hermitian_sizes =
       GridHermitian<Real, trait::dimension>::hermitianDimensions(sizes);
 
   // Copy boundary sizes
   auto sizes_it = sizes.begin();
   std::array<UInt, trait::boundary_dimension> boundary_hermitian_sizes;
   // Ignoring first dimension if model is volumetric
   if (trait::dimension > trait::boundary_dimension)
     ++sizes_it;
   std::copy(sizes_it, sizes.end(), boundary_hermitian_sizes.begin());
   boundary_hermitian_sizes =
       GridHermitian<Real, trait::boundary_dimension>::hermitianDimensions(
-          boundary_hermitian_sizes);
+	  boundary_hermitian_sizes);
 
   constexpr UInt nb_components = trait::components;
 
   buffer.setNbComponents(nb_components);
   buffer.resize(boundary_hermitian_sizes);
   influence.setNbComponents(nb_components * nb_components);
   influence.resize(hermitian_sizes);
 }
 
 /* -------------------------------------------------------------------------- */
 
 template <model_type type>
 void Westergaard<type>::initInfluence() {
   TAMAAS_EXCEPTION("Not implemented");
 }
 
 /* -------------------------------------------------------------------------- */
 
 #define INFLUENCE_BASIC_MACRO                                                  \
   do {                                                                         \
     constexpr UInt bdim = trait::boundary_dimension;                           \
     auto wavevectors =                                                         \
-        FFTransform<Real, bdim>::computeFrequencies<true>(influence.sizes());  \
+	FFTransform<Real, bdim>::computeFrequencies<true>(influence.sizes());  \
     Real E_star = model->getHertzModulus();                                    \
     StaticVector<const Real, bdim> domain(model->getSystemSize()[0]);          \
     Loop::stridedLoop(influence::Basic<bdim>(E_star, domain), wavevectors,     \
-                      influence);                                              \
+		      influence);                                              \
     influence(0) = 0;                                                          \
   } while (0)
 
 template <>
 void Westergaard<model_type::basic_2d>::initInfluence() {
   INFLUENCE_BASIC_MACRO;
 }
 
 template <>
 void Westergaard<model_type::basic_1d>::initInfluence() {
   INFLUENCE_BASIC_MACRO;
 }
 
 #undef INFLUENCE_BASIC_MACRO
 
 #define INFLUENCE_SURFACE_MACRO                                                \
   do {                                                                         \
     constexpr UInt bdim = trait::boundary_dimension;                           \
     auto wavevectors =                                                         \
-        FFTransform<Real, bdim>::computeFrequencies<true>(influence.sizes());  \
+	FFTransform<Real, bdim>::computeFrequencies<true>(influence.sizes());  \
     Real E = model->getYoungModulus();                                         \
     Real nu = model->getPoissonRatio();                                        \
     StaticVector<const Real, bdim> domain(model->getSystemSize()[0]);          \
     Loop::stridedLoop(influence::Surface<bdim>(E, nu, domain), wavevectors,    \
-                      influence);                                              \
+		      influence);                                              \
     StaticMatrix<Complex, bdim + 1, bdim + 1> inf_0(influence(0));             \
     inf_0 = 0;                                                                 \
   } while (0)
 
 template <>
 void Westergaard<model_type::surface_1d>::initInfluence() {
   INFLUENCE_SURFACE_MACRO;
 }
 
 template <>
 void Westergaard<model_type::surface_2d>::initInfluence() {
   INFLUENCE_SURFACE_MACRO;
 }
 
 #undef INFLUENCE_SURFACE_MACRO
 
   /* ------------------------------------------------------------------------ */
 
 #define NEUMANN_BASIC(type)                                                    \
   template <>                                                                  \
   void Westergaard<type>::solveNeumann(GridBase<Real>& in,                     \
-                                       GridBase<Real>& out) const {            \
+				       GridBase<Real>& out) const {            \
     auto apply = [](decltype(buffer)& buffer,                                  \
-                    const decltype(influence)& influence) {                    \
+		    const decltype(influence)& influence) {                    \
       buffer *= influence;                                                     \
       buffer(0) = 0;                                                           \
     };                                                                         \
-                                                                               \
+									       \
     fourierApply(apply, in, out);                                              \
   }
 
 #define DIRICHLET_BASIC(type)                                                  \
   template <>                                                                  \
   void Westergaard<type>::solveDirichlet(GridBase<Real>& in,                   \
-                                         GridBase<Real>& out) const {          \
+					 GridBase<Real>& out) const {          \
     auto apply = [](decltype(buffer)& buffer,                                  \
-                    const decltype(influence)& influence) {                    \
+		    const decltype(influence)& influence) {                    \
       buffer /= influence;                                                     \
       buffer(0) = 0;                                                           \
     };                                                                         \
-                                                                               \
+									       \
     fourierApply(apply, in, out);                                              \
   }
 
 NEUMANN_BASIC(model_type::basic_1d);
 NEUMANN_BASIC(model_type::basic_2d);
 DIRICHLET_BASIC(model_type::basic_1d);
 DIRICHLET_BASIC(model_type::basic_2d);
 
 #undef NEUMANN_BASIC
 #undef DIRICHLET_BASIC
 
 /* -------------------------------------------------------------------------- */
 
+namespace detail {
+/// Class needed for cuda
+template <UInt d>
+struct SurfaceApply {
+void operator()(GridHermitian<Real, d>& buffer_,
+		const GridHermitian<Real, d>& influence) {
+    Loop::stridedLoop(
+	[] CUDA_LAMBDA(StaticVector<Complex, d>&& buf,
+			StaticMatrix<const Complex, d, d>&& inf) {
+	Complex copy_[d];
+	StaticVector<Complex, d> copy(copy_[0]);
+	copy = buf;
+	buf.mul(inf, copy);
+	},
+	buffer_, influence);
+
+    StaticVector<Complex, d> buff_0(buffer_(0));
+    buff_0 = 0;
+}
+};
+} // namespace detail
+
 #define NEUMANN_SURFACE(type)                                                  \
   template <>                                                                  \
   void Westergaard<type>::solveNeumann(GridBase<Real>& in,                     \
-                                       GridBase<Real>& out) const {            \
-    auto apply = [](decltype(buffer)& buffer_,                                 \
-                    const decltype(influence)& influence) {                    \
-      constexpr UInt d = decltype(buffer)::dimension;                          \
-                                                                               \
-      Loop::stridedLoop(                                                       \
-          [] CUDA_LAMBDA(StaticVector<Complex, d>&& buf,                       \
-                         StaticMatrix<const Complex, d, d>&& inf) {            \
-            Complex copy_[d];                                                  \
-            StaticVector<Complex, d> copy(copy_[0]);                           \
-            copy = buf;                                                        \
-            buf.mul(inf, copy);                                                \
-          },                                                                   \
-          buffer_, influence);                                                 \
-                                                                               \
-      StaticVector<Complex, d> buff_0(buffer_(0));                             \
-      buff_0 = 0;                                                              \
-    };                                                                         \
-                                                                               \
-    fourierApply(apply, in, out);                                              \
+				       GridBase<Real>& out) const {            \
+    fourierApply(detail::SurfaceApply<trait::boundary_dimension>(), in, out);  \
   }
 
 /// \TODO Find a way to make it work with cuda
-// NEUMANN_SURFACE(model_type::surface_1d);
-// NEUMANN_SURFACE(model_type::surface_2d);
+NEUMANN_SURFACE(model_type::surface_1d);
+NEUMANN_SURFACE(model_type::surface_2d);
 
 #undef NEUMANN_SURFACE
 
 /* -------------------------------------------------------------------------- */
 
 template <model_type type>
 void Westergaard<type>::solveNeumann(GridBase<Real>& /*neumann*/,
-                                     GridBase<Real>& /*dirichlet*/) const {
+				     GridBase<Real>& /*dirichlet*/) const {
   TAMAAS_EXCEPTION("Not yet implemented");
 }
 
 /* -------------------------------------------------------------------------- */
 
 template <model_type type>
 void Westergaard<type>::solveDirichlet(GridBase<Real>& /*neumann*/,
-                                       GridBase<Real>& /*dirichlet*/) const {
+				       GridBase<Real>& /*dirichlet*/) const {
   TAMAAS_EXCEPTION("Not yet implemented");
 }
 
 /* -------------------------------------------------------------------------- */
 /* Template instanciation                                                     */
 /* -------------------------------------------------------------------------- */
 template class Westergaard<model_type::basic_1d>;
 template class Westergaard<model_type::basic_2d>;
 template class Westergaard<model_type::surface_1d>;
 template class Westergaard<model_type::surface_2d>;
 template class Westergaard<model_type::volume_1d>;
 template class Westergaard<model_type::volume_2d>;
 
 /* -------------------------------------------------------------------------- */
 __END_TAMAAS__
 /* -------------------------------------------------------------------------- */