!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! ! ice_type_mod - maintains the sea ice data, reads/writes restarts, reads the ! ! namelist and initializes diagnostics. - Mike Winton ! !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! module ice_type_mod use mpp_mod, only: mpp_pe, mpp_root_pe, mpp_sum, mpp_clock_id, CLOCK_COMPONENT, & CLOCK_LOOP, CLOCK_ROUTINE, stdout,input_nml_file use mpp_domains_mod, only: domain2D, mpp_update_domains, CORNER, BGRID_NE use fms_mod, only: file_exist, open_namelist_file, check_nml_error, write_version_number,& read_data, close_file, field_exist, & stderr, stdlog, error_mesg, FATAL, WARNING, NOTE, clock_flag_default use fms_io_mod, only: save_restart, restore_state, query_initialized, & register_restart_field, restart_file_type, set_domain, nullify_domain, & parse_mask_table use diag_manager_mod, only: diag_axis_init, register_diag_field, & register_static_field, send_data use time_manager_mod, only: time_type, get_time use coupler_types_mod,only: coupler_2d_bc_type, coupler_3d_bc_type use constants_mod, only: Tfreeze, radius, pi use ice_grid_mod, only: set_ice_grid, t_to_uv, ice_grid_end use ice_grid_mod, only: Domain, isc, iec, jsc, jec, isd, ied, jsd, jed, im, jm, km use ice_grid_mod, only: geo_lon, geo_lat, cell_area, sin_rot, cos_rot, wett, xb1d, yb1d use ice_grid_mod, only: geo_lonv_ib, geo_latv_ib use ice_grid_mod, only: x_cyclic, tripolar_grid, dtn, dte, wetv use ice_grid_mod, only: reproduce_siena_201303 use ice_thm_mod, only: ice_thm_param, DI, DS, e_to_melt use ice_dyn_mod, only: ice_dyn_param use constants_mod, only: LI => hlf ! latent heat of fusion - 334e3 J/(kg-ice) use ice_bergs, only: icebergs_init, icebergs_end, icebergs, icebergs_stock_pe use ice_bergs, only: icebergs_save_restart use astronomy_mod, only: astronomy_init, astronomy_end implicit none private public :: ice_data_type, ice_model_init, ice_model_end, ice_stock_pe, kmelt, & mom_rough_ice, heat_rough_ice, atmos_winds, hlim, slab_ice, & spec_ice, verbose, ice_bulk_salin, do_ice_restore, do_ice_limit, & max_ice_limit, ice_restore_timescale, do_init, h2o, heat, slp2ocean,& cm2_bugs, conservation_check, do_icebergs, ice_model_restart, & add_diurnal_sw, channel_viscosity, smag_ocn, ssh_gravity, & chan_cfl_limit, ice_data_type_chksum public :: do_sun_angle_for_alb public :: id_cn, id_hi, id_hs, id_t1, id_t2, id_ts public :: id_mi, id_sh, id_lh, id_sw, id_lw, id_snofl, id_rain, id_runoff, & id_calving, id_runoff_hflx, id_calving_hflx, & id_evap, id_saltf, id_tmelt, id_bmelt, id_bheat, id_e2m, & id_frazil, id_alb, id_xprt, id_lsrc, id_lsnk, id_bsnk, id_strna, & id_sigi, id_sigii, id_stren, id_ui, id_vi, id_fax, id_fay, id_fix, & id_fiy, id_fcx, id_fcy, id_fwx, id_fwy, id_swdn, id_lwdn, id_sn2ic, & id_slp, id_ext, id_sst, id_sss, id_ssh, id_uo, id_vo, id_ta, id_obi,& id_qfres, id_qflim, id_ix_trans, id_iy_trans, & id_sw_vis, id_sw_dir, id_sw_dif, id_sw_vis_dir, id_sw_vis_dif, & id_sw_nir_dir, id_sw_nir_dif, id_mib, id_ustar, id_vstar, id_vocean,& id_uocean, id_vchan, id_uchan, id_wnd public :: id_alb_vis_dir, id_alb_vis_dif,id_alb_nir_dir, id_alb_nir_dif public :: iceClock,iceClock1,iceClock2,iceClock3,iceClock4,iceClock5,iceClock6,iceClock7,iceClock8,iceClock9 public :: iceClocka,iceClockb,iceClockc character(len=128) :: version = '$Id: ice_type.F90,v 20.0 2013/12/13 23:28:32 fms Exp $' character(len=128) :: tagname = '$Name: tikal $' !---- id for diagnositics ------------------- integer :: id_xb, id_xt, id_yb, id_yt, id_ct, id_xv, id_yv integer :: id_cn, id_hi, id_hs, id_t1, id_t2, id_ts integer :: id_mi, id_sh, id_lh, id_sw, id_lw, id_snofl, id_rain integer :: id_runoff, id_calving, id_runoff_hflx, id_calving_hflx integer :: id_evap, id_saltf, id_tmelt, id_bmelt, id_bheat, id_e2m integer :: id_frazil, id_alb, id_xprt, id_lsrc, id_lsnk, id_bsnk, id_strna integer :: id_sigi, id_sigii, id_stren, id_ui, id_vi, id_fax, id_fay, id_fix integer :: id_fiy, id_fcx, id_fcy, id_fwx, id_fwy, id_swdn, id_lwdn, id_sn2ic integer :: id_slp, id_ext, id_sst, id_sss, id_ssh, id_uo, id_vo, id_ta, id_obi integer :: id_qfres, id_qflim, id_ix_trans, id_iy_trans integer :: id_sw_vis, id_sw_dir, id_sw_dif, id_sw_vis_dir, id_sw_vis_dif integer :: id_sw_nir_dir, id_sw_nir_dif, id_mib, id_ustar, id_vstar integer :: id_vocean, id_uocean, id_vchan, id_uchan integer :: id_alb_vis_dir, id_alb_vis_dif, id_alb_nir_dir, id_alb_nir_dif integer :: id_wnd !--- namelist interface -------------- real :: mom_rough_ice = 1.0e-4 ! momentum same, cd10=(von_k/ln(10/z0))^2 real :: heat_rough_ice = 1.0e-4 ! heat roughness length real :: kmelt = 6e-5*4e6 ! ocean/ice heat flux constant real :: ks = 0.31 ! snow conductivity (W/mK) real :: alb_sno = 0.85 ! snow albedo (less if melting) real :: alb_ice = 0.5826 ! ice albedo (less if melting) real :: pen_ice = 0.3 ! part unreflected solar penetrates ice real :: opt_dep_ice = 0.67 ! ice optical depth real :: t_range_melt = 1.0 ! melt albedos scaled in over T range real :: ice_bulk_salin = 0.0 ! ice bulk salinity (for ocean salt flux) real :: p0 = 2.75e4 ! ice strength parameter real :: c0 = 20.0 ! another ice strength parameter real :: cdw = 3.24e-3 ! water/ice drag coefficient real :: wd_turn = 25.0 ! water/ice drag turning angle real :: h_lo_lim = 0.0 ! min ice thickness for temp. calc. integer :: nsteps_dyn = 432 ! dynamics steps per slow timestep integer :: nsteps_adv = 8 ! advection steps per slow timestep integer :: num_part = 6 ! number of ice grid partitions ! partition 1 is open water ! partitions 2 to num_part-1 are ! thickness limited ice categories ! partition num_part is unlimited ice logical :: atmos_winds = .true. ! wind stress from atmosphere model over t points and has wrong sign logical :: slab_ice = .false. ! do old-style GFDL slab ice? logical :: spec_ice = .false. ! old-style GFDL slab ice with SST, ice thickness and conc. from data logical :: do_ice_restore = .false. ! restore sea-ice toward climatology logical :: do_ice_limit = .false. ! limit sea ice to max_ice_limit real :: max_ice_limit = 4.0 ! maximum sea ice height(m), ! if do_ice_limit is true ! TK: default chosen based on observed ! ice thickness data used by climate ! group, which range up to 7.31 m real :: ice_restore_timescale = 5.0 ! time scale for restoring ice (days) logical :: conservation_check = .true. ! check for heat and h2o conservation logical :: slp2ocean = .false.! apply sea level pressure to ocean surface logical :: cm2_bugs = .false.! keep cm2 bugs for reproducibility logical :: verbose = .false.! control printing message, will slow model down when turn true logical :: do_icebergs = .false.! call iceberg code to modify calving field logical :: add_diurnal_sw = .false.! apply an additional diurnal cycle to shortwave radiation real :: channel_viscosity = 0. ! viscosity used in one-cell wide channels to parameterize transport (m^2/s) real :: smag_ocn = 0.15 ! Smagorinksy coefficient for viscosity (dimensionless) real :: ssh_gravity = 9.81 ! Gravity parameter used in channel viscosity parameterization (m/s^2) real :: chan_cfl_limit = 0.25 ! CFL limit for channel viscosity parameterization (dimensionless) logical :: do_sun_angle_for_alb = .false.! find the sun angle for ocean albed in the frame of the ice model integer :: layout(2) = (/0, 0/) integer :: io_layout(2) = (/0, 0/) ! mask_table contains information for masking domain ( n_mask, layout and mask_list). ! A text file to specify n_mask, layout and mask_list to reduce number of processor ! usage by masking out some domain regions which contain all land points. ! The default file name of mask_table is "INPUT/ice_mask_table". Please note that ! the file name must begin with "INPUT/". The first ! line of mask_table will be number of region to be masked out. The second line ! of the mask_table will be the layout of the model. User need to set ice_model_nml ! variable layout to be the same as the second line of the mask table. ! The following n_mask line will be the position of the processor to be masked out. ! The mask_table could be created by tools check_mask. ! For example the mask_table will be as following if n_mask=2, layout=4,6 and ! the processor (1,2) and (3,6) will be masked out. ! 2 ! 4,6 ! 1,2 ! 3,6 character(len=128) :: mask_table = "INPUT/ice_mask_table" namelist /ice_model_nml/ mom_rough_ice, heat_rough_ice, p0, c0, cdw, wd_turn, & kmelt, alb_sno, alb_ice, pen_ice, opt_dep_ice, & nsteps_dyn, nsteps_adv, num_part, atmos_winds, & slab_ice, spec_ice, ice_bulk_salin, layout, & do_ice_restore, do_ice_limit, max_ice_limit, & ice_restore_timescale, slp2ocean, conservation_check, & t_range_melt, cm2_bugs, ks, h_lo_lim, verbose, & do_icebergs, add_diurnal_sw, io_layout, channel_viscosity,& smag_ocn, ssh_gravity, chan_cfl_limit, do_sun_angle_for_alb, & mask_table, reproduce_siena_201303 logical :: do_init = .false. real :: hlim(8) = (/ 0.0, 0.1, 0.3, 0.7, 1.1, 1.5, 2.0, 2.5 /) ! thickness limits 1...num_part-1 real :: h2o(4), heat(4) ! for conservation analysis ! 1 - initial ice h2o/heat content ! 2 - h2o/heat flux down at top of ice ! 3 - h2o/heat flux down at bottom of ice ! 4 - final ice h2o/heat content !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! ! This structure contains the ice model data (some used by calling routines); ! ! the third index is partition (1 is open water; 2 is ice cover) ! !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! type ice_data_type type(domain2D) :: Domain type (time_type) :: Time_Init, Time type (time_type) :: Time_step_fast, Time_step_slow integer :: avg_kount logical :: pe integer, pointer, dimension(:) :: pelist =>NULL() ! Used for flux-exchange. logical, pointer, dimension(:,:) :: mask =>NULL() ! where ice can be logical, pointer, dimension(:,:,:) :: ice_mask =>NULL() ! where ice actually is real, pointer, dimension(:,:,:) :: part_size =>NULL() real, pointer, dimension(:,:,:) :: part_size_uv =>NULL() real, pointer, dimension(:,:,:) :: albedo =>NULL() real, pointer, dimension(:,:,:) :: albedo_vis_dir =>NULL() real, pointer, dimension(:,:,:) :: albedo_nir_dir =>NULL() real, pointer, dimension(:,:,:) :: albedo_vis_dif =>NULL() real, pointer, dimension(:,:,:) :: albedo_nir_dif =>NULL() real, pointer, dimension(:,:,:) :: rough_mom =>NULL() real, pointer, dimension(:,:,:) :: rough_heat =>NULL() real, pointer, dimension(:,:,:) :: rough_moist =>NULL() real, pointer, dimension(:,:,:) :: t_surf =>NULL() real, pointer, dimension(:,:,:) :: u_surf =>NULL() real, pointer, dimension(:,:,:) :: v_surf =>NULL() real, pointer, dimension(:,:) :: sea_lev =>NULL() real, pointer, dimension(:,:) :: s_surf =>NULL() real, pointer, dimension(:,:) :: u_ocn =>NULL() real, pointer, dimension(:,:) :: v_ocn =>NULL() real, pointer, dimension(:,:,:) :: flux_u_top =>NULL() real, pointer, dimension(:,:,:) :: flux_v_top =>NULL() real, pointer, dimension(:,:,:) :: flux_u_top_bgrid =>NULL() real, pointer, dimension(:,:,:) :: flux_v_top_bgrid =>NULL() real, pointer, dimension(:,:,:) :: flux_t_top =>NULL() real, pointer, dimension(:,:,:) :: flux_q_top =>NULL() real, pointer, dimension(:,:,:) :: flux_lw_top =>NULL() real, pointer, dimension(:,:,:) :: flux_sw_vis_dir_top =>NULL() real, pointer, dimension(:,:,:) :: flux_sw_vis_dif_top =>NULL() real, pointer, dimension(:,:,:) :: flux_sw_nir_dir_top =>NULL() real, pointer, dimension(:,:,:) :: flux_sw_nir_dif_top =>NULL() real, pointer, dimension(:,:,:) :: flux_lh_top =>NULL() real, pointer, dimension(:,:,:) :: lprec_top =>NULL() real, pointer, dimension(:,:,:) :: fprec_top =>NULL() real, pointer, dimension(:,: ) :: flux_u =>NULL() real, pointer, dimension(:,: ) :: flux_v =>NULL() real, pointer, dimension(:,: ) :: flux_t =>NULL() real, pointer, dimension(:,: ) :: flux_q =>NULL() real, pointer, dimension(:,: ) :: flux_lw =>NULL() real, pointer, dimension(:,: ) :: flux_sw_vis_dir =>NULL() real, pointer, dimension(:,: ) :: flux_sw_vis_dif =>NULL() real, pointer, dimension(:,: ) :: flux_sw_nir_dir =>NULL() real, pointer, dimension(:,: ) :: flux_sw_nir_dif =>NULL() real, pointer, dimension(:,: ) :: flux_lh =>NULL() real, pointer, dimension(:,: ) :: lprec =>NULL() real, pointer, dimension(:,: ) :: fprec =>NULL() real, pointer, dimension(:,: ) :: p_surf =>NULL() real, pointer, dimension(:,: ) :: runoff =>NULL() real, pointer, dimension(:,: ) :: calving =>NULL() real, pointer, dimension(:,: ) :: runoff_hflx =>NULL() real, pointer, dimension(:,: ) :: calving_hflx =>NULL() real, pointer, dimension(:,: ) :: flux_salt =>NULL() real, pointer, dimension(:,:) :: lwdn =>NULL() real, pointer, dimension(:,: ) :: swdn =>NULL() ! downward long/shortwave real, pointer, dimension(:,:,:) :: pen =>NULL() real, pointer, dimension(:,:,:) :: trn =>NULL() ! ice optical parameters real, pointer, dimension(:,:,:) :: tmelt =>NULL() real, pointer, dimension(:,:,:) :: bmelt =>NULL() real, pointer, dimension(:,:,:) :: h_snow =>NULL() real, pointer, dimension(:,:,:) :: h_ice =>NULL() real, pointer, dimension(:,:,:) :: t_ice1 =>NULL() real, pointer, dimension(:,:,:) :: t_ice2 =>NULL() real, pointer, dimension(:,:) :: u_ice =>NULL() real, pointer, dimension(:,:) :: v_ice =>NULL() real, pointer, dimension(:,:) :: sig11 =>NULL() real, pointer, dimension(:,:) :: sig22 =>NULL() real, pointer, dimension(:,:) :: sig12 =>NULL() real, pointer, dimension(:,:) :: frazil =>NULL() real, pointer, dimension(:,:) :: bheat =>NULL() real, pointer, dimension(:,:) :: qflx_lim_ice =>NULL() real, pointer, dimension(:,:) :: qflx_res_ice =>NULL() real, pointer, dimension(:,:) :: area =>NULL() real, pointer, dimension(:,:) :: vmask =>NULL() ! where ice vels can be non-zero real, pointer, dimension(:,:) :: mi =>NULL() ! This is needed for the wave model. It is introduced here, ! because flux_ice_to_ocean cannot handle 3D fields. This may be ! removed, if the information on ice thickness can be derived from ! eventually from h_ice outside the ice module. real, pointer, dimension(:,:,:) :: wnd =>NULL() ! need 10m wind speed to deduce Stokes drift etc in QL_KF_17 logical, pointer, dimension(:,:) :: maskmap =>NULL() ! A pointer to an array indicating which ! logical processors are actually used for ! the ocean code. The other logical ! processors would be all land points and ! are not assigned to actual processors. ! This need not be assigned if all logical ! processors are used integer, dimension(3) :: axes type(coupler_3d_bc_type) :: ocean_fields ! array of fields used for additional tracers type(coupler_2d_bc_type) :: ocean_fluxes ! array of fluxes used for additional tracers type(coupler_3d_bc_type) :: ocean_fluxes_top ! array of fluxes for averaging type(icebergs), pointer :: icebergs end type ice_data_type integer :: iceClock, iceClock1, iceCLock2, iceCLock3, iceClock4, iceClock5, & iceClock6, iceClock7, iceClock8, iceClock9, iceClocka, iceClockb, iceClockc type(restart_file_type), save :: Ice_restart contains !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! ! ice_stock_pe - returns stocks of heat, water, etc. for conservation checks ! !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! subroutine ice_stock_pe(Ice, index, value) use stock_constants_mod, only : ISTOCK_WATER, ISTOCK_HEAT, ISTOCK_SALT use ice_grid_mod, only : all_avg type(ice_data_type) :: Ice integer, intent(in) :: index real, intent(out) :: value integer :: i, j, k real :: icebergs_value value = 0.0 if(.not.Ice%pe) return select case (index) case (ISTOCK_WATER) value = sum(cell_area*all_avg(DI*Ice%h_ice(isc:iec,jsc:jec,:)+ & DS*Ice%h_snow(isc:iec,jsc:jec,:), & Ice%part_size(isc:iec,jsc:jec,:))) & *4*pi*radius*radius case (ISTOCK_HEAT) value = 0.0 do k=2,km do j=jsc, jec do i=isc, iec if ((Ice%part_size(i,j,k)>0.0.and.Ice%h_ice(i,j,k)>0.0)) then if (slab_ice) then value = value - cell_area(i,j) * Ice%part_size(i,j,k)*Ice%h_ice(i,j,2)*DI*LI else value = value - cell_area(i,j) * Ice%part_size(i,j,k)*e_to_melt(Ice%h_snow(i,j,k), & Ice%h_ice(i,j,k)/2, Ice%t_ice1(i,j,k), Ice%h_ice(i,j,k)/2, Ice%t_ice2(i,j,k)) end if end if end do end do end do value = value*4*pi*radius*radius case (ISTOCK_SALT) !No salt in the h_snow component. value = sum(cell_area*all_avg(DI*Ice%h_ice(isc:iec,jsc:jec,:),Ice%part_size(isc:iec,jsc:jec,:))) & *ice_bulk_salin*4*pi*radius*radius case default value = 0.0 end select if (do_icebergs) then call icebergs_stock_pe(Ice%icebergs, index, icebergs_value) value = value + icebergs_value endif end subroutine ice_stock_pe !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! ! ice_model_init - initializes ice model data, parameters and diagnostics ! !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! subroutine ice_model_init (Ice, Time_Init, Time, Time_step_fast, Time_step_slow ) type (ice_data_type), intent(inout) :: Ice type (time_type) , intent(in) :: Time_Init ! starting time of model integration type (time_type) , intent(in) :: Time ! current time type (time_type) , intent(in) :: Time_step_fast ! time step for the ice_model_fast type (time_type) , intent(in) :: Time_step_slow ! time step for the ice_model_slow integer :: io, ierr, nlon, nlat, npart, unit, log_unit, k integer :: sc, dy, i, j integer :: id_restart, id_restart_albedo, id_restart_flux_sw real :: dt_slow character(len=64) :: restart_file integer :: stdlogunit, stdoutunit stdlogunit=stdlog() stdoutunit = stdout() ! ! read namelist and write to logfile ! #ifdef INTERNAL_FILE_NML read (input_nml_file, nml=ice_model_nml, iostat=io) #else unit = open_namelist_file() read (unit, ice_model_nml,iostat=io) call close_file (unit) #endif ierr = check_nml_error(io,'ice_model_nml') write (stdoutunit,'(/)') write (stdoutunit, ice_model_nml) write (stdlogunit, ice_model_nml) call write_version_number(version, tagname) if (spec_ice) then slab_ice = .true. nsteps_dyn = 0 nsteps_adv = 0 end if if (slab_ice) num_part = 2 ! open water and ice ... but never in same place if (num_part>size(hlim(:))+1) & call error_mesg ('ice_model_init', 'not enough thickness limits', FATAL) call get_time(Time_step_slow, sc, dy); dt_slow=864e2*dy+sc if(file_exist(mask_table)) then write(stdoutunit, *) '==> NOTE from ice_model_init: reading maskmap information from '//trim(mask_table) if(layout(1) == 0 .OR. layout(2) == 0 ) call error_mesg('ice_model_init', & 'ice_model_nml layout should be set when file '//trim(mask_table)//' exists', FATAL) allocate(Ice%maskmap(layout(1), layout(2))) call parse_mask_table(mask_table, Ice%maskmap, "Ice model") endif if( ASSOCIATED(Ice%maskmap) ) then call set_ice_grid(Ice%domain, dt_slow, nsteps_dyn, nsteps_adv, num_part, layout, io_layout, Ice%maskmap ) else call set_ice_grid(Ice%domain, dt_slow, nsteps_dyn, nsteps_adv, num_part, layout, io_layout ) end if call set_domain(domain) call ice_dyn_param(p0, c0, cdw, wd_turn, slab_ice) call ice_thm_param(alb_sno, alb_ice, pen_ice, opt_dep_ice, slab_ice, & t_range_melt, cm2_bugs, ks, h_lo_lim) allocate ( Ice % mask (isc:iec, jsc:jec) , & Ice % ice_mask (isc:iec, jsc:jec, km) , & Ice % t_surf (isc:iec, jsc:jec, km) , & Ice % s_surf (isc:iec, jsc:jec) , & Ice % vmask (isd:ied, jsd:jed) , & Ice % sea_lev (isd:ied, jsd:jed) , & Ice % part_size (isd:ied, jsd:jed, km) , & Ice % part_size_uv (isc:iec, jsc:jec, km) , & Ice % wnd (isc:iec, jsc:jec, km) , & Ice % u_surf (isc:iec, jsc:jec, km) , & Ice % v_surf (isc:iec, jsc:jec, km) , & Ice % u_ocn (isd:ied, jsd:jed) , & Ice % v_ocn (isd:ied, jsd:jed) , & Ice % rough_mom (isc:iec, jsc:jec, km) , & Ice % rough_heat (isc:iec, jsc:jec, km) , & Ice % rough_moist (isc:iec, jsc:jec, km) , & Ice % albedo (isc:iec, jsc:jec, km) , & Ice % albedo_vis_dir (isc:iec, jsc:jec, km) , & Ice % albedo_nir_dir (isc:iec, jsc:jec, km) , & Ice % albedo_vis_dif (isc:iec, jsc:jec, km) , & Ice % albedo_nir_dif (isc:iec, jsc:jec, km) ) allocate ( Ice % flux_u_top (isd:ied, jsd:jed, km) , & Ice % flux_v_top (isd:ied, jsd:jed, km) , & Ice % flux_t_top (isc:iec, jsc:jec, km) , & Ice % flux_q_top (isc:iec, jsc:jec, km) , & Ice % flux_sw_vis_dir_top(isc:iec, jsc:jec, km) , & Ice % flux_sw_vis_dif_top(isc:iec, jsc:jec, km) , & Ice % flux_sw_nir_dir_top(isc:iec, jsc:jec, km) , & Ice % flux_sw_nir_dif_top(isc:iec, jsc:jec, km) , & Ice % flux_lw_top (isc:iec, jsc:jec, km) , & Ice % flux_lh_top (isc:iec, jsc:jec, km) , & Ice % lprec_top (isc:iec, jsc:jec, km) , & Ice % fprec_top (isc:iec, jsc:jec, km) ) allocate ( Ice % flux_u_top_bgrid (isc:iec, jsc:jec, km) , & Ice % flux_v_top_bgrid (isc:iec, jsc:jec, km) ) allocate ( Ice % flux_u (isc:iec, jsc:jec ) , & Ice % flux_v (isc:iec, jsc:jec ) , & Ice % flux_t (isc:iec, jsc:jec ) , & Ice % flux_q (isc:iec, jsc:jec ) , & Ice % flux_sw_vis_dir (isc:iec, jsc:jec) , & Ice % flux_sw_vis_dif (isc:iec, jsc:jec) , & Ice % flux_sw_nir_dir (isc:iec, jsc:jec) , & Ice % flux_sw_nir_dif (isc:iec, jsc:jec) , & Ice % flux_lw (isc:iec, jsc:jec ) , & Ice % flux_lh (isc:iec, jsc:jec ) , & Ice % lprec (isc:iec, jsc:jec ) , & Ice % fprec (isc:iec, jsc:jec ) , & Ice % p_surf (isc:iec, jsc:jec ) , & Ice % runoff (isc:iec, jsc:jec ) , & Ice % calving (isc:iec, jsc:jec ) , & Ice % runoff_hflx (isc:iec, jsc:jec ) , & Ice % calving_hflx (isc:iec, jsc:jec ) , & Ice % flux_salt (isc:iec, jsc:jec ) , & Ice % lwdn (isc:iec, jsc:jec ) , & Ice % swdn (isc:iec, jsc:jec ) ) allocate ( Ice % frazil (isc:iec, jsc:jec), Ice % bheat (isc:iec, jsc:jec), & Ice % u_ice (isd:ied, jsd:jed), Ice % v_ice (isd:ied, jsd:jed), & Ice % sig11 (isd:ied, jsd:jed), Ice % sig22 (isd:ied, jsd:jed), & Ice % sig12 (isd:ied, jsd:jed) ) allocate ( Ice % tmelt (isc:iec, jsc:jec, 2:km), Ice % bmelt (isc:iec, jsc:jec, 2:km) , & Ice % pen (isc:iec, jsc:jec, 2:km), Ice % trn (isc:iec, jsc:jec, 2:km) , & Ice % h_snow (isd:ied, jsd:jed, 2:km), Ice % h_ice (isd:ied, jsd:jed, 2:km) , & Ice % t_ice1 (isd:ied, jsd:jed, 2:km), Ice % t_ice2 (isd:ied, jsd:jed, 2:km) ) allocate ( Ice % qflx_lim_ice (isc:iec, jsc:jec) , Ice % qflx_res_ice (isc:iec, jsc:jec) ) allocate ( Ice % area (isc:iec, jsc:jec) ) allocate ( Ice % mi (isc:iec, jsc:jec) ) Ice % flux_sw_vis_dir =0. Ice % flux_sw_vis_dif =0. Ice % flux_sw_nir_dir =0. Ice % flux_sw_nir_dif =0. Ice % flux_lh =0. Ice % lwdn =0. Ice % swdn =0. Ice % flux_u_top =0. Ice % flux_v_top =0. Ice % flux_u_top_bgrid=0. Ice % flux_v_top_bgrid=0. Ice % sea_lev =0. Ice % part_size =0. Ice % wnd =0. Ice % u_ocn =0. Ice % v_ocn =0. Ice % u_ice =0. Ice % v_ice =0. Ice % sig11 =0. Ice % sig12 =0. Ice % sig22 =0. Ice % h_snow =0. Ice % h_ice =0. Ice % t_ice1 =0. Ice % t_ice2 =0. Ice % area = cell_area * 4*PI*RADIUS*RADIUS Ice % mi =0. Ice % u_surf =0. Ice % v_surf =0. Ice % s_surf =0. Ice % flux_t_top =0. Ice % flux_q_top =0. Ice % flux_lw_top =0. Ice % flux_sw_vis_dir_top =0. Ice % flux_sw_vis_dif_top =0. Ice % flux_sw_nir_dir_top =0. Ice % flux_sw_nir_dif_top =0. Ice % flux_lh_top =0. Ice % lprec_top =0. Ice % fprec_top =0. Ice % flux_salt =0. Ice % pen =0. Ice % trn =0. Ice % bheat =0. do j = jsc, jec do i = isc, iec if( wett(i,j) > 0.5 ) then Ice % mask(i,j) = .true. else Ice % mask(i,j) = .false. end if if( wetv(i,j) > 0.5 ) then Ice % vmask(i,j) = 1. else Ice % vmask(i,j) = 0. end if enddo enddo if(reproduce_siena_201303) then call mpp_update_domains(Ice%vmask, domain=domain ) else call mpp_update_domains(Ice%vmask, domain=domain, position=CORNER ) endif Ice % Time = Time Ice % Time_Init = Time_Init Ice % Time_step_fast = Time_step_fast Ice % Time_step_slow = Time_step_slow Ice % avg_kount = 0 ! ! read restart ! !! Need to create new restart version including the needed albedos ! that have been added ?? !--- restart_file = 'ice_model.res.nc' id_restart = register_restart_field(Ice_restart, restart_file, 'part_size', Ice%part_size, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'albedo', Ice%albedo, domain=domain) id_restart_albedo = register_restart_field(Ice_restart, restart_file, 'albedo_vis_dir', Ice%albedo_vis_dir, & domain=domain, mandatory=.false.) id_restart = register_restart_field(Ice_restart, restart_file, 'albedo_nir_dir', Ice%albedo_nir_dir, & domain=domain, mandatory=.false.) id_restart = register_restart_field(Ice_restart, restart_file, 'albedo_vis_dif', Ice%albedo_vis_dif, & domain=domain, mandatory=.false.) id_restart = register_restart_field(Ice_restart, restart_file, 'albedo_nir_dif', Ice%albedo_nir_dif, & domain=domain, mandatory=.false.) id_restart = register_restart_field(Ice_restart, restart_file, 'rough_mom', Ice%rough_mom, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'rough_heat', Ice%rough_heat, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'rough_moist', Ice%rough_moist, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 't_surf', Ice%t_surf, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'h_snow', Ice%h_snow(:,:,2:km), domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'h_ice', Ice%h_ice(:,:,2:km), domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 't_ice1', Ice%t_ice1(:,:,2:km), domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 't_ice2', Ice%t_ice2(:,:,2:km), domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'u_ice', Ice%u_ice, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'v_ice', Ice%v_ice, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'sig11', Ice%sig11, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'sig22', Ice%sig22, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'sig12', Ice%sig12, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'flux_u', Ice%flux_u, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'flux_v', Ice%flux_v, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'flux_t', Ice%flux_t, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'flux_q', Ice%flux_q, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'flux_salt', Ice%flux_salt, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'flux_lw', Ice%flux_lw, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'lprec', Ice%lprec, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'fprec', Ice%fprec, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'runoff', Ice%runoff, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'calving', Ice%calving, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'runoff_hflx', Ice%runoff_hflx, domain=domain, mandatory=.false.) id_restart = register_restart_field(Ice_restart, restart_file, 'calving_hflx',Ice%calving_hflx, domain=domain, mandatory=.false.) id_restart = register_restart_field(Ice_restart, restart_file, 'p_surf', Ice%p_surf, domain=domain) id_restart = register_restart_field(Ice_restart, restart_file, 'flux_sw_vis_dir', Ice%flux_sw_vis_dir, & domain=domain, mandatory=.false.) id_restart = register_restart_field(Ice_restart, restart_file, 'flux_sw_vis_dif', Ice%flux_sw_vis_dif, & domain=domain, mandatory=.false.) id_restart_flux_sw = register_restart_field(Ice_restart, restart_file, 'flux_sw_nir_dir', Ice%flux_sw_nir_dir, & domain=domain, mandatory=.false.) id_restart = register_restart_field(Ice_restart, restart_file, 'flux_sw_nir_dif', Ice%flux_sw_nir_dif, & domain=domain, mandatory=.false.) ! ! Total SW flux is broken into 4 components in Nalanda, ! it was a single component preNalanda. ! restart_file = 'INPUT/ice_model.res.nc' if (file_exist(restart_file)) then call restore_state(Ice_restart) if( .NOT. query_initialized(Ice_restart, id_restart_albedo)) then Ice % albedo_vis_dir = 0.0 Ice % albedo_nir_dir = 0.0 Ice % albedo_vis_dif = 0.0 Ice % albedo_nir_dif = 0.0 endif if(.NOT. query_initialized(Ice_restart, id_restart_flux_sw) ) then call error_mesg ('ice_model_init', & 'Restart file does not contain flux_sw_* subcomponents!', WARNING) ! for compatibility with preN restarts, we should check for total SW flux if(field_exist(restart_file,'flux_sw')) then call read_data( restart_file, 'flux_sw',Ice%flux_sw_vis_dir , domain ) ! simplest way to brake the total flux to 4 components Ice % flux_sw_vis_dir = Ice%flux_sw_vis_dir / 4 Ice % flux_sw_vis_dif = Ice%flux_sw_vis_dir Ice % flux_sw_nir_dir = Ice%flux_sw_vis_dir Ice % flux_sw_nir_dif = Ice%flux_sw_vis_dir else call error_mesg ('ice_model_init', & 'Restart file does not contain flux_sw total or its components!', FATAL) endif endif !--- update to data domain call mpp_update_domains(Ice%part_size, Domain) call mpp_update_domains(Ice%h_snow(:,:,2:km), Domain ) call mpp_update_domains(Ice%h_ice (:,:,2:km), Domain ) call mpp_update_domains(Ice%t_ice1(:,:,2:km), Domain ) call mpp_update_domains(Ice%t_ice2(:,:,2:km), Domain ) if(reproduce_siena_201303) then call mpp_update_domains(Ice%u_ice, Domain ) call mpp_update_domains(Ice%v_ice, Domain ) else call mpp_update_domains(Ice%u_ice, Ice%v_ice, Domain, gridtype=BGRID_NE ) endif call mpp_update_domains(Ice%sig11, Domain ) call mpp_update_domains(Ice%sig22, Domain ) call mpp_update_domains(Ice%sig12, Domain ) else ! no restart => no ice Ice % part_size = 0.0 ! where (Ice%mask) Ice % part_size (:,:,1) = 1.0 - flux_exchange won't allow Ice % part_size (:,:,1) = 1.0 Ice % albedo = 0.0 Ice % albedo_vis_dir = 0.0 Ice % albedo_nir_dir = 0.0 Ice % albedo_vis_dif = 0.0 Ice % albedo_nir_dif = 0.0 Ice % rough_mom = mom_rough_ice Ice % rough_heat = heat_rough_ice Ice % rough_moist = heat_rough_ice Ice % t_surf = Tfreeze-5.0 Ice % h_snow = 0.0 Ice % h_ice = 0.0 Ice % t_ice1 = -5.0 Ice % t_ice2 = -5.0 Ice % u_ice = 0.0 Ice % v_ice = 0.0 Ice % sig11 = 0.0 Ice % sig22 = 0.0 Ice % sig12 = 0.0 Ice % flux_u = 0.0 Ice % flux_v = 0.0 Ice % flux_t = 0.0 Ice % flux_q = 0.0 Ice % flux_lw = 0.0 Ice % flux_salt = 0.0 Ice % lprec = 0.0 Ice % fprec = 0.0 Ice % p_surf = 0.0 Ice % runoff = 0.0 Ice % calving = 0.0 Ice % runoff_hflx = 0.0 Ice % calving_hflx = 0.0 Ice % frazil = 0.0 Ice % flux_sw_vis_dir = 0.0 Ice % flux_sw_vis_dif = 0.0 Ice % flux_sw_nir_dir = 0.0 Ice % flux_sw_nir_dif = 0.0 do_init = .true. ! done in ice_model end if Ice % tmelt = 0.0 Ice % bmelt = 0.0 Ice % qflx_lim_ice = 0.0 Ice % qflx_res_ice = 0.0 Ice%part_size_uv(:,:,1) = 1.0 do k=2,km call t_to_uv(Ice%part_size(:,:,k),Ice%part_size_uv(:,:,k)) Ice%part_size_uv (:,:,1) = Ice%part_size_uv(:,:,1)-Ice%part_size_uv (:,:,k) end do call ice_diagnostics_init(Ice) !Balaji iceClock = mpp_clock_id( 'Ice', flags=clock_flag_default, grain=CLOCK_COMPONENT ) iceClock1 = mpp_clock_id( 'Ice: bot to top', flags=clock_flag_default, grain=CLOCK_ROUTINE ) iceClock2 = mpp_clock_id( 'Ice: update slow (dn)', flags=clock_flag_default, grain=CLOCK_ROUTINE ) iceClock7 = mpp_clock_id( ' Ice: slow: conservation check', flags=clock_flag_default, grain=CLOCK_LOOP ) 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 ) iceClock5 = mpp_clock_id( ' Ice: slow: thermodynamics', flags=clock_flag_default, grain=CLOCK_LOOP ) iceClock6 = mpp_clock_id( ' Ice: slow: restore/limit', flags=clock_flag_default, grain=CLOCK_LOOP ) iceClock8 = mpp_clock_id( ' Ice: slow: salt to ocean', flags=clock_flag_default, grain=CLOCK_LOOP ) iceClock9 = mpp_clock_id( ' Ice: slow: thermodyn diags', flags=clock_flag_default, grain=CLOCK_LOOP ) iceClock3 = mpp_clock_id( 'Ice: update fast', flags=clock_flag_default, grain=CLOCK_ROUTINE ) ! Initialize icebergs if (do_icebergs) call icebergs_init(Ice%icebergs, & im, jm, layout, io_layout, Ice%axes(1:2), Ice%maskmap, x_cyclic, tripolar_grid, & dt_slow, Time, geo_lonv_ib, geo_latv_ib, wett, dtn, dte, cell_area, cos_rot, sin_rot ) if (add_diurnal_sw .or. do_sun_angle_for_alb) call astronomy_init call nullify_domain() end subroutine ice_model_init !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! ! ice_model_end - writes the restart file ! !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! subroutine ice_model_end (Ice) type (ice_data_type), intent(inout) :: Ice integer :: k integer :: unit character(len=22) :: restart='RESTART/ice_model.res' if (conservation_check) then do k=1,4 call mpp_sum(h2o(k)) call mpp_sum(heat(k)) end do if (mpp_pe()==mpp_root_pe()) then print * print '(a10,5a13)', 'ICE MODEL ',' AT START ', & ' TOP FLUX DN.', & ' BOT FLUX DN.', & ' AT END ', & ' ERROR ' print '(a10,5es13.5)','WATER ', h2o , h2o (4)-(h2o (1)+h2o (2)-h2o (3)) print '(a10,5es13.5)','HEAT ', heat, heat(4)-(heat(1)+heat(2)-heat(3)) print * end if end if call ice_model_restart() !--- release memory ------------------------------------------------ call ice_grid_end() deallocate(Ice % mask, Ice % ice_mask, Ice % t_surf, Ice % s_surf, Ice % sea_lev ) deallocate(Ice % vmask) deallocate(Ice % part_size, Ice % part_size_uv, Ice % u_surf, Ice % v_surf ) deallocate(Ice % u_ocn, Ice % v_ocn , Ice % rough_mom, Ice % rough_heat ) deallocate(Ice % rough_moist, Ice % albedo, Ice % flux_u_top, Ice % flux_v_top ) deallocate(Ice % flux_u_top_bgrid, Ice % flux_v_top_bgrid ) deallocate(Ice % flux_t_top, Ice % flux_q_top, Ice % flux_lw_top ) deallocate(Ice % flux_lh_top, Ice % lprec_top, Ice % fprec_top, Ice % flux_u ) deallocate(Ice % flux_v, Ice % flux_t, Ice % flux_q, Ice % flux_lw ) deallocate(Ice % flux_lh, Ice % lprec, Ice % fprec, Ice % p_surf, Ice % runoff ) deallocate(Ice % calving, Ice % runoff_hflx, Ice % calving_hflx ) deallocate(Ice % flux_salt) deallocate( Ice % lwdn) deallocate( Ice % swdn) deallocate( Ice % frazil ) deallocate(Ice % bheat, Ice % u_ice, Ice % v_ice, Ice % sig11, Ice % sig22 ) deallocate(Ice % sig12, Ice % tmelt, Ice % bmelt, Ice % pen, Ice % trn ) deallocate(Ice % h_snow, Ice % h_ice, Ice % t_ice1, Ice % t_ice2 ) deallocate(Ice % qflx_lim_ice, Ice % qflx_res_ice ) deallocate(Ice % flux_sw_vis_dir, Ice % flux_sw_vis_dif ) deallocate(Ice % flux_sw_nir_dir, Ice % flux_sw_nir_dif ) deallocate(Ice % wnd ) ! End icebergs if (do_icebergs) call icebergs_end(Ice%icebergs) if (add_diurnal_sw .or. do_sun_angle_for_alb) call astronomy_end end subroutine ice_model_end !####################################################################### ! ! ! Write out restart files registered through register_restart_file ! subroutine ice_model_restart(Ice, time_stamp) type (ice_data_type), intent(inout), optional :: Ice character(len=*), intent(in), optional :: time_stamp call save_restart(Ice_restart, time_stamp) !call icebergs_save_restart(Ice%icebergs) ! This should go here but since "Ice" is not available we have to ! rely on the restart written via ice_model_end() -AJA end subroutine ice_model_restart ! !####################################################################### subroutine ice_diagnostics_init(Ice) type (ice_data_type), intent(inout) :: Ice real, parameter :: missing = -1e34 integer, dimension(2) :: axt, axv, axtv, axvt integer, dimension(3) :: axt2 integer :: id_geo_lon, id_geo_lat, id_sin_rot, id_cos_rot, id_cell_area logical :: sent ! ! diagnostics MUST use a domain without halos otherwise same as the ! regular domain: Domain (see ice_grid.f90) ! id_xv = diag_axis_init('xv', xb1d(2:im+1), 'degrees_E', 'X','longitude', set_name='ice', Domain2=Domain ) id_yv = diag_axis_init('yv', yb1d(2:jm+1), 'degrees_N', 'Y','latitude', set_name='ice', Domain2=Domain ) id_xb = diag_axis_init('xb', xb1d, 'degrees_E', 'X', 'longitude', set_name='ice', Domain2=Domain ) id_yb = diag_axis_init('yb', yb1d, 'degrees_N', 'Y', 'latitude', set_name='ice', Domain2=Domain ) id_xt = diag_axis_init('xt', (xb1d(1:im)+xb1d(2:im+1))/2, 'degrees_E', 'X', & 'longitude',set_name='ice',edges=id_xb,Domain2=Domain) id_yt = diag_axis_init('yt', (yb1d(1:jm)+yb1d(2:jm+1))/2, 'degrees_N', 'Y', & 'latitude',set_name='ice', edges=id_yb,Domain2=Domain) id_ct = diag_axis_init('ct', hlim(1:num_part-1), 'meters','Z', 'thickness') axv = (/ id_xv, id_yv /) axt = (/ id_xt, id_yt /) axt2 = (/ id_xt, id_yt, id_ct/) Ice%axes(:) = axt2(:) axtv = (/ id_xt, id_yv /); ! for north faces of t-cells axvt = (/ id_xv, id_yt /); ! for east faces of t-cells id_sin_rot = register_static_field('ice_model', 'SINROT', axt, & '-SINROT,COSROT points north', 'none') id_cos_rot = register_static_field('ice_model', 'COSROT', axt, & 'COSROT,SINROT points east','none') id_geo_lon = register_static_field('ice_model', 'GEOLON', axt, 'longitude', & 'degrees') id_geo_lat = register_static_field('ice_model', 'GEOLAT', axt, 'latitude', & 'degrees') id_cell_area = register_static_field('ice_model', 'CELL_AREA', axt, & 'cell area', 'sphere') id_ext = register_diag_field('ice_model', 'MOI', axt, Ice%Time, & 'ice modeled', '0 or 1', missing_value=missing) if (id_ext > 0 ) then call error_mesg ('ice_model_init', & 'Diagnostic MOI has been renamed EXT. Change your diag_table.', WARNING) else id_ext = register_diag_field('ice_model', 'EXT', axt, Ice%Time, & 'ice modeled', '0 or 1', missing_value=missing) end if id_wnd = register_diag_field('ice_model', 'wnd ', axt, Ice%Time, & 'wind speed', 'm/s', missing_value=missing) id_mi = register_diag_field('ice_model', 'MI', axt, Ice%Time, & 'ice mass', 'kg/m^2', missing_value=missing) id_mib = register_diag_field('ice_model', 'MIB', axt, Ice%Time, & 'ice + bergs mass', 'kg/m^2', missing_value=missing) id_cn = register_diag_field('ice_model', 'CN', axt2, Ice%Time, & 'ice concentration', '0-1', missing_value=missing) id_hs = register_diag_field('ice_model', 'HS', axt, Ice%Time, & 'snow thickness', 'm-snow', missing_value=missing) id_hi = register_diag_field('ice_model', 'HI', axt, Ice%Time, & 'ice thickness', 'm-ice', missing_value=missing) id_t1 = register_diag_field('ice_model', 'T1', axt, Ice%Time, & 'upper ice layer temperature', 'C', missing_value=missing) id_t2 = register_diag_field('ice_model', 'T2', axt, Ice%Time, & 'lower ice layer temperature', 'C', missing_value=missing) id_ts = register_diag_field('ice_model', 'TS', axt, Ice%Time, & 'surface temperature', 'C', missing_value=missing) id_sh = register_diag_field('ice_model','SH' ,axt, Ice%Time, & 'sensible heat flux', 'W/m^2', missing_value=missing) id_lh = register_diag_field('ice_model','LH' ,axt, Ice%Time, & 'latent heat flux', 'W/m^2', missing_value=missing) id_sw = register_diag_field('ice_model','SW' ,axt, Ice%Time, & 'short wave heat flux', 'W/m^2', missing_value=missing) id_lw = register_diag_field('ice_model','LW' ,axt, Ice%Time, & 'long wave heat flux over ice', 'W/m^2', missing_value=missing) id_snofl = register_diag_field('ice_model','SNOWFL' ,axt, Ice%Time, & 'rate of snow fall', 'kg/(m^2*s)', missing_value=missing) id_rain = register_diag_field('ice_model','RAIN' ,axt, Ice%Time, & 'rate of rain fall', 'kg/(m^2*s)', missing_value=missing) id_runoff = register_diag_field('ice_model','RUNOFF' ,axt, Ice%Time, & 'liquid runoff', 'kg/(m^2*s)', missing_value=missing) id_calving = register_diag_field('ice_model','CALVING',axt, Ice%Time, & 'frozen runoff', 'kg/(m^2*s)', missing_value=missing) id_runoff_hflx = register_diag_field('ice_model','RUNOFF_HFLX' ,axt, Ice%Time, & 'liquid runoff sensible heat flux', 'W/m^2', missing_value=missing) id_calving_hflx = register_diag_field('ice_model','CALVING_HFLX',axt, Ice%Time, & 'frozen runoff sensible heat flux', 'W/m^2', missing_value=missing) id_evap = register_diag_field('ice_model','EVAP',axt, Ice%Time, & 'evaporation', 'kg/(m^2*s)', missing_value=missing) id_saltf = register_diag_field('ice_model','SALTF' ,axt, Ice%Time, & 'ice to ocean salt flux', 'kg/(m^2*s)', missing_value=missing) id_sn2ic = register_diag_field('ice_model','SN2IC' ,axt,Ice%Time, & 'rate of snow to ice conversion', 'kg/(m^2*s)', missing_value=missing) id_tmelt = register_diag_field('ice_model','TMELT' ,axt, Ice%Time, & 'upper surface melting energy flux', 'W/m^2', missing_value=missing) id_bmelt = register_diag_field('ice_model','BMELT' ,axt, Ice%Time, & 'bottom surface melting energy flux', 'W/m^2', missing_value=missing) id_bheat = register_diag_field('ice_model','BHEAT' ,axt, Ice%Time, & 'ocean to ice heat flux', 'W/m^2', missing_value=missing) id_e2m = register_diag_field('ice_model','E2MELT' ,axt, Ice%Time, & 'heat needed to melt ice', 'J/m^2', missing_value=missing) id_frazil = register_diag_field('ice_model','FRAZIL' ,axt, Ice%Time, & 'energy flux of frazil formation', 'W/m^2', missing_value=missing) id_alb = register_diag_field('ice_model','ALB',axt, Ice%Time, & 'surface albedo','0-1', missing_value=missing ) id_alb_vis_dir = register_diag_field('ice_model','alb_vis_dir',axt, Ice%Time, & 'ice surface albedo vis_dir','0-1', missing_value=missing ) id_alb_vis_dif = register_diag_field('ice_model','alb_vis_dif',axt, Ice%Time, & 'ice surface albedo vis_dif','0-1', missing_value=missing ) id_alb_nir_dir = register_diag_field('ice_model','alb_nir_dir',axt, Ice%Time, & 'ice surface albedo nir_dir','0-1', missing_value=missing ) id_alb_nir_dif = register_diag_field('ice_model','alb_nir_dif',axt, Ice%Time, & 'ice surface albedo nir_dif','0-1', missing_value=missing ) id_xprt = register_diag_field('ice_model','XPRT',axt, Ice%Time, & 'frozen water transport convergence', 'kg/(m^2*yr)', missing_value=missing) id_lsrc = register_diag_field('ice_model','LSRC', axt, Ice%Time, & 'frozen water local source', 'kg/(m^2*yr)', missing_value=missing) id_lsnk = register_diag_field('ice_model','LSNK',axt, Ice%Time, & 'frozen water local sink', 'kg/(m^2*yr)', missing_value=missing) id_bsnk = register_diag_field('ice_model','BSNK',axt, Ice%Time, & 'frozen water local bottom sink', 'kg/(m^2*yr)', missing_value=missing) id_qfres = register_diag_field('ice_model', 'QFLX_RESTORE_ICE', axt, Ice%Time, & 'Ice Restoring heat flux', 'W/m^2', missing_value=missing) id_qflim = register_diag_field('ice_model', 'QFLX_LIMIT_ICE', axt, Ice%Time, & 'Ice Limit heat flux', 'W/m^2', missing_value=missing) id_strna = register_diag_field('ice_model','STRAIN_ANGLE', axt,Ice%Time, & 'strain angle', 'none', missing_value=missing) id_sigi = register_diag_field('ice_model','SIGI' ,axt, Ice%Time, & 'first stress invariant', 'none', missing_value=missing) id_sigii = register_diag_field('ice_model','SIGII' ,axt, Ice%Time, & 'second stress invariant', 'none', missing_value=missing) id_stren = register_diag_field('ice_model','STRENGTH' ,axt, Ice%Time, & 'ice strength', 'Pa*m', missing_value=missing) id_ui = register_diag_field('ice_model', 'UI', axv, Ice%Time, & 'ice velocity - x component', 'm/s', missing_value=missing) id_vi = register_diag_field('ice_model', 'VI', axv, Ice%Time, & 'ice velocity - y component', 'm/s', missing_value=missing) id_ix_trans = register_diag_field('ice_model', 'IX_TRANS', axvt, Ice%Time, & 'x-direction ice transport', 'kg/s', missing_value=missing) id_iy_trans = register_diag_field('ice_model', 'IY_TRANS', axtv, Ice%Time, & 'y-direction ice transport', 'kg/s', missing_value=missing) id_fax = register_diag_field('ice_model', 'FA_X', axv, Ice%Time, & 'air stress on ice - x component', 'Pa', missing_value=missing) id_fay = register_diag_field('ice_model', 'FA_Y', axv, Ice%Time, & 'air stress on ice - y component', 'Pa', missing_value=missing) id_fix = register_diag_field('ice_model', 'FI_X', axv, Ice%Time, & 'ice internal stress - x component', 'Pa', missing_value=missing) id_fiy = register_diag_field('ice_model', 'FI_Y', axv, Ice%Time, & 'ice internal stress - y component', 'Pa', missing_value=missing) id_fcx = register_diag_field('ice_model', 'FC_X', axv, Ice%Time, & 'coriolis force - x component', 'Pa', missing_value=missing) id_fcy = register_diag_field('ice_model', 'FC_Y', axv, Ice%Time, & 'coriolis force - y component', 'Pa', missing_value=missing) id_fwx = register_diag_field('ice_model', 'FW_X', axv, Ice%Time, & 'water stress on ice - x component', 'Pa', missing_value=missing) id_fwy = register_diag_field('ice_model', 'FW_Y', axv, Ice%Time, & 'water stress on ice - y component', 'Pa', missing_value=missing) id_uo = register_diag_field('ice_model', 'UO', axv, Ice%Time, & 'surface current - x component', 'm/s', missing_value=missing) id_vo = register_diag_field('ice_model', 'VO', axv, Ice%Time, & 'surface current - y component', 'm/s', missing_value=missing) id_sw_vis = register_diag_field('ice_model','SW_VIS' ,axt, Ice%Time, & 'visible short wave heat flux', 'W/m^2', missing_value=missing) id_sw_dir = register_diag_field('ice_model','SW_DIR' ,axt, Ice%Time, & 'direct short wave heat flux', 'W/m^2', missing_value=missing) id_sw_dif = register_diag_field('ice_model','SW_DIF' ,axt, Ice%Time, & 'diffuse short wave heat flux', 'W/m^2', missing_value=missing) id_sw_vis_dir = register_diag_field('ice_model','SW_VIS_DIR' ,axt, Ice%Time, & 'visible direct short wave heat flux', 'W/m^2', missing_value=missing) id_sw_vis_dif = register_diag_field('ice_model','SW_VIS_DIF' ,axt, Ice%Time, & 'visible diffuse short wave heat flux', 'W/m^2', missing_value=missing) id_sw_nir_dir = register_diag_field('ice_model','SW_NIR_DIR' ,axt, Ice%Time, & 'near IR direct short wave heat flux', 'W/m^2', missing_value=missing) id_sw_nir_dif = register_diag_field('ice_model','SW_NIR_DIF' ,axt, Ice%Time, & 'near IR diffuse short wave heat flux', 'W/m^2', missing_value=missing) id_ustar = register_diag_field('ice_model', 'U_STAR', axvt, Ice%Time, & 'channel transport velocity - x component', 'm/s', missing_value=missing) id_vstar = register_diag_field('ice_model', 'V_STAR', axtv, Ice%Time, & 'channel transport velocity - y component', 'm/s', missing_value=missing) id_uocean = register_diag_field('ice_model', 'U_CHAN_OCN', axvt, Ice%Time, & 'ocean component of channel transport - x', 'm/s', missing_value=missing) id_vocean = register_diag_field('ice_model', 'V_CHAN_OCN', axtv, Ice%Time, & 'ocean component of channel transport - y', 'm/s', missing_value=missing) id_uchan = register_diag_field('ice_model', 'U_CHAN_VISC', axvt, Ice%Time, & 'viscous component of channel transport - x', 'm/s', missing_value=missing) id_vchan = register_diag_field('ice_model', 'V_CHAN_VISC', axtv, Ice%Time, & 'viscous component of channel transport - y', 'm/s', missing_value=missing) ! ! diagnostics for quantities produced outside the ice model ! id_swdn = register_diag_field('ice_model','SWDN' ,axt, Ice%Time, & 'downward shortwave flux', 'W/m^2', missing_value=missing) id_lwdn = register_diag_field('ice_model','LWDN' ,axt, Ice%Time, & 'downward longwave flux', 'W/m^2', missing_value=missing) id_ta = register_diag_field('ice_model', 'TA', axt, Ice%Time, & 'surface air temperature', 'C', missing_value=missing) id_slp = register_diag_field('ice_model', 'SLP', axt, Ice%Time, & 'sea level pressure', 'Pa', missing_value=missing) id_sst = register_diag_field('ice_model', 'SST', axt, Ice%Time, & 'sea surface temperature', 'deg-C', missing_value=missing) id_sss = register_diag_field('ice_model', 'SSS', axt, Ice%Time, & 'sea surface salinity', 'psu', missing_value=missing) id_ssh = register_diag_field('ice_model', 'SSH', axt, Ice%Time, & 'sea surface height', 'm', missing_value=missing) id_obi = register_diag_field('ice_model', 'OBI', axt, Ice%Time, & 'ice observed', '0 or 1', missing_value=missing) if (id_sin_rot>0) sent=send_data(id_sin_rot, sin_rot(isc:iec,jsc:jec), Ice%Time); if (id_cos_rot>0) sent=send_data(id_cos_rot, cos_rot(isc:iec,jsc:jec), Ice%Time); if (id_geo_lon>0) sent=send_data(id_geo_lon, geo_lon, Ice%Time); if (id_geo_lat>0) sent=send_data(id_geo_lat, geo_lat, Ice%Time); if (id_cell_area>0) sent=send_data(id_cell_area, cell_area, Ice%Time); end subroutine ice_diagnostics_init subroutine ice_data_type_chksum(id, timestep, data_type) use fms_mod, only: stdout use mpp_mod, only: mpp_chksum character(len=*), intent(in) :: id integer , intent(in) :: timestep type(ice_data_type), intent(in) :: data_type integer :: n, m, outunit outunit = stdout() 100 FORMAT(" CHECKSUM::",A32," = ",Z20) write(outunit,*) "BEGIN CHECKSUM(ice_data_type):: ", id, timestep write(outunit,100) 'ice_data_type%part_size ',mpp_chksum(data_type%part_size(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%part_size_uv ',mpp_chksum(data_type%part_size_uv(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%albedo ',mpp_chksum(data_type%albedo(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%albedo_vis_dir ',mpp_chksum(data_type%albedo_vis_dir(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%albedo_nir_dir ',mpp_chksum(data_type%albedo_nir_dir(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%albedo_vis_dif ',mpp_chksum(data_type%albedo_vis_dif(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%albedo_nir_dif ',mpp_chksum(data_type%albedo_nir_dif(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%rough_mom ',mpp_chksum(data_type%rough_mom(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%rough_heat ',mpp_chksum(data_type%rough_heat(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%rough_moist ',mpp_chksum(data_type%rough_moist(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%flux_u ',mpp_chksum(data_type%flux_u(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%flux_v ',mpp_chksum(data_type%flux_v(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%flux_t ',mpp_chksum(data_type%flux_t(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%flux_q ',mpp_chksum(data_type%flux_q(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%flux_lw ',mpp_chksum(data_type%flux_lw(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%flux_sw_vis_dir ',mpp_chksum(data_type%flux_sw_vis_dir(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%flux_sw_vis_dif ',mpp_chksum(data_type%flux_sw_vis_dif(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%flux_sw_nir_dir ',mpp_chksum(data_type%flux_sw_nir_dir(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%flux_sw_nir_dif ',mpp_chksum(data_type%flux_sw_nir_dif(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%lprec ',mpp_chksum(data_type%lprec(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%fprec ',mpp_chksum(data_type%fprec(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%p_surf ',mpp_chksum(data_type%p_surf(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%runoff ',mpp_chksum(data_type%runoff(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%calving ',mpp_chksum(data_type%calving(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%flux_salt ',mpp_chksum(data_type%flux_salt(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%h_snow ',mpp_chksum(data_type%h_snow(isc:iec,jsc:jec,:)) write(outunit,100) 'ice_data_type%h_ice ',mpp_chksum(data_type%h_ice(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%t_ice1 ',mpp_chksum(data_type%t_ice1(isc:iec,jsc:jec,:)) write(outunit,100) 'ice_data_type%t_ice2 ',mpp_chksum(data_type%t_ice2(isc:iec,jsc:jec,:)) write(outunit,100) 'ice_data_type%u_ice ',mpp_chksum(data_type%u_ice(isc:iec,jsc:jec)) write(outunit,100) 'ice_data_type%v_ice ',mpp_chksum(data_type%v_ice(isc:iec,jsc:jec)) write(outunit,100) 'ice_data_type%sig11 ',mpp_chksum(data_type%sig11(isc:iec,jsc:jec)) write(outunit,100) 'ice_data_type%sig22 ',mpp_chksum(data_type%sig22(isc:iec,jsc:jec)) write(outunit,100) 'ice_data_type%sig12 ',mpp_chksum(data_type%sig12(isc:iec,jsc:jec)) write(outunit,100) 'ice_data_type%frazil ',mpp_chksum(data_type%frazil(isc:iec,jsc:jec)) write(outunit,100) 'ice_data_type%qflx_lim_ice ',mpp_chksum(data_type%qflx_lim_ice(isc:iec,jsc:jec)) write(outunit,100) 'ice_data_type%qflx_res_ice ',mpp_chksum(data_type%qflx_res_ice(isc:iec,jsc:jec)) write(outunit,*) ' ======The following are not restart variables======' write(outunit,100) 'ice_data_type%u_surf ',mpp_chksum(data_type%u_surf(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%v_surf ',mpp_chksum(data_type%v_surf(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%sea_lev ',mpp_chksum(data_type%sea_lev(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%s_surf ',mpp_chksum(data_type%s_surf(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%u_ocn ',mpp_chksum(data_type%u_ocn(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%v_ocn ',mpp_chksum(data_type%v_ocn(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%flux_u_top ',mpp_chksum(data_type%flux_u_top(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%flux_v_top ',mpp_chksum(data_type%flux_v_top(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%flux_t_top ',mpp_chksum(data_type%flux_t_top(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%flux_q_top ',mpp_chksum(data_type%flux_q_top(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%flux_lw_top ',mpp_chksum(data_type%flux_lw_top(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%flux_sw_vis_dir_top',mpp_chksum(data_type%flux_sw_vis_dir_top(isc:iec,jsc:jec,:)) write(outunit,100) 'ice_data_type%flux_sw_vis_dif_top',mpp_chksum(data_type%flux_sw_vis_dif_top(isc:iec,jsc:jec,:)) write(outunit,100) 'ice_data_type%flux_sw_nir_dir_top',mpp_chksum(data_type%flux_sw_nir_dir_top(isc:iec,jsc:jec,:)) write(outunit,100) 'ice_data_type%flux_sw_nir_dif_top',mpp_chksum(data_type%flux_sw_nir_dif_top(isc:iec,jsc:jec,:)) write(outunit,100) 'ice_data_type%flux_lh_top ',mpp_chksum(data_type%flux_lh_top(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%lprec_top ',mpp_chksum(data_type%lprec_top(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%fprec_top ',mpp_chksum(data_type%fprec_top(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%flux_lh ',mpp_chksum(data_type%flux_lh(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%lwdn ',mpp_chksum(data_type%lwdn(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%swdn ',mpp_chksum(data_type%swdn(isc:iec,jsc:jec) ) write(outunit,100) 'ice_data_type%pen ',mpp_chksum(data_type%pen(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%trn ',mpp_chksum(data_type%trn(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%tmelt ',mpp_chksum(data_type%tmelt(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%bmelt ',mpp_chksum(data_type%bmelt(isc:iec,jsc:jec,:) ) write(outunit,100) 'ice_data_type%bheat ',mpp_chksum(data_type%bheat(isc:iec,jsc:jec)) do n = 1, data_type%ocean_fields%num_bcs !{ do m = 1, data_type%ocean_fields%bc(n)%num_fields !{ !write(outunit,101) 'ice%', m, n, mpp_chksum(Ice%ocean_fields%bc(n)%field(m)%values) write(outunit,101) 'ice%',trim(data_type%ocean_fields%bc(n)%name), & trim(data_type%ocean_fields%bc(n)%field(m)%name), & mpp_chksum(data_type%ocean_fields%bc(n)%field(m)%values) enddo !} m enddo !} n 101 FORMAT(" CHECKSUM::",A16,a,'%',a," = ",Z20) end subroutine ice_data_type_chksum end module ice_type_mod