From 3209acc6fb6a022b3e252ae87df1e0df04a6971d Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Thu, 20 Jul 2017 14:09:33 -0400 Subject: [PATCH 1/4] +Use new coupler_type routines in the SIS2 code Replaced extensive code blocks going through the detailed elements of the coupler_types with equivalent calls to the new subroutines that have been introduced to more succinctly perform tasks on the coupler types. All answers are bitwise identical, but this requires the new interfaces that have been introduced on the coupler_type_reform_rwh branch of the FMS code. --- src/combined_ice_ocean_driver.F90 | 22 +++++++--------------- src/ice_boundary_types.F90 | 11 ++--------- src/ice_type.F90 | 9 ++------- 3 files changed, 11 insertions(+), 31 deletions(-) diff --git a/src/combined_ice_ocean_driver.F90 b/src/combined_ice_ocean_driver.F90 index 512255e9..c9ff864f 100644 --- a/src/combined_ice_ocean_driver.F90 +++ b/src/combined_ice_ocean_driver.F90 @@ -38,13 +38,15 @@ module combined_ice_ocean_driver use ocean_model_mod, only : update_ocean_model, ocean_model_end! , ocean_model_init use ocean_model_mod, only : ocean_public_type, ocean_state_type, ice_ocean_boundary_type +use coupler_types_mod, only: coupler_type_send_data, coupler_type_data_override +use coupler_types_mod, only: coupler_type_copy_data, coupler_type_redistribute_data use data_override_mod, only : data_override use diag_manager_mod, only : send_data use mpp_domains_mod, only : domain2D, mpp_get_layout, mpp_get_compute_domain implicit none ; private -public update_slow_ice_and_ocean, ice_ocean_driver_init, ice_ocean_driver_end +public :: update_slow_ice_and_ocean, ice_ocean_driver_init, ice_ocean_driver_end type, public :: ice_ocean_driver_type ; private logical :: CS_is_initialized = .false. @@ -247,12 +249,7 @@ subroutine direct_flux_ice_to_IOB( Time, Ice, IOB ) if (ASSOCIATED(IOB%q_flux)) IOB%q_flux(:,:) = Ice%flux_q(:,:) ! Extra fluxes - do n=1,IOB%fluxes%num_bcs ; do m=1,IOB%fluxes%bc(n)%num_fields - if ( associated(IOB%fluxes%bc(n)%field(m)%values) ) then - IOB%fluxes%bc(n)%field(m)%values(:,:) = Ice%ocean_fluxes%bc(n)%field(m)%values(:,:) - endif - enddo ; enddo - + call coupler_type_copy_data(Ice%ocean_fluxes, IOB%fluxes) ! These lines allow the data override code to reset the fluxes to the ocean. call data_override('OCN', 'u_flux', IOB%u_flux , Time ) @@ -282,14 +279,9 @@ subroutine direct_flux_ice_to_IOB( Time, Ice, IOB ) call data_override('OCN', 'mass_berg', IOB%mass_berg , Time) ! Override and output extra fluxes of tracers or gasses - do n=1,IOB%fluxes%num_bcs ; do m=1,IOB%fluxes%bc(n)%num_fields - call data_override('OCN', IOB%fluxes%bc(n)%field(m)%name, & - IOB%fluxes%bc(n)%field(m)%values, Time) - - ! Perform diagnostic output for the extra fluxes - used = send_data(IOB%fluxes%bc(n)%field(m)%id_diag, & - IOB%fluxes%bc(n)%field(m)%values, Time) - enddo ; enddo + call coupler_type_data_override('OCN', IOB%fluxes, Time ) + + call coupler_type_send_data(IOB%fluxes, Time ) end subroutine direct_flux_ice_to_IOB diff --git a/src/ice_boundary_types.F90 b/src/ice_boundary_types.F90 index ae754177..8d3bcdfe 100644 --- a/src/ice_boundary_types.F90 +++ b/src/ice_boundary_types.F90 @@ -6,7 +6,7 @@ !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! module ice_boundary_types -use coupler_types_mod, only : coupler_2d_bc_type, coupler_3d_bc_type +use coupler_types_mod, only : coupler_2d_bc_type, coupler_3d_bc_type, coupler_type_write_chksums use fms_mod, only : stdout use mpp_mod, only : mpp_chksum use mpp_parameter_mod, only : CGRID_NE, BGRID_NE, AGRID @@ -125,14 +125,7 @@ subroutine ocn_ice_bnd_type_chksum(id, timestep, bnd_type) ! write(outunit,100) 'ocn_ice_bnd_type%data ',mpp_chksum(bnd_type%data ) 100 FORMAT("CHECKSUM::",A32," = ",Z20) - do n = 1, bnd_type%fields%num_bcs !{ - do m = 1, bnd_type%fields%bc(n)%num_fields !{ - write(outunit,101) 'oibt%',trim(bnd_type%fields%bc(n)%name), & - trim(bnd_type%fields%bc(n)%field(m)%name), & - mpp_chksum(bnd_type%fields%bc(n)%field(m)%values) - enddo !} m - enddo !} n - 101 FORMAT("CHECKSUM::",A16,a,'%',a," = ",Z20) + call coupler_type_write_chksums(bnd_type%fields, outunit, 'oibt%') end subroutine ocn_ice_bnd_type_chksum diff --git a/src/ice_type.F90 b/src/ice_type.F90 index efd61507..0ce7e2b3 100644 --- a/src/ice_type.F90 +++ b/src/ice_type.F90 @@ -11,7 +11,7 @@ module ice_type_mod use fms_io_mod, only: save_restart, restore_state, query_initialized use fms_io_mod, only: register_restart_field, restart_file_type use time_manager_mod, only: time_type, time_type_to_real -use coupler_types_mod,only: coupler_2d_bc_type, coupler_3d_bc_type +use coupler_types_mod,only: coupler_2d_bc_type, coupler_3d_bc_type, coupler_type_write_chksums use SIS_hor_grid, only : SIS_hor_grid_type use ice_grid, only : ice_grid_type @@ -598,11 +598,7 @@ subroutine ice_data_type_chksum(id, timestep, Ice) write(outunit,100) 'ice_data_type%u_surf ',mpp_chksum(Ice%u_surf ) write(outunit,100) 'ice_data_type%v_surf ',mpp_chksum(Ice%v_surf ) - do n=1,Ice%ocean_fields%num_bcs ; do m=1,Ice%ocean_fields%bc(n)%num_fields - write(outunit,101) 'ice%', trim(Ice%ocean_fields%bc(n)%name), & - trim(Ice%ocean_fields%bc(n)%field(m)%name), & - mpp_chksum(Ice%ocean_fields%bc(n)%field(m)%values) - enddo ; enddo + call coupler_type_write_chksums(Ice%ocean_fields, outunit, 'ice%') endif if (Ice%slow_ice_PE) then @@ -633,7 +629,6 @@ subroutine ice_data_type_chksum(id, timestep, Ice) endif 100 FORMAT(" CHECKSUM::",A32," = ",Z20) -101 FORMAT(" CHECKSUM::",A16,a,'%',a," = ",Z20) end subroutine ice_data_type_chksum From e2497e6207fd10ae75de16083e201526feb9a322 Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Thu, 27 Jul 2017 16:09:15 -0400 Subject: [PATCH 2/4] +Complete reformulation of SIS2 tracer fluxes Completely revised how SIS2 handles the additional tracer-related fluxes and ocean surface fields, so that they now use coupler_type variables throughout the algorithm and quasi-object-oriented routines to operate on these types. All answers are bitwise identical, and SIS2 no longer has any code that depends on the details of the implementation of the coupler_types. Some elements of publicly visible type have been changed or eliminated. These changes require an updated version of the coupler_type_reform_rwh branch of the FMS code. --- src/SIS_fast_thermo.F90 | 107 +++--------- src/SIS_slow_thermo.F90 | 33 +--- src/SIS_types.F90 | 262 ++++++------------------------ src/combined_ice_ocean_driver.F90 | 2 +- src/ice_model.F90 | 114 +++++++------ 5 files changed, 146 insertions(+), 372 deletions(-) diff --git a/src/SIS_fast_thermo.F90 b/src/SIS_fast_thermo.F90 index 724846ad..b31d4be1 100644 --- a/src/SIS_fast_thermo.F90 +++ b/src/SIS_fast_thermo.F90 @@ -47,7 +47,8 @@ module SIS_fast_thermo use MOM_time_manager, only : set_date, set_time, operator(+), operator(-) use MOM_time_manager, only : operator(>), operator(*), operator(/), operator(/=) -use coupler_types_mod, only : coupler_3d_bc_type +use coupler_types_mod, only : coupler_3d_bc_type, coupler_type_spawn +use coupler_types_mod, only : coupler_type_increment_data, coupler_type_rescale_data use SIS_optics, only : ice_optics_SIS2, bright_ice_temp, SIS_optics_CS use SIS_optics, only : VIS_DIR, VIS_DIF, NIR_DIR, NIR_DIF use SIS_types, only : ice_state_type, IST_chksum, IST_bounds_check @@ -151,30 +152,17 @@ subroutine sum_top_quantities (FIA, ABT, flux_u, flux_v, flux_sh, evap, & real :: t_sfc integer :: i, j, k, m, n, b, nb, i2, j2, k2, isc, iec, jsc, jec, i_off, j_off, ncat - integer :: ind + integer :: isd, ied, jsd, jed isc = G%isc ; iec = G%iec ; jsc = G%jsc ; jec = G%jec ; ncat = IG%CatIce + isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed nb = size(FIA%flux_sw_top,4) i_off = LBOUND(ABT%t_flux,1) - G%isc j_off = LBOUND(ABT%t_flux,2) - G%jsc - if (FIA%num_tr_fluxes < 0) then - ! Determine how many atmospheric boundary fluxes have been passed in, and - ! set up both an indexing array for these and a space to take their average. - ! This code is only exercised the first time that sum_top_quantities is called. - FIA%num_tr_fluxes = 0 - if (ABT%fluxes%num_bcs > 0) then - do n=1,ABT%fluxes%num_bcs - FIA%num_tr_fluxes = FIA%num_tr_fluxes + ABT%fluxes%bc(n)%num_fields - enddo - - if (FIA%num_tr_fluxes > 0) then - allocate(FIA%tr_flux_top(G%isd:G%ied, G%jsd:G%jed, 0:IG%CatIce, FIA%num_tr_fluxes)) - FIA%tr_flux_top(:,:,:,:) = 0.0 - endif - endif - endif + call coupler_type_spawn(ABT%fluxes, FIA%tr_flux, (/isd, isc, iec, ied/), & + (/jsd, jsc, jec, jed/), (/0, ncat/), as_needed=.true.) if (FIA%avg_count == 0) then ! zero_top_quantities - zero fluxes to begin summing in ice fast physics. @@ -191,7 +179,7 @@ subroutine sum_top_quantities (FIA, ABT, flux_u, flux_v, flux_sh, evap, & FIA%flux_lw0(:,:,:) = 0.0 ; FIA%Tskin_cat(:,:,:) = 0.0 endif - if (FIA%num_tr_fluxes > 0) FIA%tr_flux_top(:,:,:,:) = 0.0 + call coupler_type_rescale_data(FIA%tr_flux, 0.0) endif !$OMP parallel do default(shared) @@ -209,18 +197,8 @@ subroutine sum_top_quantities (FIA, ABT, flux_u, flux_v, flux_sh, evap, & ! FIA%flux_sw_dn is accumulated where the fast radiation diagnostics are output ! because it depends on arrays that are stored in the public ice_data_type. - ind = 0 - do n=1,ABT%fluxes%num_bcs ; do m=1,ABT%fluxes%bc(n)%num_fields - ind = ind + 1 - !Do not handle air_sea_deposition fluxes here, they need to be handled after atmos_down - if(ABT%fluxes%bc(n)%flux_type /= 'air_sea_deposition') then - do k=0,ncat ; do j=jsc,jec ; do i=isc,iec - i2 = i+i_off ; j2 = j+j_off ; k2 = k+1 - FIA%tr_flux_top(i,j,k,ind) = FIA%tr_flux_top(i,j,k,ind) + & - ABT%fluxes%bc(n)%field(m)%values(i2,j2,k2) - enddo ; enddo ; enddo - endif - enddo ; enddo + !Do not handle air_sea_deposition fluxes here, they need to be handled after atmos_down + call coupler_type_increment_data( ABT%fluxes, FIA%tr_flux, exclude_flux_type='air_sea_deposition') if (allocated(FIA%flux_sh0)) then !$OMP parallel do default(shared) private(t_sfc) @@ -293,10 +271,8 @@ subroutine avg_top_quantities(FIA, Rad, IST, G, IG) ! FIA%fprec_top(i,j,k) = FIA%fprec_top(i,j,k) - FIA%evap_top(i,j,k) ! FIA%evap_top(i,j,k) = 0.0 ! endif - do n=1,FIA%num_tr_fluxes - FIA%tr_flux_top(i,j,k,n) = FIA%tr_flux_top(i,j,k,n) * I_avc - enddo enddo ; enddo + do b=1,nb ; do i=isc,iec FIA%flux_sw_dn(i,j,b) = FIA%flux_sw_dn(i,j,b)*I_avc enddo ; enddo @@ -304,6 +280,7 @@ subroutine avg_top_quantities(FIA, Rad, IST, G, IG) FIA%Tskin_avg(i,j) = FIA%Tskin_avg(i,j) * I_avc enddo enddo + call coupler_type_rescale_data(FIA%tr_flux, I_avc) ! Determine the fractional ice coverage and the wind stresses averaged ! across all the ice thickness categories on an A-grid. @@ -389,20 +366,14 @@ subroutine total_top_quantities(FIA, TSF, part_size, G, IG) isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed nb = size(FIA%flux_sw_top,4) - if (TSF%num_tr_fluxes < 0) then - ! Allocate the arrays to hold the tracer fluxes. This code is only exercised - ! the first time that total_top_quantities is called. - TSF%num_tr_fluxes = FIA%num_tr_fluxes - if (TSF%num_tr_fluxes > 0) then - allocate(TSF%tr_flux(G%isd:G%ied, G%jsd:G%jed, TSF%num_tr_fluxes)) - endif - endif + call coupler_type_spawn(FIA%tr_flux, TSF%tr_flux, (/isd, isc, iec, ied/), & + (/jsd, jsc, jec, jed/), as_needed=.true. ) TSF%flux_u(:,:) = 0.0 ; TSF%flux_v(:,:) = 0.0 TSF%flux_sh(:,:) = 0.0 ; TSF%flux_lw(:,:) = 0.0 ; TSF%flux_lh(:,:) = 0.0 TSF%flux_sw(:,:,:) = 0.0 TSF%evap(:,:) = 0.0 ; TSF%fprec(:,:) = 0.0 ; TSF%lprec(:,:) = 0.0 - if (TSF%num_tr_fluxes > 0) TSF%tr_flux(:,:,:) = 0.0 + call coupler_type_rescale_data(TSF%tr_flux, 0.0) do k=0,ncat ; do j=jsc,jec ; do i=isc,iec TSF%flux_u(i,j) = TSF%flux_u(i,j) + part_size(i,j,k) * FIA%flux_u_top(i,j,k) @@ -417,11 +388,8 @@ subroutine total_top_quantities(FIA, TSF, part_size, G, IG) TSF%evap(i,j) = TSF%evap(i,j) + part_size(i,j,k) * FIA%evap_top(i,j,k) TSF%fprec(i,j) = TSF%fprec(i,j) + part_size(i,j,k) * FIA%fprec_top(i,j,k) TSF%lprec(i,j) = TSF%lprec(i,j) + part_size(i,j,k) * FIA%lprec_top(i,j,k) - - do n=1,TSF%num_tr_fluxes - TSF%tr_flux(i,j,n) = TSF%tr_flux(i,j,n) + part_size(i,j,k) * FIA%tr_flux_top(i,j,k,n) - enddo enddo ; enddo ; enddo + call coupler_type_increment_data(FIA%tr_flux, part_size, TSF%tr_flux) ! If the sum of part_size across all the ice and ocean categories is not ! exactly 1, rescaling might be advisable, but for now it is assumed that @@ -447,15 +415,11 @@ subroutine find_excess_fluxes(FIA, TSF, XSF, part_size, G, IG) isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed nb = size(FIA%flux_sw_top,4) - if ((FIA%num_tr_fluxes > 0) .and. (XSF%num_tr_fluxes < 0)) then - ! This is the first call when the number of tracer fluxes are known, and - ! the XSF tracer flux arrays need to be allocated now. This code is only - ! exercised the first or second time that total_top_quantities is called. - XSF%num_tr_fluxes = FIA%num_tr_fluxes - if (XSF%num_tr_fluxes > 0) then - allocate(XSF%tr_flux(G%isd:G%ied, G%jsd:G%jed, XSF%num_tr_fluxes)) - endif - endif + ! Check whether this is the first call when the number of tracer fluxes are + ! known, and the XSF tracer flux arrays need to be allocated now. This call + ! only does anything the first or second time that total_top_quantities is called. + call coupler_type_spawn(FIA%tr_flux, XSF%tr_flux, (/isd, isc, iec, ied/), & + (/jsd, jsc, jec, jed/), as_needed=.true. ) ! Note that XSF%flux_u and XSF%flux_v are not necessary as the changing ice cover is ! already being taken into account. @@ -465,7 +429,7 @@ subroutine find_excess_fluxes(FIA, TSF, XSF, part_size, G, IG) XSF%flux_sh(:,:) = 0.0 ; XSF%flux_lw(:,:) = 0.0 ; XSF%flux_lh(:,:) = 0.0 XSF%flux_sw(:,:,:) = 0.0 XSF%evap(:,:) = 0.0 ; XSF%fprec(:,:) = 0.0 ; XSF%lprec(:,:) = 0.0 - if (XSF%num_tr_fluxes > 0) XSF%tr_flux(:,:,:) = 0.0 + call coupler_type_rescale_data(XSF%tr_flux, 0.0) do k=0,ncat ; do j=jsc,jec ; do i=isc,iec XSF%flux_sh(i,j) = XSF%flux_sh(i,j) + part_size(i,j,k) * FIA%flux_sh_top(i,j,k) @@ -478,11 +442,8 @@ subroutine find_excess_fluxes(FIA, TSF, XSF, part_size, G, IG) XSF%evap(i,j) = XSF%evap(i,j) + part_size(i,j,k) * FIA%evap_top(i,j,k) XSF%fprec(i,j) = XSF%fprec(i,j) + part_size(i,j,k) * FIA%fprec_top(i,j,k) XSF%lprec(i,j) = XSF%lprec(i,j) + part_size(i,j,k) * FIA%lprec_top(i,j,k) - - do n=1,XSF%num_tr_fluxes - XSF%tr_flux(i,j,n) = XSF%tr_flux(i,j,n) + part_size(i,j,k) * FIA%tr_flux_top(i,j,k,n) - enddo enddo ; enddo ; enddo + call coupler_type_increment_data(FIA%tr_flux, part_size, XSF%tr_flux) do j=jsc,jec ; do i=isc,iec XSF%flux_sh(i,j) = XSF%flux_sh(i,j) - TSF%flux_sh(i,j) @@ -492,11 +453,8 @@ subroutine find_excess_fluxes(FIA, TSF, XSF, part_size, G, IG) XSF%evap(i,j) = XSF%evap(i,j) - TSF%evap(i,j) XSF%fprec(i,j) = XSF%fprec(i,j) - TSF%fprec(i,j) XSF%lprec(i,j) = XSF%lprec(i,j) - TSF%lprec(i,j) - - do n=1,XSF%num_tr_fluxes - XSF%tr_flux(i,j,n) = XSF%tr_flux(i,j,n) - TSF%tr_flux(i,j,n) - enddo enddo ; enddo + call coupler_type_increment_data(TSF%tr_flux, XSF%tr_flux, scale_factor=-1.0) end subroutine find_excess_fluxes @@ -872,24 +830,7 @@ subroutine accumulate_deposition_fluxes(ABT, FIA, G, IG) type(SIS_hor_grid_type), intent(in) :: G type(ice_grid_type), intent(in) :: IG - integer :: i, j, k, m, n, i_off, j_off, i2, j2, k2, isc, iec, jsc, jec, ncat, ind - - isc = G%isc ; iec = G%iec ; jsc = G%jsc ; jec = G%jec ; ncat = IG%CatIce - - i_off = LBOUND(ABT%t_flux,1) - G%isc - j_off = LBOUND(ABT%t_flux,2) - G%jsc - - ind = 0 - do n=1,ABT%fluxes%num_bcs ; do m=1,ABT%fluxes%bc(n)%num_fields - ind = ind + 1 - if (ABT%fluxes%bc(n)%flux_type == 'air_sea_deposition') then - do k=0,ncat ; do j=jsc,jec ; do i=isc,iec - i2 = i+i_off ; j2 = j+j_off ; k2 = k+1 - FIA%tr_flux_top(i,j,k,ind) = FIA%tr_flux_top(i,j,k,ind) + & - ABT%fluxes%bc(n)%field(m)%values(i2,j2,k2) - enddo ; enddo ; enddo - endif - enddo; enddo + call coupler_type_increment_data( ABT%fluxes, FIA%tr_flux, only_flux_type='air_sea_deposition') end subroutine accumulate_deposition_fluxes diff --git a/src/SIS_slow_thermo.F90 b/src/SIS_slow_thermo.F90 index 8e75dee2..f2de2b42 100644 --- a/src/SIS_slow_thermo.F90 +++ b/src/SIS_slow_thermo.F90 @@ -48,6 +48,8 @@ module SIS_slow_thermo use MOM_hor_index, only : hor_index_type use MOM_EOS, only : EOS_type, calculate_density_derivs +use coupler_types_mod, only : coupler_type_spawn +use coupler_types_mod, only : coupler_type_increment_data, coupler_type_rescale_data use fms_mod, only : clock_flag_default use mpp_mod, only : mpp_clock_id, mpp_clock_begin, mpp_clock_end use mpp_mod, only : CLOCK_COMPONENT, CLOCK_LOOP, CLOCK_ROUTINE @@ -396,29 +398,10 @@ subroutine slow_thermodynamics(IST, dt_slow, CS, OSS, FIA, XSF, IOF, G, IG) ! IOF must be updated regardless of whether the ice is specified or the ! prognostic model is being used. - if (FIA%num_tr_fluxes>0) then - ! Only one OMP thread executes the following block because IOF is shared. - !$OMP SINGLE - if (IOF%num_tr_fluxes < 0) then - ! This is the first call when the number of tracer fluxes are known, and - ! the IOF tracer flux arrays need to be allocated now. - IOF%num_tr_fluxes = FIA%num_tr_fluxes - - allocate(IOF%tr_flux_ocn_top(SZI_(G), SZJ_(G), IOF%num_tr_fluxes)) - IOF%tr_flux_ocn_top(:,:,:) = 0.0 - endif - !$OMP END SINGLE - !$OMP parallel do default(shared) - do m=1,FIA%num_tr_fluxes - do j=jsc,jec ; do i=isc,iec - IOF%tr_flux_ocn_top(i,j,m) = IST%part_size(i,j,0) * FIA%tr_flux_top(i,j,0,m) - enddo ; enddo - do k=1,ncat ; do j=jsc,jec ; do i=isc,iec - IOF%tr_flux_ocn_top(i,j,m) = IOF%tr_flux_ocn_top(i,j,m) + & - IST%part_size(i,j,k) * FIA%tr_flux_top(i,j,k,m) - enddo ; enddo ; enddo - enddo - endif + call coupler_type_spawn(FIA%tr_flux, IOF%tr_flux_ocn_top, (/isd, isc, iec, ied/), & + (/jsd, jsc, jec, jed/), as_needed=.true. ) + call coupler_type_rescale_data(IOF%tr_flux_ocn_top, 0.0) + call coupler_type_increment_data(FIA%tr_flux, IST%part_size, IOF%tr_flux_ocn_top) ! No other thermodynamics need to be done for ice that is specified, if (CS%specified_ice) then @@ -522,9 +505,7 @@ subroutine add_excess_fluxes(IOF, XSF, G) IOF%lprec_ocn_top(i,j) = IOF%lprec_ocn_top(i,j) - XSF%lprec(i,j) IOF%fprec_ocn_top(i,j) = IOF%fprec_ocn_top(i,j) - XSF%fprec(i,j) - do n=1,XSF%num_tr_fluxes - IOF%tr_flux_ocn_top(i,j,n) = IOF%tr_flux_ocn_top(i,j,n) - XSF%tr_flux(i,j,n) - enddo + call coupler_type_increment_data(XSF%tr_flux, IOF%tr_flux_ocn_top, scale_factor=-1.0) ! The shortwave fluxes are more complicated because there are multiple bands ! and none of these should have negative fluxes if it can be avoided. diff --git a/src/SIS_types.F90 b/src/SIS_types.F90 index 45ffaa6b..26ac2d63 100644 --- a/src/SIS_types.F90 +++ b/src/SIS_types.F90 @@ -15,7 +15,8 @@ module SIS_types use fms_io_mod, only : restore_state, query_initialized use time_manager_mod, only : time_type, time_type_to_real use coupler_types_mod, only : coupler_2d_bc_type, coupler_3d_bc_type - +use coupler_types_mod, only : coupler_type_spawn, coupler_type_copy_data +use coupler_types_mod, only : coupler_type_redistribute_data use SIS_hor_grid, only : SIS_hor_grid_type use ice_grid, only : ice_grid_type @@ -131,10 +132,9 @@ module SIS_types ! (partially) converted back to its equivalent by the ! ocean. - real, allocatable, dimension(:,:,:) :: & - tr_array ! An array of fields related to properties for additional tracers. + type (coupler_2d_bc_type) :: & + tr_fields ! A structure of fields related to properties for additional tracers. - integer :: num_tr = -1 ! The number of additional tracer-related arrays. ! type(coupler_3d_bc_type) :: ocean_fields ! array of fields used for additional tracers real :: kmelt ! A constant that is used in the calculation of the ocean/ice @@ -166,11 +166,8 @@ module SIS_types bheat ! The upward diffusive heat flux from the ocean ! to the ice at the base of the ice, in W m-2. - real, allocatable, dimension(:,:,:) :: & - tr_array ! An array of fields related to properties for additional tracers. - integer :: num_tr = -1 ! The number of additional tracer-related arrays. - logical :: first_copy = .true. - + type (coupler_2d_bc_type) :: & + tr_fields ! A structure of fields related to properties for additional tracers. end type simple_OSS_type @@ -257,10 +254,9 @@ module SIS_types integer :: copy_calls = 0 ! The number of times this structure has been ! copied from the fast ice to the slow ice. - integer :: num_tr_fluxes = -1 ! The number of tracer flux fields - real, allocatable, dimension(:,:,:,:) :: & - tr_flux_top ! An array of tracer fluxes at the top of the - ! sea ice. + type (coupler_3d_bc_type) :: & + tr_flux ! A structure of additional tracer fluxes at the top + ! of the sea-ice ! These are the arrays that are averaged over the fast thermodynamics and ! then interpolated into unoccupied categories for the purpose of redoing @@ -329,10 +325,9 @@ module SIS_types ! from this module helping to distinguish them. integer :: copy_calls = 0 ! The number of times this structure has been ! copied from the fast ice to the slow ice. - integer :: num_tr_fluxes = -1 ! The number of tracer flux fields - real, allocatable, dimension(:,:,:) :: & - tr_flux ! An array of tracer fluxes at the top of the - ! sea ice. + type (coupler_2d_bc_type) :: & + tr_flux ! A structure of additional tracer fluxes at the top + ! of the sea-ice end type total_sfc_flux_type @@ -455,11 +450,9 @@ module SIS_types ! to determine its value. logical :: slp2ocean ! If true, apply sea level pressure to ocean surface. -! type(coupler_2d_bc_type) :: ocean_fluxes ! array of fluxes used for additional tracers - - integer :: num_tr_fluxes = -1 ! The number of tracer flux fields - real, allocatable, dimension(:,:,:) :: & - tr_flux_ocn_top ! An array of tracer fluxes at the ocean's surface. + type (coupler_2d_bc_type) :: & + tr_flux_ocn_top ! A structure of additional tracer fluxes at the top + ! of the ocean ! diagnostic IDs for ice-to-ocean fluxes. integer :: id_saltf=-1 @@ -1127,8 +1120,10 @@ subroutine translate_OSS_to_sOSS(OSS, IST, sOSS, G) type(SIS_hor_grid_type), intent(in) :: G integer :: i, j, k, m, n, i2, j2, k2, isc, iec, jsc, jec, i_off, j_off + integer :: isd, ied, jsd, jed isc = G%isc ; iec = G%iec ; jsc = G%jsc ; jec = G%jec + isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed !$OMP parallel do default(none) shared(isc,iec,jsc,jec,G,sOSS,OSS,IST) do j=jsc,jec ; do i=isc,iec @@ -1164,15 +1159,10 @@ subroutine translate_OSS_to_sOSS(OSS, IST, sOSS, G) endif enddo ; enddo - if (sOSS%num_tr<0) then - sOSS%num_tr = OSS%num_tr - if (sOSS%num_tr > 0) then - allocate(sOSS%tr_array(G%isd:G%ied,G%jsd:G%jed,sOSS%num_tr)) ; sOSS%tr_array(:,:,:) = 0.0 - endif - endif - do m=1,OSS%num_tr ; do j=jsc,jec ; do i=isc,iec - sOSS%tr_array(i,j,m) = OSS%tr_array(i,j,m) - enddo ; enddo ; enddo + call coupler_type_spawn(OSS%tr_fields, sOSS%tr_fields, (/isd, isc, iec, ied/), & + (/jsd, jsc, jec, jed/), as_needed=.true. ) + + call coupler_type_copy_data(OSS%tr_fields, sOSS%tr_fields) end subroutine translate_OSS_to_sOSS @@ -1186,7 +1176,7 @@ subroutine copy_sOSS_to_sOSS(OSS_in, OSS_out, HI_in, HI_out) type(simple_OSS_type), intent(inout) :: OSS_out type(hor_index_type), intent(in) :: HI_in, HI_out - integer :: i, j, m, isc, iec, jsc, jec + integer :: i, j, isc, iec, jsc, jec integer :: i2, j2, i_off, j_off isc = HI_in%isc ; iec = HI_in%iec ; jsc = HI_in%jsc ; jec = HI_in%jec @@ -1209,18 +1199,7 @@ subroutine copy_sOSS_to_sOSS(OSS_in, OSS_out, HI_in, HI_out) OSS_out%v_ice_A(i2,j2) = OSS_in%v_ice_A(i,j) enddo ; enddo - if (OSS_out%first_copy) then - OSS_in%first_copy = .false. ; OSS_out%first_copy = .false. - OSS_out%num_tr = OSS_in%num_tr - if (OSS_out%num_tr > 0) then - allocate(OSS_out%tr_array(HI_out%isd:HI_out%ied,HI_out%jsd:HI_out%jed,OSS_out%num_tr)) - OSS_out%tr_array(:,:,:) = 0.0 - endif - endif - - do m=1,OSS_in%num_tr ; do j=jsc,jec ; do i=isc,iec ; i2=i+i_off ; j2=j+j_off - OSS_out%tr_array(i2,j2,m) = OSS_in%tr_array(i,j,m) - enddo ; enddo ; enddo + call coupler_type_copy_data(OSS_in%tr_fields, OSS_out%tr_fields) end subroutine copy_sOSS_to_sOSS @@ -1235,46 +1214,17 @@ subroutine redistribute_sOSS_to_sOSS(OSS_in, OSS_out, domain_in, domain_out, HI_ type(domain2d), intent(in) :: domain_out !< The target data domain. type(hor_index_type), optional, intent(in) :: HI_out !< The hor_index_type on the target domain; HI_out !! may be omitted if this is not a target PE. - real, pointer, dimension(:,:) :: null_ptr => NULL() - logical :: first_copy - integer :: m, num_tr + type(coupler_2d_bc_type) :: null_bc if (.not. (associated(OSS_out) .or. associated(OSS_in))) & call SIS_error(FATAL, "redistribute_sOSS_to_sOSS called with "//& "neither OSS_in nor OSS_out associated.") - first_copy = .false. - if (associated(OSS_out)) first_copy = OSS_out%first_copy - if (associated(OSS_in)) first_copy = first_copy .or. OSS_in%first_copy - - if (first_copy) then - ! Determine the number of fluxes. - num_tr = 0 ; if (associated(OSS_in)) num_tr = OSS_in%num_tr - call max_across_PEs(num_tr) - - if (associated(OSS_out)) then - if (.not. present(HI_out)) & - call SIS_error(FATAL, "redistribute_sOSS_to_sOSS called with an "//& - "associated OSS_out but without HI_out.") - OSS_out%num_tr = num_tr - if ((num_tr > 0) .and. .not.allocated(OSS_out%tr_array)) then - allocate(OSS_out%tr_array(HI_out%isd:HI_out%ied,HI_out%jsd:HI_out%jed,num_tr)) - OSS_out%tr_array(:,:,:) = 0.0 - endif - OSS_out%first_copy = .false. - endif - - if (associated(OSS_in)) OSS_in%first_copy = .false. - endif if (associated(OSS_out) .and. associated(OSS_in)) then - ! The extra tracer arrays are copied first so that they can all have - ! complete=.false. - do m=1,OSS_in%num_tr - call mpp_redistribute(domain_in, OSS_in%tr_array(:,:,m), domain_out, & - OSS_out%tr_array(:,:,m), complete=.false.) - enddo - + ! This could have complete set to .false. if the halo sizes matched. + call coupler_type_redistribute_data(OSS_in%tr_fields, domain_in, & + OSS_out%tr_fields, domain_out, complete=.false.) call mpp_redistribute(domain_in, OSS_in%SST_C, domain_out, & OSS_out%SST_C, complete=.false.) call mpp_redistribute(domain_in, OSS_in%s_surf, domain_out, & @@ -1293,11 +1243,9 @@ subroutine redistribute_sOSS_to_sOSS(OSS_in, OSS_out, domain_in, domain_out, HI_ OSS_out%v_ice_A, complete=.true.) elseif (associated(OSS_out)) then ! Use the null pointer in place of the unneeded input arrays. - do m=1,OSS_out%num_tr - call mpp_redistribute(domain_in, null_ptr, domain_out, & - OSS_out%tr_array(:,:,m), complete=.false.) - enddo + call coupler_type_redistribute_data(null_bc, domain_in, & + OSS_out%tr_fields, domain_out, complete=.false.) call mpp_redistribute(domain_in, null_ptr, domain_out, & OSS_out%SST_C, complete=.false.) call mpp_redistribute(domain_in, null_ptr, domain_out, & @@ -1316,11 +1264,8 @@ subroutine redistribute_sOSS_to_sOSS(OSS_in, OSS_out, domain_in, domain_out, HI_ OSS_out%v_ice_A, complete=.true.) elseif (associated(OSS_in)) then ! Use the null pointer in place of the unneeded output arrays. - do m=1,OSS_in%num_tr - call mpp_redistribute(domain_in, OSS_in%tr_array(:,:,m), domain_out, & - null_ptr, complete=.false.) - enddo - + call coupler_type_redistribute_data(OSS_in%tr_fields, domain_in, & + null_bc, domain_out, complete=.false.) call mpp_redistribute(domain_in, OSS_in%SST_C, domain_out, & null_ptr, complete=.false.) call mpp_redistribute(domain_in, OSS_in%s_surf, domain_out, & @@ -1426,29 +1371,8 @@ subroutine copy_FIA_to_FIA(FIA_in, FIA_out, HI_in, HI_out, IG) "copy_FIA_to_FIA called an inconsistent number of time for the input and output types.") FIA_in%copy_calls = FIA_in%copy_calls + 1 ; FIA_out%copy_calls = FIA_out%copy_calls + 1 - if (FIA_in%copy_calls <= 2) then - if ((FIA_in%num_tr_fluxes >= 0) .and. (FIA_out%num_tr_fluxes < 0)) then - ! Allocate the tr_flux_top arrays to accommodate the size of the input - ! fluxes. This only occurs the first time FIA_out is copied from a fully - ! initialized FIA_in. - FIA_out%num_tr_fluxes = FIA_in%num_tr_fluxes - if (FIA_out%num_tr_fluxes > 0) then - isd = HI_out%isd ; ied = HI_out%ied ; jsd = HI_out%jsd ; jed = HI_out%jed - allocate(FIA_out%tr_flux_top(isd:ied, jsd:jed, 0:ncat, FIA_out%num_tr_fluxes)) - FIA_out%tr_flux_top(:,:,:,:) = 0.0 - endif - endif - endif - - if (FIA_in%num_tr_fluxes >= 0) then - if (FIA_in%num_tr_fluxes /= FIA_out%num_tr_fluxes) & - call SIS_error(FATAL, "copy_FIA_to_FIA called with different num_tr_fluxes.") - do n=1,FIA_in%num_tr_fluxes ; do k=0,ncat ; do j=jsc,jec ; do i=isc,iec - i2 = i+i_off ; j2 = j+j_off - FIA_out%tr_flux_top(i2,j2,k,n) = FIA_in%tr_flux_top(i,j,k,n) - enddo ; enddo ; enddo ; enddo - endif + call coupler_type_copy_data(FIA_in%tr_flux, FIA_out%tr_flux) ! avg_count, atmos_winds, and the IDs are deliberately not being copied. end subroutine copy_FIA_to_FIA @@ -1467,9 +1391,9 @@ subroutine redistribute_FIA_to_FIA(FIA_in, FIA_out, domain_in, domain_out, G_out real, pointer, dimension(:,:) :: null_ptr2D => NULL() real, pointer, dimension(:,:,:) :: null_ptr3D => NULL() real, pointer, dimension(:,:,:,:) :: null_ptr4D => NULL() + type(coupler_3d_bc_type) :: null_bc integer :: copy_calls ! The number of times these FIA_types have been copied. integer :: i, j, b, isd, ied, jsd, jed, ncat - integer :: num_tr copy_calls = 0 if (associated(FIA_out)) then @@ -1481,28 +1405,6 @@ subroutine redistribute_FIA_to_FIA(FIA_in, FIA_out, domain_in, domain_out, G_out copy_calls = max(copy_calls,FIA_in%copy_calls) endif - if (copy_calls <= 2) then - ! Determine the number of fluxes. - num_tr = 0 ; if (associated(FIA_in)) num_tr = FIA_in%num_tr_fluxes - call max_across_PEs(num_tr) - - if (associated(FIA_out)) then - if (.not. present(G_out)) & - call SIS_error(FATAL, "redistribute_sFIA_to_sFIA called with an "//& - "associated FIA_out but without G_out.") - if (.not. present(IG)) & - call SIS_error(FATAL, "redistribute_sFIA_to_sFIA called with an "//& - "associated FIA_out but without IG.") - FIA_out%num_tr_fluxes = num_tr - if ((num_tr > 0) .and. .not.allocated(FIA_out%tr_flux_top)) then - isd = G_out%isd ; ied = G_out%ied ; jsd = G_out%jsd ; jed = G_out%jed - ncat = IG%CatIce - allocate(FIA_out%tr_flux_top(isd:ied, jsd:jed, 0:ncat, num_tr)) - FIA_out%tr_flux_top(:,:,:,:) = 0.0 - endif - endif - endif - ! FIA%flux_u_top and flux_v_top are deliberately not being copied, as they ! are only needed on the fast_ice_PEs ! FIA%frazil_left is deliberately not being copied, as it is only valid on @@ -1512,6 +1414,8 @@ subroutine redistribute_FIA_to_FIA(FIA_in, FIA_out, domain_in, domain_out, G_out ! avg_count, atmos_winds, and the IDs are deliberately not being copied. if (associated(FIA_out) .and. associated(FIA_in)) then + call coupler_type_redistribute_data(FIA_in%tr_flux, domain_in, & + FIA_out%tr_flux, domain_out, complete=.false.) do b=1,size(FIA_in%flux_sw_top,4) call mpp_redistribute(domain_in, FIA_in%flux_sw_top(:,:,:,b), domain_out, & FIA_out%flux_sw_top(:,:,:,b), complete=.false.) @@ -1582,12 +1486,10 @@ subroutine redistribute_FIA_to_FIA(FIA_in, FIA_out, domain_in, domain_out, G_out FIA_out%Tskin_cat, complete=.true.) endif - if (FIA_in%num_tr_fluxes > 0) then - call mpp_redistribute(domain_in, FIA_in%tr_flux_top, domain_out, & - FIA_out%tr_flux_top) - endif elseif (associated(FIA_out)) then ! Use the null pointers in place of the unneeded input arrays. + call coupler_type_redistribute_data(null_bc, domain_in, & + FIA_out%tr_flux, domain_out, complete=.false.) do b=1,size(FIA_out%flux_sw_top,4) call mpp_redistribute(domain_in, null_ptr3D, domain_out, & FIA_out%flux_sw_top(:,:,:,b), complete=.false.) @@ -1659,12 +1561,10 @@ subroutine redistribute_FIA_to_FIA(FIA_in, FIA_out, domain_in, domain_out, G_out endif - if (FIA_out%num_tr_fluxes > 0) then - call mpp_redistribute(domain_in, null_ptr4D, domain_out, & - FIA_out%tr_flux_top) - endif elseif (associated(FIA_in)) then ! Use the null pointers in place of the unneeded output arrays. + call coupler_type_redistribute_data(FIA_in%tr_flux, domain_in, & + null_bc, domain_out, complete=.false.) do b=1,size(FIA_in%flux_sw_top,4) call mpp_redistribute(domain_in, FIA_in%flux_sw_top(:,:,:,b), domain_out, & null_ptr3D, complete=.false.) @@ -1735,10 +1635,6 @@ subroutine redistribute_FIA_to_FIA(FIA_in, FIA_out, domain_in, domain_out, G_out null_ptr3D, complete=.true.) endif - if (FIA_in%num_tr_fluxes > 0) then - call mpp_redistribute(domain_in, FIA_in%tr_flux_top, domain_out, & - null_ptr4D) - endif else call SIS_error(FATAL, "redistribute_FIA_to_FIA called with "//& "neither FIA_in nor FIA_out associated.") @@ -1783,31 +1679,9 @@ subroutine copy_TSF_to_TSF(TSF_in, TSF_out, HI_in, HI_out) if (TSF_in%copy_calls /= TSF_out%copy_calls) call SIS_error(WARNING, & "copy_TSF_to_TSF called an inconsistent number of time for the input and output types.") - TSF_in%copy_calls = TSF_in%copy_calls + 1 ; TSF_out%copy_calls = TSF_out%copy_calls + 1 - if (TSF_in%copy_calls <= 2) then - if ((TSF_in%num_tr_fluxes >= 0) .and. (TSF_out%num_tr_fluxes < 0)) then - ! Allocate the tr_flux_top arrays to accommodate the size of the input - ! fluxes. This only occurs the first time TSF_out is copied from a fully - ! initialized TSF_in. - TSF_out%num_tr_fluxes = TSF_in%num_tr_fluxes - if (TSF_out%num_tr_fluxes > 0) then - isd = HI_out%isd ; ied = HI_out%ied ; jsd = HI_out%jsd ; jed = HI_out%jed - allocate(TSF_out%tr_flux(isd:ied, jsd:jed, TSF_out%num_tr_fluxes)) - TSF_out%tr_flux(:,:,:) = 0.0 - endif - endif - endif - if (TSF_in%num_tr_fluxes >= 0) then - if (TSF_in%num_tr_fluxes /= TSF_out%num_tr_fluxes) & - call SIS_error(FATAL, "copy_TSF_to_TSF called with different num_tr_fluxes.") - - do n=1,TSF_in%num_tr_fluxes ; do j=jsc,jec ; do i=isc,iec - i2 = i+i_off ; j2 = j+j_off - TSF_out%tr_flux(i2,j2,n) = TSF_in%tr_flux(i,j,n) - enddo ; enddo ; enddo - endif + call coupler_type_copy_data(TSF_in%tr_flux, TSF_out%tr_flux) end subroutine copy_TSF_to_TSF @@ -1823,6 +1697,7 @@ subroutine redistribute_TSF_to_TSF(TSF_in, TSF_out, domain_in, domain_out, HI_ou !! may be omitted if this is not a target PE. real, pointer, dimension(:,:) :: null_ptr2D => NULL() + type(coupler_2d_bc_type) :: null_bc integer :: copy_calls ! The number of times these TSF_types have been copied. integer :: b, m, num_tr @@ -1839,30 +1714,11 @@ subroutine redistribute_TSF_to_TSF(TSF_in, TSF_out, domain_in, domain_out, HI_ou copy_calls = max(copy_calls,TSF_in%copy_calls) endif - if (copy_calls <= 2) then - ! Determine the number of fluxes. - num_tr = 0 ; if (associated(TSF_in)) num_tr = TSF_in%num_tr_fluxes - call max_across_PEs(num_tr) - - if (associated(TSF_out)) then - if (.not. present(HI_out)) & - call SIS_error(FATAL, "redistribute_TSF_to_TSF called with an "//& - "associated TSF_out but without HI_out.") - TSF_out%num_tr_fluxes = num_tr - if ((num_tr > 0) .and. .not.allocated(TSF_out%tr_flux)) then - allocate(TSF_out%tr_flux(HI_out%isd:HI_out%ied,HI_out%jsd:HI_out%jed,num_tr)) - TSF_out%tr_flux(:,:,:) = 0.0 - endif - endif - endif - if (associated(TSF_out) .and. associated(TSF_in)) then ! The extra tracer arrays are copied first so that they can all have ! complete=.false. - do m=1,TSF_in%num_tr_fluxes - call mpp_redistribute(domain_in, TSF_in%tr_flux(:,:,m), domain_out, & - TSF_out%tr_flux(:,:,m), complete=.false.) - enddo + call coupler_type_redistribute_data(TSF_in%tr_flux, domain_in, & + TSF_out%tr_flux, domain_out, complete=.false.) do b=1,size(TSF_in%flux_sw,3) call mpp_redistribute(domain_in, TSF_in%flux_sw(:,:,b), domain_out, & TSF_out%flux_sw(:,:,b), complete=.false.) @@ -1881,10 +1737,8 @@ subroutine redistribute_TSF_to_TSF(TSF_in, TSF_out, domain_in, domain_out, HI_ou TSF_out%fprec, complete=.true.) elseif (associated(TSF_out)) then ! Use the null pointer in place of the unneeded input arrays. - do m=1,TSF_out%num_tr_fluxes - call mpp_redistribute(domain_in, null_ptr2D, domain_out, & - TSF_out%tr_flux(:,:,m), complete=.false.) - enddo + call coupler_type_redistribute_data(null_bc, domain_in, & + TSF_out%tr_flux, domain_out, complete=.false.) do b=1,size(TSF_out%flux_sw,3) call mpp_redistribute(domain_in, null_ptr2D, domain_out, & TSF_out%flux_sw(:,:,b), complete=.false.) @@ -1903,10 +1757,8 @@ subroutine redistribute_TSF_to_TSF(TSF_in, TSF_out, domain_in, domain_out, HI_ou TSF_out%fprec, complete=.true.) elseif (associated(TSF_in)) then ! Use the null pointer in place of the unneeded output arrays. - do m=1,TSF_in%num_tr_fluxes - call mpp_redistribute(domain_in, TSF_in%tr_flux(:,:,m), domain_out, & - null_ptr2D, complete=.false.) - enddo + call coupler_type_redistribute_data(TSF_in%tr_flux, domain_in, & + null_bc, domain_out, complete=.false.) do b=1,size(TSF_in%flux_sw,3) call mpp_redistribute(domain_in, TSF_in%flux_sw(:,:,b), domain_out, & null_ptr2D, complete=.false.) @@ -2108,19 +1960,13 @@ subroutine register_fast_to_slow_restarts(FIA, Rad, TSF, mpp_domain, Ice_restart !### These tracer fluxes will need to be dealt with, but because tracer packages can !### be turned on or off at a model restart, these will have to be kept as - !### coupler_2d_bc_type and coupler_3d_bc_type structures, and be dealt with outside - !### at the coupler level, using the whole coupler types package. - -! if (FIA_in%num_tr_fluxes >= 0) then -! do n=1,TSF_in%num_tr_fluxes ; do j=jsc,jec ; do i=isc,iec -! i2 = i+i_off ; j2 = j+j_off -! TSF_out%tr_flux(i2,j2,n) = TSF_in%tr_flux(i,j,n) - -! do n=1,FIA_in%num_tr_fluxes ; do k=0,ncat ; do j=jsc,jec ; do i=isc,iec -! i2 = i+i_off ; j2 = j+j_off -! FIA_out%tr_flux_top(i2,j2,k,n) = FIA_in%tr_flux_top(i,j,k,n) -! enddo ; enddo ; enddo ; enddo -! endif + !### coupler_2d_bc_type and coupler_3d_bc_type structures. To enable these, the + !### coupler types will have to be read before the ice model is initialized and + !### then passed in. + +!### call coupler_type_register_restarts(TSF%tr_flux, restart_file, Ice_restart, mpp_domain) +!### call coupler_type_register_restarts(FIA%tr_flux, restart_file, Ice_restart, mpp_domain) + end subroutine register_fast_to_slow_restarts !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! diff --git a/src/combined_ice_ocean_driver.F90 b/src/combined_ice_ocean_driver.F90 index c9ff864f..17774e75 100644 --- a/src/combined_ice_ocean_driver.F90 +++ b/src/combined_ice_ocean_driver.F90 @@ -39,7 +39,7 @@ module combined_ice_ocean_driver use ocean_model_mod, only : ocean_public_type, ocean_state_type, ice_ocean_boundary_type use coupler_types_mod, only: coupler_type_send_data, coupler_type_data_override -use coupler_types_mod, only: coupler_type_copy_data, coupler_type_redistribute_data +use coupler_types_mod, only: coupler_type_copy_data use data_override_mod, only : data_override use diag_manager_mod, only : send_data use mpp_domains_mod, only : domain2D, mpp_get_layout, mpp_get_compute_domain diff --git a/src/ice_model.F90 b/src/ice_model.F90 index 96688d25..9964872f 100644 --- a/src/ice_model.F90 +++ b/src/ice_model.F90 @@ -63,6 +63,9 @@ module ice_model_mod use MOM_time_manager, only : set_date, set_time, operator(+), operator(-) use MOM_time_manager, only : operator(>), operator(*), operator(/), operator(/=) +use coupler_types_mod, only : coupler_2d_bc_type, coupler_3d_bc_type +use coupler_types_mod, only : coupler_type_spawn, coupler_type_initialized +use coupler_types_mod, only : coupler_type_rescale_data, coupler_type_copy_data use fms_mod, only : file_exist, clock_flag_default use fms_io_mod, only : set_domain, nullify_domain, restore_state, query_initialized use fms_io_mod, only : restore_state, query_initialized @@ -73,7 +76,6 @@ module ice_model_mod use astronomy_mod, only : astronomy_init, astronomy_end use astronomy_mod, only : universal_time, orbital_time, diurnal_solar, daily_mean_solar -use coupler_types_mod, only : coupler_3d_bc_type use ocean_albedo_mod, only : compute_ocean_albedo ! ice sets ocean surface use ocean_rough_mod, only : compute_ocean_roughness ! properties over water @@ -442,12 +444,44 @@ subroutine exchange_fast_to_slow_ice(Ice) type(ice_state_type), pointer :: IST_null => NULL() type(ice_rad_type), pointer :: Rad_null => NULL() type(total_sfc_flux_type), pointer :: TSF_null => NULL() + + integer :: isc, iec, jsc, jec, isd, ied, jsd, jed logical :: redo_fast_update redo_fast_update = .false. if (associated(Ice%fCS)) redo_fast_update = Ice%fCS%redo_fast_update if (associated(Ice%sCS)) redo_fast_update = Ice%sCS%redo_fast_update + if (associated(Ice%fCS)) then + isc = Ice%fCS%G%isc ; iec = Ice%fCS%G%iec ; jsc = Ice%fCS%G%jsc ; jec = Ice%fCS%G%jec + isd = Ice%fCS%G%isd ; ied = Ice%fCS%G%ied ; jsd = Ice%fCS%G%jsd ; jed = Ice%fCS%G%jed + + ! Propagate the coupler_type info to Ice%fCS%FIA%tr_flux and allocate its arrays. + call coupler_type_spawn(Ice%ocean_fluxes, Ice%fCS%FIA%tr_flux, & + (/isd, isc, iec, ied/), (/jsd, jsc, jec, jed/), & + (/0, Ice%fCS%IG%CatIce/), as_needed=.true.) + + if (redo_fast_update) & + ! Propagate the coupler_type info to Ice%fCS%TSF%tr_flux and allocate its arrays. + call coupler_type_spawn(Ice%ocean_fluxes, Ice%fCS%TSF%tr_flux, & + (/isd, isc, iec, ied/), (/jsd, jsc, jec, jed/), as_needed=.true.) + endif + + if (associated(Ice%sCS)) then + isc = Ice%sCS%G%isc ; iec = Ice%sCS%G%iec ; jsc = Ice%sCS%G%jsc ; jec = Ice%sCS%G%jec + isd = Ice%sCS%G%isd ; ied = Ice%sCS%G%ied ; jsd = Ice%sCS%G%jsd ; jed = Ice%sCS%G%jed + + ! Propagate the coupler_type info to Ice%sCS%FIA%tr_flux and allocate its arrays. + call coupler_type_spawn(Ice%ocean_fluxes, Ice%sCS%FIA%tr_flux, & + (/isd, isc, iec, ied/), (/jsd, jsc, jec, jed/), & + (/0, Ice%sCS%IG%CatIce/), as_needed=.true.) + + if (redo_fast_update) & + ! Propagate the coupler_type info to Ice%sCS%TSF%tr_flux and allocate its arrays. + call coupler_type_spawn(Ice%ocean_fluxes, Ice%sCS%TSF%tr_flux, & + (/isd, isc, iec, ied/), (/jsd, jsc, jec, jed/), as_needed=.true.) + endif + if(Ice%xtype == DIRECT) then if (.not.associated(Ice%fCS) .or. .not.associated(Ice%sCS)) call SIS_error(FATAL, & "With xtype=DIRECT, both the pointer to Ice%sCS and the pointer to Ice%fCS must be "//& @@ -540,7 +574,7 @@ subroutine set_ocean_top_fluxes(Ice, IST, IOF, FIA, OSS, G, IG, sCS) real :: I_count integer :: i, j, k, isc, iec, jsc, jec, m, n - integer :: i2, j2, i_off, j_off, ind, ncat, NkIce + integer :: i2, j2, i_off, j_off, ncat, NkIce isc = G%isc ; iec = G%iec ; jsc = G%jsc ; jec = G%jec ncat = IG%CatIce ; NkIce = IG%NkIce @@ -568,9 +602,7 @@ subroutine set_ocean_top_fluxes(Ice, IST, IOF, FIA, OSS, G, IG, sCS) Ice%flux_sw_vis_dir(:,:) = 0.0 ; Ice%flux_sw_vis_dif(:,:) = 0.0 Ice%flux_lw(:,:) = 0.0 ; Ice%flux_lh(:,:) = 0.0 Ice%fprec(:,:) = 0.0 ; Ice%lprec(:,:) = 0.0 - do n=1,Ice%ocean_fluxes%num_bcs ; do m=1,Ice%ocean_fluxes%bc(n)%num_fields - Ice%ocean_fluxes%bc(n)%field(m)%values(:,:) = 0.0 - enddo ; enddo + call coupler_type_rescale_data(Ice%ocean_fluxes, 0.0) i_off = LBOUND(Ice%flux_t,1) - G%isc ; j_off = LBOUND(Ice%flux_t,2) - G%jsc !$OMP parallel do default(none) shared(isc,iec,jsc,jec,Ice,IST,IOF,FIA,i_off,j_off,G,OSS) & @@ -610,22 +642,10 @@ subroutine set_ocean_top_fluxes(Ice, IST, IOF, FIA, OSS, G, IG, sCS) enddo ; enddo endif - ind = 0 - do n=1,Ice%ocean_fluxes%num_bcs ; do m=1,Ice%ocean_fluxes%bc(n)%num_fields - ind = ind + 1 - if (ind <= IOF%num_tr_fluxes) then - do j=jsc,jec ; do i=isc,iec - i2 = i+i_off ; j2 = j+j_off ! Use these to correct for indexing differences. - Ice%ocean_fluxes%bc(n)%field(m)%values(i2,j2) = IOF%tr_flux_ocn_top(i,j,ind) - enddo ; enddo - else - ! This can occur the first step of a cold-start run with lagged ice - ! coupling, but otherwise it may indicate a problem that should be trapped. - do j2=jsc+j_off,jec+j_off ; do i2=isc+i_off,iec+i_off - Ice%ocean_fluxes%bc(n)%field(m)%values(i2,j2) = 0.0 - enddo ; enddo - endif - enddo ; enddo + ! This copy may need to be skipped in the first step of a cold-start run with lagged ice + ! coupling, but otherwise if it is skipped may indicate a problem that should be trapped. + if (coupler_type_initialized(IOF%tr_flux_ocn_top)) & + call coupler_type_copy_data(IOF%tr_flux_ocn_top, Ice%ocean_fluxes) if (sCS%debug) then call Ice_public_type_chksum("End set_ocean_top_fluxes", Ice, check_slow=.true.) @@ -753,8 +773,18 @@ subroutine exchange_slow_to_fast_ice(Ice) type(simple_OSS_type), pointer :: sOSS_null => NULL() type(ice_state_type), pointer :: IST_null => NULL() + integer :: isc, iec, jsc, jec, isd, ied, jsd, jed + call mpp_clock_begin(iceClock) ; call mpp_clock_begin(ice_clock_exchange) + if (associated(Ice%fCS)) then + isc = Ice%fCS%G%isc ; iec = Ice%fCS%G%iec ; jsc = Ice%fCS%G%jsc ; jec = Ice%fCS%G%jec + isd = Ice%fCS%G%isd ; ied = Ice%fCS%G%ied ; jsd = Ice%fCS%G%jsd ; jed = Ice%fCS%G%jed + + ! Propagate the coupler_type info to Ice%fCS%sOSS%tr_fields and allocate its arrays. + call coupler_type_spawn(Ice%ocean_fields, Ice%fCS%sOSS%tr_fields, & + (/isd, isc, iec, ied/), (/jsd, jsc, jec, jed/), as_needed=.true. ) + endif if (Ice%xtype == DIRECT) then if (.not.associated(Ice%fCS) .or. .not.associated(Ice%sCS)) call SIS_error(FATAL, & @@ -841,7 +871,10 @@ subroutine unpack_ocn_ice_bdry(OIB, OSS, ITV, G, specified_ice, ocean_fields) real, parameter :: T_0degC = 273.15 ! 0 degrees C in Kelvin logical :: Cgrid_ocn integer :: i, j, k, m, n, i2, j2, k2, isc, iec, jsc, jec, i_off, j_off, index + integer :: isd, ied, jsd, jed + isc = G%isc ; iec = G%iec ; jsc = G%jsc ; jec = G%jec + isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed i_off = LBOUND(OIB%t,1) - G%isc ; j_off = LBOUND(OIB%t,2) - G%jsc call mpp_clock_begin(iceClock) ; call mpp_clock_begin(ice_clock_slow) @@ -948,27 +981,9 @@ subroutine unpack_ocn_ice_bdry(OIB, OSS, ITV, G, specified_ice, ocean_fields) call pass_var(OSS%sea_lev, G%Domain) ! Transfer the ocean state for extra tracer fluxes. - if (OSS%num_tr<0) then - ! Determine the number of tracer arrays and allocate them. - OSS%num_tr = 0 - do n=1,OIB%fields%num_bcs - OSS%num_tr = OSS%num_tr + OIB%fields%bc(n)%num_fields - enddo - if (OSS%num_tr > 0) then - allocate(OSS%tr_array(G%isd:G%ied,G%jsd:G%jed,OSS%num_tr)) ; OSS%tr_array(:,:,:) = 0.0 - endif - endif - index = 0 - do n=1,OIB%fields%num_bcs ; do m=1,OIB%fields%bc(n)%num_fields - index = index+1 - if (index == 1) then - i_off = LBOUND(OIB%fields%bc(n)%field(m)%values,1) - G%isc - j_off = LBOUND(OIB%fields%bc(n)%field(m)%values,2) - G%jsc - endif - do j=jsc,jec ; do i=isc,iec ; i2 = i+i_off ; j2 = j+j_off - OSS%tr_array(i,j,index) = OIB%fields%bc(n)%field(m)%values(i2,j2) - enddo ; enddo - enddo ; enddo + call coupler_type_spawn(OIB%fields, OSS%tr_fields, (/isd, isc, iec, ied/), & + (/jsd, jsc, jec, jed/), as_needed=.true. ) + call coupler_type_copy_data(OIB%fields, OSS%tr_fields) call mpp_clock_end(ice_clock_slow) ; call mpp_clock_end(iceClock) @@ -1169,18 +1184,9 @@ subroutine set_ice_surface_state(Ice, IST, OSS, Rad, FIA, G, IG, fCS) enddo endif - ! Copy over the additional tracer fields into the ocean_fields structure. - index = 0 - do n=1,Ice%ocean_fields%num_bcs ; do m=1,Ice%ocean_fields%bc(n)%num_fields - index = index+1 - if (index == 1) then - i_off = LBOUND(Ice%ocean_fields%bc(n)%field(m)%values,1) - G%isc - j_off = LBOUND(Ice%ocean_fields%bc(n)%field(m)%values,2) - G%jsc - endif - do j=jsc,jec ; do i=isc,iec ; i2 = i+i_off ; j2 = j+j_off - Ice%ocean_fields%bc(n)%field(m)%values(i2,j2,1) = OSS%tr_array(i,j,index) - enddo ; enddo - enddo ; enddo + ! Copy over the additional tracer fields into the the open-ocean category of + ! the Ice%ocean_fields structure. + call coupler_type_copy_data(OSS%tr_fields, Ice%ocean_fields, ind3_start=1, ind3_end=1) if (fCS%debug) then call IST_chksum("End set_ice_surface_state", IST, G, IG) From 9020bf13e1462c46a0ed7156bdb5effe39aa654a Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Mon, 31 Jul 2017 16:02:33 -0400 Subject: [PATCH 3/4] +Added optional gas_flux args to ice_model_init Added optional arguments gas_flux and gas_fields_ocn arguments of type coupler_1d_bc_type to several SIS2 initialization routines, including ice_model_init. If these are present, they are used to spawn the types used to describe additional gas or tracer fluxes and fields at the same time as other SIS2 arrays are allocated. There is a new version of the coupler that makes use of these new arguments, but because they are optional, older versions of the coupler will work also. All answers are bitwise identical. --- src/SIS_types.F90 | 77 ++++++++++++++++++++++++++++++++++++----------- src/ice_model.F90 | 31 +++++++++++++------ src/ice_type.F90 | 33 +++++++++++++++----- 3 files changed, 105 insertions(+), 36 deletions(-) diff --git a/src/SIS_types.F90 b/src/SIS_types.F90 index 26ac2d63..7e8e5e0a 100644 --- a/src/SIS_types.F90 +++ b/src/SIS_types.F90 @@ -14,7 +14,7 @@ module SIS_types use fms_io_mod, only : register_restart_field, restart_file_type use fms_io_mod, only : restore_state, query_initialized use time_manager_mod, only : time_type, time_type_to_real -use coupler_types_mod, only : coupler_2d_bc_type, coupler_3d_bc_type +use coupler_types_mod, only : coupler_1d_bc_type, coupler_2d_bc_type, coupler_3d_bc_type use coupler_types_mod, only : coupler_type_spawn, coupler_type_copy_data use coupler_types_mod, only : coupler_type_redistribute_data use SIS_hor_grid, only : SIS_hor_grid_type @@ -729,16 +729,20 @@ end subroutine ice_state_read_alt_restarts !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! !> alloc_fast_ice_avg allocates and zeros out the arrays in a fast_ice_avg_type. -subroutine alloc_fast_ice_avg(FIA, HI, IG, interp_fluxes) +subroutine alloc_fast_ice_avg(FIA, HI, IG, interp_fluxes, gas_fluxes) type(fast_ice_avg_type), pointer :: FIA type(hor_index_type), intent(in) :: HI type(ice_grid_type), intent(in) :: IG logical, intent(in) :: interp_fluxes - - integer :: isd, ied, jsd, jed, CatIce + type(coupler_1d_bc_type), & + optional, intent(in) :: gas_fluxes !< If present, this type describes the + !! additional gas or other tracer fluxes between the + !! ocean, ice, and atmosphere. + integer :: isc, iec, jsc, jec, isd, ied, jsd, jed, CatIce if (.not.associated(FIA)) allocate(FIA) CatIce = IG%CatIce + isc = HI%isc ; iec = HI%iec ; jsc = HI%jsc ; jec = HI%jec isd = HI%isd ; ied = HI%ied ; jsd = HI%jsd ; jed = HI%jed FIA%avg_count = 0 @@ -783,17 +787,27 @@ subroutine alloc_fast_ice_avg(FIA, HI, IG, interp_fluxes) allocate(FIA%flux_sw_dn(isd:ied, jsd:jed, NBANDS)) ; FIA%flux_sw_dn(:,:,:) = 0.0 allocate(FIA%sw_abs_ocn(isd:ied, jsd:jed, CatIce)) ; FIA%sw_abs_ocn(:,:,:) = 0.0 + if (present(gas_fluxes)) & + call coupler_type_spawn(gas_fluxes, FIA%tr_flux, (/isd, isc, iec, ied/), & + (/jsd, jsc, jec, jed/), (/0, CatIce/)) + end subroutine alloc_fast_ice_avg !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! !> alloc_total_sfc_flux allocates and zeros out the arrays in a total_sfc_flux_type. -subroutine alloc_total_sfc_flux(TSF, HI) - type(total_sfc_flux_type), pointer :: TSF - type(hor_index_type), intent(in) :: HI +subroutine alloc_total_sfc_flux(TSF, HI, gas_fluxes) + type(total_sfc_flux_type), pointer :: TSF !< The total surface flux type being allocated + type(hor_index_type), intent(in) :: HI !< The hor_index_type with information about the + !! array extents to be allocated. + type(coupler_1d_bc_type), & + optional, intent(in) :: gas_fluxes !< If present, this type describes the + !! additional gas or other tracer fluxes between the + !! ocean, ice, and atmosphere. - integer :: isd, ied, jsd, jed + integer :: isc, iec, jsc, jec, isd, ied, jsd, jed if (.not.associated(TSF)) allocate(TSF) + isc = HI%isc ; iec = HI%iec ; jsc = HI%jsc ; jec = HI%jec isd = HI%isd ; ied = HI%ied ; jsd = HI%jsd ; jed = HI%jed allocate(TSF%flux_u(isd:ied, jsd:jed)) ; TSF%flux_u(:,:) = 0.0 @@ -805,6 +819,9 @@ subroutine alloc_total_sfc_flux(TSF, HI) allocate(TSF%evap(isd:ied, jsd:jed)) ; TSF%evap(:,:) = 0.0 allocate(TSF%lprec(isd:ied, jsd:jed)) ; TSF%lprec(:,:) = 0.0 allocate(TSF%fprec(isd:ied, jsd:jed)) ; TSF%fprec(:,:) = 0.0 + if (present(gas_fluxes)) & + call coupler_type_spawn(gas_fluxes, TSF%tr_flux, (/isd, isc, iec, ied/), & + (/jsd, jsc, jec, jed/)) end subroutine alloc_total_sfc_flux @@ -917,10 +934,21 @@ end subroutine alloc_ice_ocean_flux !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! !> alloc_ocean_sfc_state allocates and zeros out the arrays in an ocean_sfc_state_type. -subroutine alloc_ocean_sfc_state(OSS, HI, Cgrid_dyn) - type(ocean_sfc_state_type), pointer :: OSS - type(hor_index_type), intent(in) :: HI - logical, intent(in) :: Cgrid_dyn +subroutine alloc_ocean_sfc_state(OSS, HI, Cgrid_dyn, gas_fields_ocn) + type(ocean_sfc_state_type), pointer :: OSS !< The ocean_sfc_state_type being allocated + type(hor_index_type), intent(in) :: HI !< The hor_index_type with information about the + !! array extents to be allocated. + logical, intent(in) :: Cgrid_dyn !< A variable indicating whether the ice + !! ice dynamics are calculated on a C-grid (true) + !! or on a B-grid (false). + type(coupler_1d_bc_type), & + optional, intent(in) :: gas_fields_ocn !< If present, this type describes the + !! ocean and surface-ice fields that will participate + !! in the calculation of additional gas or other + !! tracer fluxes. + integer :: isc, iec, jsc, jec, isd, ied, jsd, jed + isc = HI%isc ; iec = HI%iec ; jsc = HI%jsc ; jec = HI%jec + isd = HI%isd ; ied = HI%ied ; jsd = HI%jsd ; jed = HI%jed if (.not.associated(OSS)) allocate(OSS) @@ -942,19 +970,30 @@ subroutine alloc_ocean_sfc_state(OSS, HI, Cgrid_dyn) OSS%Cgrid_dyn = Cgrid_dyn + if (present(gas_fields_ocn)) & + call coupler_type_spawn(gas_fields_ocn, OSS%tr_fields, (/isd, isc, iec, ied/), & + (/jsd, jsc, jec, jed/)) + end subroutine alloc_ocean_sfc_state !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! !> alloc_simple_ocean_sfc_state allocates and zeros out the arrays in a !! simple_OSS_type. -subroutine alloc_simple_OSS(OSS, HI) - type(simple_OSS_type), pointer :: OSS - type(hor_index_type), intent(in) :: HI - - integer :: isd, ied, jsd, jed +subroutine alloc_simple_OSS(OSS, HI, gas_fields_ocn) + type(simple_OSS_type), pointer :: OSS !< The simple_OSS_type being allocated + type(hor_index_type), intent(in) :: HI !< The hor_index_type with information about the + !! array extents to be allocated. + type(coupler_1d_bc_type), & + optional, intent(in) :: gas_fields_ocn !< If present, this type describes the + !! ocean and surface-ice fields that will participate + !! in the calculation of additional gas or other + !! tracer fluxes. + + integer :: isc, iec, jsc, jec, isd, ied, jsd, jed if (.not.associated(OSS)) allocate(OSS) + isc = HI%isc ; iec = HI%iec ; jsc = HI%jsc ; jec = HI%jec isd = HI%isd ; ied = HI%ied ; jsd = HI%jsd ; jed = HI%jed allocate(OSS%s_surf(isd:ied, jsd:jed)) ; OSS%s_surf(:,:) = 0.0 @@ -965,10 +1004,12 @@ subroutine alloc_simple_OSS(OSS, HI) allocate(OSS%v_ocn_A(isd:ied, jsd:jed)) ; OSS%v_ocn_A(:,:) = 0.0 allocate(OSS%u_ice_A(isd:ied, jsd:jed)) ; OSS%u_ice_A(:,:) = 0.0 allocate(OSS%v_ice_A(isd:ied, jsd:jed)) ; OSS%v_ice_A(:,:) = 0.0 + if (present(gas_fields_ocn)) & + call coupler_type_spawn(gas_fields_ocn, OSS%tr_fields, (/isd, isc, iec, ied/), & + (/jsd, jsc, jec, jed/)) end subroutine alloc_simple_OSS - !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! !> copy_IST_to_IST copies the computational domain of one ice state type into !! the computational domain of another ice_state_type. Both must use the same diff --git a/src/ice_model.F90 b/src/ice_model.F90 index 9964872f..94ac172d 100644 --- a/src/ice_model.F90 +++ b/src/ice_model.F90 @@ -63,7 +63,7 @@ module ice_model_mod use MOM_time_manager, only : set_date, set_time, operator(+), operator(-) use MOM_time_manager, only : operator(>), operator(*), operator(/), operator(/=) -use coupler_types_mod, only : coupler_2d_bc_type, coupler_3d_bc_type +use coupler_types_mod, only : coupler_1d_bc_type, coupler_2d_bc_type, coupler_3d_bc_type use coupler_types_mod, only : coupler_type_spawn, coupler_type_initialized use coupler_types_mod, only : coupler_type_rescale_data, coupler_type_copy_data use fms_mod, only : file_exist, clock_flag_default @@ -1607,7 +1607,7 @@ end subroutine add_diurnal_sw !> ice_model_init - initializes ice model data, parameters and diagnostics. It !! might operate on the fast ice processors, the slow ice processors or both. subroutine ice_model_init(Ice, Time_Init, Time, Time_step_fast, Time_step_slow, & - Verona_coupler, Concurrent_ice ) + Verona_coupler, Concurrent_ice, gas_fluxes, gas_fields_ocn ) type(ice_data_type), intent(inout) :: Ice !< The ice data type that is being initialized. type(time_type) , intent(in) :: Time_Init !< The starting time of the model integration @@ -1625,6 +1625,17 @@ subroutine ice_model_init(Ice, Time_Init, Time, Time_step_fast, Time_step_slow, !! settings appropriate for running the atmosphere and !! slow ice simultaneously, including embedding the !! slow sea-ice time stepping in the ocean model. + type(coupler_1d_bc_type), & + optional, intent(in) :: gas_fluxes !< If present, this type describes the + !! additional gas or other tracer fluxes between the + !! ocean, ice, and atmosphere, and can be used to + !! spawn related internal variables in the ice model. + type(coupler_1d_bc_type), & + optional, intent(in) :: gas_fields_ocn !< If present, this type describes the + !! ocean and surface-ice fields that will participate + !! in the calculation of additional gas or other + !! tracer fluxes, and can be used to spawn related + !! internal variables in the ice model. ! This include declares and sets the variable "version". #include "version_variable.h" @@ -2076,19 +2087,19 @@ subroutine ice_model_init(Ice, Time_Init, Time, Time_step_fast, Time_step_slow, call alloc_IST_arrays(sHI, sIG, sIST, omit_tsurf=Eulerian_tsurf) call ice_state_register_restarts(sIST, sG, sIG, Ice%Ice_restart, restart_file) - call alloc_ocean_sfc_state(Ice%sCS%OSS, sHI, sIST%Cgrid_dyn) + call alloc_ocean_sfc_state(Ice%sCS%OSS, sHI, sIST%Cgrid_dyn, gas_fields_ocn) Ice%sCS%OSS%kmelt = kmelt - call alloc_simple_OSS(Ice%sCS%sOSS, sHI) + call alloc_simple_OSS(Ice%sCS%sOSS, sHI, gas_fields_ocn) call alloc_ice_ocean_flux(Ice%sCS%IOF, sHI, do_iceberg_fields=Ice%sCS%do_icebergs) Ice%sCS%IOF%slp2ocean = slp2ocean Ice%sCS%IOF%flux_uv_stagger = Ice%flux_uv_stagger - call alloc_fast_ice_avg(Ice%sCS%FIA, sHI, sIG, interp_fluxes) + call alloc_fast_ice_avg(Ice%sCS%FIA, sHI, sIG, interp_fluxes, gas_fluxes) if (Ice%sCS%redo_fast_update) then - call alloc_total_sfc_flux(Ice%sCS%TSF, sHI) - call alloc_total_sfc_flux(Ice%sCS%XSF, sHI) + call alloc_total_sfc_flux(Ice%sCS%TSF, sHI, gas_fluxes) + call alloc_total_sfc_flux(Ice%sCS%XSF, sHI, gas_fluxes) call alloc_ice_rad(Ice%sCS%Rad, sHI, sIG) endif @@ -2223,11 +2234,11 @@ subroutine ice_model_init(Ice, Time_Init, Time, Time_step_fast, Time_step_slow, omit_velocities=.true., omit_tsurf=Eulerian_tsurf) endif if (.not.single_IST) then - call alloc_fast_ice_avg(Ice%fCS%FIA, fHI, Ice%fCS%IG, interp_fluxes) + call alloc_fast_ice_avg(Ice%fCS%FIA, fHI, Ice%fCS%IG, interp_fluxes, gas_fluxes) - call alloc_simple_OSS(Ice%fCS%sOSS, fHI) + call alloc_simple_OSS(Ice%fCS%sOSS, fHI, gas_fields_ocn) endif - call alloc_total_sfc_flux(Ice%fCS%TSF, fHI) + call alloc_total_sfc_flux(Ice%fCS%TSF, fHI, gas_fluxes) Ice%fCS%FIA%atmos_winds = atmos_winds call ice_rad_register_restarts(fGD%mpp_domain, fHI, Ice%fCS%IG, param_file, & diff --git a/src/ice_type.F90 b/src/ice_type.F90 index 0ce7e2b3..95266cda 100644 --- a/src/ice_type.F90 +++ b/src/ice_type.F90 @@ -11,7 +11,8 @@ module ice_type_mod use fms_io_mod, only: save_restart, restore_state, query_initialized use fms_io_mod, only: register_restart_field, restart_file_type use time_manager_mod, only: time_type, time_type_to_real -use coupler_types_mod,only: coupler_2d_bc_type, coupler_3d_bc_type, coupler_type_write_chksums +use coupler_types_mod,only: coupler_1d_bc_type, coupler_2d_bc_type, coupler_3d_bc_type +use coupler_types_mod,only: coupler_type_spawn, coupler_type_write_chksums use SIS_hor_grid, only : SIS_hor_grid_type use ice_grid, only : ice_grid_type @@ -156,13 +157,17 @@ module ice_type_mod ! ice restart files. ! !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! subroutine ice_type_slow_reg_restarts(domain, CatIce, param_file, Ice, & - Ice_restart, restart_file) + Ice_restart, restart_file, gas_fluxes) type(domain2d), intent(in) :: domain integer, intent(in) :: CatIce type(param_file_type), intent(in) :: param_file type(ice_data_type), intent(inout) :: Ice type(restart_file_type), pointer :: Ice_restart character(len=*), intent(in) :: restart_file + type(coupler_1d_bc_type), & + optional, intent(in) :: gas_fluxes !< If present, this type describes the + !! additional gas or other tracer fluxes between the + !! ocean, ice, and atmosphere. ! This subroutine allocates the externally visible ice_data_type's arrays and ! registers the appopriate ones for inclusion in the restart file. @@ -202,6 +207,10 @@ subroutine ice_type_slow_reg_restarts(domain, CatIce, param_file, Ice, & allocate(Ice%mass_berg(isc:iec, jsc:jec)) ; Ice%mass_berg(:,:) = 0.0 endif ; endif + if (present(gas_fluxes)) & + call coupler_type_spawn(gas_fluxes, Ice%ocean_fluxes, (/isc,isc,iec,iec/), & + (/jsc,jsc,jec,jec/), suffix = '_ice') + ! These are used by the ocean model, and need to be in the slow PE restarts. if (associated(Ice_restart)) then idr = register_restart_field(Ice_restart, restart_file, 'flux_u', Ice%flux_u, domain=domain) @@ -229,19 +238,23 @@ subroutine ice_type_slow_reg_restarts(domain, CatIce, param_file, Ice, & end subroutine ice_type_slow_reg_restarts !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! -! ice_type_slow_reg_restarts - allocate the arrays in the ice_data_type ! -! that are predominantly associated with the fast processors, and register ! -! any variables in the ice data type that need to be included in the fast ! -! ice restart files. ! -!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! +!> ice_type_fast_reg_restarts allocates the arrays in the ice_data_type that are +!! predominantly associated with the fast processors, and registers any +!! variables in the ice data type that need to be included in the fast +!! ice restart files. subroutine ice_type_fast_reg_restarts(domain, CatIce, param_file, Ice, & - Ice_restart, restart_file) + Ice_restart, restart_file, gas_fields_ocn) type(domain2d), intent(in) :: domain integer, intent(in) :: CatIce type(param_file_type), intent(in) :: param_file type(ice_data_type), intent(inout) :: Ice type(restart_file_type), pointer :: Ice_restart character(len=*), intent(in) :: restart_file + type(coupler_1d_bc_type), & + optional, intent(in) :: gas_fields_ocn !< If present, this type describes the + !! ocean and surface-ice fields that will participate + !! in the calculation of additional gas or other + !! tracer fluxes. ! This subroutine allocates the externally visible ice_data_type's arrays and ! registers the appopriate ones for inclusion in the restart file. @@ -268,6 +281,10 @@ subroutine ice_type_fast_reg_restarts(domain, CatIce, param_file, Ice, & allocate(Ice%albedo_vis_dif(isc:iec, jsc:jec, km)) ; Ice%albedo_vis_dif(:,:,:) = 0.0 allocate(Ice%albedo_nir_dif(isc:iec, jsc:jec, km)) ; Ice%albedo_nir_dif(:,:,:) = 0.0 + if (present(gas_fields_ocn)) & + call coupler_type_spawn(gas_fields_ocn, Ice%ocean_fields, (/isc,isc,iec,iec/), & + (/jsc,jsc,jec,jec/), (/1, km/), suffix = '_ice') + ! Now register some of these arrays to be read from the restart files. ! These are used by the atmospheric model, and need to be in the fast PE restarts. if (associated(Ice_restart)) then From 725d709142938e045e089ad78297f14b62ff1007 Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Mon, 14 Aug 2017 16:22:53 -0400 Subject: [PATCH 4/4] (*)Save gas fluxes to restarts with concurrent ice Added calls to save any gas fluxes in the fast-ice-average and total_surface_flux types to a restart file with concurrent ice coupling, if the newer version of the coupler, with the optional gas_flux arguments to ice_model_init, is used. All answers are bitwise identical, but there is now a good chance that coupled configurations with gas fluxes and concurrent ice coupling will reproduce across restarts (this is not well tested by the current suite of MOM6-examples test cases). These changes require an updated version of the coupler_type_reform_rwh branch of the FMS code, which should be made available more generally as a patch to FMS in or about September 2017. --- src/SIS_types.F90 | 31 +++++++++++++------------------ 1 file changed, 13 insertions(+), 18 deletions(-) diff --git a/src/SIS_types.F90 b/src/SIS_types.F90 index 7e8e5e0a..5f4b1390 100644 --- a/src/SIS_types.F90 +++ b/src/SIS_types.F90 @@ -5,18 +5,14 @@ !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! module SIS_types -! use mpp_mod, only: mpp_sum, stdout, input_nml_file, PE_here => mpp_pe -! use mpp_domains_mod, only: domain2D, mpp_get_compute_domain, CORNER, EAST, NORTH use mpp_domains_mod, only : domain2D, CORNER, EAST, NORTH, mpp_redistribute -! use mpp_parameter_mod, only: CGRID_NE, BGRID_NE, AGRID -! use fms_mod, only: open_namelist_file, check_nml_error, close_file -! use fms_io_mod, only: save_restart, restore_state, query_initialized use fms_io_mod, only : register_restart_field, restart_file_type use fms_io_mod, only : restore_state, query_initialized use time_manager_mod, only : time_type, time_type_to_real use coupler_types_mod, only : coupler_1d_bc_type, coupler_2d_bc_type, coupler_3d_bc_type -use coupler_types_mod, only : coupler_type_spawn, coupler_type_copy_data -use coupler_types_mod, only : coupler_type_redistribute_data +use coupler_types_mod, only : coupler_type_spawn, coupler_type_initialized +use coupler_types_mod, only : coupler_type_redistribute_data, coupler_type_copy_data +use coupler_types_mod, only : coupler_type_register_restarts use SIS_hor_grid, only : SIS_hor_grid_type use ice_grid, only : ice_grid_type @@ -1906,9 +1902,9 @@ subroutine register_fast_to_slow_restarts(FIA, Rad, TSF, mpp_domain, Ice_restart type(fast_ice_avg_type), pointer :: FIA !< The fast ice model's fast_ice_avg_type type(ice_rad_type), pointer :: Rad !< The fast ice model's ice_rad_type type(total_sfc_flux_type), pointer :: TSF !< The fast ice model's total_sfc_flux_type - type(domain2d), intent(in) :: mpp_domain !< The mpp domain descriptor - type(restart_file_type), intent(inout) :: Ice_restart !< The restart_file_type for these restarts - character(len=*), intent(in) :: restart_file !< The name and path to the restart file + type(domain2d), intent(in) :: mpp_domain !< The mpp domain descriptor + type(restart_file_type), pointer :: Ice_restart !< The restart_file_type for these restarts + character(len=*), intent(in) :: restart_file !< The name and path to the restart file integer :: idr @@ -1999,14 +1995,13 @@ subroutine register_fast_to_slow_restarts(FIA, Rad, TSF, mpp_domain, Ice_restart longname="Total shortwave flux by frequency and angular band", & domain=mpp_domain, mandatory=.false., units="W m-2") - !### These tracer fluxes will need to be dealt with, but because tracer packages can - !### be turned on or off at a model restart, these will have to be kept as - !### coupler_2d_bc_type and coupler_3d_bc_type structures. To enable these, the - !### coupler types will have to be read before the ice model is initialized and - !### then passed in. - -!### call coupler_type_register_restarts(TSF%tr_flux, restart_file, Ice_restart, mpp_domain) -!### call coupler_type_register_restarts(FIA%tr_flux, restart_file, Ice_restart, mpp_domain) + if (coupler_type_initialized(TSF%tr_flux) .and. & + coupler_type_initialized(FIA%tr_flux)) then + call coupler_type_register_restarts(TSF%tr_flux, restart_file, & + Ice_restart, mpp_domain, "TSF_") + call coupler_type_register_restarts(FIA%tr_flux, restart_file, & + Ice_restart, mpp_domain, "FIA_") + endif end subroutine register_fast_to_slow_restarts