diff --git a/src/schemes/fluid/surface_normal_sph.jl b/src/schemes/fluid/surface_normal_sph.jl
index b07c57b1f..6eebb650f 100644
--- a/src/schemes/fluid/surface_normal_sph.jl
+++ b/src/schemes/fluid/surface_normal_sph.jl
@@ -66,16 +66,18 @@ function calc_normal!(system::FluidSystem, neighbor_system::FluidSystem, u_syste
 end
 
 # Section 2.2 in Akinci et al. 2013 "Versatile Surface Tension and Adhesion for SPH Fluids"
-# and Section 5 in Morris 2000 "Simulating surface tension with smoothed particle hydrodynamics"
-function calc_normal!(system::FluidSystem, neighbor_system::BoundarySystem, u_system, v,
-                      v_neighbor_system, u_neighbor_system, semi, surfn, nsurfn)
+# Note: This is the simplest form of normal approximation commonly used in SPH and comes
+# with serious deficits in accuracy especially at corners, small neighborhoods and boundaries
+function calc_normal_akinci!(system::FluidSystem, neighbor_system::BoundarySystem, u_system,
+                             v,
+                             v_neighbor_system, u_neighbor_system, semi, surfn)
     (; cache) = system
     (; colorfield, colorfield_bnd) = neighbor_system.boundary_model.cache
     (; boundary_contact_threshold) = surfn
 
     system_coords = current_coordinates(u_system, system)
     neighbor_system_coords = current_coordinates(u_neighbor_system, neighbor_system)
-    nhs = get_neighborhood_search(system, neighbor_system, semi)
+    nhs = get_neighborhood_search(neighbor_system, system, semi)
 
     # First we need to calculate the smoothed colorfield values of the boundary
     # TODO: move colorfield to extra step