Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Katetc/dglc negfluxes #303

Merged
merged 10 commits into from
Sep 6, 2024
Merged

Katetc/dglc negfluxes #303

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
8e1f68d
Fixed usrf bug and first draft simple hole filling for neg qice
Katetc Aug 29, 2024
da714ce
Passing logunit to run comp subroutine
Katetc Aug 29, 2024
aedd722
Fix step 2, scaled removal by fraction of local positives
Katetc Aug 30, 2024
b73abc3
Third step complete, added global sums to spread calculation out over…
Katetc Aug 30, 2024
9a8afbc
Fixed area weighting on application
Katetc Aug 30, 2024
5251638
Removed second area weighting and second if clause to make sure the 0…
Katetc Sep 3, 2024
2e0a332
Merge up new changes from main
Katetc Sep 4, 2024
34838d2
One more review comment addressed
Katetc Sep 5, 2024
9b22a50
Final code cleanup and review comments
Katetc Sep 5, 2024
5b604f9
Merge branch 'main' into katetc/dglc_negfluxes
Katetc Sep 6, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
103 changes: 89 additions & 14 deletions dglc/dglc_datamode_noevolve_mod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,9 @@ module dglc_datamode_noevolve_mod
use ESMF , only : ESMF_Mesh, ESMF_DistGrid, ESMF_FieldBundle, ESMF_Field
use ESMF , only : ESMF_FieldBundleCreate, ESMF_FieldCreate, ESMF_MeshLoc_Element
use ESMF , only : ESMF_FieldBundleAdd, ESMF_MeshGet, ESMF_DistGridGet, ESMF_Typekind_R8
use ESMF , only : ESMF_VMGetCurrent, ESMF_VMBroadCast, ESMF_VM
use ESMF , only : ESMF_GridComp, ESMF_GridCompGet
use ESMF , only : ESMF_VM, ESMF_VMAllreduce, ESMF_REDUCE_SUM
use ESMF , only : ESMF_VMGetCurrent, ESMF_VMBroadCast
use NUOPC , only : NUOPC_Advertise, NUOPC_IsConnected
use shr_kind_mod , only : r8=>shr_kind_r8, i8=>shr_kind_i8, cl=>shr_kind_cl, cs=>shr_kind_cs
use shr_sys_mod , only : shr_sys_abort
Expand Down Expand Up @@ -217,15 +219,17 @@ subroutine dglc_datamode_noevolve_init_pointers(NStateExp, NstateImp, rc)
end subroutine dglc_datamode_noevolve_init_pointers

!===============================================================================
subroutine dglc_datamode_noevolve_advance(pio_subsystem, io_type, io_format, &
model_meshes, model_internal_gridsize, model_datafiles, rc)
subroutine dglc_datamode_noevolve_advance(gcomp, pio_subsystem, io_type, io_format, &
logunit, model_meshes, model_internal_gridsize, model_datafiles, rc)

! Assume that the model mesh is the same as the input data mesh

