diff --git a/src/regolith/module_regolith.f90 b/src/regolith/module_regolith.f90 index 3f6fc9c1..dd0288f4 100644 --- a/src/regolith/module_regolith.f90 +++ b/src/regolith/module_regolith.f90 @@ -95,14 +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,meltinejecta) + age_collector,xmints,xsfints,rsh,depthb,meltinejecta,totvol) 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 - real(DP),intent(out) :: meltinejecta + real(DP),intent(out) :: meltinejecta,totvol real(DP),intent(out) :: vmare,totseb real(SP),dimension(:),intent(inout) :: age_collector real(DP),intent(in) :: xmints @@ -112,14 +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,meltinejecta) + age_collector,xmints,xsfints,vol,meltinejecta,totvol) 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 - real(DP),intent(out) :: meltinejecta + real(DP),intent(out) :: meltinejecta,totvol real(DP),intent(out) :: vmare,totseb real(SP),dimension(:),intent(inout) :: age_collector real(DP),intent(in) :: xmints @@ -130,14 +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,meltinejecta) + totseb,age_collector,xmints,xsfints,depthb,meltinejecta,totvol) 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 - real(DP),intent(inout) :: meltinejecta + real(DP),intent(inout) :: meltinejecta,totvol real(DP),intent(inout) :: vmare,totseb real(SP),dimension(:),intent(inout) :: age_collector real(DP),intent(in) :: xmints @@ -146,12 +146,12 @@ end subroutine regolith_streamtube_lineseg end interface interface - subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector,meltinejecta) + subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector,meltinejecta,totvol) use module_globals implicit none type(usertype),intent(in) :: user type(surftype),intent(in) :: surfi - real(DP),intent(inout) :: meltinejecta + real(DP),intent(inout) :: meltinejecta,totvol 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 35228df0..d3110dbe 100644 --- a/src/regolith/regolith_streamtube.f90 +++ b/src/regolith/regolith_streamtube.f90 @@ -89,7 +89,7 @@ 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 - real(DP) :: meltinejecta + real(DP) :: meltinejecta, totvol type(regodatatype) :: newlayer ! Constrain the tangital tube's volume with CTEM result @@ -256,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,meltinejecta) + age_collector,xmints,xsfints,vol,meltinejecta,totvol) newlayer%age(:) = newlayer%age(:) + age_collector(:) @@ -269,7 +269,7 @@ subroutine regolith_streamtube(user,surf,crater,domain,ejb,ejtble,xp,yp,xpi,ypi, 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,& - newlayer%meltvolume, newlayer%totvolume, newlayer%ejm, newlayer%ejmf, newlayer%ejm/newlayer%ejmf + newlayer%meltvolume, newlayer%totvolume, newlayer%ejm, newlayer%ejmf, totvol end if else @@ -283,7 +283,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,meltinejecta) + totseb,age_collector,xmints,xsfints,rsh,depthb,meltinejecta,totvol) totmare = totmare + vmare tots = tots + totseb end if @@ -300,7 +300,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,meltinejecta) + totseb,age_collector,xmints,xsfints,rsh,depthb,meltinejecta,totvol) totmare = totmare + vmare tots = tots + totseb end if @@ -309,7 +309,7 @@ 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,meltinejecta) + call regolith_streamtube_head(user,surf(xstpi,ystpi),deltar,totmare,tots,age_collector,meltinejecta,totvol) newlayer%thickness = ebh newlayer%comp = min(totmare/tots, 1.0_DP) @@ -319,7 +319,7 @@ subroutine regolith_streamtube(user,surf,crater,domain,ejb,ejtble,xp,yp,xpi,ypi, 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,& - newlayer%meltvolume, newlayer%totvolume, newlayer%ejm, newlayer%ejmf, newlayer%ejm/newlayer%ejmf + newlayer%meltvolume, newlayer%totvolume, newlayer%ejm, newlayer%ejmf, totvol end if end if diff --git a/src/regolith/regolith_streamtube_head.f90 b/src/regolith/regolith_streamtube_head.f90 index bf8d886b..6d16bd00 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,meltinejecta) +subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector,meltinejecta,totvol) 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 - real(DP),intent(inout) :: meltinejecta + real(DP),intent(inout) :: meltinejecta,totvol real(DP),intent(in) :: deltar real(DP),intent(inout) :: totmare,tots real(SP),dimension(:),intent(inout) :: age_collector @@ -67,6 +67,7 @@ subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector age_collector(:) = age_collector(:) + current(N)%age(:) * recyratio headmeltvol = current(N)%meltfrac * recyratio * vsgly meltinejecta = meltinejecta + headmeltvol + totvol = totvol + vsgly else ! head is not intersected with layers. do @@ -80,6 +81,7 @@ subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector age_collector(:) = age_collector(:) + current(N)%age(:) * recyratio headmeltvol = current(N)%meltfrac * recyratio * tothead meltinejecta = meltinejecta + headmeltvol + totvol = totvol + tothead !current => current%next N = N - 1 z = z + current(N)%thickness @@ -92,6 +94,7 @@ subroutine regolith_streamtube_head(user,surfi,deltar,totmare,tots,age_collector age_collector(:) = age_collector(:) + current(N)%age(:) * recyratio headmeltvol = current(N)%meltfrac * recyratio * tothead meltinejecta = meltinejecta + headmeltvol + totvol = totvol + tothead exit end if end do diff --git a/src/regolith/regolith_streamtube_lineseg.f90 b/src/regolith/regolith_streamtube_lineseg.f90 index ddd2f104..70e648aa 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,meltinejecta) + age_collector,xmints,xsfints,depthb,meltinejecta,totvol) use module_globals use module_regolith, EXCEPT_THIS_ONE => regolith_streamtube_lineseg implicit none @@ -28,7 +28,7 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad type(surftype),intent(in) :: surfi real(DP),intent(in) :: thetast,ri,rip1,zmin,zmax,erad,eradi,deltar type(regodatatype),intent(inout) :: newlayer - real(DP),intent(inout) :: meltinejecta + real(DP),intent(inout) :: meltinejecta,totvol real(DP),intent(inout) :: vmare,totseb real(SP),dimension(:),intent(inout) :: age_collector real(DP),intent(in) :: xmints @@ -85,6 +85,7 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad vol = vol + sum(current(N-1)%age(:)) * recyratio linmelt = current(N-1)%meltfrac * vsgly * recyratio meltinejecta = meltinejecta + linmelt + totvol = totvol + vsgly 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) @@ -93,6 +94,7 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad vol = vol + sum(current(N-1)%age(:)) * recyratio linmelt = current(N-1)%meltfrac * vseg * recyratio meltinejecta = meltinejecta + linmelt + totvol = totvol + vseg end if exit else @@ -121,6 +123,7 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad vol = vol + sum(current(N)%age(:)) * recyratio linmelt = current(N)%meltfrac * vseg * recyratio meltinejecta = meltinejecta + linmelt + totvol = totvol + vseg end if ! A segment coming from the side of emerging location of a streamtube rip1 @@ -132,6 +135,7 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad vol = vol + sum(current(N)%age(:)) * recyratio linmelt = current(N)%meltfrac * vseg * recyratio meltinejecta = meltinejecta + linmelt + totvol = totvol + vseg end if !current => current%next N = N - 1 @@ -155,6 +159,7 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad vol = vol + sum(current(N)%age(:)) * recyratio linmelt = current(N)%meltfrac * vseg * recyratio meltinejecta = meltinejecta + linmelt + totvol = totvol + vseg end if exit @@ -166,6 +171,7 @@ subroutine regolith_streamtube_lineseg(user,surfi,thetast,ri,rip1,zmin,zmax,erad totseb = vsgly linmelt = current(N)%meltfrac * vsgly * recyratio meltinejecta = meltinejecta + linmelt + totvol = totvol + vsgly exit end if diff --git a/src/regolith/regolith_subpixel_streamtube.f90 b/src/regolith/regolith_subpixel_streamtube.f90 index 1fe3d24d..f9566a1d 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,meltinejecta) + age_collector,xmints,xsfints,vol,meltinejecta,totvol) use module_globals use module_regolith, EXCEPT_THIS_ONE => regolith_subpixel_streamtube implicit none @@ -61,7 +61,7 @@ subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer type(surftype),intent(in) :: surfi real(DP),intent(in) :: deltar,ri,rip1,eradi type(regodatatype),intent(inout) :: newlayer - real(DP),intent(out) :: meltinejecta + real(DP),intent(out) :: meltinejecta, totvol real(DP),intent(out) :: vmare,totseb real(SP),dimension(:),intent(inout) :: age_collector real(DP),intent(in) :: xmints @@ -99,6 +99,7 @@ subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer vsgly1 = 0.0_DP vsgly2 = 0.0_DP meltinejecta = 0.0_DP + totvol = 0.0_DP ! Two cases: subpixel is inside the first layer, and its volume is simply the landing ejecta blanket. if (zend>=zmax) then @@ -109,6 +110,7 @@ subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer vsh = regolith_shock_damage(eradi,deltar,xmints,xsfints,0.0_DP,eradi) recyratio = max(vseg-vsh,0.0 )/ (user%pix**2) / (surfi%regolayer(N)%thickness) meltinejecta = surfi%regolayer(N)%meltfrac * vseg * recyratio + totvol = vseg age_collector(:) = age_collector(:) + surfi%regolayer(N)%age(:) * recyratio vol = vol + sum(age_collector(:)) ! write(*,*) '1',eradi, xmints, xsfints, & @@ -197,6 +199,7 @@ subroutine regolith_subpixel_streamtube(user,surfi,deltar,ri,rip1,eradi,newlayer !current => current%next meltinejecta = meltinejecta + mvl + mvr + totvol = totvol + vsgly1 + vsgly2 N = N - 1 z = z + current(N)%thickness zstart = zend @@ -213,7 +216,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,meltinejecta) + call regolith_streamtube_head(user,surfi,deltar,vmare,totseb,age_collector,meltinejecta,totvol) return end subroutine regolith_subpixel_streamtube diff --git a/src/regolith/regolith_traverse_streamtube.f90 b/src/regolith/regolith_traverse_streamtube.f90 index ea083f9e..549780ea 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,meltinejecta) + age_collector,xmints,xsfints,depthb,meltinejecta,totvol) use module_globals use module_regolith, EXCEPT_THIS_ONE => regolith_traverse_streamtube implicit none @@ -29,7 +29,7 @@ subroutine regolith_traverse_streamtube(user,surfi,deltar,ri,rip1,eradi,erado,ne type(surftype),intent(inout) :: surfi real(DP),intent(in) :: deltar,ri,rip1,eradi,erado type(regodatatype),intent(inout) :: newlayer - real(DP),intent(out) :: meltinejecta + real(DP),intent(out) :: meltinejecta,totvol real(DP),intent(out) :: vmare,totseb real(SP),dimension(:),intent(inout) :: age_collector real(DP),intent(in) :: xmints @@ -53,6 +53,7 @@ subroutine regolith_traverse_streamtube(user,surfi,deltar,ri,rip1,eradi,erado,ne !executable code meltinejecta = 0.0_DP + totvol = 0.0_DP erad = (eradi + erado)/2.0 rzmax = erad * sqrt(3.0)/4.0 @@ -90,12 +91,13 @@ subroutine regolith_traverse_streamtube(user,surfi,deltar,ri,rip1,eradi,erado,ne recyratio = max(vseg-vsh,0.0_DP) / (user%pix**2) / surfi%regolayer(N)%thickness age_collector(:) = age_collector(:) + surfi%regolayer(N)%age(:) * recyratio meltinejecta = surfi%regolayer(N)%meltfrac * vseg * recyratio + totvol = vseg 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,meltinejecta) + newlayer,vmare,totseb,age_collector,xmints,xsfints,depthb,meltinejecta,totvol) end if