From 8dc16287c04a405c6c766c0b8e70713e82533262 Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Sun, 2 Aug 2020 15:54:44 -0400 Subject: [PATCH] Use MOM6 cpu_clock timers in place of mpp_clock Replace calls to mpp_clock routines from mpp_mod with calls to cpu_clock from the module MOM_cpu_clock, to limit the direct use of fms routines and expand the standardized use of the MOM6 infrastructure. Also corrected the granularity of two of the sea-ice clocks, to avoid double counting with clocks at the same level of granularity, and added clocks around calls that to dierctly exchange information between the ice and ocean. All answers are bitwise identical. --- src/SIS_dyn_trans.F90 | 139 ++++++++++++++---------------- src/combined_ice_ocean_driver.F90 | 33 +++---- 2 files changed, 84 insertions(+), 88 deletions(-) diff --git a/src/SIS_dyn_trans.F90 b/src/SIS_dyn_trans.F90 index 13c46a2e..7235f932 100644 --- a/src/SIS_dyn_trans.F90 +++ b/src/SIS_dyn_trans.F90 @@ -16,6 +16,8 @@ module SIS_dyn_trans ! and lateral transport. ! !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! +use MOM_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end +use MOM_cpu_clock, only : CLOCK_COMPONENT, CLOCK_SUBCOMPONENT, CLOCK_LOOP, CLOCK_ROUTINE use MOM_domains, only : pass_var, pass_vector, AGRID, BGRID_NE, CGRID_NE use MOM_error_handler, only : SIS_error=>MOM_error, FATAL, WARNING, SIS_mesg=>MOM_mesg use MOM_error_handler, only : callTree_enter, callTree_leave, callTree_waypoint @@ -28,11 +30,8 @@ module SIS_dyn_trans use MOM_EOS, only : EOS_type, calculate_density_derivs use coupler_types_mod, only: coupler_type_initialized, coupler_type_send_data -use fms_mod, only : clock_flag_default use fms_io_mod, only : restart_file_type use mpp_domains_mod, only : domain2D -use mpp_mod, only : mpp_clock_id, mpp_clock_begin, mpp_clock_end -use mpp_mod, only : CLOCK_COMPONENT, CLOCK_LOOP, CLOCK_ROUTINE use SIS_continuity, only : SIS_continuity_CS, summed_continuity, ice_cover_transport use SIS_debugging, only : chksum, Bchksum, hchksum @@ -373,12 +372,6 @@ subroutine SIS_dynamics_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, U 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 - !### if (CS%merged_cont) then !### This is here for debugging only. Delete it later. - ! call SIS_multi_dyn_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, IG, tracer_CSp, & - ! .true., .true., dt_slow) - ! return - !endif - CS%n_calls = CS%n_calls + 1 IOF%stress_count = 0 @@ -412,7 +405,7 @@ subroutine SIS_dynamics_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, U do nds=1,ndyn_steps - call mpp_clock_begin(iceClock4) + call cpu_clock_begin(iceClock4) ! The code timed by iceClock4 is the non-merged-cont equivalent of convert_IST_to_simple_state. ! Convert the category-resolved ice state into the simplified 2-d ice state. @@ -437,7 +430,7 @@ subroutine SIS_dynamics_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, U if (.not.CS%Warsaw_sum_order) then do j=jsd,jed ; do i=isd,ied ; ice_free(i,j) = max(1.0 - ice_cover(i,j), 0.0) ; enddo ; enddo endif - call mpp_clock_end(iceClock4) + call cpu_clock_end(iceClock4) ! ! Dynamics - update ice velocities. @@ -450,7 +443,7 @@ subroutine SIS_dynamics_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, U ! equation) are not included in the dynamics. All of the thickness categories ! are merged together. - call mpp_clock_begin(iceClock4) + call cpu_clock_begin(iceClock4) ! The code timed by iceClock4 is the non-merged-cont equivalent of SIS_merged_dyn_cont. if (CS%Cgrid_dyn) then @@ -476,7 +469,7 @@ subroutine SIS_dynamics_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, U endif !### Ridging needs to be added with C-grid dynamics. - call mpp_clock_begin(iceClocka) + call cpu_clock_begin(iceClocka) if (CS%do_ridging) rdg_rate(:,:) = 0.0 if (CS%Warsaw_sum_order) then call SIS_C_dynamics(1.0-ice_free(:,:), misp_sum, mi_sum, IST%u_ice_C, IST%v_ice_C, & @@ -487,17 +480,17 @@ subroutine SIS_dynamics_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, U OSS%u_ocn_C, OSS%v_ocn_C, WindStr_x_Cu, WindStr_y_Cv, OSS%sea_lev, & str_x_ice_ocn_Cu, str_y_ice_ocn_Cv, dt_slow_dyn, G, US, CS%SIS_C_dyn_CSp) endif - call mpp_clock_end(iceClocka) + call cpu_clock_end(iceClocka) if (CS%debug) call uvchksum("After ice_dynamics [uv]_ice_C", IST%u_ice_C, IST%v_ice_C, G, scale=US%L_T_to_m_s) - call mpp_clock_begin(iceClockb) + call cpu_clock_begin(iceClockb) call pass_vector(IST%u_ice_C, IST%v_ice_C, G%Domain, stagger=CGRID_NE) call pass_vector(str_x_ice_ocn_Cu, str_y_ice_ocn_Cv, G%Domain, stagger=CGRID_NE) - call mpp_clock_end(iceClockb) + call cpu_clock_end(iceClockb) ! Dynamics diagnostics - call mpp_clock_begin(iceClockc) + call cpu_clock_begin(iceClockc) if (CS%id_fax>0) call post_data(CS%id_fax, WindStr_x_Cu, CS%diag) if (CS%id_fay>0) call post_data(CS%id_fay, WindStr_y_Cv, CS%diag) @@ -513,7 +506,7 @@ subroutine SIS_dynamics_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, U str_x_ice_ocn_Cu, str_y_ice_ocn_Cv, ice_free, ice_cover, G, US) endif - call mpp_clock_end(iceClockc) + call cpu_clock_end(iceClockc) else ! B-grid dynamics. @@ -535,7 +528,7 @@ subroutine SIS_dynamics_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, U scale=US%RZ_T_to_kg_m2s*US%L_T_to_m_s) endif - call mpp_clock_begin(iceClocka) + call cpu_clock_begin(iceClocka) if (CS%do_ridging) rdg_rate(:,:) = 0.0 if (CS%Warsaw_sum_order) then call SIS_B_dynamics(1.0-ice_free(:,:), misp_sum, mi_sum, IST%u_ice_B, IST%v_ice_B, & @@ -548,16 +541,16 @@ subroutine SIS_dynamics_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, U str_x_ice_ocn_B, str_y_ice_ocn_B, CS%do_ridging, & rdg_rate, dt_slow_dyn, G, US, CS%SIS_B_dyn_CSp) endif - call mpp_clock_end(iceClocka) + call cpu_clock_end(iceClocka) if (CS%debug) call Bchksum_pair("After dynamics [uv]_ice_B", IST%u_ice_B, IST%v_ice_B, G, scale=US%L_T_to_m_s) - call mpp_clock_begin(iceClockb) + call cpu_clock_begin(iceClockb) call pass_vector(IST%u_ice_B, IST%v_ice_B, G%Domain, stagger=BGRID_NE) - call mpp_clock_end(iceClockb) + call cpu_clock_end(iceClockb) ! Dynamics diagnostics - call mpp_clock_begin(iceClockc) + call cpu_clock_begin(iceClockc) if ((CS%id_fax>0) .or. (CS%id_fay>0)) then !$OMP parallel do default(shared) private(ps_vel) do J=jsc-1,jec ; do I=isc-1,iec @@ -582,7 +575,7 @@ subroutine SIS_dynamics_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, U call set_ocean_top_stress_B2(IOF, WindStr_x_ocn_B, WindStr_y_ocn_B, & str_x_ice_ocn_B, str_y_ice_ocn_B, ice_free, ice_cover, G, US) endif - call mpp_clock_end(iceClockc) + call cpu_clock_end(iceClockc) ! Convert the velocities to C-grid points for use in transport. do j=jsc,jec ; do I=isc-1,iec @@ -601,12 +594,12 @@ subroutine SIS_dynamics_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, U enddo ; enddo endif - call mpp_clock_begin(iceClock4) + call cpu_clock_end(iceClock4) enddo ! nds=1,ndyn_steps ! Do ice mass transport and related tracer transport. This updates the category-decomposed ice state. - call mpp_clock_begin(iceClock8) + call cpu_clock_begin(iceClock8) ! The code timed by iceClock8 is the non-merged_cont equivalent to complete_IST_transport. if (CS%debug) call uvchksum("Before ice_transport [uv]_ice_C", IST%u_ice_C, IST%v_ice_C, G, scale=US%L_T_to_m_s) call enable_SIS_averaging(dt_slow_dyn_sec, Time_cycle_start + real_to_time(nds*dt_slow_dyn_sec), CS%diag) @@ -623,7 +616,7 @@ subroutine SIS_dynamics_trans(IST, OSS, FIA, IOF, dt_slow, CS, icebergs_CS, G, U call finish_ice_transport(CS%CAS, IST, IST%TrReg, G, US, IG, CS%SIS_transport_CSp) endif endif - call mpp_clock_end(iceClock8) + call cpu_clock_end(iceClock8) endif ! (.not.CS%merged_cont) @@ -745,7 +738,7 @@ subroutine complete_IST_transport(DS2d, CAS, IST, dt_adv_cycle, G, US, IG, CS) isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB - call mpp_clock_begin(iceClock8) + call cpu_clock_begin(iceClock8) ! Do the transport of mass and tracers by category and vertical layer. call ice_cat_transport(CS%CAS, IST%TrReg, dt_adv_cycle, DS2d%nts, G, US, IG, & CS%SIS_transport_CSp, mca_tot=DS2d%mca_step(:,:,0:DS2d%nts), & @@ -772,7 +765,7 @@ subroutine complete_IST_transport(DS2d, CAS, IST, dt_adv_cycle, G, US, IG, CS) endif IST%valid_IST = .true. - call mpp_clock_end(iceClock8) + call cpu_clock_end(iceClock8) end subroutine complete_IST_transport @@ -806,7 +799,7 @@ subroutine ice_state_cleanup(IST, OSS, IOF, dt_slow, G, US, IG, CS, tracer_CSp) endif ! Calculate and output various diagnostics of the ice state. - call mpp_clock_begin(iceClock9) + call cpu_clock_begin(iceClock9) call enable_SIS_averaging(US%T_to_s*dt_slow, CS%Time, CS%diag) call post_ice_state_diagnostics(CS%IDs, IST, OSS, IOF, dt_slow, CS%Time, G, US, IG, CS%diag) @@ -824,7 +817,7 @@ subroutine ice_state_cleanup(IST, OSS, IOF, dt_slow, G, US, IG, CS, tracer_CSp) call write_ice_statistics(IST, CS%Time, CS%n_calls, G, US, IG, CS%sum_output_CSp) endif - call mpp_clock_end(iceClock9) + call cpu_clock_end(iceClock9) end subroutine ice_state_cleanup @@ -851,7 +844,7 @@ subroutine convert_IST_to_simple_state(IST, DS2d, CAS, G, US, IG, CS) isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB - call mpp_clock_begin(iceClock4) + call cpu_clock_begin(iceClock4) DS2d%mca_step(:,:,0) = 0.0 ; DS2d%mi_sum(:,:) = 0.0 ; DS2d%ice_cover(:,:) = 0.0 !$OMP parallel do default(shared) @@ -883,7 +876,7 @@ subroutine convert_IST_to_simple_state(IST, DS2d, CAS, G, US, IG, CS) IST%valid_IST = .false. - call mpp_clock_end(iceClock4) + call cpu_clock_end(iceClock4) end subroutine convert_IST_to_simple_state @@ -961,7 +954,7 @@ subroutine SIS_merged_dyn_cont(OSS, FIA, IOF, DS2d, dt_cycle, Time_start, G, US, continuing_call = .false. ; if (present(end_call)) continuing_call = .not.end_call do nds=1,ndyn_steps - call mpp_clock_begin(iceClock4) + call cpu_clock_begin(iceClock4) call enable_SIS_averaging(dt_slow_dyn_sec, Time_start + real_to_time(nds*dt_slow_dyn_sec), CS%diag) do j=jsd,jed ; do i=isd,ied ; ice_free(i,j) = max(1.0 - DS2d%ice_cover(i,j), 0.0) ; enddo ; enddo @@ -991,7 +984,7 @@ subroutine SIS_merged_dyn_cont(OSS, FIA, IOF, DS2d, dt_cycle, Time_start, G, US, ! call hchksum_pair("WindStr_[xy]_A before SIS_C_dynamics", WindStr_x_A, WindStr_y_A, G, halos=1) endif - call mpp_clock_begin(iceClocka) + call cpu_clock_begin(iceClocka) !### Ridging needs to be added with C-grid dynamics. if (CS%do_ridging) rdg_rate(:,:) = 0.0 call SIS_C_dynamics(DS2d%ice_cover, DS2d%mca_step(:,:,DS2d%nts), DS2d%mi_sum, & @@ -999,17 +992,17 @@ subroutine SIS_merged_dyn_cont(OSS, FIA, IOF, DS2d, dt_cycle, Time_start, G, US, OSS%u_ocn_C, OSS%v_ocn_C, WindStr_x_Cu, WindStr_y_Cv, OSS%sea_lev, & str_x_ice_ocn_Cu, str_y_ice_ocn_Cv, dt_slow_dyn, G, US, CS%SIS_C_dyn_CSp) - call mpp_clock_end(iceClocka) + call cpu_clock_end(iceClocka) if (CS%debug) call uvchksum("After ice_dynamics [uv]_ice_C", DS2d%u_ice_C, DS2d%v_ice_C, G, scale=US%L_T_to_m_s) - call mpp_clock_begin(iceClockb) + call cpu_clock_begin(iceClockb) call pass_vector(DS2d%u_ice_C, DS2d%v_ice_C, G%Domain, stagger=CGRID_NE) call pass_vector(str_x_ice_ocn_Cu, str_y_ice_ocn_Cv, G%Domain, stagger=CGRID_NE) - call mpp_clock_end(iceClockb) + call cpu_clock_end(iceClockb) ! Dynamics diagnostics - call mpp_clock_begin(iceClockc) + call cpu_clock_begin(iceClockc) if (CS%id_fax>0) call post_data(CS%id_fax, WindStr_x_Cu, CS%diag) if (CS%id_fay>0) call post_data(CS%id_fay, WindStr_y_Cv, CS%diag) if (CS%debug) call uvchksum("Before set_ocean_top_stress_Cgrid [uv]_ice_C", & @@ -1018,7 +1011,7 @@ subroutine SIS_merged_dyn_cont(OSS, FIA, IOF, DS2d, dt_cycle, Time_start, G, US, ! Store all mechanical ocean forcing. call set_ocean_top_stress_C2(IOF, WindStr_x_ocn_Cu, WindStr_y_ocn_Cv, & str_x_ice_ocn_Cu, str_y_ice_ocn_Cv, ice_free, DS2d%ice_cover, G, US) - call mpp_clock_end(iceClockc) + call cpu_clock_end(iceClockc) else ! B-grid dynamics. @@ -1041,23 +1034,23 @@ subroutine SIS_merged_dyn_cont(OSS, FIA, IOF, DS2d, dt_cycle, Time_start, G, US, scale=US%RZ_T_to_kg_m2s*US%L_T_to_m_s) endif - call mpp_clock_begin(iceClocka) + call cpu_clock_begin(iceClocka) if (CS%do_ridging) rdg_rate(:,:) = 0.0 call SIS_B_dynamics(DS2d%ice_cover, DS2d%mca_step(:,:,DS2d%nts), DS2d%mi_sum, & DS2d%u_ice_B, DS2d%v_ice_B, & OSS%u_ocn_B, OSS%v_ocn_B, WindStr_x_B, WindStr_y_B, OSS%sea_lev, & str_x_ice_ocn_B, str_y_ice_ocn_B, CS%do_ridging, & rdg_rate, dt_slow_dyn, G, US, CS%SIS_B_dyn_CSp) - call mpp_clock_end(iceClocka) + call cpu_clock_end(iceClocka) if (CS%debug) call Bchksum_pair("After dynamics [uv]_ice_B", DS2d%u_ice_B, DS2d%v_ice_B, G, scale=US%L_T_to_m_s) - call mpp_clock_begin(iceClockb) + call cpu_clock_begin(iceClockb) call pass_vector(DS2d%u_ice_B, DS2d%v_ice_B, G%Domain, stagger=BGRID_NE) - call mpp_clock_end(iceClockb) + call cpu_clock_end(iceClockb) ! Dynamics diagnostics - call mpp_clock_begin(iceClockc) + call cpu_clock_begin(iceClockc) if ((CS%id_fax>0) .or. (CS%id_fay>0)) then !$OMP parallel do default(shared) private(ps_vel) do J=jsc-1,jec ; do I=isc-1,iec @@ -1077,7 +1070,7 @@ subroutine SIS_merged_dyn_cont(OSS, FIA, IOF, DS2d, dt_cycle, Time_start, G, US, ! Store all mechanical ocean forcing. call set_ocean_top_stress_B2(IOF, WindStr_x_ocn_B, WindStr_y_ocn_B, & str_x_ice_ocn_B, str_y_ice_ocn_B, ice_free, DS2d%ice_cover, G, US) - call mpp_clock_end(iceClockc) + call cpu_clock_end(iceClockc) ! Convert the velocities to C-grid points for use in transport. do j=jsc,jec ; do I=isc-1,iec @@ -1120,7 +1113,7 @@ subroutine SIS_merged_dyn_cont(OSS, FIA, IOF, DS2d, dt_cycle, Time_start, G, US, endif enddo DS2d%nts = DS2d%nts + CS%adv_substeps - call mpp_clock_end(iceClock4) + call cpu_clock_end(iceClock4) enddo ! nds=1,ndyn_steps @@ -1190,14 +1183,14 @@ subroutine slab_ice_dyn_trans(IST, OSS, FIA, IOF, dt_slow, CS, G, US, IG, tracer call enable_SIS_averaging(dt_slow_dyn_sec, CS%Time - real_to_time((ndyn_steps-nds)*dt_slow_dyn_sec), CS%diag) - call mpp_clock_begin(iceClock4) + call cpu_clock_begin(iceClock4) !$OMP parallel do default(shared) do j=jsd,jed ; do i=isd,ied mi_sum(i,j) = IST%mH_ice(i,j,1) * IST%part_size(i,j,1) misp_sum(i,j) = mi_sum(i,j) + IST%part_size(i,j,1) * & (IST%mH_snow(i,j,1) + IST%mH_pond(i,j,1)) enddo ; enddo - call mpp_clock_end(iceClock4) + call cpu_clock_end(iceClock4) ! ! Dynamics - update ice velocities. @@ -1209,7 +1202,7 @@ subroutine slab_ice_dyn_trans(IST, OSS, FIA, IOF, dt_slow, CS, G, US, IG, tracer ! are merged together. if (CS%Cgrid_dyn) then - call mpp_clock_begin(iceClock4) + call cpu_clock_begin(iceClock4) ! Correct the wind stresses for changes in the fractional ice-coverage and set ! the wind stresses on the ice and the open ocean for a C-grid staggering. ! This block of code must be executed if ice_cover and ice_free or the various wind @@ -1230,20 +1223,20 @@ subroutine slab_ice_dyn_trans(IST, OSS, FIA, IOF, dt_slow, CS, G, US, IG, tracer ! call hchksum_pair("WindStr_[xy]_A before SIS_C_dynamics", WindStr_x_A, WindStr_y_A, G, halos=1) endif - call mpp_clock_begin(iceClocka) + call cpu_clock_begin(iceClocka) call slab_ice_dynamics(IST%u_ice_C, IST%v_ice_C, OSS%u_ocn_C, OSS%v_ocn_C, & WindStr_x_Cu, WindStr_y_Cv, str_x_ice_ocn_Cu, str_y_ice_ocn_Cv) - call mpp_clock_end(iceClocka) + call cpu_clock_end(iceClocka) if (CS%debug) call uvchksum("After ice_dynamics [uv]_ice_C", IST%u_ice_C, IST%v_ice_C, G, scale=US%L_T_to_m_s) - call mpp_clock_begin(iceClockb) + call cpu_clock_begin(iceClockb) call pass_vector(IST%u_ice_C, IST%v_ice_C, G%Domain, stagger=CGRID_NE) call pass_vector(str_x_ice_ocn_Cu, str_y_ice_ocn_Cv, G%Domain, stagger=CGRID_NE) - call mpp_clock_end(iceClockb) + call cpu_clock_end(iceClockb) ! Dynamics diagnostics - call mpp_clock_begin(iceClockc) + call cpu_clock_begin(iceClockc) if (CS%id_fax>0) call post_data(CS%id_fax, WindStr_x_Cu, CS%diag) if (CS%id_fay>0) call post_data(CS%id_fay, WindStr_y_Cv, CS%diag) @@ -1253,13 +1246,13 @@ subroutine slab_ice_dyn_trans(IST, OSS, FIA, IOF, dt_slow, CS, G, US, IG, tracer ! Store all mechanical ocean forcing. call set_ocean_top_stress_C2(IOF, WindStr_x_ocn_Cu, WindStr_y_ocn_Cv, str_x_ice_ocn_Cu, str_y_ice_ocn_Cv, & IST%part_size(:,:,0), IST%part_size(:,:,1), G, US) - call mpp_clock_end(iceClockc) + call cpu_clock_end(iceClockc) - call mpp_clock_end(iceClock4) + call cpu_clock_end(iceClock4) else ! B-grid dynamics. - call mpp_clock_begin(iceClock4) + call cpu_clock_begin(iceClock4) ! Correct the wind stresses for changes in the fractional ice-coverage and set ! the wind stresses on the ice and the open ocean for a C-grid staggering. ! This block of code must be executed if ice_cover and ice_free or the various wind @@ -1279,19 +1272,19 @@ subroutine slab_ice_dyn_trans(IST, OSS, FIA, IOF, dt_slow, CS, G, US, IG, tracer scale=US%RZ_T_to_kg_m2s*US%L_T_to_m_s) endif - call mpp_clock_begin(iceClocka) + call cpu_clock_begin(iceClocka) call slab_ice_dynamics(IST%u_ice_B, IST%v_ice_B, OSS%u_ocn_B, OSS%v_ocn_B, & WindStr_x_B, WindStr_y_B, str_x_ice_ocn_B, str_y_ice_ocn_B) - call mpp_clock_end(iceClocka) + call cpu_clock_end(iceClocka) if (CS%debug) call Bchksum_pair("After dynamics [uv]_ice_B", IST%u_ice_B, IST%v_ice_B, G, scale=US%L_T_to_m_s) - call mpp_clock_begin(iceClockb) + call cpu_clock_begin(iceClockb) call pass_vector(IST%u_ice_B, IST%v_ice_B, G%Domain, stagger=BGRID_NE) - call mpp_clock_end(iceClockb) + call cpu_clock_end(iceClockb) ! Dynamics diagnostics - call mpp_clock_begin(iceClockc) + call cpu_clock_begin(iceClockc) if ((CS%id_fax>0) .or. (CS%id_fay>0)) then !$OMP parallel do default(shared) private(ps_vel) do J=jsc-1,jec ; do I=isc-1,iec @@ -1311,7 +1304,7 @@ subroutine slab_ice_dyn_trans(IST, OSS, FIA, IOF, dt_slow, CS, G, US, IG, tracer ! Store all mechanical ocean forcing. call set_ocean_top_stress_B2(IOF, WindStr_x_ocn_B, WindStr_y_ocn_B, str_x_ice_ocn_B, str_y_ice_ocn_B, & IST%part_size(:,:,0), IST%part_size(:,:,1), G, US) - call mpp_clock_end(iceClockc) + call cpu_clock_end(iceClockc) ! Convert the B-grid velocities to C-grid points for transport. if (CS%debug) call Bchksum_pair("Before ice_transport [uv]_ice_B", IST%u_ice_B, IST%v_ice_B, G, & @@ -1323,18 +1316,18 @@ subroutine slab_ice_dyn_trans(IST, OSS, FIA, IOF, dt_slow, CS, G, US, IG, tracer IST%v_ice_C(i,J) = 0.5 * ( IST%v_ice_B(I-1,J) + IST%v_ice_B(I,J) ) enddo ; enddo - call mpp_clock_end(iceClock4) + call cpu_clock_end(iceClock4) endif ! End of B-grid dynamics ! Do ice mass transport and related tracer transport. This updates the category-decomposed ice state. - call mpp_clock_begin(iceClock8) + call cpu_clock_begin(iceClock8) if (CS%debug) call uvchksum("Before ice_transport [uv]_ice_C", IST%u_ice_C, IST%v_ice_C, G, scale=US%L_T_to_m_s) call enable_SIS_averaging(dt_slow_dyn_sec, CS%Time - real_to_time((ndyn_steps-nds)*dt_slow_dyn_sec), CS%diag) call slab_ice_advect(IST%u_ice_C, IST%v_ice_C, IST%mH_ice(:,:,1), 4.0*US%kg_m3_to_R*US%m_to_Z, & dt_slow_dyn, G, US, IST%part_size(:,:,1), nsteps=CS%adv_substeps) - call mpp_clock_end(iceClock8) + call cpu_clock_end(iceClock8) if (CS%column_check) & call write_ice_statistics(IST, CS%Time, CS%n_calls, G, US, IG, CS%sum_output_CSp, & @@ -2345,12 +2338,12 @@ subroutine SIS_dyn_trans_init(Time, G, US, IG, param_file, diag, CS, output_dir, call register_ice_state_diagnostics(Time, IG, US, param_file, diag, CS%IDs) - iceClock4 = mpp_clock_id( ' Ice: slow: dynamics', flags=clock_flag_default, grain=CLOCK_LOOP ) - iceClocka = mpp_clock_id( ' slow: ice_dynamics', flags=clock_flag_default, grain=CLOCK_LOOP ) - iceClockb = mpp_clock_id( ' slow: comm/cut check ', flags=clock_flag_default, grain=CLOCK_LOOP ) - iceClockc = mpp_clock_id( ' slow: diags', flags=clock_flag_default, grain=CLOCK_LOOP ) - iceClock8 = mpp_clock_id( ' Ice: slow: transport', flags=clock_flag_default, grain=CLOCK_LOOP ) - iceClock9 = mpp_clock_id( ' Ice: slow: thermodyn diags', flags=clock_flag_default, grain=CLOCK_LOOP ) + iceClock4 = cpu_clock_id( ' Ice: slow: dynamics', grain=CLOCK_ROUTINE ) + iceClocka = cpu_clock_id( ' slow: ice_dynamics', grain=CLOCK_LOOP ) + iceClockb = cpu_clock_id( ' slow: comm/cut check ', grain=CLOCK_LOOP ) + iceClockc = cpu_clock_id( ' slow: diags', grain=CLOCK_LOOP ) + iceClock8 = cpu_clock_id( ' Ice: slow: transport', grain=CLOCK_ROUTINE ) + iceClock9 = cpu_clock_id( ' Ice: slow: thermodyn diags', grain=CLOCK_ROUTINE ) call callTree_leave("SIS_dyn_trans_init()") diff --git a/src/combined_ice_ocean_driver.F90 b/src/combined_ice_ocean_driver.F90 index dfccdac2..53801a07 100644 --- a/src/combined_ice_ocean_driver.F90 +++ b/src/combined_ice_ocean_driver.F90 @@ -6,12 +6,11 @@ module combined_ice_ocean_driver !----------------------------------------------------------------------- ! -! This module provides a common interface for jointly stepping SIS2 and -! MOM6, and will evolve as a platform for tightly integrating the ocean -! and sea ice models. +! This module provides a common interface for jointly stepping SIS2 and MOM6, and +! will evolve as a platform for tightly integrating the ocean and sea ice models. -use MOM_error_handler, only : MOM_error, FATAL, WARNING, is_root_pe -use MOM_error_handler, only : callTree_enter, callTree_leave +use MOM_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end, CLOCK_COMPONENT +use MOM_error_handler, only : MOM_error, FATAL, WARNING, callTree_enter, callTree_leave use MOM_file_parser, only : param_file_type, open_param_file, close_param_file use MOM_file_parser, only : read_param, get_param, log_param, log_version use MOM_io, only : file_exists, close_file, slasher, ensembler @@ -19,10 +18,10 @@ module combined_ice_ocean_driver use MOM_time_manager, only : time_type, time_type_to_real, real_to_time_type use MOM_time_manager, only : operator(+), operator(-), operator(>) -use ice_model_mod, only : ice_data_type, ice_model_end -use ice_model_mod, only : update_ice_slow_thermo, update_ice_dynamics_trans -use ocean_model_mod, only : update_ocean_model, ocean_model_end -use ocean_model_mod, only : ocean_public_type, ocean_state_type, ice_ocean_boundary_type +use ice_model_mod, only : ice_data_type, ice_model_end +use ice_model_mod, only : update_ice_slow_thermo, update_ice_dynamics_trans +use ocean_model_mod, only : update_ocean_model, ocean_model_end +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 @@ -48,6 +47,10 @@ module combined_ice_ocean_driver !! The default is -1. end type ice_ocean_driver_type +!>@{ CPU time clock IDs +integer :: fluxIceOceanClock +!!@} + contains !======================================================================= @@ -128,6 +131,8 @@ subroutine ice_ocean_driver_init(CS, Time_init, Time_in) call close_param_file(param_file, component="CIOD") CS%CS_is_initialized = .true. + fluxIceOceanClock = cpu_clock_id('Direct Ice Ocean Exchange', grain=CLOCK_COMPONENT ) + call callTree_leave("ice_ocean_driver_init(") end subroutine ice_ocean_driver_init @@ -183,7 +188,6 @@ subroutine update_slow_ice_and_ocean(CS, Ice, Ocn, Ocean_sfc, IOB, & call MOM_error(FATAL, "update_slow_ice_and_ocean can only be used if the "//& "ocean and slow ice layouts and domain sizes are identical.") - !### Add clocks of the various calls. if (CS%intersperse_ice_ocn) then ! First step the ice, then ocean thermodynamics. call update_ice_slow_thermo(Ice) @@ -221,15 +225,11 @@ subroutine update_slow_ice_and_ocean(CS, Ice, Ocn, Ocean_sfc, IOB, & call update_ice_dynamics_trans(Ice) -! call mpp_clock_begin(fluxIceOceanClock) call direct_flux_ice_to_IOB(time_start_update, Ice, IOB) -! call mpp_clock_end(fluxIceOceanClock) if (CS%single_MOM_call) then call update_ocean_model(IOB, Ocn, Ocean_sfc, time_start_update, coupling_time_step ) else - !### This line could be a temporary measure to avoid using newer arguments to update_ocean_model. - ! call update_ocean_model(IOB, Ocn, Ocean_sfc, time_start_update, coupling_time_step ) call update_ocean_model(IOB, Ocn, Ocean_sfc, time_start_update, coupling_time_step, & update_dyn=.true., update_thermo=.false., & start_cycle=.true., end_cycle=.false., cycle_length=dt_coupling) @@ -279,6 +279,8 @@ subroutine direct_flux_ice_to_IOB(Time, Ice, IOB, do_thermo) integer :: i, j, is, ie, js, je, i_off, j_off, n, m logical :: used, do_therm + call cpu_clock_begin(fluxIceOceanClock) + do_therm = .true. ; if (present(do_thermo)) do_therm = do_thermo ! Do a direct copy of the ice surface fluxes into the Ice-ocean-boundary type. @@ -336,7 +338,6 @@ subroutine direct_flux_ice_to_IOB(Time, Ice, IOB, do_thermo) call data_override('OCN', 'mi', IOB%mi , Time) if (ASSOCIATED(IOB%stress_mag) ) & call data_override('OCN', 'stress_mag', IOB%stress_mag, Time ) - !Are these if statements needed, or does data_override routine check if variable is associated? if (ASSOCIATED(IOB%ustar_berg)) & call data_override('OCN', 'ustar_berg', IOB%ustar_berg, Time) if (ASSOCIATED(IOB%area_berg) ) & @@ -351,6 +352,8 @@ subroutine direct_flux_ice_to_IOB(Time, Ice, IOB, do_thermo) call coupler_type_send_data(IOB%fluxes, Time ) endif + call cpu_clock_end(fluxIceOceanClock) + end subroutine direct_flux_ice_to_IOB !=======================================================================