fixed a few bugs

src/crater/crater_populate.f90

@@ -268,7 +268,8 @@ subroutine crater_populate(user,surf,crater,domain,prod,production_list,vdist,nt
                call ejecta_table_define(user,crater,domain,ejb,ejtble)
                !call ejecta_interpolate(crater,domain,crater%frad,ejb(1:ejtble),ejtble,crater%ejrim)
             end if
-            call ejecta_emplace(user,surf,crater,domain,ejb(1:ejtble),ejtble,ejbmass,age,age_resolution,ejecta_dem,vmeltsheet)
+            call ejecta_emplace(user,surf,crater,domain,ejb(1:ejtble),ejtble,ejbmass,age,age_resolution,&
+               ejecta_dem,nmeltsheet,vmeltsheet)
          else
             ejtble = 0
          end if

src/ejecta/ejecta_emplace.f90

@@ -75,7 +75,8 @@
 !                The cutoff of ejecta thickness is still buggy.  
 !
 !**********************************************************************************************************************************
-subroutine ejecta_emplace(user,surf,crater,domain,ejb,ejtble,deltaMtot,age,age_resolution,cumulative_elchange,vmeltsheet)
+subroutine ejecta_emplace(user,surf,crater,domain,ejb,ejtble,deltaMtot,age,age_resolution,cumulative_elchange,&
+   nmeltsheet,vmeltsheet)
    use module_globals
    use module_util
    use module_io
@@ -95,6 +96,7 @@ subroutine ejecta_emplace(user,surf,crater,domain,ejb,ejtble,deltaMtot,age,age_r
    real(DP),intent(in)  :: age
    real(DP),intent(in)  :: age_resolution
    real(DP),dimension(:,:),allocatable,intent(out) :: cumulative_elchange
+   integer(I4B),intent(in) :: nmeltsheet
    real(DP),intent(out) :: vmeltsheet
 
    ! Internal variables
@@ -342,7 +344,11 @@ subroutine ejecta_emplace(user,surf,crater,domain,ejb,ejtble,deltaMtot,age,age_r
       end do
    end if
 
-   vmeltsheet = totmelt - vmelt
+   if (totmelt > vmelt) then
+      vmeltsheet = totmelt - vmelt
+   else !give the craters a melt sheet of 1mm
+      vmeltsheet = 1.0_DP * user%pix * user%pix * nmeltsheet
+   end if
 
    ! Create box for soften calculation (will be no bigger than the grid itself)
    if (2 * inc + 1 < user%gridsize) then

src/ejecta/module_ejecta.f90

@@ -25,7 +25,8 @@ module module_ejecta
 save
 
    interface
-      subroutine ejecta_emplace(user,surf,crater,domain,ejb,ejtble,deltaMtot,age,age_resolution,cumulative_elchange,vmeltsheet)
+      subroutine ejecta_emplace(user,surf,crater,domain,ejb,ejtble,deltaMtot,age,age_resolution,cumulative_elchange,&
+         nmeltsheet,vmeltsheet)
       use module_globals
       implicit none
       type(usertype),intent(in) :: user
@@ -38,6 +39,7 @@ subroutine ejecta_emplace(user,surf,crater,domain,ejb,ejtble,deltaMtot,age,age_r
       real(DP),intent(in)  :: age
       real(DP),intent(in)  :: age_resolution
       real(DP),dimension(:,:),allocatable,intent(inout) :: cumulative_elchange
+      integer(I4B),intent(in) :: nmeltsheet
       real(DP),intent(out) :: vmeltsheet
       end subroutine ejecta_emplace
    end interface

src/regolith/regolith_subpixel_streamtube.f90

@@ -113,9 +113,9 @@ subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer
          vseg             = regolith_streamtube_volume_func(eradi,xmints,eradi,deltar)
          vsh              = regolith_shock_damage(eradi,deltar,xmints,xsfints,0.0_DP,eradi)
          recyratio        = max(vseg-vsh,0.0 )/ (user%pix**2) / (surfi%regolayer(M)%thickness)
-         meltinejecta     = surfi%regolayer(M)%meltfrac * (vseg-vsh)! * recyratio
-         distvol(:)       = distvol(:) + (surfi%regolayer(M)%meltdist(:)*(vseg-vsh))!*recyratio)
-         totvol           = vseg - vsh
+         meltinejecta     = surfi%regolayer(M)%meltfrac * max((vseg-vsh),0.0_DP)! * recyratio
+         distvol(:)       = distvol(:) + (surfi%regolayer(M)%meltdist(:)*max((vseg-vsh),0.0_DP))!*recyratio)
+         totvol           = vseg
          age_collector(:) = age_collector(:) + surfi%regolayer(M)%age(:) * recyratio
          vol              = vol + sum(age_collector(:))
 !         write(*,*) '1',eradi, xmints, xsfints, &
@@ -207,7 +207,7 @@ subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer
          !current => current%next
          meltinejecta = meltinejecta + mvl + mvr
          !distvol(:) = distvol(:) + (current(N)%meltdist(:)*((vsgly1-vsh)*recyratio)+((vsgly2-vsh2)*recyratio2))
-         distvol(:) = distvol(:) + (current(N)%meltdist(:)*(max((vsgly1-vsh),0.0_DP)+max((vsgly2-vsh),0.0_DP)))
+         distvol(:) = distvol(:) + (current(N)%meltdist(:)*(max((vsgly1-vsh),0.0_DP)+max((vsgly2-vsh2),0.0_DP)))
          totvol = totvol + vsgly1 + vsgly2
          !N = N - 1
          z = z + current(N-1)%thickness