diff --git a/src/modules/mesh/moduleMesh@boundaryEM.f90 b/src/modules/mesh/moduleMesh@boundaryEM.f90
index e7567c8..17b0036 100644
--- a/src/modules/mesh/moduleMesh@boundaryEM.f90
+++ b/src/modules/mesh/moduleMesh@boundaryEM.f90
@@ -450,10 +450,10 @@ submodule(moduleMesh) boundaryEM
       integer:: e
 
       write(fileID, '(A)') self%name
+      write(fileID, '(A,",",A)') 'Edge id', 'Electric Field (V m^-1)'
       select type(self)
       type is(boundaryEMFreeCurrent)
         do e = 1, self%nEdges
-          write(fileID, '(A,",",A)') 'Edge id', 'Electric Field (V m^-1)'
           write(fileID, '('//fmtColInt//','//fmtColReal//')') self%edges(e)%obj%n, self%electricField(e)*EF_ref
 
         end do
diff --git a/src/modules/mesh/moduleMesh@boundaryParticle.f90 b/src/modules/mesh/moduleMesh@boundaryParticle.f90
index 514e613..1f2d3f0 100644
--- a/src/modules/mesh/moduleMesh@boundaryParticle.f90
+++ b/src/modules/mesh/moduleMesh@boundaryParticle.f90
@@ -576,13 +576,13 @@ submodule(moduleMesh) boundaryParticle
               do n = 1, edge%nNodes
                 node => mesh%nodes(nodes(n))%obj
 
-                den_nodes(n) = node%output(s)%den/(node%v*Vol_ref*n_ref)
-                mom_nodes(n) = dot_product(node%output(s)%mom, edge%normal)/(node%v*Vol_ref*n_ref)
+                den_nodes(n) =  node%output(s)%den/(node%v*Vol_ref*n_ref)
+                mom_nodes(n) = -dot_product(node%output(s)%mom, edge%normal)/(node%v*Vol_ref*n_ref)
 
               end do
 
               den_center = edge%gatherF(edge%centerXi(), edge%nNodes, den_nodes)
-              if (den_center > 1.0d-10) then
+              if (den_center > 1.0d-12) then
                 mom_center = edge%gatherF(edge%centerXi(), edge%nNodes, mom_nodes)
                 vel_center = mom_center / den_center
 
@@ -604,7 +604,7 @@ submodule(moduleMesh) boundaryParticle
           alpha = vel_rest - vel_incident
 
           ! Apply correction with a factor of 0.1 to avoid fast changes
-          self%velocity_shift(edge%n) = self%velocity_shift(edge%n) + 1.0d-2 * alpha
+          self%velocity_shift(edge%n) = self%velocity_shift(edge%n) + 1.0d-3 * alpha
 
           deallocate(den_nodes)
           deallocate(mom_nodes)