!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !! !! !! GNU General Public License !! !! !! !! This file is part of the Flexible Modeling System (FMS). !! !! !! !! FMS is free software; you can redistribute it and/or modify !! !! it and are expected to follow the terms of the GNU General Public !! !! License as published by the Free Software Foundation. !! !! !! !! FMS is distributed in the hope that it will be useful, !! !! but WITHOUT ANY WARRANTY; without even the implied warranty of !! !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the !! !! GNU General Public License for more details. !! !! !! !! You should have received a copy of the GNU General Public License !! !! along with FMS; if not, write to: !! !! Free Software Foundation, Inc. !! !! 59 Temple Place, Suite 330 !! !! Boston, MA 02111-1307 USA !! !! or see: !! !! http://www.gnu.org/licenses/gpl.txt !! !! !! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! ! ice_spec_mod - sea ice and SST specified from data as per GFDL climate group ! !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! module ice_spec_mod use fms_mod, only: open_namelist_file, check_nml_error, close_file, & stdlog, mpp_pe, mpp_root_pe, write_version_number use time_manager_mod, only: time_type, get_date, set_date use data_override_mod,only: data_override use constants_mod, only: Tfreeze implicit none include 'netcdf.inc' private public :: get_sea_surface character(len=128), parameter :: version = '$Id: ice_spec.F90,v 17.0 2009/07/21 03:01:45 fms Exp $' character(len=128), parameter :: tagname = '$Name: mom4p1_pubrel_dec2009_nnz $' logical :: module_is_initialized = .false. logical :: mcm_ice = .false. ! When mcm_ice=.true., ice is handled as in supersource real :: sst_pert = 0. ! global temperature perturbation used for sensitivity experiments real :: minimum_ice_concentration = 0.2 real :: minimum_ice_thickness = 1.0 logical :: do_leads = .true. ! when do_leads=false there is no fractional ice concentration ! also you should set the minimum_ice_concentration = 0.5 logical :: sst_degk = .false. ! when sst_degk=true the input sst data is in degrees Kelvin ! otherwise it is assumed to be in degrees Celsius !amip date for repeating single day (rsd) option integer :: repeat_date(3)=(/-1,-1,-1/) namelist / ice_spec_nml / mcm_ice, do_leads, minimum_ice_concentration, & minimum_ice_thickness, sst_degk, sst_pert, repeat_date contains !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! ! get_sea_surface - get SST, ice concentration and thickness from data ! !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~! subroutine get_sea_surface(Time, ts, cn, iceh) type (time_type), intent(in) :: Time real, dimension(:, :), intent(out) :: ts real, dimension(size(ts,1),size(ts,2),2), intent(out) :: cn real, dimension(size(ts,1),size(ts,2)), intent(out) :: iceh real, dimension(size(ts,1),size(ts,2)) :: sst, icec real :: t_sw_freeze0 = -1.8 real :: t_sw_freeze integer :: ierr, io, unit type(time_type) :: Spec_Time integer :: tod(3),dum if(.not.module_is_initialized) then unit = open_namelist_file() ierr=1 do while (ierr /= 0) read(unit, nml=ice_spec_nml, iostat=io, end=20) ierr = check_nml_error (io, 'ice_spec_nml') enddo 20 call close_file (unit) call write_version_number(version, tagname) unit = stdlog() if(mpp_pe() == mpp_root_pe()) write (unit, nml=ice_spec_nml) module_is_initialized = .true. endif ! modify time repeating single day option if (all(repeat_date>0)) then call get_date(Time,dum,dum,dum,tod(1),tod(2),tod(3)) Spec_Time = set_date(repeat_date(1),repeat_date(2),repeat_date(3),tod(1),tod(2),tod(3)) else Spec_Time = Time endif t_sw_freeze = t_sw_freeze0 if (sst_degk) then t_sw_freeze = t_sw_freeze0 + Tfreeze ! convert sea water freeze point to degK endif icec = 0.0; iceh = 0.0; sst = t_sw_freeze; call data_override('ICE', 'sic_obs', icec, Spec_Time) call data_override('ICE', 'sit_obs', iceh, Spec_Time) call data_override('ICE', 'sst_obs', sst, Spec_Time) if(mcm_ice) then icec = 0.0 ! TK Mod: Limit minimum non-zero sea ice thickness to 0.01m. ! This is to eliminate some very thin but non-zero ! sea ice thickness values, where they really should be zero ! but have become nonzero due to spatial interpolation ! where the input grid and model grid are not ! EXACTLY the same. 0.01 was obtained by trial and ! error to roughly match supersource behavior. ! 5/22/01; 8/23/01 where (iceh < 0.01) iceh=0.0 where (iceh>0.0) icec = 1.0 sst = t_sw_freeze end where else where (icec >= minimum_ice_concentration) iceh = max(iceh, minimum_ice_thickness) sst = t_sw_freeze elsewhere icec = 0.0 iceh = 0.0 end where if (.not.do_leads) then where (icec >= minimum_ice_concentration) icec = 1.0 endif endif where (icec==0.0 .and. sst<=t_sw_freeze) sst = t_sw_freeze+1e-10 ! add on non-zero sea surface temperature perturbation (namelist option) ! this perturbation may be useful in accessing model sensitivities if ( abs(sst_pert) > 0.0001 ) then sst = sst + sst_pert endif cn(:,:,2) = icec cn(:,:,1) = 1-cn(:,:,2) if (sst_degk) then ts = sst else ts = sst+Tfreeze endif return end subroutine get_sea_surface end module ice_spec_mod