From 0687dbcfcdad96fe8b9e6665994674a803b5073b Mon Sep 17 00:00:00 2001 From: David Minton Date: Tue, 10 Sep 2019 11:02:45 -0400 Subject: [PATCH] Added code to test a variable fe hypothesis --- src/crater/crater_populate.f90 | 6 ------ src/crater/crater_soften.f90 | 14 ++++++++++-- src/globals/module_globals.f90 | 7 ++++++ src/io/io_input.f90 | 39 ++++++++++++++++++++++++++++++++++ 4 files changed, 58 insertions(+), 8 deletions(-) diff --git a/src/crater/crater_populate.f90 b/src/crater/crater_populate.f90 index 21c96455..3ec949bd 100644 --- a/src/crater/crater_populate.f90 +++ b/src/crater/crater_populate.f90 @@ -190,12 +190,6 @@ subroutine crater_populate(user,surf,crater,domain,prod,production_list,vdist,nt ! Do seismic shaking if (user%doseismic) call seismic_shake(user,surf,crater,domain) - ! Generate interior anomalous diffusion - !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) - ! find the average height and slope at crater location call crater_averages(user,surf,crater) diff --git a/src/crater/crater_soften.f90 b/src/crater/crater_soften.f90 index e20d4f1a..f9d81235 100644 --- a/src/crater/crater_soften.f90 +++ b/src/crater/crater_soften.f90 @@ -40,8 +40,18 @@ subroutine crater_soften(user,surf,crater,domain) integer(I4B) :: inc,incsq,N,xpi,ypi,iradsq,i,j,ss,ioerr,bigi,bigj real(DP) :: xp,yp,fradsq,xbar,ybar,areafrac,kdiffmax,sf + real(DP) :: mfe,bfe + + ! TEST VARIABLE FE MODEL + 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 = min(user%femax,crater%fe) + !^^^^^^^^^^^^^^^^^^^^^^^^^^ + kdiffmax = user%Kd1 * crater%frad**user%psi - inc = int(min(user%fe * crater%frad / user%pix,SQRT2 * user%gridsize)) + 3 + inc = int(min(crater%fe * crater%frad / user%pix,SQRT2 * user%gridsize)) + 3 maxhits = (1 + (inc / (user%gridsize / 2)))**2 crater%maxinc = max(crater%maxinc,inc) incsq = (inc - 2)**2 @@ -78,7 +88,7 @@ subroutine crater_soften(user,surf,crater,domain) xbar = xp - crater%xl ybar = yp - crater%yl - areafrac = util_area_intersection(user%fe * crater%frad,xbar,ybar,user%pix) + areafrac = util_area_intersection(crater%fe * crater%frad,xbar,ybar,user%pix) kdiff(i,j) = kdiffmax * areafrac ! Apply user-defined diffusive degradation to degradation region !end if diff --git a/src/globals/module_globals.f90 b/src/globals/module_globals.f90 index f3f271ae..95a8b191 100644 --- a/src/globals/module_globals.f90 +++ b/src/globals/module_globals.f90 @@ -187,6 +187,13 @@ module module_globals real(DP) :: Kd1 ! Degradation function coefficient (from Minton et al. (2019)) real(DP) :: psi ! Degradation function exponent (from Minton et al. (2019)) 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 diff --git a/src/io/io_input.f90 b/src/io/io_input.f90 index 258a83ae..ab880bfb 100644 --- a/src/io/io_input.f90 +++ b/src/io/io_input.f90 @@ -89,7 +89,12 @@ subroutine io_input(infile,user) user%Kd1 = 0.00_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 open(unit=LUN,file=infile,status="old",iostat=ierr) if (ierr /= 0) then @@ -452,6 +457,40 @@ subroutine io_input(infile,user) 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") + ifirst = ilast + 1 + call io_get_token(line, ilength, ifirst, ilast, ierr) + token = line(ifirst:ilast) + read(token, *) user%dovariablefe + case ("FEMIN") + 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") + ifirst = ilast + 1 + call io_get_token(line, ilength, ifirst, ilast, ierr) + token = line(ifirst:ilast) + read(token, *) user%rminfe + + + + + + + !************************************************************************** ! The following is for backwards compatibility with older style input files !**************************************************************************