From 0e240abbc93ba8269c98f091d60f37ca4375c8a7 Mon Sep 17 00:00:00 2001 From: Spencer Clark Date: Fri, 17 Apr 2020 16:59:33 -0400 Subject: [PATCH] Some initial infrastructure for online coarse-graining (#26) * Add basic example of a coarse-grained diagnostic * Slight refactor; added error checking * Minor cleanup * Remove placeholder implementation of weighted_block_average; add static diagnostic fields * Make sure fine array bounds are properly initialized * Reorder subroutines * Clean up static diagnostics * Remove unneeded block_sum_3d * Reduce things down to a minimal example * Simplify interface some * Shorten names somewhat * Address Jeremy's commemnts * Add more comments * Address Oli's comment; a few other minor cleanups * Move do_coarse_graining flag up a level --- .../driver/fvGFS/atmosphere.F90 | 17 ++ FV3/atmos_cubed_sphere/makefile | 2 + FV3/atmos_cubed_sphere/model/fv_arrays.F90 | 36 ++- FV3/atmos_cubed_sphere/model/fv_control.F90 | 4 +- .../tools/coarse_grained_diagnostics.F90 | 135 +++++++++ .../tools/coarse_graining.F90 | 271 ++++++++++++++++++ 6 files changed, 462 insertions(+), 3 deletions(-) create mode 100644 FV3/atmos_cubed_sphere/tools/coarse_grained_diagnostics.F90 create mode 100644 FV3/atmos_cubed_sphere/tools/coarse_graining.F90 diff --git a/FV3/atmos_cubed_sphere/driver/fvGFS/atmosphere.F90 b/FV3/atmos_cubed_sphere/driver/fvGFS/atmosphere.F90 index 6930442b7..9767ffafa 100644 --- a/FV3/atmos_cubed_sphere/driver/fvGFS/atmosphere.F90 +++ b/FV3/atmos_cubed_sphere/driver/fvGFS/atmosphere.F90 @@ -195,6 +195,8 @@ module atmosphere_mod a_step, p_step, current_time_in_seconds use mpp_domains_mod, only: mpp_get_data_domain, mpp_get_compute_domain +use coarse_graining_mod, only: coarse_graining_init +use coarse_grained_diagnostics_mod, only: fv_coarse_diag_init, fv_coarse_diag !$ser verbatim use k_checkpoint, only: set_nz implicit none @@ -335,6 +337,11 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Grid_box, area) if (grids_on_this_pe(n)) mytile = n enddo + if (Atm(mytile)%flagstruct%do_coarse_graining) then + call coarse_graining_init(Atm(mytile)%flagstruct%npx, Atm(mytile)%npz, & + Atm(mytile)%layout, Atm(mytile)%bd, Atm(mytile)%coarse_graining) + endif + Atm(mytile)%Time_init = Time_init a_step = 0 @@ -426,6 +433,10 @@ subroutine atmosphere_init (Time_init, Time, Time_step, Grid_box, area) !I've had trouble getting this to work with multiple grids at a time; worth revisiting? call fv_diag_init(Atm(mytile:mytile), Atm(mytile)%atmos_axes, Time, npx, npy, npz, Atm(mytile)%flagstruct%p_ref) + if (Atm(mytile)%coarse_graining%do_coarse_graining) then + call fv_coarse_diag_init(Atm(mytile)%bd, Time, Atm(mytile)%atmos_axes(3), & + Atm(mytile)%atmos_axes(4), Atm(mytile)%coarse_graining) + endif !---------- reference profile ----------- ps1 = 101325. ps2 = 81060. @@ -779,6 +790,9 @@ subroutine atmosphere_end (Time, Grid_box) call timing_on('FV_DIAG') call fv_diag(Atm(mytile:mytile), zvir, fv_time, Atm(mytile)%flagstruct%print_freq) call fv_nggps_diag(Atm(mytile:mytile), zvir, fv_time) + if (Atm(mytile)%coarse_graining%do_coarse_graining) then + call fv_coarse_diag(Atm(mytile:mytile), fv_time) + endif first_diag = .false. call timing_off('FV_DIAG') endif @@ -1627,6 +1641,9 @@ subroutine atmosphere_state_update (Time, IPD_Data, IAU_Data, Atm_block, flip_vc call nullify_domain() call timing_on('FV_DIAG') call fv_diag(Atm(mytile:mytile), zvir, fv_time, Atm(mytile)%flagstruct%print_freq) + if (Atm(mytile)%coarse_graining%do_coarse_graining) then + call fv_coarse_diag(Atm(mytile:mytile), fv_time) + endif first_diag = .false. call timing_off('FV_DIAG') diff --git a/FV3/atmos_cubed_sphere/makefile b/FV3/atmos_cubed_sphere/makefile index 58efdcbdf..b981f758a 100644 --- a/FV3/atmos_cubed_sphere/makefile +++ b/FV3/atmos_cubed_sphere/makefile @@ -55,6 +55,7 @@ SRCS_F90 = \ ./model/nh_utils.F90 \ ./tools/external_ic.F90 \ ./tools/external_sst.F90 \ + ./tools/coarse_grained_diagnostics.F90 \ ./tools/fv_diagnostics.F90 \ ./tools/fv_eta.F90 \ ./tools/fv_grid_tools.F90 \ @@ -67,6 +68,7 @@ SRCS_F90 = \ ./tools/fv_surf_map.F90 \ ./tools/fv_timing.F90 \ ./tools/init_hydro.F90 \ + ./tools/coarse_graining.F90 \ ./tools/sim_nc_mod.F90 \ ./tools/sorted_index.F90 \ ./tools/test_cases.F90 \ diff --git a/FV3/atmos_cubed_sphere/model/fv_arrays.F90 b/FV3/atmos_cubed_sphere/model/fv_arrays.F90 index dbfff7a0f..bc7741278 100644 --- a/FV3/atmos_cubed_sphere/model/fv_arrays.F90 +++ b/FV3/atmos_cubed_sphere/model/fv_arrays.F90 @@ -116,6 +116,12 @@ module fv_arrays_mod end type fv_diag_type + type fv_coarse_diag_type + + integer :: id_omega_coarse + integer :: n_3d_diagnostics = 1 + + end type fv_coarse_diag_type !>@brief The type 'fv_grid_type' is made up of grid-dependent information from fv_grid_tools and fv_grid_utils. !>@details It should not contain any user options (that goes in a different structure) nor data which @@ -1015,10 +1021,10 @@ module fv_arrays_mod !f1p logical :: adj_mass_vmr = .false. !TER: This is to reproduce answers for verona patch. This default can be changed ! to .true. in the next city release if desired - + + logical :: do_coarse_graining = .false. ! Whether to enable online coarse-graining of restart files and diagnostics !integer, pointer :: test_case !real, pointer :: alpha - end type fv_flags_type type fv_nest_BC_type_3D @@ -1158,6 +1164,30 @@ module fv_arrays_mod end type fv_grid_bounds_type + type fv_coarse_grid_bounds_type + + integer :: is_coarse, ie_coarse, js_coarse, je_coarse + + end type fv_coarse_grid_bounds_type + + type fv_coarse_graining_type + + type(fv_coarse_grid_bounds_type) :: bd + type(domain2d) :: domain + integer :: factor + integer :: nx_coarse + integer :: id_x_coarse ! diagnostic x-axis id for data on x-edges + integer :: id_y_coarse ! diagnostic y-axis id for data on y-edges + integer :: id_xt_coarse ! diagnostic x-axis id for data on x-centers + integer :: id_yt_coarse ! diagnostic y-axis id for data on y-centers + integer :: id_pfull ! diagnostic vertical axis id for data on z-centers + integer :: id_phalf ! diagnostic vertical axis id for data on z-edges + character(len=64) :: strategy ! Current valid values are: 'model_level' + logical :: do_coarse_graining = .false. + type(fv_coarse_diag_type) :: idiag ! container for coarse diagnostic ids + + end type fv_coarse_graining_type + type fv_regional_bc_bounds_type integer :: is_north ,ie_north ,js_north ,je_north & @@ -1328,6 +1358,8 @@ module fv_arrays_mod type(nudge_diag_type) :: nudge_diag + type(fv_coarse_graining_type) :: coarse_graining + end type fv_atmos_type contains diff --git a/FV3/atmos_cubed_sphere/model/fv_control.F90 b/FV3/atmos_cubed_sphere/model/fv_control.F90 index e87b01963..1b9f2715d 100644 --- a/FV3/atmos_cubed_sphere/model/fv_control.F90 +++ b/FV3/atmos_cubed_sphere/model/fv_control.F90 @@ -335,6 +335,7 @@ module fv_control_mod real, pointer :: s_weight, update_blend integer, pointer :: layout(:), io_layout(:) + logical, pointer :: do_coarse_graining integer :: ntilesMe ! Number of tiles on this process =1 for now @@ -669,7 +670,7 @@ subroutine run_setup(Atm, dt_atmos, grids_on_this_pe, p_split) nested, twowaynest, parent_grid_num, parent_tile, nudge_qv, & refinement, nestbctype, nestupdate, nsponge, s_weight, & ioffset, joffset, check_negative, nudge_ic, halo_update_type, gfs_phil, agrid_vel_rst, & - do_uni_zfull, adj_mass_vmr, fac_n_spl, fhouri, regional, bc_update_interval + do_uni_zfull, adj_mass_vmr, fac_n_spl, fhouri, regional, bc_update_interval, do_coarse_graining namelist /test_case_nml/test_case, bubble_do, alpha, nsolitons, soliton_Umax, soliton_size #ifdef MULTI_GASES @@ -1340,6 +1341,7 @@ subroutine setup_pointers(Atm) layout => Atm%layout io_layout => Atm%io_layout + do_coarse_graining => Atm%flagstruct%do_coarse_graining end subroutine setup_pointers diff --git a/FV3/atmos_cubed_sphere/tools/coarse_grained_diagnostics.F90 b/FV3/atmos_cubed_sphere/tools/coarse_grained_diagnostics.F90 new file mode 100644 index 000000000..f1b9470d4 --- /dev/null +++ b/FV3/atmos_cubed_sphere/tools/coarse_grained_diagnostics.F90 @@ -0,0 +1,135 @@ +module coarse_grained_diagnostics_mod + + use diag_manager_mod, only: diag_axis_init, register_diag_field, register_static_field, send_data + use fv_arrays_mod, only: fv_atmos_type, fv_coarse_diag_type, fv_coarse_graining_type, fv_grid_bounds_type + use mpp_domains_mod, only: domain2d + use mpp_mod, only: FATAL, mpp_error + use coarse_graining_mod, only: block_sum, get_fine_array_bounds, get_coarse_array_bounds, MODEL_LEVEL, weighted_block_average + use time_manager_mod, only: time_type + + implicit none + private + + public :: fv_coarse_diag_init, fv_coarse_diag + + integer :: tile_count = 1 ! Following fv_diagnostics.F90 + +contains + + subroutine fv_coarse_diag_init(bd, Time, id_pfull, id_phalf, coarse_graining) + type(fv_grid_bounds_type), intent(in) :: bd + type(time_type), intent(in) :: Time + integer, intent(in) :: id_pfull, id_phalf + type(fv_coarse_graining_type), intent(inout) :: coarse_graining + + integer :: is, ie, js, je, is_coarse, ie_coarse, js_coarse, je_coarse + + call get_fine_array_bounds(bd, is, ie, js, je) + call get_coarse_array_bounds(coarse_graining%bd, is_coarse, ie_coarse, js_coarse, je_coarse) + call initialize_coarse_diagnostic_axes(coarse_graining%domain, coarse_graining%nx_coarse, & + coarse_graining%id_x_coarse, coarse_graining%id_y_coarse, coarse_graining%id_xt_coarse, & + coarse_graining%id_yt_coarse) + + coarse_graining%id_pfull = id_pfull + coarse_graining%id_phalf = id_phalf + + call register_coarse_diagnostics(coarse_graining%idiag, Time, & + coarse_graining%id_xt_coarse, coarse_graining%id_yt_coarse, id_pfull) + end subroutine fv_coarse_diag_init + + subroutine initialize_coarse_diagnostic_axes(coarse_domain, & + nx_coarse, id_x_coarse, id_y_coarse, id_xt_coarse, id_yt_coarse) + type(domain2d), intent(in) :: coarse_domain + integer, intent(in) :: nx_coarse + integer, intent(inout) :: id_x_coarse, id_y_coarse, id_xt_coarse, id_yt_coarse + + integer :: i, j + real, allocatable :: grid_x_coarse(:), grid_y_coarse(:), grid_xt_coarse(:), grid_yt_coarse(:) + + allocate(grid_x_coarse(nx_coarse + 1)) + allocate(grid_y_coarse(nx_coarse + 1)) + allocate(grid_xt_coarse(nx_coarse)) + allocate(grid_yt_coarse(nx_coarse)) + + grid_x_coarse = (/ (i, i=1, nx_coarse + 1) /) + grid_y_coarse = (/ (j, j=1, nx_coarse + 1) /) + grid_xt_coarse = (/ (i, i=1, nx_coarse) /) + grid_yt_coarse = (/ (j, j=1, nx_coarse) /) + + id_x_coarse = diag_axis_init('grid_x_coarse', grid_x_coarse, & + 'index', 'x', 'x-index of cell corner points', set_name='coarse_grid', & + Domain2=coarse_domain, tile_count=tile_count) + id_y_coarse = diag_axis_init('grid_y_coarse', grid_y_coarse, & + 'index', 'y', 'y-index of cell corner points', set_name='coarse_grid', & + Domain2=coarse_domain, tile_count=tile_count) + + id_xt_coarse = diag_axis_init('grid_xt_coarse', grid_xt_coarse, & + 'index', 'x', 'x-index of cell center points', set_name='coarse_grid', & + Domain2=coarse_domain, tile_count=tile_count) + id_yt_coarse = diag_axis_init('grid_yt_coarse', grid_yt_coarse, & + 'index', 'y', 'y-index of cell center points', set_name='coarse_grid', & + Domain2=coarse_domain, tile_count=tile_count) + end subroutine initialize_coarse_diagnostic_axes + + subroutine register_coarse_diagnostics(idiag_coarse, Time, id_xt_coarse,& + id_yt_coarse, id_pfull) + type(fv_coarse_diag_type), intent(inout) :: idiag_coarse + type(time_type), intent(in) :: Time + integer, intent(in) :: id_xt_coarse, id_yt_coarse, id_pfull + + integer :: coarse_axes_t(3) + real :: missing_value = -1.0e10 ! Following fv_diagnostics.F90 + + coarse_axes_t = (/ id_xt_coarse, id_yt_coarse, id_pfull /) + idiag_coarse%id_omega_coarse = register_diag_field('dynamics', & + 'omega_coarse', coarse_axes_t(1:3), Time, & + 'coarse-grained omega', & + 'Pa/s', missing_value=missing_value) + end subroutine register_coarse_diagnostics + + subroutine fv_coarse_diag(Atm, Time) + type(fv_atmos_type), intent(in), target :: Atm(:) + type(time_type), intent(in) :: Time + + character(len=256) :: error_message + + if (trim(Atm(tile_count)%coarse_graining%strategy) .eq. MODEL_LEVEL) then + call fv_coarse_diag_model_levels(Atm, Time) + endif + end subroutine fv_coarse_diag + + subroutine fv_coarse_diag_model_levels(Atm, Time) + type(fv_atmos_type), intent(in), target :: Atm(:) + type(time_type), intent(in) :: Time + + real, allocatable :: work_3d_coarse(:,:,:) + integer :: diagnostic_ids_3d(Atm(tile_count)%coarse_graining%idiag%n_3d_diagnostics) + integer :: is, ie, js, je, is_coarse, ie_coarse, js_coarse, je_coarse, npz + logical :: used + + call get_diagnostic_ids_3d(Atm(tile_count)%coarse_graining%idiag, diagnostic_ids_3d) + call get_fine_array_bounds(Atm(tile_count)%bd, is, ie, js, je) + call get_coarse_array_bounds(Atm(tile_count)%coarse_graining%bd, is_coarse, ie_coarse, js_coarse, je_coarse) + npz = Atm(tile_count)%npz + + if (any(diagnostic_ids_3d > 0)) then + allocate(work_3d_coarse(is_coarse:ie_coarse,js_coarse:je_coarse,1:npz)) + endif + + if (Atm(tile_count)%coarse_graining%idiag%id_omega_coarse > 0) then + call weighted_block_average( & + Atm(tile_count)%gridstruct%area(is:ie,js:je), & + Atm(tile_count)%omga(is:ie,js:je,1:npz), & + work_3d_coarse) + used = send_data(Atm(tile_count)%coarse_graining%idiag%id_omega_coarse, work_3d_coarse, Time) + endif + end subroutine fv_coarse_diag_model_levels + + subroutine get_diagnostic_ids_3d(idiag_coarse, diagnostic_ids_3d) + type(fv_coarse_diag_type), intent(in) :: idiag_coarse + integer, intent(out) :: diagnostic_ids_3d(idiag_coarse%n_3d_diagnostics) + + diagnostic_ids_3d = (/ idiag_coarse%id_omega_coarse /) + end subroutine + +end module coarse_grained_diagnostics_mod diff --git a/FV3/atmos_cubed_sphere/tools/coarse_graining.F90 b/FV3/atmos_cubed_sphere/tools/coarse_graining.F90 new file mode 100644 index 000000000..73791e65f --- /dev/null +++ b/FV3/atmos_cubed_sphere/tools/coarse_graining.F90 @@ -0,0 +1,271 @@ +module coarse_graining_mod + + use fv_arrays_mod, only: fv_coarse_grid_bounds_type, fv_grid_bounds_type, fv_coarse_graining_type + use fms_mod, only: check_nml_error, close_file, open_namelist_file + use mpp_domains_mod, only: domain2d, mpp_define_io_domain, mpp_define_mosaic, mpp_get_compute_domain + use mpp_mod, only: FATAL, input_nml_file, mpp_error, mpp_npes + + implicit none + private + + public :: block_sum, get_fine_array_bounds, get_coarse_array_bounds, coarse_graining_init, weighted_block_average, MODEL_LEVEL + + interface block_sum + module procedure block_sum_2d + end interface block_sum + + interface weighted_block_average + module procedure weighted_block_average_2d + module procedure weighted_block_average_3d_field_2d_weights + end interface weighted_block_average + + ! Global variables for the module, initialized in coarse_graining_init + type(domain2d) :: coarse_domain + integer :: is, ie, js, je, npz + integer :: is_coarse, ie_coarse, js_coarse, je_coarse + integer :: nx_coarse + character(len=11) :: MODEL_LEVEL = 'model_level' + + ! Namelist parameters initialized with default values + integer :: coarsening_factor = 8 + integer :: coarse_io_layout(2) = (/1, 1/) + character(len=64) :: strategy = 'model_level' ! Valid values are 'model_level' + + namelist /coarse_graining_nml/ coarsening_factor, coarse_io_layout, strategy + +contains + + subroutine coarse_graining_init(npx, atm_npz, layout, bd, coarse_graining) + integer, intent(in) :: npx + integer, intent(in) :: atm_npz + integer, intent(in) :: layout(2) + type(fv_grid_bounds_type), intent(in) :: bd + type(fv_coarse_graining_type), intent(inout) :: coarse_graining + + character(len=256) :: error_message + logical :: exists + integer :: error_code, iostat + + read(input_nml_file, coarse_graining_nml, iostat=iostat) + error_code = check_nml_error(iostat, 'coarse_graining_nml') + coarse_graining%do_coarse_graining = .true. + + call assert_valid_strategy(strategy) + call compute_nx_coarse(npx, coarsening_factor, nx_coarse) + call assert_valid_domain_layout(nx_coarse, layout) + call define_cubic_mosaic(coarse_domain, nx_coarse, nx_coarse, layout) + call mpp_define_io_domain(coarse_domain, coarse_io_layout) + call mpp_get_compute_domain(coarse_domain, is_coarse, ie_coarse, js_coarse, je_coarse) + call get_fine_array_bounds(bd, is, ie, js, je) + npz = atm_npz + + coarse_graining%bd%is_coarse = is_coarse + coarse_graining%bd%ie_coarse = ie_coarse + coarse_graining%bd%js_coarse = js_coarse + coarse_graining%bd%je_coarse = je_coarse + coarse_graining%nx_coarse = nx_coarse + coarse_graining%domain = coarse_domain + coarse_graining%strategy = strategy + end subroutine coarse_graining_init + + subroutine compute_nx_coarse(npx, coarsening_factor, nx_coarse) + integer, intent(in) :: npx + integer, intent(in) :: coarsening_factor + integer, intent(out) :: nx_coarse + + character(len=256) :: error_message + integer :: nx + + nx = npx - 1 + if (mod(nx, coarsening_factor) > 0) then + write(error_message, *) 'coarse_graining_init: coarsening_factor does not evenly divide the native resolution.' + call mpp_error(FATAL, error_message) + endif + nx_coarse = nx / coarsening_factor + end subroutine compute_nx_coarse + + subroutine assert_valid_domain_layout(nx_coarse, layout) + integer, intent(in) :: nx_coarse + integer, intent(in) :: layout(2) + + character(len=256) :: error_message + integer :: layout_x, layout_y + layout_x = layout(1) + layout_y = layout(2) + + if (mod(nx_coarse, layout_x) > 0 .or. mod(nx_coarse, layout_y) > 0) then + write(error_message, *) 'coarse_graining_init: domain decomposition layout does not evenly divide the coarse grid.' + call mpp_error(FATAL, error_message) + endif + end subroutine assert_valid_domain_layout + + subroutine assert_valid_strategy(strategy) + character(len=64), intent(in) :: strategy + + character(len=256) :: error_message + + if (trim(strategy) .ne. MODEL_LEVEL) then + write(error_message, *) 'Invalid coarse graining strategy provided.' + call mpp_error(FATAL, error_message) + endif + end subroutine assert_valid_strategy + + subroutine get_fine_array_bounds(bd, is, ie, js, je) + type(fv_grid_bounds_type), intent(in) :: bd + integer, intent(out) :: is, ie, js, je + + is = bd%is + ie = bd%ie + js = bd%js + je = bd%je + end subroutine get_fine_array_bounds + + subroutine get_coarse_array_bounds(coarse_bd, is_coarse, ie_coarse, js_coarse, je_coarse) + type(fv_coarse_grid_bounds_type), intent(in) :: coarse_bd + integer, intent(out) :: is_coarse, ie_coarse, js_coarse, je_coarse + + is_coarse = coarse_bd%is_coarse + ie_coarse = coarse_bd%ie_coarse + js_coarse = coarse_bd%js_coarse + je_coarse = coarse_bd%je_coarse + end subroutine get_coarse_array_bounds + + subroutine block_sum_2d(fine, coarse) + real, intent(in) :: fine(is:ie,js:je) + real, intent(out) :: coarse(is_coarse:ie_coarse,js_coarse:je_coarse) + + integer :: i, j, i_coarse, j_coarse, offset + + offset = coarsening_factor - 1 + do i = is, ie, coarsening_factor + i_coarse = (i - 1) / coarsening_factor + 1 + do j = js, je, coarsening_factor + j_coarse = (j - 1) / coarsening_factor + 1 + coarse(i_coarse,j_coarse) = sum(fine(i:i+offset,j:j+offset)) + enddo + enddo + end subroutine + + subroutine weighted_block_average_2d(weights, fine, coarse) + real, intent(in) :: weights(is:ie,js:je), fine(is:ie,js:je) + real, intent(out) :: coarse(is_coarse:ie_coarse,js_coarse:je_coarse) + + real, allocatable :: weighted_fine(:,:), weighted_block_sum(:,:), block_sum_weights(:,:) + + allocate(weighted_fine(is:ie,js:je)) + allocate(weighted_block_sum(is_coarse:ie_coarse,js_coarse:je_coarse)) + allocate(block_sum_weights(is_coarse:ie_coarse,js_coarse:je_coarse)) + + weighted_fine = weights * fine + call block_sum_2d(weighted_fine, weighted_block_sum) + call block_sum_2d(weights, block_sum_weights) + coarse = weighted_block_sum / block_sum_weights + end subroutine weighted_block_average_2d + + subroutine weighted_block_average_3d_field_2d_weights(weights, fine, coarse) + real, intent(in) :: weights(is:ie,js:je), fine(is:ie,js:je,1:npz) + real, intent(out) :: coarse(is_coarse:ie_coarse,js_coarse:je_coarse,1:npz) + + integer :: k + + do k = 1, npz + call weighted_block_average_2d(weights, fine(is:ie,js:je,k), coarse(is_coarse:ie_coarse,js_coarse:je_coarse,k)) + enddo + end subroutine weighted_block_average_3d_field_2d_weights + + ! This subroutine is copied from FMS/test_fms/horiz_interp/test2_horiz_interp.F90. + ! domain_decomp in fv_mp_mod.F90 does something similar, but it does a + ! few other unnecessary things (and requires more arguments). + subroutine define_cubic_mosaic(domain, ni, nj, layout) + type(domain2d), intent(inout) :: domain + integer, intent(in) :: layout(2) + integer, intent(in) :: ni, nj + integer :: pe_start(6), pe_end(6) + integer :: global_indices(4,6), layout2d(2,6) + integer, dimension(12) :: istart1, iend1, jstart1, jend1, tile1 + integer, dimension(12) :: istart2, iend2, jstart2, jend2, tile2 + integer :: ntiles, num_contact + integer :: p, npes_per_tile, i + + ntiles = 6 + num_contact = 12 + p = 0 + npes_per_tile = mpp_npes()/ntiles + do i = 1, 6 + layout2d(:,i) = layout(:) + global_indices(1,i) = 1 + global_indices(2,i) = ni + global_indices(3,i) = 1 + global_indices(4,i) = nj + pe_start(i) = p + p = p + npes_per_tile + pe_end(i) = p-1 + enddo + + !--- Contact line 1, between tile 1 (EAST) and tile 2 (WEST) + tile1(1) = 1; tile2(1) = 2 + istart1(1) = ni; iend1(1) = ni; jstart1(1) = 1; jend1(1) = nj + istart2(1) = 1; iend2(1) = 1; jstart2(1) = 1; jend2(1) = nj + + !--- Contact line 2, between tile 1 (NORTH) and tile 3 (WEST) + tile1(2) = 1; tile2(2) = 3 + istart1(2) = 1; iend1(2) = ni; jstart1(2) = nj; jend1(2) = nj + istart2(2) = 1; iend2(2) = 1; jstart2(2) = nj; jend2(2) = 1 + + !--- Contact line 3, between tile 1 (WEST) and tile 5 (NORTH) + tile1(3) = 1; tile2(3) = 5 + istart1(3) = 1; iend1(3) = 1; jstart1(3) = 1; jend1(3) = nj + istart2(3) = ni; iend2(3) = 1; jstart2(3) = nj; jend2(3) = nj + + !--- Contact line 4, between tile 1 (SOUTH) and tile 6 (NORTH) + tile1(4) = 1; tile2(4) = 6 + istart1(4) = 1; iend1(4) = ni; jstart1(4) = 1; jend1(4) = 1 + istart2(4) = 1; iend2(4) = ni; jstart2(4) = nj; jend2(4) = nj + + !--- Contact line 5, between tile 2 (NORTH) and tile 3 (SOUTH) + tile1(5) = 2; tile2(5) = 3 + istart1(5) = 1; iend1(5) = ni; jstart1(5) = nj; jend1(5) = nj + istart2(5) = 1; iend2(5) = ni; jstart2(5) = 1; jend2(5) = 1 + + !--- Contact line 6, between tile 2 (EAST) and tile 4 (SOUTH) + tile1(6) = 2; tile2(6) = 4 + istart1(6) = ni; iend1(6) = ni; jstart1(6) = 1; jend1(6) = nj + istart2(6) = ni; iend2(6) = 1; jstart2(6) = 1; jend2(6) = 1 + + !--- Contact line 7, between tile 2 (SOUTH) and tile 6 (EAST) + tile1(7) = 2; tile2(7) = 6 + istart1(7) = 1; iend1(7) = ni; jstart1(7) = 1; jend1(7) = 1 + istart2(7) = ni; iend2(7) = ni; jstart2(7) = nj; jend2(7) = 1 + + !--- Contact line 8, between tile 3 (EAST) and tile 4 (WEST) + tile1(8) = 3; tile2(8) = 4 + istart1(8) = ni; iend1(8) = ni; jstart1(8) = 1; jend1(8) = nj + istart2(8) = 1; iend2(8) = 1; jstart2(8) = 1; jend2(8) = nj + + !--- Contact line 9, between tile 3 (NORTH) and tile 5 (WEST) + tile1(9) = 3; tile2(9) = 5 + istart1(9) = 1; iend1(9) = ni; jstart1(9) = nj; jend1(9) = nj + istart2(9) = 1; iend2(9) = 1; jstart2(9) = nj; jend2(9) = 1 + + !--- Contact line 10, between tile 4 (NORTH) and tile 5 (SOUTH) + tile1(10) = 4; tile2(10) = 5 + istart1(10) = 1; iend1(10) = ni; jstart1(10) = nj; jend1(10) = nj + istart2(10) = 1; iend2(10) = ni; jstart2(10) = 1; jend2(10) = 1 + + !--- Contact line 11, between tile 4 (EAST) and tile 6 (SOUTH) + tile1(11) = 4; tile2(11) = 6 + istart1(11) = ni; iend1(11) = ni; jstart1(11) = 1; jend1(11) = nj + istart2(11) = ni; iend2(11) = 1; jstart2(11) = 1; jend2(11) = 1 + + !--- Contact line 12, between tile 5 (EAST) and tile 6 (WEST) + tile1(12) = 5; tile2(12) = 6 + istart1(12) = ni; iend1(12) = ni; jstart1(12) = 1; jend1(12) = nj + istart2(12) = 1; iend2(12) = 1; jstart2(12) = 1; jend2(12) = nj + + call mpp_define_mosaic(global_indices, layout2d, domain, ntiles, & + num_contact, tile1, tile2, istart1, iend1, jstart1, jend1, & + istart2, iend2, jstart2, jend2, pe_start, pe_end, & + symmetry=.true., name='coarse cubic mosaic') + end subroutine define_cubic_mosaic + +end module coarse_graining_mod