changes to crater subroutines since the SVN-GitHub switchover

MintonGroup · May 10, 2017 · 47e43b2 · 47e43b2
1 parent f237874
commit 47e43b2
Show file tree

Hide file tree

Showing 4 changed files with 15 additions and 9 deletions.
diff --git a/src/crater/crater_populate.f90 b/src/crater/crater_populate.f90
@@ -184,10 +184,10 @@ subroutine crater_populate(user,surf,crater,domain,prod,production_list,vdist,nt
          if (user%doseismic) call seismic_shake(user,surf,crater,domain)
 
          ! Generate interior anomalous diffusion
-         if (user%dosoftening) call crater_soften(user,surf,crater,domain)
+         !if (user%dosoftening) call crater_soften(user,surf,crater,domain)
 
          ! Generate distal anomalous diffusion
-         if (user%dosoftening) call crater_soften_accumulate(user,surf,crater,domain,kdiff)
+         !if (user%dosoftening) call crater_soften_accumulate(user,surf,crater,domain,kdiff)
 
          ! find the average height and slope at crater location
          call crater_averages(user,surf,crater)

diff --git a/src/crater/crater_soften_accumulate.f90 b/src/crater/crater_soften_accumulate.f90
@@ -22,6 +22,7 @@
 subroutine crater_soften_accumulate(user,surf,crater,domain,kdiff)
    use module_globals
    use module_util
+   use module_ejecta
    use module_crater, EXCEPT_THIS_ONE => crater_soften_accumulate
    implicit none
 
@@ -96,7 +97,8 @@ subroutine crater_soften_accumulate(user,surf,crater,domain,kdiff)
             areafrac = areafrac * craterhole(xpi,ypi)
 
             if (.not.hit(xpi,ypi)) then 
-               kdiff(xpi,ypi) = kdiff(xpi,ypi) + kappatmax * areafrac 
+               lrad = sqrt(lradsq)
+               kdiff(xpi,ypi) = kdiff(xpi,ypi) * lrad**2 + kappatmax * areafrac 
                hit(xpi,ypi) = .true.
             end if
 

diff --git a/src/crater/crater_tally_observed.f90 b/src/crater/crater_tally_observed.f90
@@ -58,9 +58,9 @@ subroutine crater_tally_observed(user,surf,domain,nkilled,onum,obsdist,obslist,o
    integer(I4B) :: inc,xpi,ypi
    logical :: killable
    integer(I4B) :: nrim,nbowl,nouter
-   real(DP) :: rim,bowl,outer,rad,baseline
+   real(DP) :: rim,bowl,outer,rad,baseline,dd
    ! Counting parameters from Howl study
-   real(DP),parameter :: DDCUTOFF = 5.e-2_DP
+   real(DP),parameter :: DDCUTOFF = 5.0e-2_DP
    real(DP),parameter :: OCUTOFF = 5.5e-2_DP
    real(DP),parameter :: RIMDI = 1.0_DP
    real(DP),parameter :: RIMDO = 1.2_DP
@@ -200,8 +200,12 @@ subroutine crater_tally_observed(user,surf,domain,nkilled,onum,obsdist,obslist,o
       bowl = bowl / nbowl 
       outer = outer / nouter
       tmp_depthdiam(craternum) = (rim - bowl) / crater%fcrat
-
-      if (((tmp_depthdiam(craternum) > DDCUTOFF).and.((outer - rim) / crater%fcrat < OCUTOFF)).and.&
+      if (crater%fcrat < 20e3_DP) then
+         dd = DDCUTOFF
+      else
+         dd = DDCUTOFF - (crater%fcrat - 20e3_DP) * 2e-7
+      end if
+      if (((tmp_depthdiam(craternum) > dd).and.((outer - rim) / crater%fcrat < OCUTOFF)).and.&
           (nrim /= 0).and.(nbowl /= 0).and.(nouter /= 0)) then
          countable(craternum) = .true.
          killable = .false.

diff --git a/src/crater/crater_tally_true.f90 b/src/crater/crater_tally_true.f90
@@ -51,14 +51,14 @@ subroutine crater_tally_true(domain,truelist,ntrue,truedist)
 
    ! Bin the true crater distribution
    do craternum = 1,ntrue
-      ! Find out what bin this crate belongs in
+      ! Find out what bin this crater belongs in
       i = ceiling(log(truelist(1,craternum)/1e3_DP)/LOGSQRT2) - domain%plo
       if (i < 1) then
          write(*,*) 'fcrat = ',truelist(1,craternum)
          write(*,*) 'smallest bin: ',1e3_DP*SQRT2**(domain%plo)
          write(*,*) 'domain%smallest_crater = ',domain%smallest_crater
          write(*,*) 'domain%subcrater_limit = ',domain%subcrater_limit
-         read(*,*)
+         cycle
       end if
       truedist(3,i) = truedist(3,i) + log(truelist(1,craternum)) ! Geometric mean (intermediate step)
       truedist(4,i) = truedist(4,i) + 1 ! Differential number