Stripped away variable fe model and now allow for variable psi

src/crater/crater_degradation_function.f90

@@ -32,9 +32,9 @@ function crater_degradation_function(user,crater) result(Kd)
 
    ! Testing values that match both mare scale and highlands scale
    alpha_1 = user%psi
-   alpha_2 = 1.1_DP
+   alpha_2 = user%psi2
 
-   r_break = 0.10e3_DP
+   r_break = user%rbreak
    delta = 1.0_DP
 
    A = user%Kd1 * r_break**(alpha_1) / (1_DP + (alpha_1 - alpha_2) / alpha_1)**2 

src/crater/crater_generate.f90

@@ -194,20 +194,7 @@ subroutine crater_generate(user,crater,domain,prod,production_list,vdist,surf)
 
    ! Get pixel space values
    crater%fcratpx = nint(crater%fcrat / user%pix)
-
-   ! TEST VARIABLE FE MODEL
-   if (user%dovariablefe) then
-      !mfe = (user%femax - user%femin) / (log10(user%rmaxfe) - log10(user%rminfe))
-      !bfe = 0.5_DP * ((user%femax + user%femin) - mfe * (log10(user%rmaxfe) + log10(user%rminfe)))
-
-      !crater%fe = mfe * log10(crater%frad) + bfe
-      !crater%fe = max(min(user%femax,crater%fe),0.d0)
-      crater%fe = max(user%femin * (crater%frad*1e-3_DP)**(0.27_DP),user%femin)
-      crater%fe = min(crater%fe * crater%frad, 3.08e6_DP) / crater%frad
-   else
-      crater%fe = user%fe
-   end if
-   !^^^^^^^^^^^^^^^^^^^^^^^^^^
+   crater%fe = user%fe
    return
 end subroutine crater_generate
 

src/crater/crater_subpixel_diffusion.f90

@@ -160,19 +160,12 @@ subroutine crater_subpixel_diffusion(user,surf,nflux,domain,finterval,kdiffin)
             crater%ylpx = nint(crater%yl / user%pix)
 
             crater%frad = 0.5_DP * (diam + dD * rn(3))
+
             crater%continuous = RCONT * crater%frad**(EXPCONT) 
-            krad = fd * crater%frad
-
-            if (user%dovariablefe) then
-               !mfe = (user%femax - user%femin) / (log10(user%rmaxfe) - log10(user%rminfe))
-               !bfe = 0.5_DP * ((user%femax + user%femin) - mfe * (log10(user%rmaxfe) + log10(user%rminfe)))
-               !crater%fe = mfe * log10(crater%frad) + bfe
-               !crater%fe = max(min(user%femax,crater%fe),0.0_DP) 
-               crater%fe = max(user%femin * (crater%frad*1e-3_DP)**(0.27_DP),user%femin)
-               crater%fe = min(crater%fe * crater%frad, 3.08e6_DP) / crater%frad
-            else
-               crater%fe = user%fe
-            end if
+            crater%fe = user%fe
+            krad = max(fd * crater%frad,crater%fe * crater%frad)
+
+            crater%fe = user%fe
 
 
             !dKdN = user%Kd1 * crater%frad**(user%psi)

src/globals/module_globals.f90

@@ -185,16 +185,11 @@ module module_globals
 
    ! Crater diffusion input parameters
    real(DP) :: Kd1 ! Degradation function coefficient (from Minton et al. (2019))
-   real(DP) :: psi ! Degradation function exponent (from Minton et al. (2019))
+   real(DP) :: psi ! Degradation function exponent (from Minton et al. (2029))
+   real(DP) :: psi2 ! Degradation function large size exponent (from Minton et al. (2020))
+   real(DP) :: rbreak ! Degradation function break in exponent (from Minton et al. (2020))
    real(DP) :: fe  ! Scale factor for size of degradation region (from Minton et al. (2019))
 
-   ! Testing variable fe model
-   logical  :: dovariablefe ! Testing a variable fe model
-   real(DP) :: femin  ! fe at largest crater in linear interpolation
-   real(DP) :: femax  ! fe at smallest crate in linear interpolation
-   real(DP) :: rminfe   ! crater size with lowest fe 
-   real(DP) :: rmaxfe   ! crater size with highest fe 
-
    ! Ejecta softening variables
    logical           :: dosoftening  ! Set T to use the extra crater softening model
    real(DP)          :: ejecta_truncation ! Set the number of crater diameters to truncate the ejecta

src/io/io_input.f90

@@ -86,15 +86,12 @@ subroutine io_input(infile,user)
    user%testtally = .false.
    user%killatmaxcrater = .false.
    user%doporosity = .false.
-   user%Kd1 = 0.00_DP
+   user%Kd1 = 0.0315151258095091_DP
    user%psi = 2.0_DP
-   user%fe  = 1.0_DP
-   user%dovariablefe = .false.
-   user%ejecta_truncation = 10.0_DP
-   user%femax = 10.0_DP
-   user%femin = 2.0_DP
-   user%rmaxfe = 100_DP
-   user%rminfe = 400e3_DP
+   user%psi2 = 1.25_DP
+   user%rbreak = 0.1e3_DP
+   user%fe  = 5.0_DP
+   user%ejecta_truncation = 5.0_DP
 
    open(unit=LUN,file=infile,status="old",iostat=ierr)
    if (ierr /= 0) then
@@ -452,44 +449,22 @@ subroutine io_input(infile,user)
             call io_get_token(line, ilength, ifirst, ilast, ierr)
             token = line(ifirst:ilast)
             read(token, *) user%psi
-         case ("FE")
-            ifirst = ilast + 1
-            call io_get_token(line, ilength, ifirst, ilast, ierr)
-            token = line(ifirst:ilast)
-            read(token, *) user%fe
 
-!Used to test variable fe model
-         case ("DOVARIABLEFE")
+        case ("PSI2")
             ifirst = ilast + 1
             call io_get_token(line, ilength, ifirst, ilast, ierr)
             token = line(ifirst:ilast)
-            read(token, *) user%dovariablefe
-         case ("FEMIN")
+            read(token, *) user%psi2
+        case ("RBREAK")
             ifirst = ilast + 1
             call io_get_token(line, ilength, ifirst, ilast, ierr)
             token = line(ifirst:ilast)
-            read(token, *) user%femin
-         case ("FEMAX")
-            ifirst = ilast + 1
-            call io_get_token(line, ilength, ifirst, ilast, ierr)
-            token = line(ifirst:ilast)
-            read(token, *) user%femax
-         case ("RMAXFE")
-            ifirst = ilast + 1
-            call io_get_token(line, ilength, ifirst, ilast, ierr)
-            token = line(ifirst:ilast)
-            read(token, *) user%rmaxfe
-         case ("RMINFE")
+            read(token, *) user%rbreak
+         case ("FE")
             ifirst = ilast + 1
             call io_get_token(line, ilength, ifirst, ilast, ierr)
             token = line(ifirst:ilast)
-            read(token, *) user%rminfe
-
-
-
-
-
-
+            read(token, *) user%fe
 
          !**************************************************************************
          ! The following is for backwards compatibility with older style input files