diff --git a/cicecore/cicedynB/dynamics/ice_dyn_evp.F90 b/cicecore/cicedynB/dynamics/ice_dyn_evp.F90 index d8ce42681..6c1fb48c7 100644 --- a/cicecore/cicedynB/dynamics/ice_dyn_evp.F90 +++ b/cicecore/cicedynB/dynamics/ice_dyn_evp.F90 @@ -88,7 +88,7 @@ subroutine evp (dt) stress12_1, stress12_2, stress12_3, stress12_4 use ice_grid, only: tmask, umask, dxt, dyt, dxhy, dyhx, cxp, cyp, cxm, cym, & tarear, uarear, tinyarea, to_ugrid, t2ugrid_vector, u2tgrid_vector, & - grid_type, HTE, HTN + grid_type use ice_state, only: aice, vice, vsno, uvel, vvel, divu, shear, & aice_init, aice0, aicen, vicen, strength use ice_timers, only: timer_dynamics, timer_bound, & @@ -366,7 +366,6 @@ subroutine evp (dt) endif call ice_dyn_evp_1d_copyin( & nx_block,ny_block,nblocks,nx_global+2*nghost,ny_global+2*nghost, & - HTE,HTN, & !v1 dxhy,dyhx,cyp,cxp,cym,cxm,tinyarea, & !v1 waterx,watery, & icetmask, iceumask, & diff --git a/cicecore/cicedynB/dynamics/ice_dyn_evp_1d.F90 b/cicecore/cicedynB/dynamics/ice_dyn_evp_1d.F90 index 78469cc86..351bbaed4 100644 --- a/cicecore/cicedynB/dynamics/ice_dyn_evp_1d.F90 +++ b/cicecore/cicedynB/dynamics/ice_dyn_evp_1d.F90 @@ -1135,7 +1135,6 @@ end subroutine dealloc1d !---------------------------------------------------------------------------- subroutine evp_copyin_v2(nx,ny,nblk,nx_glob,ny_glob, & - I_HTE,I_HTN, & !v1 I_dxhy,I_dyhx,I_cyp,I_cxp,I_cym,I_cxm,I_tinyarea, & !v1 I_waterx,I_watery, & I_icetmask,I_iceumask, & @@ -1149,6 +1148,8 @@ subroutine evp_copyin_v2(nx,ny,nblk,nx_glob,ny_glob, use ice_gather_scatter, only: gather_global_ext use ice_domain, only: distrb_info use ice_communicate, only: my_task, master_task + use ice_grid, only: G_HTE, G_HTN + use ice_constants, only: c0 implicit none @@ -1156,7 +1157,6 @@ subroutine evp_copyin_v2(nx,ny,nblk,nx_glob,ny_glob, integer (kind=int_kind),dimension (nx,ny,nblk), intent(in) :: I_icetmask logical (kind=log_kind),dimension (nx,ny,nblk), intent(in) :: I_iceumask real (kind=dbl_kind), dimension(nx,ny,nblk), intent(in) :: & - I_HTE,I_HTN, & !v1 I_dxhy,I_dyhx,I_cyp,I_cxp,I_cym,I_cxm,I_tinyarea, & !v1 I_waterx,I_watery, & I_cdn_ocn,I_aiu,I_uocn,I_vocn,I_forcex,I_forcey,I_Tbu, & @@ -1171,7 +1171,6 @@ subroutine evp_copyin_v2(nx,ny,nblk,nx_glob,ny_glob, integer (kind=int_kind),dimension (nx_glob,ny_glob) :: G_icetmask logical (kind=log_kind),dimension (nx_glob,ny_glob) :: G_iceumask real (kind=dbl_kind), dimension(nx_glob,ny_glob) :: & - G_HTE,G_HTN, & !v1 G_dxhy,G_dyhx,G_cyp,G_cxp,G_cym,G_cxm,G_tinyarea, & !v1 G_waterx,G_watery, & G_cdn_ocn,G_aiu,G_uocn,G_vocn,G_forcex,G_forcey,G_Tbu, & @@ -1186,51 +1185,49 @@ subroutine evp_copyin_v2(nx,ny,nblk,nx_glob,ny_glob, !--------------------------------------- !-- Gather data into one single block -- - call gather_global_ext(G_icetmask, I_icetmask, master_task, distrb_info) - call gather_global_ext(G_iceumask, I_iceumask, master_task, distrb_info) - call gather_global_ext(G_HTE, I_HTE, master_task, distrb_info) - call gather_global_ext(G_HTN, I_HTN, master_task, distrb_info) -!v1 call gather_global_ext(G_dxhy, I_dxhy, master_task, distrb_info) -!v1 call gather_global_ext(G_dyhx, I_dyhx, master_task, distrb_info) -!v1 call gather_global_ext(G_cyp, I_cyp, master_task, distrb_info) -!v1 call gather_global_ext(G_cxp, I_cxp, master_task, distrb_info) -!v1 call gather_global_ext(G_cym, I_cym, master_task, distrb_info) -!v1 call gather_global_ext(G_cxm, I_cxm, master_task, distrb_info) -!v1 call gather_global_ext(G_tinyarea, I_tinyarea, master_task, distrb_info) -!v1 call gather_global_ext(G_waterx, I_waterx, master_task, distrb_info) -!v1 call gather_global_ext(G_watery, I_watery, master_task, distrb_info) - call gather_global_ext(G_cdn_ocn, I_cdn_ocn, master_task, distrb_info) - call gather_global_ext(G_aiu, I_aiu, master_task, distrb_info) - call gather_global_ext(G_uocn, I_uocn, master_task, distrb_info) - call gather_global_ext(G_vocn, I_vocn, master_task, distrb_info) - call gather_global_ext(G_forcex, I_forcex, master_task, distrb_info) - call gather_global_ext(G_forcey, I_forcey, master_task, distrb_info) - call gather_global_ext(G_Tbu, I_Tbu, master_task, distrb_info) - call gather_global_ext(G_umassdti, I_umassdti, master_task, distrb_info) - call gather_global_ext(G_fm, I_fm, master_task, distrb_info) - call gather_global_ext(G_uarear, I_uarear, master_task, distrb_info) - call gather_global_ext(G_tarear, I_tarear, master_task, distrb_info) - call gather_global_ext(G_strintx, I_strintx, master_task, distrb_info) - call gather_global_ext(G_strinty, I_strinty, master_task, distrb_info) - call gather_global_ext(G_uvel_init, I_uvel_init, master_task, distrb_info) - call gather_global_ext(G_vvel_init, I_vvel_init, master_task, distrb_info) - call gather_global_ext(G_strength, I_strength, master_task, distrb_info) - call gather_global_ext(G_uvel, I_uvel, master_task, distrb_info) - call gather_global_ext(G_vvel, I_vvel, master_task, distrb_info) - call gather_global_ext(G_dxt, I_dxt, master_task, distrb_info) - call gather_global_ext(G_dyt, I_dyt, master_task, distrb_info) - call gather_global_ext(G_stressp_1, I_stressp_1, master_task, distrb_info) - call gather_global_ext(G_stressp_2, I_stressp_2, master_task, distrb_info) - call gather_global_ext(G_stressp_3, I_stressp_3, master_task, distrb_info) - call gather_global_ext(G_stressp_4, I_stressp_4, master_task, distrb_info) - call gather_global_ext(G_stressm_1, I_stressm_1, master_task, distrb_info) - call gather_global_ext(G_stressm_2, I_stressm_2, master_task, distrb_info) - call gather_global_ext(G_stressm_3, I_stressm_3, master_task, distrb_info) - call gather_global_ext(G_stressm_4, I_stressm_4, master_task, distrb_info) - call gather_global_ext(G_stress12_1, I_stress12_1, master_task, distrb_info) - call gather_global_ext(G_stress12_2, I_stress12_2, master_task, distrb_info) - call gather_global_ext(G_stress12_3, I_stress12_3, master_task, distrb_info) - call gather_global_ext(G_stress12_4, I_stress12_4, master_task, distrb_info) + call gather_global_ext(G_icetmask, I_icetmask, master_task, distrb_info, 0 ) + call gather_global_ext(G_iceumask, I_iceumask, master_task, distrb_info, .false.) +!v1 call gather_global_ext(G_dxhy, I_dxhy, master_task, distrb_info ) +!v1 call gather_global_ext(G_dyhx, I_dyhx, master_task, distrb_info ) +!v1 call gather_global_ext(G_cyp, I_cyp, master_task, distrb_info ) +!v1 call gather_global_ext(G_cxp, I_cxp, master_task, distrb_info ) +!v1 call gather_global_ext(G_cym, I_cym, master_task, distrb_info ) +!v1 call gather_global_ext(G_cxm, I_cxm, master_task, distrb_info ) +!v1 call gather_global_ext(G_tinyarea, I_tinyarea, master_task, distrb_info ) +!v1 call gather_global_ext(G_waterx, I_waterx, master_task, distrb_info ) +!v1 call gather_global_ext(G_watery, I_watery, master_task, distrb_info ) + call gather_global_ext(G_cdn_ocn, I_cdn_ocn, master_task, distrb_info ) + call gather_global_ext(G_aiu, I_aiu, master_task, distrb_info ) + call gather_global_ext(G_uocn, I_uocn, master_task, distrb_info ) + call gather_global_ext(G_vocn, I_vocn, master_task, distrb_info ) + call gather_global_ext(G_forcex, I_forcex, master_task, distrb_info ) + call gather_global_ext(G_forcey, I_forcey, master_task, distrb_info ) + call gather_global_ext(G_Tbu, I_Tbu, master_task, distrb_info ) + call gather_global_ext(G_umassdti, I_umassdti, master_task, distrb_info ) + call gather_global_ext(G_fm, I_fm, master_task, distrb_info ) + call gather_global_ext(G_uarear, I_uarear, master_task, distrb_info ) + call gather_global_ext(G_tarear, I_tarear, master_task, distrb_info ) + call gather_global_ext(G_strintx, I_strintx, master_task, distrb_info ) + call gather_global_ext(G_strinty, I_strinty, master_task, distrb_info ) + call gather_global_ext(G_uvel_init, I_uvel_init, master_task, distrb_info ) + call gather_global_ext(G_vvel_init, I_vvel_init, master_task, distrb_info ) + call gather_global_ext(G_strength, I_strength, master_task, distrb_info ) + call gather_global_ext(G_uvel, I_uvel, master_task, distrb_info, c0 ) + call gather_global_ext(G_vvel, I_vvel, master_task, distrb_info, c0 ) + call gather_global_ext(G_dxt, I_dxt, master_task, distrb_info ) + call gather_global_ext(G_dyt, I_dyt, master_task, distrb_info ) + call gather_global_ext(G_stressp_1, I_stressp_1, master_task, distrb_info ) + call gather_global_ext(G_stressp_2, I_stressp_2, master_task, distrb_info ) + call gather_global_ext(G_stressp_3, I_stressp_3, master_task, distrb_info ) + call gather_global_ext(G_stressp_4, I_stressp_4, master_task, distrb_info ) + call gather_global_ext(G_stressm_1, I_stressm_1, master_task, distrb_info ) + call gather_global_ext(G_stressm_2, I_stressm_2, master_task, distrb_info ) + call gather_global_ext(G_stressm_3, I_stressm_3, master_task, distrb_info ) + call gather_global_ext(G_stressm_4, I_stressm_4, master_task, distrb_info ) + call gather_global_ext(G_stress12_1, I_stress12_1, master_task, distrb_info ) + call gather_global_ext(G_stress12_2, I_stress12_2, master_task, distrb_info ) + call gather_global_ext(G_stress12_3, I_stress12_3, master_task, distrb_info ) + call gather_global_ext(G_stress12_4, I_stress12_4, master_task, distrb_info ) !-- All calculations has to be done on the master-task -- diff --git a/cicecore/cicedynB/general/ice_init.F90 b/cicecore/cicedynB/general/ice_init.F90 index 91c0d10b1..f98c16b10 100644 --- a/cicecore/cicedynB/general/ice_init.F90 +++ b/cicecore/cicedynB/general/ice_init.F90 @@ -94,7 +94,8 @@ subroutine input_data bathymetry_file, use_bathymetry, & bathymetry_format, & grid_type, grid_format, & - dxrect, dyrect + dxrect, dyrect, & + pgl_global_ext use ice_dyn_shared, only: ndte, kdyn, revised_evp, yield_curve, & kevp_kernel, & seabed_stress, seabed_stress_method, & @@ -314,6 +315,7 @@ subroutine input_data ndtd = 1 ! dynamic time steps per thermodynamic time step ndte = 120 ! subcycles per dynamics timestep: ndte=dt_dyn/dte kevp_kernel = 0 ! EVP kernel (0 = 2D, >0: 1D. Only ver. 2 is implemented yet) + pgl_global_ext = .false. ! if true, init primary grid lebgths (global ext.) brlx = 300.0_dbl_kind ! revised_evp values. Otherwise overwritten in ice_dyn_shared arlx = 300.0_dbl_kind ! revised_evp values. Otherwise overwritten in ice_dyn_shared revised_evp = .false. ! if true, use revised procedure for evp dynamics @@ -637,6 +639,7 @@ subroutine input_data call broadcast_scalar(ndtd, master_task) call broadcast_scalar(ndte, master_task) call broadcast_scalar(kevp_kernel, master_task) + call broadcast_scalar(pgl_global_ext, master_task) call broadcast_scalar(brlx, master_task) call broadcast_scalar(arlx, master_task) call broadcast_scalar(revised_evp, master_task) @@ -1749,6 +1752,7 @@ subroutine input_data if (kevp_kernel /= 0) then if (kevp_kernel == 102) then kevp_kernel = 2 + if (my_task == master_task) pgl_global_ext = .true. else if (my_task == master_task) write(nu_diag,*) subname//' ERROR: kevp_kernel = ',kevp_kernel if (kevp_kernel == 2) then diff --git a/cicecore/cicedynB/infrastructure/ice_grid.F90 b/cicecore/cicedynB/infrastructure/ice_grid.F90 index a354efb6b..874902d1f 100644 --- a/cicecore/cicedynB/infrastructure/ice_grid.F90 +++ b/cicecore/cicedynB/infrastructure/ice_grid.F90 @@ -77,6 +77,10 @@ module ice_grid ocn_gridcell_frac ! only relevant for lat-lon grids ! gridcell value of [1 - (land fraction)] (T-cell) + real (kind=dbl_kind), dimension (:,:), allocatable, public :: & + G_HTE , & ! length of eastern edge of T-cell (global ext.) + G_HTN ! length of northern edge of T-cell (global ext.) + real (kind=dbl_kind), dimension (:,:,:), allocatable, public :: & cyp , & ! 1.5*HTE - 0.5*HTE cxp , & ! 1.5*HTN - 0.5*HTN @@ -125,7 +129,8 @@ module ice_grid kmt ! ocean topography mask for bathymetry (T-cell) logical (kind=log_kind), public :: & - use_bathymetry ! flag for reading in bathymetry_file + use_bathymetry, & ! flag for reading in bathymetry_file + pgl_global_ext ! flag for init primary grid lengths (global ext.) logical (kind=log_kind), & dimension (:,:,:), allocatable, public :: & @@ -153,6 +158,8 @@ subroutine alloc_grid integer (int_kind) :: ierr + character(len=*), parameter :: subname = '(alloc_grid)' + allocate( & dxt (nx_block,ny_block,max_blocks), & ! width of T-cell through the middle (m) dyt (nx_block,ny_block,max_blocks), & ! height of T-cell through the middle (m) @@ -203,7 +210,15 @@ subroutine alloc_grid mse (2,2,nx_block,ny_block,max_blocks), & msw (2,2,nx_block,ny_block,max_blocks), & stat=ierr) - if (ierr/=0) call abort_ice('(alloc_grid): Out of memory') + if (ierr/=0) call abort_ice(subname//'ERROR: Out of memory') + + if (pgl_global_ext) then + allocate( & + G_HTE(nx_global+2*nghost, ny_global+2*nghost), & ! length of eastern edge of T-cell (global ext.) + G_HTN(nx_global+2*nghost, ny_global+2*nghost), & ! length of northern edge of T-cell (global ext.) + stat=ierr) + if (ierr/=0) call abort_ice(subname//'ERROR: Out of memory') + endif end subroutine alloc_grid @@ -1498,6 +1513,9 @@ subroutine primary_grid_lengths_HTN(work_g) enddo enddo endif + if (pgl_global_ext) then + call primary_grid_lengths_global_ext(G_HTN, work_g) + endif call scatter_global(HTN, work_g, master_task, distrb_info, & field_loc_Nface, field_type_scalar) call scatter_global(dxu, work_g2, master_task, distrb_info, & @@ -1572,6 +1590,9 @@ subroutine primary_grid_lengths_HTE(work_g) enddo endif endif + if (pgl_global_ext) then + call primary_grid_lengths_global_ext(G_HTE, work_g) + endif call scatter_global(HTE, work_g, master_task, distrb_info, & field_loc_Eface, field_type_scalar) call scatter_global(dyu, work_g2, master_task, distrb_info, & @@ -2555,6 +2576,115 @@ subroutine read_seabedstress_bathy end subroutine read_seabedstress_bathy +!======================================================================= +! Initialize global primary grid lengths array with ghost cells from +! global primary grid lengths array + + subroutine primary_grid_lengths_global_ext(ARRAY_O, ARRAY_I) + + use ice_constants, only: c0 + + real (kind=dbl_kind), dimension(:,:), intent(in) :: & + ARRAY_I + + real (kind=dbl_kind), dimension(:,:), intent(out) :: & + ARRAY_O + + ! Local variables + + integer (kind=int_kind) :: & + ii, io, ji, jo + + character(len=*), parameter :: & + subname = '(primary_grid_lengths_global_ext)' + + if ((ns_boundary_type == 'tripole' ) .or. & + (ns_boundary_type == 'tripoleT')) then + call abort_ice(subname // 'ERROR: ' // & + ns_boundary_type // ' bndy type not impl for cfg') + endif + + do jo = 1, (ny_global + 2 * nghost) + ji = -nghost + jo + + ! Southern ghost cells + + if (ji < 1) then + select case (ns_boundary_type) + case ('cyclic') + ji = ji + ny_global + case ('open') + ji = nghost - jo + 1 + case ('closed') + ji = 0 + case default + call abort_ice( & + subname // 'ERROR: unknown n-s bndy type') + end select + endif + + ! Northern ghost cells + + if (ji > ny_global) then + select case (ns_boundary_type) + case ('cyclic') + ji = ji - ny_global + case ('open') + ji = 2 * ny_global - ji + 1 + case ('closed') + ji = 0 + case default + call abort_ice( & + subname // 'ERROR: unknown n-s bndy type') + end select + endif + + do io = 1, (nx_global + 2 * nghost) + ii = -nghost + io + + ! Western ghost cells + + if (ii < 1) then + select case (ew_boundary_type) + case ('cyclic') + ii = ii + nx_global + case ('open') + ii = nghost - io + 1 + case ('closed') + ii = 0 + case default + call abort_ice( & + subname//'ERROR: unknown e-w bndy type') + end select + endif + + ! Eastern ghost cells + + if (ii > nx_global) then + select case (ew_boundary_type) + case ('cyclic') + ii = ii - nx_global + case ('open') + ii = 2 * nx_global - ii + 1 + case ('closed') + ii = 0 + case default + call abort_ice( & + subname//'ERROR: unknown e-w bndy type') + end select + endif + + if ((ii == 0) .or. (ji == 0)) then + ARRAY_O(io, jo) = c0 + else + ARRAY_O(io, jo) = ARRAY_I(ii, ji) + endif + + enddo + enddo + + end subroutine primary_grid_lengths_global_ext + !======================================================================= end module ice_grid