@@ -371,38 +371,45 @@ subroutine li_advection_thickness_tracers(&
371371 ! update masks
372372 call li_calculate_mask(meshPool, velocityPool, geometryPool, domain, err_tmp)
373373 err = ior (err, err_tmp)
374+ ! Update halo on edgeMask for calculating flux across grounding line.
375+ call mpas_timer_start(" halo updates" )
376+ call mpas_dmpar_field_halo_exch(domain, ' edgeMask' )
377+ call mpas_timer_stop(" halo updates" )
374378
375- layerThicknessEdgeFlux(:,:) = 0.0_RKIND
376379 !-----------------------------------------------------------------
377380 ! Horizontal transport of thickness and tracers
378381 !-----------------------------------------------------------------
379382
380383 if ( (trim (config_thickness_advection) .eq. ' fo' ) .or. &
381384 (trim (config_thickness_advection) .eq. ' fct' ) ) then
382-
383- ! Calculate flux across grounding line before applying
384- ! sfc/ basal mass balance and advection
385+ ! Calculate flux across grounding line
386+ ! Do this after new thickness & mask have been calculated, including SMB/ BMB
385387 fluxAcrossGroundingLine(:) = 0.0_RKIND
386388 fluxAcrossGroundingLineOnCells(:) = 0.0_RKIND
387389 do iEdge = 1 , nEdges
388- if (li_mask_is_grounding_line(edgeMask(iEdge))) then
390+ if (li_mask_is_grounding_line(edgeMask(iEdge))) then
389391 iCell1 = cellsOnEdge(1 ,iEdge)
390392 iCell2 = cellsOnEdge(2 ,iEdge)
391- if (li_mask_is_grounded_ice(cellMask(iCell1))) then
393+ if (li_mask_is_grounded_ice(cellMask(iCell1))) then
392394 GLfluxSign = 1.0_RKIND ! edge sign convention is positive from iCell1 to iCell2 on an edge
393395 theGroundedCell = iCell1
394- else
396+ else
395397 GLfluxSign = - 1.0_RKIND
396398 theGroundedCell = iCell2
397399 endif
398- do k = 1 , nVertLevels
399- thicknessFluxEdge = layerNormalVelocity(k,iEdge) * dvEdge(iEdge) * layerThicknessEdge(k,iEdge)
400- fluxAcrossGroundingLine(iEdge) = fluxAcrossGroundingLine(iEdge) + GLfluxSign * thicknessFluxEdge / dvEdge(iEdge)
401- enddo
400+ if (trim (config_thickness_advection) == ' fct' ) then
401+ do k = 1 , nVertLevels
402+ fluxAcrossGroundingLine(iEdge) = fluxAcrossGroundingLine(iEdge) + GLfluxSign * layerThicknessEdgeFlux(k,iEdge)
403+ enddo
404+ else
405+ do k = 1 , nVertLevels
406+ layerThicknessEdgeFlux(k,iEdge) = layerNormalVelocity(k,iEdge) * layerThicknessEdge(k,iEdge)
407+ fluxAcrossGroundingLine(iEdge) = fluxAcrossGroundingLine(iEdge) + GLfluxSign * layerThicknessEdgeFlux(k,iEdge)
408+ enddo
409+ endif
402410 ! assign to grounded cell in fluxAcrossGroundingLineOnCells
403- if (thickness(theGroundedCell) <= 0.0_RKIND ) then
404- ! This should never be the case, but checking to avoid
405- ! possible divide by zero
411+ if (thickness(theGroundedCell) <= 0.0_RKIND ) then
412+ ! This should never be the case, but checking to avoid possible divide by zero
406413 call mpas_log_write(" thickness at a grounding line is unexepectedly <=0" , MPAS_LOG_ERR)
407414 err = ior (err, 1 )
408415 return
@@ -476,7 +483,7 @@ subroutine li_advection_thickness_tracers(&
476483 ! copy old (input) values of layer thickness and tracers to scratch arrays
477484 layerThicknessOld(:,:) = layerThickness(:,:)
478485 advectedTracersOld(:,:,:) = advectedTracers(:,:,:)
479-
486+ layerThicknessEdgeFlux(:,:) = 0.0_RKIND
480487 ! compute new values of layer thickness and tracers
481488 if ((trim (config_thickness_advection) .eq. ' fo' ) .and. &
482489 ( (trim (config_tracer_advection) .eq. ' fo' ) .or. &
0 commit comments