! input/output variables
! input/output variables
type(ESMF_GridComp) :: gcomp
type(iosystem_desc_t) , pointer :: pio_subsystem ! pio info
integer , intent(in) :: io_type ! pio info
integer , intent(in) :: io_format ! pio info
integer , intent(in) :: logunit ! For writing logs
type(ESMF_Mesh) , intent(in) :: model_meshes(:) ! ice sheets meshes
real(r8) , intent(in) :: model_internal_gridsize(:) ! internal model gridsizes (m)
character(len=*) , intent(in) :: model_datafiles(:) ! input file names
Expand All @@ -235,6 +239,7 @@ subroutine dglc_datamode_noevolve_advance(pio_subsystem, io_type, io_format, &
type(ESMF_FieldBundle) :: fldbun_noevolve
type(ESMF_DistGrid) :: distgrid
type(ESMF_Field) :: field_noevolve
type(ESMF_VM) :: vm
type(file_desc_t) :: pioid
type(io_desc_t) :: pio_iodesc
integer :: ns ! ice sheet index
Expand All @@ -254,6 +259,10 @@ subroutine dglc_datamode_noevolve_advance(pio_subsystem, io_type, io_format, &
real(r8) :: eus ! eustatic sea level
real(r8) :: lsrf ! lower surface elevation (m) on ice grid
real(r8) :: usrf ! upper surface elevation (m) on ice grid
real(r8) :: loc_pos_smb(1), Tot_pos_smb(1) ! Sum of positive smb values on each ice sheet for hole-filling
real(r8) :: loc_neg_smb(1), Tot_neg_smb(1) ! Sum of negative smb values on each ice sheet for hole-filling
real(r8) :: rat ! Ratio of hole-filling flux to apply

character(len=*), parameter :: subname='(dglc_datamode_noevolve_advance): '
!-------------------------------------------------------------------------------

Expand Down Expand Up @@ -380,25 +389,91 @@ subroutine dglc_datamode_noevolve_advance(pio_subsystem, io_type, io_format, &
end if

if (initialized_noevolve) then

! Compute Fgrg_rofi
do ns = 1,num_icesheets
call ESMF_MeshGet(model_meshes(ns), elementdistGrid=distGrid, rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call ESMF_DistGridGet(distGrid, localDe=0, elementCount=lsize, rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return

! Get number of grid cells per ice sheet
lsize = size(Fgrg_rofi(ns)%ptr)
billsacks marked this conversation as resolved.
Show resolved Hide resolved

! reset output variables to zero
Fgrg_rofi(ns)%ptr(:) = 0.d0
loc_pos_smb(1) = 0.d0
Tot_pos_smb(1) = 0.d0
loc_neg_smb(1) = 0.d0
Tot_neg_smb(1) = 0.d0
rat = 0.d0

! For No Evolve to reduce negative ice fluxes from DGLC, we will
! Calculate the total positive and total negative fluxes on each
! processor first (local totals).
do ng = 1,lsize
if (Sg_icemask_coupled_fluxes(ns)%ptr(ng) > 0.d0) then
Fgrg_rofi(ns)%ptr(ng) = Flgl_qice(ns)%ptr(ng)
else
Fgrg_rofi(ns)%ptr(ng) = 0._r8
if(Flgl_qice(ns)%ptr(ng) > 0.d0) then
loc_pos_smb(1) = loc_pos_smb(1)+Flgl_qice(ns)%ptr(ng)*Sg_area(ns)%ptr(ng)
end if
! Ignore places that are exactly 0.d0
if(Flgl_qice(ns)%ptr(ng) < 0.d0) then
loc_neg_smb(1) = loc_neg_smb(1)+Flgl_qice(ns)%ptr(ng)*Sg_area(ns)%ptr(ng)
end if
end if
end do
end do
end if
! Now do two global sums to get the ice sheet total positive
! and negative ice fluxes
call ESMF_GridCompGet(gcomp, vm=vm, rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call ESMF_VMAllreduce(vm, senddata=loc_pos_smb, recvdata=Tot_pos_smb, count=1, &
reduceflag=ESMF_REDUCE_SUM, rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
call ESMF_VMAllreduce(vm, senddata=loc_neg_smb, recvdata=Tot_neg_smb, count=1, &
reduceflag=ESMF_REDUCE_SUM, rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return

! If there's more positive than negative, then set all
! negative to zero and destribute the negative flux amount
! across the positive values, scaled by the size of the
! positive value. This section also applies to any chunks
! where there is no negative smb. In that case the ice
! runoff is exactly equal to the input smb.
if(abs(Tot_pos_smb(1)) >= abs(Tot_neg_smb(1))) then
do ng = 1,lsize
if (Sg_icemask_coupled_fluxes(ns)%ptr(ng) > 0.d0) then
if(Flgl_qice(ns)%ptr(ng) > 0.d0) then
rat = Flgl_qice(ns)%ptr(ng)/Tot_pos_smb(1)
Fgrg_rofi(ns)%ptr(ng) = Flgl_qice(ns)%ptr(ng) + rat*Tot_neg_smb(1)
else
Fgrg_rofi(ns)%ptr(ng) = 0.d0
end if
else
Fgrg_rofi(ns)%ptr(ng) = 0.d0
end if
end do
else
! If there's more negative than positive, set the positive to zero
! and distribute the positive amount to the negative spaces to
! reduce their negativity a bit. This shouldn't happen often.
! This section of code also applies if Tot_pos_smb is zero.
do ng = 1,lsize
if (Sg_icemask_coupled_fluxes(ns)%ptr(ng) > 0.d0) then
if(Flgl_qice(ns)%ptr(ng) < 0.d0) then
rat = Flgl_qice(ns)%ptr(ng)/Tot_neg_smb(1)
Fgrg_rofi(ns)%ptr(ng) = Flgl_qice(ns)%ptr(ng) + rat*Tot_pos_smb(1)
else
Fgrg_rofi(ns)%ptr(ng) = 0.d0
end if
else
Fgrg_rofi(ns)%ptr(ng) = 0.d0
end if
end do

end if ! More neg or pos smb

end do ! Each ice sheet
end if ! if initialized

! Set initialized flag
initialized_noevolve = .true.

end subroutine dglc_datamode_noevolve_advance

!===============================================================================
Expand Down
11 changes: 6 additions & 5 deletions dglc/glc_comp_nuopc.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -430,7 +430,7 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc)
end do ! end loop over ice sheets

! Run dglc
call dglc_comp_run(clock, current_ymd, current_tod, restart_write=.false., valid_inputs=.true., rc=rc)
call dglc_comp_run(gcomp, clock, current_ymd, current_tod, restart_write=.false., valid_inputs=.true., rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return

call ESMF_TraceRegionExit('dglc_strdata_init')
Expand Down Expand Up @@ -503,21 +503,22 @@ subroutine ModelAdvance(gcomp, rc)
end if

! run dglc
call dglc_comp_run(clock, next_ymd, next_tod, restart_write, valid_inputs, rc=rc)
call dglc_comp_run(gcomp, clock, next_ymd, next_tod, restart_write, valid_inputs, rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return

call ESMF_TraceRegionExit(subname)

end subroutine ModelAdvance

!===============================================================================
subroutine dglc_comp_run(clock, target_ymd, target_tod, restart_write, valid_inputs, rc)
subroutine dglc_comp_run(gcomp, clock, target_ymd, target_tod, restart_write, valid_inputs, rc)

! --------------------------
! advance dglc
! --------------------------

! input/output variables:
type(ESMF_GridComp) :: gcomp
type(ESMF_Clock) , intent(in) :: clock
integer , intent(in) :: target_ymd ! model date
integer , intent(in) :: target_tod ! model sec into model date
Expand Down Expand Up @@ -598,8 +599,8 @@ subroutine dglc_comp_run(clock, target_ymd, target_tod, restart_write, valid_inp
select case (trim(datamode))
case('noevolve')
if (valid_inputs) then
call dglc_datamode_noevolve_advance(sdat(1)%pio_subsystem, sdat(1)%io_type, sdat(1)%io_format, &
model_meshes, model_internal_gridsize, model_datafiles, rc)
call dglc_datamode_noevolve_advance(gcomp, sdat(1)%pio_subsystem, sdat(1)%io_type, sdat(1)%io_format, &
logunit, model_meshes, model_internal_gridsize, model_datafiles, rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
end if
end select
Expand Down
Loading