diff --git a/src/regolith/module_regolith.f90 b/src/regolith/module_regolith.f90
index 09fb2349..9c623b4e 100644
--- a/src/regolith/module_regolith.f90
+++ b/src/regolith/module_regolith.f90
@@ -227,10 +227,11 @@ end subroutine regolith_subcrater_mix
    end interface
 
    interface
-      subroutine regolith_mix(surf,d)
+      subroutine regolith_mix(surfi,mixing_depth,domain)
       use module_globals
-      type(surftype),intent(inout) :: surf
-      real(DP),intent(in) :: d
+      type(surftype),intent(inout) :: surfi
+      real(DP),intent(in) :: mixing_depth
+      type(domaintype),intent(in) :: domain
       end subroutine regolith_mix
    end interface
 
diff --git a/src/regolith/regolith_mix.f90 b/src/regolith/regolith_mix.f90
index bf7f4966..accf7885 100644
--- a/src/regolith/regolith_mix.f90
+++ b/src/regolith/regolith_mix.f90
@@ -16,7 +16,7 @@
 !  Notes       :  
 !
 !**********************************************************************************************************************************
-subroutine regolith_mix(surfi,mixing_depth)
+subroutine regolith_mix(surfi,mixing_depth,domain)
    use module_globals
    use module_util
    use module_regolith, EXCEPT_THIS_ONE => regolith_mix
@@ -25,6 +25,7 @@ subroutine regolith_mix(surfi,mixing_depth)
    ! Arguments
    type(surftype),intent(inout) :: surfi
    real(DP), intent(in) :: mixing_depth
+   type(domaintype),intent(in) :: domain
 
    ! Internal variables
    type(regodatatype) :: newlayer
@@ -41,6 +42,7 @@ subroutine regolith_mix(surfi,mixing_depth)
    newlayer%comp      = 0.0_DP
    newlayer%meltfrac  = 0.0_DP
    newlayer%age(:)    = 0.0_DP
+   allocate(newlayer%meltdist(domain%rcnum))
 
    !poppedlist => poppedlist_top
    !do while(associated(poppedlist%next))
@@ -50,12 +52,14 @@ subroutine regolith_mix(surfi,mixing_depth)
       newlayer%comp      = newlayer%comp + poppedarray(i)%thickness * poppedarray(i)%comp       
       newlayer%meltfrac  = newlayer%meltfrac + poppedarray(i)%thickness * poppedarray(i)%meltfrac
       newlayer%age(:)    = newlayer%age(:) + poppedarray(i)%age(:)
+      newlayer%meltdist(:) = newlayer%meltdist(:) + poppedarray(i)%thickness * poppedarray(i)%meltdist(:)
       !poppedlist => poppedlist%next
    end do
 
    ! Get average values of composition and melt fraction
    newlayer%comp = newlayer%comp / newlayer%thickness 
    newlayer%meltfrac = newlayer%meltfrac / newlayer%thickness 
+   newlayer%meltdist(:) = newlayer%meltdist(:) / newlayer%thickness
    
    call util_push_array(surfi%regolayer, newlayer)
    !call util_destroy_list(poppedlist_top)
diff --git a/src/regolith/regolith_subcrater_mix.f90 b/src/regolith/regolith_subcrater_mix.f90
index 5065b96e..ae40822f 100644
--- a/src/regolith/regolith_subcrater_mix.f90
+++ b/src/regolith/regolith_subcrater_mix.f90
@@ -66,7 +66,7 @@ subroutine regolith_subcrater_mix(user,surf,domain,nflux,finterval,p)
          N = size(surf(i,j)%regolayer)
 
          if (surf(i,j)%regolayer(N)%thickness < dd) then             
-            call regolith_mix(surf(i,j),dd)
+            call regolith_mix(surf(i,j),dd,domain)
          end if
          
        end do