From 3a3b28266d213814db0c8667467368bf060697c8 Mon Sep 17 00:00:00 2001 From: Austin Blevins Date: Fri, 3 Feb 2023 15:59:54 -0500 Subject: [PATCH] Overhauled melt tracking and made it much simpler --- src/crater/crater_populate.f90 | 2 - src/regolith/module_regolith.f90 | 23 +++++---- src/regolith/regolith_streamtube.f90 | 24 +++++----- src/regolith/regolith_streamtube_head.f90 | 25 ++-------- src/regolith/regolith_streamtube_lineseg.f90 | 47 ++++--------------- src/regolith/regolith_subpixel_streamtube.f90 | 21 ++++----- src/regolith/regolith_traverse_streamtube.f90 | 14 +++--- 7 files changed, 54 insertions(+), 102 deletions(-) diff --git a/src/crater/crater_populate.f90 b/src/crater/crater_populate.f90 index ea40c8a0..6ae9faf1 100644 --- a/src/crater/crater_populate.f90 +++ b/src/crater/crater_populate.f90 @@ -387,7 +387,5 @@ subroutine crater_populate(user,surf,crater,domain,prod,production_list,vdist,nt end if - write(*,*) 'Total times mixing was called:',nmixingtimes - return end subroutine crater_populate diff --git a/src/regolith/module_regolith.f90 b/src/regolith/module_regolith.f90 index 36150d9e..3f6fc9c1 100644 --- a/src/regolith/module_regolith.f90 +++ b/src/regolith/module_regolith.f90 @@ -35,7 +35,7 @@ module module_regolith save interface - subroutine regolith_melt_zone(user,crater,dimp,vimp,rmelt,depthb, volm) + subroutine regolith_melt_zone(user,crater,dimp,vimp,rmelt,depthb,volm) use module_globals type(usertype),intent(in) :: user type(cratertype),intent(inout) :: crater @@ -77,7 +77,7 @@ end subroutine regolith_transport interface subroutine regolith_streamtube(user,surf,crater,domain,ejb,ejtble,xp,yp,xpi,ypi,lrad,ebh,& - rm,vsq,age,age_resolution, volm) + rm,vsq,age,age_resolution,volm) use module_globals implicit none type(usertype),intent(in) :: user @@ -95,13 +95,14 @@ end subroutine regolith_streamtube interface subroutine regolith_traverse_streamtube(user,surfi,deltar,ri,rip1,eradi,erado,newlayer,vmare,totseb,& - age_collector,xmints,xsfints,rsh,depthb,mixedregodata) + age_collector,xmints,xsfints,rsh,depthb,meltinejecta) use module_globals implicit none type(usertype),intent(in) :: user type(surftype),intent(inout) :: surfi real(DP),intent(in) :: deltar,ri,rip1,eradi,erado - type(regodatatype),intent(inout) :: newlayer, mixedregodata + type(regodatatype),intent(inout) :: newlayer + real(DP),intent(out) :: meltinejecta real(DP),intent(out) :: vmare,totseb real(SP),dimension(:),intent(inout) :: age_collector real(DP),intent(in) :: xmints @@ -111,13 +112,14 @@ end subroutine regolith_traverse_streamtube interface subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer,vmare,totseb,& - age_collector,xmints,xsfints,vol,mixedregodata) + age_collector,xmints,xsfints,vol,meltinejecta) use module_globals implicit none type(usertype),intent(in) :: user type(surftype),intent(in) :: surfi real(DP),intent(in) :: deltar,ri,rip1,eradi - type(regodatatype),intent(inout) :: newlayer, mixedregodata + type(regodatatype),intent(inout) :: newlayer + real(DP),intent(out) :: meltinejecta real(DP),intent(out) :: vmare,totseb real(SP),dimension(:),intent(inout) :: age_collector real(DP),intent(in) :: xmints @@ -128,13 +130,14 @@ end subroutine regolith_subpixel_streamtube interface subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad,eradi,deltar,newlayer,vmare,& - totseb,age_collector,xmints,xsfints,depthb,mixedregodata) + totseb,age_collector,xmints,xsfints,depthb,meltinejecta) use module_globals implicit none type(usertype),intent(in) :: user type(surftype),intent(in) :: surfi real(DP),intent(in) :: thetast,ri,rip1,zmin,zmax,erad,eradi,deltar - type(regodatatype),intent(inout) :: newlayer,mixedregodata + type(regodatatype),intent(inout) :: newlayer + real(DP),intent(inout) :: meltinejecta real(DP),intent(inout) :: vmare,totseb real(SP),dimension(:),intent(inout) :: age_collector real(DP),intent(in) :: xmints @@ -143,12 +146,12 @@ end subroutine regolith_streamtube_lineseg end interface interface - subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector,mixedregodata) + subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector,meltinejecta) use module_globals implicit none type(usertype),intent(in) :: user type(surftype),intent(in) :: surfi - type(regodatatype),intent(inout) :: mixedregodata + real(DP),intent(inout) :: meltinejecta real(DP),intent(in) :: deltar real(DP),intent(inout) :: totmare,tots real(SP),dimension(:),intent(inout) :: age_collector diff --git a/src/regolith/regolith_streamtube.f90 b/src/regolith/regolith_streamtube.f90 index 69ffe66c..8bfe83af 100644 --- a/src/regolith/regolith_streamtube.f90 +++ b/src/regolith/regolith_streamtube.f90 @@ -89,7 +89,8 @@ subroutine regolith_streamtube(user,surf,crater,domain,ejb,ejtble,xp,yp,xpi,ypi, integer(I4B) :: i,j,k,toti,totj,toty,cnt,xstpi,ystpi real(DP) :: vtot,vseg,ri,rip1,xc,yc,thetast real(DP) :: vst,vbody,rbody,vmare,totmare,totseb,tots - type(regodatatype) :: newlayer, mixedregodata + real(DP) :: meltinejecta + type(regodatatype) :: newlayer ! Constrain the tangital tube's volume with CTEM result real(DP) :: k1,k2,k3,k4,c1,c2 @@ -117,10 +118,6 @@ subroutine regolith_streamtube(user,surf,crater,domain,ejb,ejtble,xp,yp,xpi,ypi, ! Executalbe code - mixedregodata%totvolume = 0 - mixedregodata%meltvolume = 0 - mixedregodata%meltfrac = 0 - ! ****** Interpolate radial distance, erad, for a given pixel ******* ! outeredge = crater%frad + domain%ejbres * (EJBTABSIZE - 0.5_DP) ! inneredge = crater%frad + 0.5_DP * domain%ejbres @@ -259,7 +256,7 @@ subroutine regolith_streamtube(user,surf,crater,domain,ejb,ejtble,xp,yp,xpi,ypi, newlayer%thickness = vseg/(user%pix**2) call util_periodic(xstpi,ystpi,user%gridsize) call regolith_subpixel_streamtube(user,surf(xstpi,ystpi),deltar,ri,rip1,eradi,newlayer,vmare,totseb,& - age_collector,xmints,xsfints,vol,mixedregodata) + age_collector,xmints,xsfints,vol,meltinejecta) newlayer%age(:) = newlayer%age(:) + age_collector(:) @@ -268,8 +265,7 @@ subroutine regolith_streamtube(user,surf,crater,domain,ejb,ejtble,xp,yp,xpi,ypi, newlayer%thickness = ebh newlayer%comp = min(totmare/tots, 1.0_DP) newlayer%age(:) = newlayer%age(:) * min( (ebh * user%pix**2) / tots, 1.0_DP) - newlayer%meltvolume = newlayer%meltvolume + mixedregodata%meltvolume - newlayer%totvolume = newlayer%totvolume + mixedregodata%totvolume + newlayer%meltvolume = newlayer%meltvolume + meltinejecta newlayer%meltfrac = newlayer%meltvolume / newlayer%totvolume else @@ -283,7 +279,7 @@ subroutine regolith_streamtube(user,surf,crater,domain,ejb,ejtble,xp,yp,xpi,ypi, newlayer%thickness = vseg/(user%pix**2) call util_periodic(xstpi,ystpi,user%gridsize) call regolith_traverse_streamtube(user,surf(xstpi,ystpi),deltar,rbody,eradi,eradi,erado,newlayer,vmare,& - totseb,age_collector,xmints,xsfints,rsh,depthb,mixedregodata) + totseb,age_collector,xmints,xsfints,rsh,depthb,meltinejecta) totmare = totmare + vmare tots = tots + totseb end if @@ -300,7 +296,7 @@ subroutine regolith_streamtube(user,surf,crater,domain,ejb,ejtble,xp,yp,xpi,ypi, newlayer%thickness = vseg/(user%pix**2) call util_periodic(xstpi,ystpi,user%gridsize) call regolith_traverse_streamtube(user,surf(xstpi,ystpi),deltar,ri,rip1,eradi,erado,newlayer,vmare,& - totseb,age_collector,xmints,xsfints,rsh,depthb,mixedregodata) + totseb,age_collector,xmints,xsfints,rsh,depthb,meltinejecta) totmare = totmare + vmare tots = tots + totseb end if @@ -309,15 +305,17 @@ subroutine regolith_streamtube(user,surf,crater,domain,ejb,ejtble,xp,yp,xpi,ypi, xstpi = crater%xlpx + nint(eradc*xl/lrad/user%pix) ystpi = crater%ylpx + nint(eradc*yl/lrad/user%pix) call util_periodic(xstpi,ystpi,user%gridsize) - call regolith_streamtube_head(user,surf(xstpi,ystpi),deltar,totmare,tots,age_collector,mixedregodata) + call regolith_streamtube_head(user,surf(xstpi,ystpi),deltar,totmare,tots,age_collector,meltinejecta) newlayer%thickness = ebh newlayer%comp = min(totmare/tots, 1.0_DP) newlayer%age(:) = newlayer%age(:) + age_collector(:) newlayer%age(:) = newlayer%age(:) * min( (ebh * user%pix**2) / tots, 1.0_DP) - newlayer%meltvolume = newlayer%meltvolume + mixedregodata%meltvolume - newlayer%totvolume = newlayer%totvolume + mixedregodata%totvolume + newlayer%meltvolume = newlayer%meltvolume + meltinejecta newlayer%meltfrac = newlayer%meltvolume / newlayer%totvolume + if (newlayer%meltfrac > 1.0_DP) then + write(*,*) "Melt fraction >1!", xpi,ypi,crater%timestamp,crater%fcrat,crater%xlpx,crater%ylpx + end if end if diff --git a/src/regolith/regolith_streamtube_head.f90 b/src/regolith/regolith_streamtube_head.f90 index a13cb3bd..bf8d886b 100644 --- a/src/regolith/regolith_streamtube_head.f90 +++ b/src/regolith/regolith_streamtube_head.f90 @@ -18,14 +18,14 @@ ! Notes : The stream tube's head is always attached to the surface. ! !********************************************************************************************************************************** -subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector,mixedregodata) +subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector,meltinejecta) use module_globals use module_regolith, EXCEPT_THIS_ONE => regolith_streamtube_head implicit none ! arguemnts type(usertype),intent(in) :: user type(surftype),intent(in) :: surfi - type(regodatatype),intent(inout) :: mixedregodata + real(DP),intent(inout) :: meltinejecta real(DP),intent(in) :: deltar real(DP),intent(inout) :: totmare,tots real(SP),dimension(:),intent(inout) :: age_collector @@ -66,12 +66,7 @@ subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector recyratio = vsgly / (user%pix**2) /current(N)%thickness age_collector(:) = age_collector(:) + current(N)%age(:) * recyratio headmeltvol = current(N)%meltfrac * recyratio * vsgly - mixedregodata%totvolume = mixedregodata%totvolume + vsgly - mixedregodata%meltvolume = mixedregodata%meltvolume + headmeltvol - mixedregodata%meltfrac = mixedregodata%meltvolume / mixedregodata%totvolume - if (mixedregodata%meltfrac > 1.0_DP) then - write(*,*) "ERROR! mixedregodata%meltfrac >1! (HEAD)" - end if + meltinejecta = meltinejecta + headmeltvol else ! head is not intersected with layers. do @@ -83,13 +78,8 @@ subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector totmarehead = totmarehead + vhead * vratio * current(N)%comp recyratio = vhead * vratio / (user%pix**2) / current(N)%thickness age_collector(:) = age_collector(:) + current(N)%age(:) * recyratio - mixedregodata%totvolume = mixedregodata%totvolume + tothead headmeltvol = current(N)%meltfrac * recyratio * tothead - mixedregodata%meltvolume = mixedregodata%meltvolume + headmeltvol - mixedregodata%meltfrac = mixedregodata%meltvolume / mixedregodata%totvolume - if (mixedregodata%meltfrac > 1.0_DP) then - write(*,*) "ERROR! mixedregodata%meltfrac >1! (HEAD)" - end if + meltinejecta = meltinejecta + headmeltvol !current => current%next N = N - 1 z = z + current(N)%thickness @@ -100,13 +90,8 @@ subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector tothead = vsgly recyratio = (vsgly - tothead) / (user%pix**2) / current(N)%thickness age_collector(:) = age_collector(:) + current(N)%age(:) * recyratio - mixedregodata%totvolume = mixedregodata%totvolume + tothead headmeltvol = current(N)%meltfrac * recyratio * tothead - mixedregodata%meltvolume = mixedregodata%meltvolume + headmeltvol - mixedregodata%meltfrac = mixedregodata%meltvolume / mixedregodata%totvolume - if (mixedregodata%meltfrac > 1.0_DP) then - write(*,*) "ERROR! mixedregodata%meltfrac >1! (HEAD)" - end if + meltinejecta = meltinejecta + headmeltvol exit end if end do diff --git a/src/regolith/regolith_streamtube_lineseg.f90 b/src/regolith/regolith_streamtube_lineseg.f90 index fe99112c..ddd2f104 100644 --- a/src/regolith/regolith_streamtube_lineseg.f90 +++ b/src/regolith/regolith_streamtube_lineseg.f90 @@ -19,7 +19,7 @@ ! !********************************************************************************************************************************** subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad,eradi,deltar,newlayer,vmare,totseb,& - age_collector,xmints,xsfints,depthb,mixedregodata) + age_collector,xmints,xsfints,depthb,meltinejecta) use module_globals use module_regolith, EXCEPT_THIS_ONE => regolith_streamtube_lineseg implicit none @@ -27,7 +27,8 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad type(usertype),intent(in) :: user type(surftype),intent(in) :: surfi real(DP),intent(in) :: thetast,ri,rip1,zmin,zmax,erad,eradi,deltar - type(regodatatype),intent(inout) :: newlayer,mixedregodata + type(regodatatype),intent(inout) :: newlayer + real(DP),intent(inout) :: meltinejecta real(DP),intent(inout) :: vmare,totseb real(SP),dimension(:),intent(inout) :: age_collector real(DP),intent(in) :: xmints @@ -82,26 +83,16 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad recyratio = max(ebh_recyl,0.0_DP) / current(N-1)%thickness age_collector(:) = age_collector(:) + current(N-1)%age(:) * recyratio vol = vol + sum(current(N-1)%age(:)) * recyratio - mixedregodata%totvolume = mixedregodata%totvolume + (vsgly-vsh) linmelt = current(N-1)%meltfrac * vsgly * recyratio - mixedregodata%meltvolume = mixedregodata%meltvolume + linmelt - mixedregodata%meltfrac = mixedregodata%meltvolume / mixedregodata%totvolume - if (mixedregodata%meltfrac > 1.0_DP) then - write(*,*) "ERROR! mixedregodata%meltfrac >1! (LINESEG)" - end if + meltinejecta = meltinejecta + linmelt else if (ri <= xmints .and. rip1 > xmints) then vseg = regolith_streamtube_volume_func(eradi,xmints,rip1,deltar) vsh = regolith_shock_damage(eradi,deltar,xmints,xsfints,ri,rip1) recyratio = max((vseg-vsh),0.0_DP) / (user%pix**2) / current(N-1)%thickness age_collector(:) = age_collector(:) + current(N-1)%age(:) * recyratio vol = vol + sum(current(N-1)%age(:)) * recyratio - mixedregodata%totvolume = mixedregodata%totvolume + (vseg-vsh) linmelt = current(N-1)%meltfrac * vseg * recyratio - mixedregodata%meltvolume = mixedregodata%meltvolume + linmelt - mixedregodata%meltfrac = mixedregodata%meltvolume / mixedregodata%totvolume - if (mixedregodata%meltfrac > 1.0_DP) then - write(*,*) "ERROR! mixedregodata%meltfrac >1! (LINESEG)" - end if + meltinejecta = meltinejecta + linmelt end if exit else @@ -128,13 +119,8 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad recyratio = max((vseg-vsh),0.0_DP) / (user%pix**2) / current(N)%thickness age_collector(:) = age_collector(:) + current(N)%age(:) * recyratio vol = vol + sum(current(N)%age(:)) * recyratio - mixedregodata%totvolume = mixedregodata%totvolume + (vseg-vsh) linmelt = current(N)%meltfrac * vseg * recyratio - mixedregodata%meltvolume = mixedregodata%meltvolume + linmelt - mixedregodata%meltfrac = mixedregodata%meltvolume / mixedregodata%totvolume - if (mixedregodata%meltfrac > 1.0_DP) then - write(*,*) "ERROR! mixedregodata%meltfrac >1! (LINESEG)" - end if + meltinejecta = meltinejecta + linmelt end if ! A segment coming from the side of emerging location of a streamtube rip1 @@ -144,13 +130,8 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad recyratio = max((vseg-vsh),0.0_DP) / (user%pix**2) / current(N)%thickness age_collector(:) = age_collector(:) + current(N)%age(:) * recyratio vol = vol + sum(current(N)%age(:)) * recyratio - mixedregodata%totvolume = mixedregodata%totvolume + (vseg-vsh) linmelt = current(N)%meltfrac * vseg * recyratio - mixedregodata%meltvolume = mixedregodata%meltvolume + linmelt - mixedregodata%meltfrac = mixedregodata%meltvolume / mixedregodata%totvolume - if (mixedregodata%meltfrac > 1.0_DP) then - write(*,*) "ERROR! mixedregodata%meltfrac >1! (LINESEG)" - end if + meltinejecta = meltinejecta + linmelt end if !current => current%next N = N - 1 @@ -172,13 +153,8 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad recyratio = max((vseg-vsh),0.0_DP) / (user%pix**2) / current(N)%thickness age_collector(:) = age_collector(:) + current(N)%age(:) * recyratio vol = vol + sum(current(N)%age(:)) * recyratio - mixedregodata%totvolume = mixedregodata%totvolume + (vseg-vsh) linmelt = current(N)%meltfrac * vseg * recyratio - mixedregodata%meltvolume = mixedregodata%meltvolume + linmelt - mixedregodata%meltfrac = mixedregodata%meltvolume / mixedregodata%totvolume - if (mixedregodata%meltfrac > 1.0_DP) then - write(*,*) "ERROR! mixedregodata%meltfrac >1! (LINESEG)" - end if + meltinejecta = meltinejecta + linmelt end if exit @@ -188,13 +164,8 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad vsgly = regolith_streamtube_volume_func(eradi,ri,rip1,deltar) vmare = vmare + (vsgly - totseb) * current(N)%comp totseb = vsgly - mixedregodata%totvolume = mixedregodata%totvolume + vsgly linmelt = current(N)%meltfrac * vsgly * recyratio - mixedregodata%meltvolume = mixedregodata%meltvolume + linmelt - mixedregodata%meltfrac = mixedregodata%meltvolume / mixedregodata%totvolume - if (mixedregodata%meltfrac > 1.0_DP) then - write(*,*) "ERROR! mixedregodata%meltfrac >1! (LINESEG)" - end if + meltinejecta = meltinejecta + linmelt exit end if diff --git a/src/regolith/regolith_subpixel_streamtube.f90 b/src/regolith/regolith_subpixel_streamtube.f90 index aba3df9f..1fe3d24d 100644 --- a/src/regolith/regolith_subpixel_streamtube.f90 +++ b/src/regolith/regolith_subpixel_streamtube.f90 @@ -51,7 +51,7 @@ ! !********************************************************************************************************************************** subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer,vmare,totseb,& - age_collector,xmints,xsfints,vol,mixedregodata) + age_collector,xmints,xsfints,vol,meltinejecta) use module_globals use module_regolith, EXCEPT_THIS_ONE => regolith_subpixel_streamtube implicit none @@ -60,7 +60,8 @@ subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer type(usertype),intent(in) :: user type(surftype),intent(in) :: surfi real(DP),intent(in) :: deltar,ri,rip1,eradi - type(regodatatype),intent(inout) :: newlayer, mixedregodata + type(regodatatype),intent(inout) :: newlayer + real(DP),intent(out) :: meltinejecta real(DP),intent(out) :: vmare,totseb real(SP),dimension(:),intent(inout) :: age_collector real(DP),intent(in) :: xmints @@ -97,6 +98,7 @@ subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer vseg = 0.0_DP vsgly1 = 0.0_DP vsgly2 = 0.0_DP + meltinejecta = 0.0_DP ! Two cases: subpixel is inside the first layer, and its volume is simply the landing ejecta blanket. if (zend>=zmax) then @@ -105,10 +107,8 @@ subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer if (eradi>xmints) then vseg = regolith_streamtube_volume_func(eradi,xmints,eradi,deltar) vsh = regolith_shock_damage(eradi,deltar,xmints,xsfints,0.0_DP,eradi) - mixedregodata%totvolume = vseg-vsh recyratio = max(vseg-vsh,0.0 )/ (user%pix**2) / (surfi%regolayer(N)%thickness) - mixedregodata%meltvolume = surfi%regolayer(N)%meltfrac * vseg * recyratio - mixedregodata%meltfrac = mixedregodata%meltvolume / mixedregodata%totvolume + meltinejecta = surfi%regolayer(N)%meltfrac * vseg * recyratio age_collector(:) = age_collector(:) + surfi%regolayer(N)%age(:) * recyratio vol = vol + sum(age_collector(:)) ! write(*,*) '1',eradi, xmints, xsfints, & @@ -196,19 +196,14 @@ subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer end if !current => current%next - mixedregodata%meltvolume = mixedregodata%meltvolume + mvl + mvr - mixedregodata%totvolume = mixedregodata%totvolume + vsgly - mixedregodata%meltfrac = mixedregodata%meltvolume / mixedregodata%totvolume - if (mixedregodata%meltfrac > 1.0_DP) then - write(*,*) "ERROR! mixedregodata%meltfrac >1! (SUBPIXEL)" - end if + meltinejecta = meltinejecta + mvl + mvr N = N - 1 z = z + current(N)%thickness zstart = zend zend = z rlefti = rleftf rrightf = rrighti - ! final part of a stream tube + ! final part of a stream tube else !I think this means it's in the melt zone like Ya-Huei said above.. if so, nothing needed here. vsgly = 0.25 * PI * deltar**2 * a**2 * eradi / b * (abs(tan(b) - b)) vmare = vmare + (vsgly - totseb) * current(N)%comp @@ -218,7 +213,7 @@ subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer end do end if - call regolith_streamtube_head(user,surfi,deltar,vmare,totseb,age_collector, mixedregodata) + call regolith_streamtube_head(user,surfi,deltar,vmare,totseb,age_collector,meltinejecta) return end subroutine regolith_subpixel_streamtube diff --git a/src/regolith/regolith_traverse_streamtube.f90 b/src/regolith/regolith_traverse_streamtube.f90 index 6c8636f5..ea083f9e 100644 --- a/src/regolith/regolith_traverse_streamtube.f90 +++ b/src/regolith/regolith_traverse_streamtube.f90 @@ -19,7 +19,7 @@ ! !********************************************************************************************************************************** subroutine regolith_traverse_streamtube(user,surfi,deltar,ri,rip1,eradi,erado,newlayer,vmare,totseb,& - age_collector,xmints,xsfints,depthb,mixedregodata) + age_collector,xmints,xsfints,depthb,meltinejecta) use module_globals use module_regolith, EXCEPT_THIS_ONE => regolith_traverse_streamtube implicit none @@ -28,7 +28,8 @@ subroutine regolith_traverse_streamtube(user,surfi,deltar,ri,rip1,eradi,erado,ne type(usertype),intent(in) :: user type(surftype),intent(inout) :: surfi real(DP),intent(in) :: deltar,ri,rip1,eradi,erado - type(regodatatype),intent(inout) :: newlayer, mixedregodata + type(regodatatype),intent(inout) :: newlayer + real(DP),intent(out) :: meltinejecta real(DP),intent(out) :: vmare,totseb real(SP),dimension(:),intent(inout) :: age_collector real(DP),intent(in) :: xmints @@ -50,6 +51,9 @@ subroutine regolith_traverse_streamtube(user,surfi,deltar,ri,rip1,eradi,erado,ne ! Shock zone real(DP) :: vsh + !executable code + meltinejecta = 0.0_DP + erad = (eradi + erado)/2.0 rzmax = erad * sqrt(3.0)/4.0 if (ri .eq. 0) then @@ -85,15 +89,13 @@ subroutine regolith_traverse_streamtube(user,surfi,deltar,ri,rip1,eradi,erado,ne vsh = regolith_shock_damage(eradi,deltar,xmints,xsfints,ri,rip1) recyratio = max(vseg-vsh,0.0_DP) / (user%pix**2) / surfi%regolayer(N)%thickness age_collector(:) = age_collector(:) + surfi%regolayer(N)%age(:) * recyratio - mixedregodata%totvolume = vseg-vsh - mixedregodata%meltvolume = surfi%regolayer(N)%meltfrac * vseg * recyratio - mixedregodata%meltfrac = mixedregodata%meltvolume / mixedregodata%totvolume + meltinejecta = surfi%regolayer(N)%meltfrac * vseg * recyratio end if else call regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad,eradi,deltar,& - newlayer,vmare,totseb,age_collector,xmints,xsfints,depthb,mixedregodata) + newlayer,vmare,totseb,age_collector,xmints,xsfints,depthb,meltinejecta) end if