diff --git a/src/geometry_plastic_nonlocal.f90 b/src/geometry_plastic_nonlocal.f90
index b042943d0..056c7cf57 100644
--- a/src/geometry_plastic_nonlocal.f90
+++ b/src/geometry_plastic_nonlocal.f90
@@ -120,26 +120,35 @@ end subroutine geometry_plastic_nonlocal_disable
 !---------------------------------------------------------------------------------------------------
 subroutine geometry_plastic_nonlocal_results
   
-  integer,     dimension(:),   allocatable :: s
+  integer,     dimension(:),   allocatable :: shp
 
 #if defined(DAMASK_HDF5)
   call results_openJobFile
 
   writeVolume: block
     real(pReal), dimension(:), allocatable :: temp
-    s = shape(geometry_plastic_nonlocal_IPvolume0)
-    temp = reshape(geometry_plastic_nonlocal_IPvolume0,[s(1)*s(2)])
+    shp = shape(geometry_plastic_nonlocal_IPvolume0)
+    temp = reshape(geometry_plastic_nonlocal_IPvolume0,[shp(1)*shp(2)])
     call results_writeDataset('geometry',temp,'v_0',&
                               'initial cell volume','m³')
   end block writeVolume
 
-  writeArea: block
+  writeAreas: block
     real(pReal), dimension(:,:), allocatable :: temp
-    s = shape(geometry_plastic_nonlocal_IParea0)
-    temp = reshape(geometry_plastic_nonlocal_IParea0,[s(1),s(2)*s(3)])
+    shp = shape(geometry_plastic_nonlocal_IParea0)
+    temp = reshape(geometry_plastic_nonlocal_IParea0,[shp(1),shp(2)*shp(3)])
     call results_writeDataset('geometry',temp,'a_0',&
                               'initial cell face area','m²')
-  end block writeArea
+  end block writeAreas
+
+  writeNormals: block
+    real(pReal), dimension(:,:,:), allocatable :: temp
+    shp = shape(geometry_plastic_nonlocal_IPareaNormal0)
+    temp = reshape(geometry_plastic_nonlocal_IPareaNormal0,[shp(1),shp(2),shp(3)*shp(4)])
+    call results_writeDataset('geometry',temp,'n_0',&
+                              'initial cell face normals','-',transposed=.false.)
+  end block writeNormals
+
 
   call results_closeJobFile
 #endif
diff --git a/src/mesh_marc.f90 b/src/mesh_marc.f90
index 4154e75d2..66906207a 100644
--- a/src/mesh_marc.f90
+++ b/src/mesh_marc.f90
@@ -77,8 +77,8 @@ subroutine mesh_init(ip,el)
     connectivity_cell                                                                                !< cell connectivity for each element,ip/cell
   integer, dimension(:,:), allocatable :: &
     connectivity_elem
-
-  real(pReal), dimension(:,:,:,:),allocatable :: x
+  real(pReal), dimension(:,:,:,:),allocatable :: &
+    unscaledNormals
 
   write(6,'(/,a)')   ' <<<+-  mesh init  -+>>>'
  
@@ -120,8 +120,9 @@ subroutine mesh_init(ip,el)
                      node0_cell,ip_reshaped)
 
   call geometry_plastic_nonlocal_setIPvolume(IPvolume(elem,node0_cell,connectivity_cell))
-  x = IPareaNormal(elem,nElems,connectivity_cell,node0_cell)
-  call geometry_plastic_nonlocal_setIParea(norm2(x,1))
+  unscaledNormals = IPareaNormal(elem,nElems,connectivity_cell,node0_cell)
+  call geometry_plastic_nonlocal_setIParea(norm2(unscaledNormals,1))
+  call geometry_plastic_nonlocal_setIPareaNormal(unscaledNormals/spread(norm2(unscaledNormals,1),1,3))
   call geometry_plastic_nonlocal_results
   
 
@@ -1034,7 +1035,7 @@ function IPareaNormal(elem,nElem,connectivity,node)
 
   do e = 1,nElem
     do i = 1,elem%nIPs
-      do f = 1,size(elem%cellFace,2)
+      do f = 1,elem%nIPneighbors
         nodePos = node(1:3,connectivity(elem%cellface(1:m,f),i,e))
         
         select case (elem%cellType)