From eee8f13a7adab47d60794fb999d99178cfb5c65b Mon Sep 17 00:00:00 2001 From: David Minton Date: Fri, 10 Apr 2020 20:28:12 -0400 Subject: [PATCH] Tweaked some of the ejecta noise parameters and increased the efficiency a bit --- src/crater/crater_realistic_topography.f90 | 17 ++++++----------- 1 file changed, 6 insertions(+), 11 deletions(-) diff --git a/src/crater/crater_realistic_topography.f90 b/src/crater/crater_realistic_topography.f90 index bc5c6239..f55c73f3 100644 --- a/src/crater/crater_realistic_topography.f90 +++ b/src/crater/crater_realistic_topography.f90 @@ -311,7 +311,7 @@ subroutine complex_wall_texture(user,surf,crater,domain,deltaMtot) real(DP), parameter :: rad = 0.20_DP ! Radial size of central peak relative to frad integer(I4B), parameter :: num_octaves = 10 ! Number of Perlin noise octaves integer(I4B), parameter :: offset = 4000 ! Scales the random xy-offset so that each crater's random noise is unique - real(DP), parameter :: xy_noise_fac = 4.0_DP ! Spatial "size" of noise features at the first octave + real(DP), parameter :: xy_noise_fac = 3.0_DP ! Spatial "size" of noise features at the first octave real(DP), parameter :: noise_height = 7.0e-3_DP ! Spatial "size" of noise features at the first octave real(DP), parameter :: freq = 2.0_DP ! Spatial size scale factor multiplier at each octave level real(DP), parameter :: pers = 1.20_DP ! The relative size scaling at each octave level @@ -341,7 +341,7 @@ subroutine complex_wall_texture(user,surf,crater,domain,deltaMtot) areafrac = 1.0 - util_area_intersection(0.3_DP * crater%floordiam,xbar,ybar,user%pix) areafrac = areafrac * util_area_intersection(2.1_DP * crater%frad,xbar,ybar,user%pix) areafrac = areafrac * min((r / flr)**12,1.0_DP) ! Smooth out interface between wall and floor - areafrac = areafrac * max(min(2._DP - 1 * r,1.0_DP),0.0_DP) ! Smooth out region outside of the rim + areafrac = areafrac * max(min(2._DP - r,1.0_DP),0.0_DP) ! Smooth out region outside of the rim ! Add some roughness to the walls noise = 0.0_DP @@ -438,7 +438,7 @@ subroutine complex_terrace(user,surf,crater,deltaMtot) noise = (dnoise**2)**noise_slope hprof = (r / router)**(-1) - tprof = crater_profile(user,crater,rinner) * (r/rinner)**2 + tprof = crater_profile(user,crater,rinner) * (r/rinner)**2 ! This is the floor profile that replaces the old one at each terrace if (1.0_DP - noise < hprof) then newdem = min(tprof,newdem) @@ -537,7 +537,6 @@ subroutine ejecta_texture(user,surf,crater,deltaMtot,inc,ejecta_dem) real(DP) :: newdem,elchange,xbar,ybar,xstretch,ystretch,ejbmass,ejbmassnew,fmasscons,r,phi,areafrac integer(I4B) :: i,j,xpi,ypi integer(I4B),dimension(2,-inc:inc,-inc:inc) :: indarray - real(DP),dimension(-inc:inc,-inc:inc) :: noise_dem real(DP), dimension(2) :: rn ! Topographic noise parameters @@ -593,7 +592,6 @@ subroutine ejecta_texture(user,surf,crater,deltaMtot,inc,ejecta_dem) ! Add the base texture to the ejecta proportional to the thickness - noise_dem = 0.0_DP do j = -inc,inc do i = -inc,inc @@ -614,7 +612,7 @@ subroutine ejecta_texture(user,surf,crater,deltaMtot,inc,ejecta_dem) areafrac = util_area_intersection(inc * user%pix,xbar,ybar,user%pix) areafrac = areafrac * (1.0_DP - util_area_intersection(crater%frad,xbar,ybar,user%pix)) - areafrac = areafrac * (1.0_DP - max(min(2._DP - 1 * r,1.0_DP),0.0_DP)) ! Blend in with the wall texture + areafrac = areafrac * (1.0_DP - max(min(2._DP - r,1.0_DP),0.0_DP)) ! Blend in with the wall texture ! Make the splat pattern @@ -628,7 +626,7 @@ subroutine ejecta_texture(user,surf,crater,deltaMtot,inc,ejecta_dem) dsplat = (dsplat**2)**(splat_slope) splatnoise = splatnoise + dsplat end do - hprof = (1.0_DP * r)**(1.5_DP) + hprof = r**(1.0_DP) insplat = 1.0_DP + splatnoise > hprof ! make base texture and then add extra layers if we are in one @@ -642,14 +640,11 @@ subroutine ejecta_texture(user,surf,crater,deltaMtot,inc,ejecta_dem) if (insplat) noise = noise + (1.0_DP + splatnoise - hprof) * ejecta_dem(i,j) * splatmag end do - noise_dem(i,j) = noise_dem(i,j) + noise * areafrac + ejecta_dem(i,j) = max(ejecta_dem(i,j) + noise * areafrac,0.0_DP) end do end do - - ejecta_dem(-inc:inc,-inc:inc) = max(ejecta_dem(-inc:inc,-inc:inc) + noise_dem(-inc:inc,-inc:inc),0.0_DP) - ! Save new total mass for mass conservation calculation ejbmassnew = sum(ejecta_dem)