module ice_comp_mct 2,32
!---------------------------------------------------------------------------
!BOP
!
! !MODULE: ice_comp_mct
!
! !DESCRIPTION:
! CICE interface routine for the ccsm cpl7 mct system
!
! !USES:
use shr_kind_mod
, only : r8 => shr_kind_r8
use shr_sys_mod, only : shr_sys_abort, shr_sys_flush
! use shr_mem_mod, only : shr_get_memusage, shr_init_memusage
use shr_file_mod, only : shr_file_getlogunit, shr_file_getloglevel, &
shr_file_setloglevel, shr_file_setlogunit
use mct_mod
#ifdef USE_ESMF_LIB
use esmf_mod
#else
use esmf_mod, only: ESMF_clock
#endif
use seq_flds_mod
use seq_flds_indices
use seq_cdata_mod, only : seq_cdata, seq_cdata_setptrs
use seq_infodata_mod,only : seq_infodata_type, seq_infodata_getdata, &
seq_infodata_putdata, seq_infodata_start_type_cont, &
seq_infodata_start_type_brnch, seq_infodata_start_type_start
use seq_timemgr_mod, only : seq_timemgr_eclockgetdata, &
seq_timemgr_restartalarmison, &
seq_timemgr_eclockdateinsync, &
seq_timemgr_stopalarmison
use perf_mod, only : t_startf, t_stopf
use ice_flux, only : strairxt, strairyt, strocnxt, strocnyt, &
alvdr, alidr, alvdf, alidf, tref, qref, flat, &
fsens, flwout, evap, fswabs, fhocn, fswthru, &
fresh, fsalt, zlvl, uatm, vatm, potT, Tair, Qa, &
rhoa, swvdr, swvdf, swidr, swidf, flw, frain, &
fsnow, uocn, vocn, sst, ss_tltx, ss_tlty, frzmlt,&
sss, tf, wind, fsw, init_flux_atm, init_flux_ocn,&
faero
use ice_state, only : vice, vsno, aice, trcr, filename_aero, filename_iage, &
filename_volpn, filename_FY, filename_lvl, &
tr_aero, tr_iage, tr_FY, tr_pond, tr_lvl
use ice_domain_size, only : nx_global, ny_global, block_size_x, block_size_y, max_blocks
use ice_domain, only : nblocks, blocks_ice, halo_info, distrb_info, profile_barrier
use ice_blocks, only : block, get_block, nx_block, ny_block
use ice_grid, only : tlon, tlat, tarea, tmask, anglet, hm, ocn_gridcell_frac, &
grid_type, t2ugrid_vector
use ice_constants, only : c0, c1, puny, tffresh, spval_dbl, rad_to_deg, radius, &
field_loc_center, field_type_scalar, field_type_vector, c100
use ice_communicate, only : my_task, master_task, lprint_stats, MPI_COMM_ICE
use ice_calendar, only : idate, mday, time, month, daycal, secday, &
sec, dt, dt_dyn, xndt_dyn, calendar, &
calendar_type, nextsw_cday, days_per_year,&
get_daycal, leap_year_count
use ice_timers
use ice_probability, only : init_numIceCells, print_numIceCells, &
write_numIceCells, accum_numIceCells2
use ice_kinds_mod, only : int_kind, dbl_kind, char_len_long, log_kind
! use ice_init
use ice_boundary, only : ice_HaloUpdate
use ice_scam, only : scmlat, scmlon, single_column
use ice_fileunits, only : nu_diag
use ice_dyn_evp, only : kdyn
use ice_prescribed_mod
use ice_step_mod
use CICE_RunMod
use ice_global_reductions
use ice_broadcast
! !PUBLIC MEMBER FUNCTIONS:
implicit none
public :: ice_init_mct
public :: ice_run_mct
public :: ice_final_mct
SAVE
private ! By default make data private
!
! ! PUBLIC DATA:
!
! !REVISION HISTORY:
! Author: Jacob Sewall, Mariana Vertenstein
!
!EOP
! !PRIVATE MEMBER FUNCTIONS:
private :: ice_export_mct
private :: ice_import_mct
private :: ice_SetGSMap_mct
private :: ice_domain_mct
!
! !PRIVATE VARIABLES
integer (kind=int_kind) :: ICEID
logical (kind=log_kind) :: atm_aero
!=======================================================================
contains
!=======================================================================
!BOP
!
! !IROUTINE: ice_init_mct
!
! !INTERFACE:
subroutine ice_init_mct( EClock, cdata_i, x2i_i, i2x_i, NLFilename ) 1,30
!
! !DESCRIPTION:
! Initialize thermodynamic ice model and obtain relevant atmospheric model
! arrays back from driver
!
! !USES:
use CICE_InitMod
use ice_restart, only: runid, runtype, restart_dir, restart_format
use ice_history, only: history_dir, history_file
!
! !ARGUMENTS:
type(ESMF_Clock) , intent(in) :: EClock
type(seq_cdata) , intent(inout) :: cdata_i
type(mct_aVect) , intent(inout) :: x2i_i, i2x_i
character(len=*), optional , intent(in) :: NLFilename ! Namelist filename
!
! !LOCAL VARIABLES:
!
type(mct_gsMap) , pointer :: gsMap_ice
type(mct_gGrid) , pointer :: dom_i
type(seq_infodata_type) , pointer :: infodata ! Input init object
integer :: lsize
character(len=256) :: drvarchdir ! driver archive directory
character(len=32) :: starttype ! infodata start type
integer :: start_ymd ! Start date (YYYYMMDD)
integer :: start_tod ! start time of day (s)
integer :: ref_ymd ! Reference date (YYYYMMDD)
integer :: ref_tod ! reference time of day (s)
integer :: iyear ! yyyy
integer :: dtime ! time step
integer :: shrlogunit,shrloglev ! old values
integer :: iam,ierr
integer :: lbnum
integer :: daycal(13) !number of cumulative days per month
integer :: nleaps ! number of leap days before current year
integer :: mpicom_loc ! temporary mpicom
real(r8) :: mrss, mrss0,msize,msize0
! !REVISION HISTORY:
! Author: Mariana Vertenstein
!EOP
!-----------------------------------------------------------------------
!---------------------------------------------------------------------------
! Set cdata pointers
!---------------------------------------------------------------------------
call seq_cdata_setptrs(cdata_i, ID=ICEID, mpicom=mpicom_loc, &
gsMap=gsMap_ice, dom=dom_i, infodata=infodata)
! Determine time of next atmospheric shortwave calculation
call seq_infodata_GetData(infodata, nextsw_cday=nextsw_cday )
! Determine if aerosols are coming from the coupler
call seq_infodata_GetData(infodata, atm_aero=atm_aero )
! call shr_init_memusage()
!---------------------------------------------------------------------------
! use infodata to determine type of run
!---------------------------------------------------------------------------
! Preset single column values
single_column = .false.
scmlat = -999.
scmlon = -999.
call seq_infodata_GetData( infodata, case_name=runid , &
single_column=single_column ,scmlat=scmlat,scmlon=scmlon)
call seq_infodata_GetData( infodata, start_type=starttype)
if ( trim(starttype) == trim(seq_infodata_start_type_start)) then
runtype = "initial"
else if (trim(starttype) == trim(seq_infodata_start_type_cont) ) then
runtype = "continue"
else if (trim(starttype) == trim(seq_infodata_start_type_brnch)) then
runtype = "branch"
else
write(nu_diag,*) 'ice_comp_mct ERROR: unknown starttype'
call shr_sys_abort()
end if
! Set nextsw_cday to -1 for continue and branch runs.
if (trim(runtype) /= 'initial') nextsw_cday = -1
!=============================================================
! Set ice dtime to ice coupling frequency
!=============================================================
call seq_timemgr_EClockGetData(EClock, dtime=dtime, calendar=calendar_type)
dt = real(dtime)
!=============================================================
! Initialize cice because grid information is needed for
! creation of GSMap_ice. cice_init also sets time manager info
!=============================================================
call t_startf ('cice_init')
call cice_init( mpicom_loc )
call t_stopf ('cice_init')
call init_numIceCells
!---------------------------------------------------------------------------
! Reset shr logging to my log file
!---------------------------------------------------------------------------
call shr_file_getLogUnit (shrlogunit)
call shr_file_getLogLevel(shrloglev)
call shr_file_setLogUnit (nu_diag)
!---------------------------------------------------------------------------
! use EClock to reset calendar information on initial start
!---------------------------------------------------------------------------
! - the following logic duplicates the logic for the concurrent system -
! cice_init is called then init_cpl is called where the start date is received
! from the flux coupler
! - in the following calculation for the variable time, iyear is used
! rather than (iyear-1) (as in init_cpl) since the sequential system permits you
! to start with year "0" not year "1"
! - on restart run
! - istep0, time and time_forc are read from restart file
! - istep1 is set to istep0
! - idate is determined from time via the call to calendar (see below)
! - on initial run
! - iyear, month and mday obtained from sync clock
! - time determined from iyear, month and mday
! - istep0 and istep1 are set to 0
if (runtype == 'initial') then
call seq_timemgr_EClockGetData(EClock, &
start_ymd=start_ymd, start_tod=start_tod, &
ref_ymd=ref_ymd, ref_tod=ref_tod)
if (ref_ymd /= start_ymd .or. ref_tod /= start_tod) then
if (my_task == master_task) then
write(nu_diag,*) 'ice_comp_mct: ref_ymd ',ref_ymd, &
' must equal start_ymd ',start_ymd
write(nu_diag,*) 'ice_comp_mct: ref_ymd ',ref_tod, &
' must equal start_ymd ',start_tod
end if
call shr_sys_abort()
end if
if (my_task == master_task) then
write(nu_diag,*) '(ice_init_mct) idate from sync clock = ', &
start_ymd
write(nu_diag,*) '(ice_init_mct) tod from sync clock = ', &
start_tod
write(nu_diag,*) &
'(ice_init_mct) resetting idate to match sync clock'
end if
idate = start_ymd
iyear = (idate/10000) ! integer year of basedate
month = (idate-iyear*10000)/100 ! integer month of basedate
mday = idate-iyear*10000-month*100-1 ! day of month of basedate
! (starts at 0)
call get_daycal(year=iyear,days_per_year_in=days_per_year, &
daycal_out=daycal)
nleaps = leap_year_count(iyear) ! this sets nleaps in ice_calendar
time = (((iyear)*days_per_year + nleaps + daycal(month)+mday)*secday) &
+ start_tod
call shr_sys_flush(nu_diag)
end if
call calendar(time) ! update calendar info
!---------------------------------------------------------------------------
! Initialize MCT attribute vectors and indices
!---------------------------------------------------------------------------
call t_startf ('cice_mct_init')
! Initialize ice gsMap
call ice_SetGSMap_mct( MPI_COMM_ICE, ICEID, GSMap_ice )
lsize = mct_gsMap_lsize(gsMap_ice, MPI_COMM_ICE)
! Initialize mct ice domain (needs ice initialization info)
call ice_domain_mct( lsize, gsMap_ice, dom_i )
! Inialize mct attribute vectors
call mct_aVect_init(x2i_i, rList=seq_flds_x2i_fields, lsize=lsize)
call mct_aVect_zero(x2i_i)
call mct_aVect_init(i2x_i, rList=seq_flds_i2x_fields, lsize=lsize)
call mct_aVect_zero(i2x_i)
!-----------------------------------------------------------------
! Prescribed ice initialization
!-----------------------------------------------------------------
call ice_prescribed_init(ICEID, gsmap_ice, dom_i)
!-----------------------------------------------------------------
! Get ready for coupling
!-----------------------------------------------------------------
call coupling_prep
!---------------------------------------------------------------------------
! Fill in export state for driver
!---------------------------------------------------------------------------
call ice_export_mct (i2x_i) !Send initial state to driver
call seq_infodata_PutData( infodata, ice_prognostic=.true., &
ice_nx = nx_global, ice_ny = ny_global )
call t_stopf ('cice_mct_init')
!---------------------------------------------------------------------------
! Reset shr logging to original values
!---------------------------------------------------------------------------
call shr_file_setLogUnit (shrlogunit)
call shr_file_setLogLevel(shrloglev)
call ice_timer_stop(timer_total) ! time entire run
! call shr_get_memusage(msize,mrss)
! call shr_mpi_max(mrss, mrss0, MPI_COMM_ICE,'ice_init_mct mrss0')
! call shr_mpi_max(msize,msize0,MPI_COMM_ICE,'ice_init_mct msize0')
! if(my_task == 0) then
! write(shrlogunit,105) 'ice_init_mct: memory_write: model date = ',start_ymd,start_tod, &
! ' memory = ',msize0,' MB (highwater) ',mrss0,' MB (usage)'
! endif
105 format( A, 2i8, A, f10.2, A, f10.2, A)
end subroutine ice_init_mct
!---------------------------------------------------------------------------
!BOP
!
! !IROUTINE: ice_run_mct
!
! !INTERFACE:
subroutine ice_run_mct( EClock, cdata_i, x2i_i, i2x_i ) 1,69
!
! !DESCRIPTION:
! Run thermodynamic CICE
!
! !USES:
use ice_history
use ice_restart
use ice_diagnostics
use ice_aerosol, only: write_restart_aero
use ice_age, only: write_restart_age
use ice_meltpond, only: write_restart_pond
use ice_FY, only: write_restart_FY
use ice_lvl, only: write_restart_lvl
use ice_restoring, only: restore_ice, ice_HaloRestore
use ice_shortwave, only: init_shortwave
! !ARGUMENTS:
type(ESMF_Clock),intent(in) :: EClock
type(seq_cdata), intent(inout) :: cdata_i
type(mct_aVect), intent(inout) :: x2i_i
type(mct_aVect), intent(inout) :: i2x_i
! !LOCAL VARIABLES:
integer :: k ! index
logical :: rstwr ! .true. ==> write a restart file
logical :: stop_now ! .true. ==> stop at the end of this run phase
integer :: ymd ! Current date (YYYYMMDD)
integer :: tod ! Current time of day (sec)
integer :: yr_sync ! Sync current year
integer :: mon_sync ! Sync current month
integer :: day_sync ! Sync current day
integer :: tod_sync ! Sync current time of day (sec)
integer :: ymd_sync ! Current year of sync clock
integer :: shrlogunit,shrloglev ! old values
integer :: lbnum
integer :: n
type(mct_gGrid) , pointer :: dom_i
type(seq_infodata_type), pointer :: infodata
type(mct_gsMap) , pointer :: gsMap_i
character(len=char_len_long) :: fname
character(len=char_len_long) :: string1, string2
character(len=*), parameter :: SubName = "ice_run_mct"
real(r8) :: mrss, mrss0,msize,msize0
logical, save :: first_time = .true.
!
! !REVISION HISTORY:
! Author: Jacob Sewall, Mariana Vertenstein
!
!EOP
!---------------------------------------------------------------------------
call ice_timer_start(timer_total) ! time entire run
!---------------------------------------------------------------------------
! Reset shr logging to my log file
!---------------------------------------------------------------------------
call shr_file_getLogUnit (shrlogunit)
call shr_file_getLogLevel(shrloglev)
call shr_file_setLogUnit (nu_diag)
! Determine time of next atmospheric shortwave calculation
call seq_cdata_setptrs(cdata_i, infodata=infodata, dom=dom_i, &
gsMap=gsMap_i)
call seq_infodata_GetData(infodata, nextsw_cday=nextsw_cday )
!-------------------------------------------------------------------
! get import state
!-------------------------------------------------------------------
call t_startf ('cice_import')
call ice_timer_start(timer_cplrecv)
call ice_import_mct( x2i_i )
call ice_timer_stop(timer_cplrecv)
call t_stopf ('cice_import')
!--------------------------------------------------------------------
! timestep update
!--------------------------------------------------------------------
call ice_timer_start(timer_step)
istep = istep + 1 ! update time step counters
istep1 = istep1 + 1
time = time + dt ! determine the time and date
call calendar(time) ! at the end of the timestep
if (resttype == 'old' .and. istep == 1) &
call init_shortwave ! initialize radiative transfer using current swdn
!-----------------------------------------------------------------
! restoring on grid boundaries
!-----------------------------------------------------------------
if (restore_ice) call ice_HaloRestore
call t_startf ('cice_initmd')
call init_mass_diags ! diagnostics per timestep
call t_stopf ('cice_initmd')
if(prescribed_ice) then ! read prescribed ice
call t_startf ('cice_presc')
call ice_prescribed_run(idate, sec)
call t_stopf ('cice_presc')
endif
call init_flux_atm ! initialize atmosphere fluxes sent to coupler
call init_flux_ocn ! initialize ocean fluxes sent to coupler
!-----------------------------------------------------------------
! Scale radiation fields
!-----------------------------------------------------------------
call t_startf ('cice_prep_radiation')
call ice_timer_start(timer_sw)
call prep_radiation(dt)
call ice_timer_stop(timer_sw)
call t_stopf ('cice_prep_radiation')
!-----------------------------------------------------------------
! thermodynamics1
!-----------------------------------------------------------------
call t_startf ('cice_therm1')
call step_therm1(dt)
call t_stopf ('cice_therm1')
!-----------------------------------------------------------------
! thermodynamics2
!-----------------------------------------------------------------
if (.not.prescribed_ice) then
call t_startf ('cice_therm2')
call step_therm2 (dt) ! post-coupler thermodynamics
call t_stopf ('cice_therm2')
end if
!-----------------------------------------------------------------
! dynamics, transport, ridging
!-----------------------------------------------------------------
if (.not.prescribed_ice .and. kdyn>0) then
if (xndt_dyn > c1) then
call t_startf ('cice_dyn')
do k = 1, nint(xndt_dyn)
call step_dynamics(dt_dyn,dt) ! dynamics, transport, ridging
enddo
call t_stopf ('cice_dyn')
else
if (mod(time, dt_dyn) == c0) then
call t_startf ('cice_dyn')
call step_dynamics(dt_dyn,dt) ! dynamics, transport, ridging
call t_stopf ('cice_dyn')
endif
endif
endif ! not prescribed_ice
!-----------------------------------------------------------------
! radiation
!-----------------------------------------------------------------
call t_startf ('cice_radiation')
call ice_timer_start(timer_sw)
call step_radiation(dt)
call ice_timer_stop(timer_sw)
call t_stopf ('cice_radiation')
!-----------------------------------------------------------------
! get ready for coupling
!-----------------------------------------------------------------
call coupling_prep
call ice_timer_stop(timer_step)
!-----------------------------------------------------------------
! write data
!-----------------------------------------------------------------
call t_startf ('cice_diag')
call ice_timer_start(timer_diags)
if (mod(istep,diagfreq) == 0) call runtime_diags(dt) ! log file
call ice_timer_stop(timer_diags)
call t_stopf ('cice_diag')
call t_startf ('cice_hist')
call ice_timer_start(timer_hist)
#if (defined _NOIO)
! Not enought memory on BGL to write a history file yet!
! call ice_write_hist (dt) ! history file
#else
call ice_write_hist (dt) ! history file
#endif
call ice_timer_stop(timer_hist)
call t_stopf ('cice_hist')
!--------------------------------------------
! Accumualate the number of active ice cells
!--------------------------------------------
call accum_numIceCells2(aice)
rstwr = seq_timemgr_RestartAlarmIsOn(EClock)
if (rstwr) then
call seq_timemgr_EClockGetData(EClock, curr_ymd=ymd_sync, curr_tod=tod_sync, &
curr_yr=yr_sync,curr_mon=mon_sync,curr_day=day_sync)
fname = restart_filename(yr_sync, mon_sync, day_sync, tod_sync)
if (my_task == master_task) then
write(nu_diag,*) &
'ice_comp_mct: calling dumpfile for restart filename= ', fname
endif
if (restart_format /= 'nc') then
if (tr_pond) then
n = index(fname,'cice.r') + 6
string1 = trim(fname(1:n-1))
string2 = trim(fname(n:lenstr(fname)))
write(filename_volpn,'(a,a,a,a)') &
string1(1:lenstr(string1)),'.volpn', &
string2(1:lenstr(string2))
endif
if (tr_aero) then
n = index(fname,'cice.r') + 6
string1 = trim(fname(1:n-1))
string2 = trim(fname(n:lenstr(fname)))
write(filename_aero,'(a,a,a,a)') &
string1(1:lenstr(string1)),'.aero', &
string2(1:lenstr(string2))
endif
if (tr_iage) then
n = index(fname,'cice.r') + 6
string1 = trim(fname(1:n-1))
string2 = trim(fname(n:lenstr(fname)))
write(filename_iage,'(a,a,a,a)') &
string1(1:lenstr(string1)),'.age', &
string2(1:lenstr(string2))
endif
if (tr_FY) then
n = index(fname,'cice.r') + 6
string1 = trim(fname(1:n-1))
string2 = trim(fname(n:lenstr(fname)))
write(filename_FY,'(a,a,a,a)') &
string1(1:lenstr(string1)),'.FY', &
string2(1:lenstr(string2))
endif
if (tr_lvl) then
n = index(fname,'cice.r') + 6
string1 = trim(fname(1:n-1))
string2 = trim(fname(n:lenstr(fname)))
write(filename_lvl,'(a,a,a,a)') &
string1(1:lenstr(string1)),'.lvl', &
string2(1:lenstr(string2))
endif
endif ! restart_format
#if (defined _NOIO)
! Not enought memory on BGL to call dumpfile file yet!
! call dumpfile(fname)
#else
call ice_timer_start(timer_readwrite)
call dumpfile(fname)
if (restart_format /= 'nc') then
if (tr_aero) call write_restart_aero(filename_aero)
if (tr_iage) call write_restart_age(filename_iage)
if (tr_FY) call write_restart_FY(filename_FY)
if (tr_lvl) call write_restart_lvl(filename_lvl)
if (tr_pond) call write_restart_pond(filename_volpn)
endif
call ice_timer_stop(timer_readwrite)
#endif
end if
!-----------------------------------------------------------------
! send export state to driver
!-----------------------------------------------------------------
call t_startf ('cice_export')
call ice_timer_start(timer_cplsend)
call ice_export_mct ( i2x_i )
call ice_timer_stop(timer_cplsend)
call t_stopf ('cice_export')
!--------------------------------------------------------------------
! check that internal clock is in sync with master clock
!--------------------------------------------------------------------
tod = sec
ymd = idate
if (.not. seq_timemgr_EClockDateInSync( EClock, ymd, tod )) then
call seq_timemgr_EClockGetData( EClock, curr_ymd=ymd_sync, &
curr_tod=tod_sync )
write(nu_diag,*)' cice ymd=',ymd ,' cice tod= ',tod
write(nu_diag,*)' sync ymd=',ymd_sync,' sync tod= ',tod_sync
call shr_sys_abort( SubName// &
":: Internal sea-ice clock not in sync with Sync Clock")
end if
! reset shr logging to my original values
call shr_file_setLogUnit (shrlogunit)
call shr_file_setLogLevel(shrloglev)
!-------------------------------------------------------------------
! stop timers and print timer info
!-------------------------------------------------------------------
! Need to have this logic here instead of in ice_final_mct since
! the ice_final_mct.F90 will still be called even in aqua-planet mode
! Could put this logic in the driver - but it seems easier here
! Need to stop this at the end of every run phase in a coupled run.
call ice_timer_stop(timer_total) ! stop timing
stop_now = seq_timemgr_StopAlarmIsOn( EClock )
if (stop_now) then
call ice_timer_print_all(stats=.true.) ! print timing information
if(lprint_stats) then
call write_numIceCells()
endif
call release_all_fileunits
end if
! if(tod == 0) then
! call shr_get_memusage(msize,mrss)
! call shr_mpi_max(mrss, mrss0, MPI_COMM_ICE,'ice_run_mct mrss0')
! call shr_mpi_max(msize,msize0,MPI_COMM_ICE,'ice_run_mct msize0')
! if(my_task == 0 ) then
! write(shrlogunit,105) 'ice_run_mct: memory_write: model date = ',ymd,tod, &
! ' memory = ',msize0,' MB (highwater) ',mrss0,' MB (usage)'
! endif
! endif
105 format( A, 2i8, A, f10.2, A, f10.2, A)
end subroutine ice_run_mct
!---------------------------------------------------------------------------
!BOP
!
! !IROUTINE: ice_final_mct
!
! !INTERFACE:
subroutine ice_final_mct( ) 1
!
! !DESCRIPTION:
! Finalize CICE
!
! !USES:
!
!------------------------------------------------------------------------------
!BOP
!
! !ARGUMENTS:
!
! !REVISION HISTORY:
!
!EOP
!---------------------------------------------------------------------------
end subroutine ice_final_mct
!===============================================================================
subroutine ice_SetGSMap_mct( mpicom, ID, gsMap_ice ) 1,2
!-------------------------------------------------------------------
!
! Arguments
!
implicit none
integer , intent(in) :: mpicom
integer , intent(in) :: ID
type(mct_gsMap), intent(inout) :: gsMap_ice
!
! Local variables
!
integer,allocatable :: gindex(:)
integer :: lat
integer :: lon
integer :: i, j, iblk, n, gi
integer :: lsize,gsize
integer :: ier
integer :: ilo, ihi, jlo, jhi ! beginning and end of physical domain
type(block) :: this_block ! block information for current block
!-------------------------------------------------------------------
! Build the CICE grid numbering for MCT
! NOTE: Numbering scheme is: West to East and South to North
! starting at south pole. Should be the same as what's used
! in SCRIP
! number the local grid
n=0
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi
do j = jlo, jhi
do i = ilo, ihi
n = n+1
enddo !i
enddo !j
enddo !iblk
lsize = n
! not valid for padded decomps
! lsize = block_size_x*block_size_y*nblocks
gsize = nx_global*ny_global
allocate(gindex(lsize),stat=ier)
n=0
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi
do j = jlo, jhi
do i = ilo, ihi
n = n+1
lon = this_block%i_glob(i)
lat = this_block%j_glob(j)
gi = (lat-1)*nx_global + lon
gindex(n) = gi
enddo !i
enddo !j
enddo !iblk
call mct_gsMap_init( gsMap_ice, gindex, mpicom, ID, lsize, gsize )
deallocate(gindex)
end subroutine ice_SetGSMap_mct
!===============================================================================
subroutine ice_export_mct( i2x_i ) 2,2
!-----------------------------------------------------
type(mct_aVect) , intent(inout) :: i2x_i
integer :: i, j, iblk, n, ij
integer :: ilo, ihi, jlo, jhi !beginning and end of physical domain
integer (kind=int_kind) :: icells ! number of ocean/ice cells
integer (kind=int_kind), dimension (nx_block*ny_block) :: indxi ! compressed indices in i
integer (kind=int_kind), dimension (nx_block*ny_block) :: indxj ! compressed indices in i
real (kind=dbl_kind), dimension(nx_block,ny_block,max_blocks) :: &
Tsrf & ! surface temperature
, tauxa & ! atmo/ice stress
, tauya &
, tauxo & ! ice/ocean stress
, tauyo &
, sicthk & ! needed for cam/som only
, ailohi ! fractional ice area
real (kind=dbl_kind) :: &
workx, worky ! tmps for converting grid
type(block) :: this_block ! block information for current block
logical :: flag
!-----------------------------------------------------
flag=.false.
!calculate ice thickness from aice and vice. Also
!create Tsrf from the first tracer (trcr) in ice_state.F
!$OMP PARALLEL DO PRIVATE(iblk,i,j,workx,worky)
do iblk = 1, nblocks
do j = 1, ny_block
do i = 1, nx_block
! sea-ice thickness needed for cam-som only
sicthk(i,j,iblk) = vice(i,j,iblk)/(aice(i,j,iblk)+puny)
! ice fraction
ailohi(i,j,iblk) = min(aice(i,j,iblk), c1)
! surface temperature
Tsrf(i,j,iblk) = Tffresh + trcr(i,j,1,iblk) !Kelvin (original ???)
! wind stress (on POP T-grid: convert to lat-lon)
workx = strairxT(i,j,iblk) ! N/m^2
worky = strairyT(i,j,iblk) ! N/m^2
tauxa(i,j,iblk) = workx*cos(ANGLET(i,j,iblk)) &
- worky*sin(ANGLET(i,j,iblk))
tauya(i,j,iblk) = worky*cos(ANGLET(i,j,iblk)) &
+ workx*sin(ANGLET(i,j,iblk))
! ice/ocean stress (on POP T-grid: convert to lat-lon)
workx = -strocnxT(i,j,iblk) ! N/m^2
worky = -strocnyT(i,j,iblk) ! N/m^2
tauxo(i,j,iblk) = workx*cos(ANGLET(i,j,iblk)) &
- worky*sin(ANGLET(i,j,iblk))
tauyo(i,j,iblk) = worky*cos(ANGLET(i,j,iblk)) &
+ workx*sin(ANGLET(i,j,iblk))
enddo
enddo
enddo
!$OMP END PARALLEL DO
do iblk = 1, nblocks
do j = 1, ny_block
do i = 1, nx_block
if (tmask(i,j,iblk) .and. ailohi(i,j,iblk) < c0 ) then
flag = .true.
endif
end do
end do
end do
if (flag) then
do iblk = 1, nblocks
do j = 1, ny_block
do i = 1, nx_block
if (tmask(i,j,iblk) .and. ailohi(i,j,iblk) < c0 ) then
write(nu_diag,*) &
' (ice) send: ERROR ailohi < 0.0 ',i,j,ailohi(i,j,iblk)
call shr_sys_flush(nu_diag)
endif
end do
end do
end do
endif
! Fill export state i2x_i
i2x_i%rAttr(:,:) = spval_dbl
n=0
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi
do j = jlo, jhi
do i = ilo, ihi
n = n+1
!-------states--------------------
i2x_i%rAttr(index_i2x_Si_sicthk,n) = sicthk(i,j,iblk) ! (needed by CAM/SOM only)
i2x_i%rAttr(index_i2x_Si_ifrac ,n) = ailohi(i,j,iblk)
if ( tmask(i,j,iblk) .and. ailohi(i,j,iblk) > c0 ) then
!-------states--------------------
i2x_i%rAttr(index_i2x_Si_t ,n) = Tsrf(i,j,iblk)
i2x_i%rAttr(index_i2x_Si_avsdr ,n) = alvdr(i,j,iblk)
i2x_i%rAttr(index_i2x_Si_anidr ,n) = alidr(i,j,iblk)
i2x_i%rAttr(index_i2x_Si_avsdf ,n) = alvdf(i,j,iblk)
i2x_i%rAttr(index_i2x_Si_anidf ,n) = alidf(i,j,iblk)
i2x_i%rAttr(index_i2x_Si_tref ,n) = Tref(i,j,iblk)
i2x_i%rAttr(index_i2x_Si_qref ,n) = Qref(i,j,iblk)
i2x_i%rAttr(index_i2x_Si_snowh ,n) = vsno(i,j,iblk) &
/ ailohi(i,j,iblk)
!--- a/i fluxes computed by ice
i2x_i%rAttr(index_i2x_Faii_taux ,n) = tauxa(i,j,iblk)
i2x_i%rAttr(index_i2x_Faii_tauy ,n) = tauya(i,j,iblk)
i2x_i%rAttr(index_i2x_Faii_lat ,n) = flat(i,j,iblk)
i2x_i%rAttr(index_i2x_Faii_sen ,n) = fsens(i,j,iblk)
i2x_i%rAttr(index_i2x_Faii_lwup ,n) = flwout(i,j,iblk)
i2x_i%rAttr(index_i2x_Faii_evap ,n) = evap(i,j,iblk)
i2x_i%rAttr(index_i2x_Faii_swnet,n) = fswabs(i,j,iblk)
!--- i/o fluxes computed by ice
i2x_i%rAttr(index_i2x_Fioi_melth,n) = fhocn(i,j,iblk)
i2x_i%rAttr(index_i2x_Fioi_swpen,n) = fswthru(i,j,iblk) ! hf from melting
i2x_i%rAttr(index_i2x_Fioi_meltw,n) = fresh(i,j,iblk) ! h2o flux from melting ???
i2x_i%rAttr(index_i2x_Fioi_salt ,n) = fsalt(i,j,iblk) ! salt flux from melting ???
i2x_i%rAttr(index_i2x_Fioi_taux ,n) = tauxo(i,j,iblk) ! stress : i/o zonal ???
i2x_i%rAttr(index_i2x_Fioi_tauy ,n) = tauyo(i,j,iblk) ! stress : i/o meridional ???
end if
enddo !i
enddo !j
enddo !iblk
end subroutine ice_export_mct
!==============================================================================
subroutine ice_import_mct( x2i_i ) 1,10
!-----------------------------------------------------
implicit none
type(mct_aVect) , intent(inout) :: x2i_i
integer :: i, j, iblk, n
integer :: ilo, ihi, jlo, jhi !beginning and end of physical domain
type(block) :: this_block ! block information for current block
integer,parameter :: nflds=15,nfldv=6
real (kind=dbl_kind),allocatable :: aflds(:,:,:,:)
real (kind=dbl_kind) :: &
workx, worky
logical (kind=log_kind) :: &
first_call = .true.
logical (kind=log_kind) :: &
prescribed_aero_tmp
!-----------------------------------------------------
! Note that the precipitation fluxes received from the coupler
! are in units of kg/s/m^2 which is what CICE requires.
! Note also that the read in below includes only values needed
! by the thermodynamic component of CICE. Variables uocn, vocn,
! ss_tltx, and ss_tlty are excluded. Also, because the SOM and
! DOM don't compute SSS. SSS is not read in and is left at
! the initilized value (see ice_flux.F init_coupler_flux) of
! 34 ppt
! Use aflds to gather the halo updates of multiple fields
! Need to separate the scalar from the vector halo updates
allocate(aflds(nx_block,ny_block,nflds,nblocks))
aflds = c0
n=0
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi
do j = jlo, jhi
do i = ilo, ihi
n = n+1
aflds(i,j, 1,iblk) = x2i_i%rAttr(index_x2i_So_t,n)
aflds(i,j, 2,iblk) = x2i_i%rAttr(index_x2i_So_s,n)
aflds(i,j, 3,iblk) = x2i_i%rAttr(index_x2i_Sa_z,n)
aflds(i,j, 4,iblk) = x2i_i%rAttr(index_x2i_Sa_ptem,n)
aflds(i,j, 5,iblk) = x2i_i%rAttr(index_x2i_Sa_tbot,n)
aflds(i,j, 6,iblk) = x2i_i%rAttr(index_x2i_Sa_shum,n)
aflds(i,j, 7,iblk) = x2i_i%rAttr(index_x2i_Sa_dens,n)
aflds(i,j, 8,iblk) = x2i_i%rAttr(index_x2i_Fioo_q,n)
aflds(i,j, 9,iblk) = x2i_i%rAttr(index_x2i_Faxa_swvdr,n)
aflds(i,j,10,iblk) = x2i_i%rAttr(index_x2i_Faxa_swndr,n)
aflds(i,j,11,iblk) = x2i_i%rAttr(index_x2i_Faxa_swvdf,n)
aflds(i,j,12,iblk) = x2i_i%rAttr(index_x2i_Faxa_swndf,n)
aflds(i,j,13,iblk) = x2i_i%rAttr(index_x2i_Faxa_lwdn,n)
aflds(i,j,14,iblk) = x2i_i%rAttr(index_x2i_Faxa_rain,n)
aflds(i,j,15,iblk) = x2i_i%rAttr(index_x2i_Faxa_snow,n)
enddo !i
enddo !j
enddo !iblk
if (.not.prescribed_ice) then
call t_startf ('cice_imp_halo')
call ice_HaloUpdate(aflds, halo_info, field_loc_center, &
field_type_scalar)
call t_stopf ('cice_imp_halo')
endif
!$OMP PARALLEL DO PRIVATE(iblk,i,j)
do iblk = 1, nblocks
do j = 1,ny_block
do i = 1,nx_block
sst (i,j,iblk) = aflds(i,j, 1,iblk)
sss (i,j,iblk) = aflds(i,j, 2,iblk)
zlvl (i,j,iblk) = aflds(i,j, 3,iblk)
potT (i,j,iblk) = aflds(i,j, 4,iblk)
Tair (i,j,iblk) = aflds(i,j, 5,iblk)
Qa (i,j,iblk) = aflds(i,j, 6,iblk)
rhoa (i,j,iblk) = aflds(i,j, 7,iblk)
frzmlt (i,j,iblk) = aflds(i,j, 8,iblk)
swvdr(i,j,iblk) = aflds(i,j, 9,iblk)
swidr(i,j,iblk) = aflds(i,j,10,iblk)
swvdf(i,j,iblk) = aflds(i,j,11,iblk)
swidf(i,j,iblk) = aflds(i,j,12,iblk)
flw (i,j,iblk) = aflds(i,j,13,iblk)
frain(i,j,iblk) = aflds(i,j,14,iblk)
fsnow(i,j,iblk) = aflds(i,j,15,iblk)
enddo !i
enddo !j
enddo !iblk
!$OMP END PARALLEL DO
deallocate(aflds)
allocate(aflds(nx_block,ny_block,nfldv,nblocks))
aflds = c0
n=0
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi
do j = jlo, jhi
do i = ilo, ihi
n = n+1
aflds(i,j, 1,iblk) = x2i_i%rAttr(index_x2i_So_u,n)
aflds(i,j, 2,iblk) = x2i_i%rAttr(index_x2i_So_v,n)
aflds(i,j, 3,iblk) = x2i_i%rAttr(index_x2i_Sa_u,n)
aflds(i,j, 4,iblk) = x2i_i%rAttr(index_x2i_Sa_v,n)
aflds(i,j, 5,iblk) = x2i_i%rAttr(index_x2i_So_dhdx,n)
aflds(i,j, 6,iblk) = x2i_i%rAttr(index_x2i_So_dhdy,n)
enddo
enddo
enddo
if (.not.prescribed_ice) then
call t_startf ('cice_imp_halo')
call ice_HaloUpdate(aflds, halo_info, field_loc_center, &
field_type_vector)
call t_stopf ('cice_imp_halo')
endif
!$OMP PARALLEL DO PRIVATE(iblk,i,j)
do iblk = 1, nblocks
do j = 1,ny_block
do i = 1,nx_block
uocn (i,j,iblk) = aflds(i,j, 1,iblk)
vocn (i,j,iblk) = aflds(i,j, 2,iblk)
uatm (i,j,iblk) = aflds(i,j, 3,iblk)
vatm (i,j,iblk) = aflds(i,j, 4,iblk)
ss_tltx(i,j,iblk) = aflds(i,j, 5,iblk)
ss_tlty(i,j,iblk) = aflds(i,j, 6,iblk)
enddo !i
enddo !j
enddo !iblk
!$OMP END PARALLEL DO
deallocate(aflds)
!-------------------------------------------------------
! Set aerosols from coupler
!-------------------------------------------------------
if (first_call) then
if (tr_aero .and. .not. atm_aero) then
write(nu_diag,*) 'ice_comp_mct ERROR: atm_aero must be set for tr_aero'
call shr_sys_abort()
end if
first_call = .false.
end if
n=0
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi
do j = jlo, jhi
do i = ilo, ihi
n = n+1
faero(i,j,1,iblk) = x2i_i%rAttr(index_x2i_Faxa_bcphodry,n)
faero(i,j,2,iblk) = x2i_i%rAttr(index_x2i_Faxa_bcphidry,n) &
+ x2i_i%rAttr(index_x2i_Faxa_bcphiwet,n)
! Combine all of the dust into one category
faero(i,j,3,iblk) = x2i_i%rAttr(index_x2i_Faxa_dstwet1,n) &
+ x2i_i%rAttr(index_x2i_Faxa_dstdry1,n) &
+ x2i_i%rAttr(index_x2i_Faxa_dstwet2,n) &
+ x2i_i%rAttr(index_x2i_Faxa_dstdry2,n) &
+ x2i_i%rAttr(index_x2i_Faxa_dstwet3,n) &
+ x2i_i%rAttr(index_x2i_Faxa_dstdry3,n) &
+ x2i_i%rAttr(index_x2i_Faxa_dstwet4,n) &
+ x2i_i%rAttr(index_x2i_Faxa_dstdry4,n)
enddo !i
enddo !j
enddo !iblk
!-----------------------------------------------------------------
! rotate zonal/meridional vectors to local coordinates
! compute data derived quantities
!-----------------------------------------------------------------
! Vector fields come in on T grid, but are oriented geographically
! need to rotate to pop-grid FIRST using ANGLET
! then interpolate to the U-cell centers (otherwise we
! interpolate across the pole)
! use ANGLET which is on the T grid !
call t_startf ('cice_imp_ocn')
!$OMP PARALLEL DO PRIVATE(iblk,i,j,workx,worky)
do iblk = 1, nblocks
do j = 1,ny_block
do i = 1,nx_block
! ocean
workx = uocn (i,j,iblk) ! currents, m/s
worky = vocn (i,j,iblk)
uocn(i,j,iblk) = workx*cos(ANGLET(i,j,iblk)) & ! convert to POP grid
+ worky*sin(ANGLET(i,j,iblk))
vocn(i,j,iblk) = worky*cos(ANGLET(i,j,iblk)) &
- workx*sin(ANGLET(i,j,iblk))
workx = ss_tltx (i,j,iblk) ! sea sfc tilt, m/m
worky = ss_tlty (i,j,iblk)
ss_tltx(i,j,iblk) = workx*cos(ANGLET(i,j,iblk)) & ! convert to POP grid
+ worky*sin(ANGLET(i,j,iblk))
ss_tlty(i,j,iblk) = worky*cos(ANGLET(i,j,iblk)) &
- workx*sin(ANGLET(i,j,iblk))
sst(i,j,iblk) = sst(i,j,iblk) - Tffresh ! sea sfc temp (C)
Tf (i,j,iblk) = -1.8_dbl_kind ! hardwired for NCOM
! Tf (i,j,iblk) = -depressT*sss(i,j,iblk) ! freezing temp (C)
! Tf (i,j,iblk) = -depressT*max(sss(i,j,iblk),ice_ref_salinity)
enddo
enddo
enddo
!$OMP END PARALLEL DO
call t_stopf ('cice_imp_ocn')
! Interpolate ocean dynamics variables from T-cell centers to
! U-cell centers.
if (.not.prescribed_ice) then
call t_startf ('cice_imp_t2u')
call t2ugrid_vector(uocn)
call t2ugrid_vector(vocn)
call t2ugrid_vector(ss_tltx)
call t2ugrid_vector(ss_tlty)
call t_stopf ('cice_imp_t2u')
end if
! Atmosphere variables are needed in T cell centers in
! subroutine stability and are interpolated to the U grid
! later as necessary.
call t_startf ('cice_imp_atm')
!$OMP PARALLEL DO PRIVATE(iblk,i,j,workx,worky)
do iblk = 1, nblocks
do j = 1, ny_block
do i = 1, nx_block
! atmosphere
workx = uatm(i,j,iblk) ! wind velocity, m/s
worky = vatm(i,j,iblk)
uatm (i,j,iblk) = workx*cos(ANGLET(i,j,iblk)) & ! convert to POP grid
+ worky*sin(ANGLET(i,j,iblk)) ! note uatm, vatm, wind
vatm (i,j,iblk) = worky*cos(ANGLET(i,j,iblk)) & ! are on the T-grid here
- workx*sin(ANGLET(i,j,iblk))
wind (i,j,iblk) = sqrt(uatm(i,j,iblk)**2 + vatm(i,j,iblk)**2)
fsw (i,j,iblk) = swvdr(i,j,iblk) + swvdf(i,j,iblk) &
+ swidr(i,j,iblk) + swidf(i,j,iblk)
enddo
enddo
enddo
!$OMP END PARALLEL DO
call t_stopf ('cice_imp_atm')
end subroutine ice_import_mct
!=======================================================================
subroutine ice_domain_mct( lsize, gsMap_i, dom_i ) 1,5
!-------------------------------------------------------------------
!
! Arguments
!
integer , intent(in) :: lsize
type(mct_gsMap), intent(in) :: gsMap_i
type(mct_ggrid), intent(inout) :: dom_i
!
! Local Variables
!
integer :: i, j, iblk, n, gi ! indices
integer :: ilo, ihi, jlo, jhi ! beginning and end of physical domain
real(dbl_kind), pointer :: work_dom(:) ! temporary
real(dbl_kind), pointer :: data(:) ! temporary
integer , pointer :: idata(:) ! temporary
type(block) :: this_block ! block information for current block
!-------------------------------------------------------------------
!
! Initialize mct domain type
! lat/lon in degrees, area in radians^2, mask is 1 (ocean), 0 (non-ocean)
!
call mct_gGrid_init( GGrid=dom_i, CoordChars=trim(seq_flds_dom_coord), &
OtherChars=trim(seq_flds_dom_other), lsize=lsize )
call mct_aVect_zero(dom_i%data)
!
allocate(data(lsize))
!
! Determine global gridpoint number attribute, GlobGridNum, which is set automatically by MCT
!
call mct_gsMap_orderedPoints(gsMap_i, my_task, idata)
call mct_gGrid_importIAttr(dom_i,'GlobGridNum',idata,lsize)
!
! Determine domain (numbering scheme is: West to East and South to North to South pole)
! Initialize attribute vector with special value
!
data(:) = -9999.0_R8
call mct_gGrid_importRAttr(dom_i,"lat" ,data,lsize)
call mct_gGrid_importRAttr(dom_i,"lon" ,data,lsize)
call mct_gGrid_importRAttr(dom_i,"area" ,data,lsize)
call mct_gGrid_importRAttr(dom_i,"aream",data,lsize)
data(:) = 0.0_R8
call mct_gGrid_importRAttr(dom_i,"mask",data,lsize)
call mct_gGrid_importRAttr(dom_i,"frac",data,lsize)
!
! Fill in correct values for domain components
!
allocate(work_dom(lsize))
work_dom(:) = 0.0_dbl_kind
n=0
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi
do j = jlo, jhi
do i = ilo, ihi
n = n+1
data(n) = TLON(i,j,iblk)*rad_to_deg
enddo !i
enddo !j
enddo !iblk
call mct_gGrid_importRattr(dom_i,"lon",data,lsize)
n=0
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi
do j = jlo, jhi
do i = ilo, ihi
n = n+1
data(n) = TLAT(i,j,iblk)*rad_to_deg
enddo !i
enddo !j
enddo !iblk
call mct_gGrid_importRattr(dom_i,"lat",data,lsize)
n=0
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi
do j = jlo, jhi
do i = ilo, ihi
n = n+1
data(n) = tarea(i,j,iblk)/(radius*radius)
enddo !i
enddo !j
enddo !iblk
call mct_gGrid_importRattr(dom_i,"area",data,lsize)
n=0
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi
do j = jlo, jhi
do i = ilo, ihi
n = n+1
data(n) = real(nint(hm(i,j,iblk)),kind=dbl_kind)
enddo !i
enddo !j
enddo !iblk
call mct_gGrid_importRattr(dom_i,"mask",data,lsize)
n=0
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi
do j = jlo, jhi
do i = ilo, ihi
n = n+1
if (trim(grid_type) == 'latlon') then
data(n) = ocn_gridcell_frac(i,j,iblk)
else
data(n) = real(nint(hm(i,j,iblk)),kind=dbl_kind)
end if
enddo !i
enddo !j
enddo !iblk
call mct_gGrid_importRattr(dom_i,"frac",data,lsize)
deallocate(data)
deallocate(idata)
deallocate(work_dom)
end subroutine ice_domain_mct
!=======================================================================
! BOP
!
! !ROUTINE: restart_filename
!
! !INTERFACE:
character(len=char_len_long) function restart_filename( yr_spec, mon_spec, day_spec, sec_spec ) 1,6
!
! !DESCRIPTION: Create a filename from a filename specifyer. Interpret filename specifyer
! string with:
! %c for case,
! %t for optional number argument sent into function
! %y for year
! %m for month
! %d for day
! %s for second
! %% for the "%" character
! If the filename specifyer has spaces " ", they will be trimmed out
! of the resulting filename.
!
! !USES:
use ice_restart, only: runid
!
! !INPUT/OUTPUT PARAMETERS:
integer , intent(in) :: yr_spec ! Simulation year
integer , intent(in) :: mon_spec ! Simulation month
integer , intent(in) :: day_spec ! Simulation day
integer , intent(in) :: sec_spec ! Seconds into current simulation day
!
! EOP
!
integer :: i, n ! Loop variables
integer :: year ! Simulation year
integer :: month ! Simulation month
integer :: day ! Simulation day
integer :: ncsec ! Seconds into current simulation day
character(len=char_len_long) :: string ! Temporary character string
character(len=char_len_long) :: format ! Format character string
character(len=char_len_long) :: filename_spec = '%c.cice.r.%y-%m-%d-%s' ! ice restarts
!-----------------------------------------------------------------
! Determine year, month, day and sec to put in filename
!-----------------------------------------------------------------
if ( len_trim(filename_spec) == 0 )then
call shr_sys_abort ('restart_filename: filename specifier is empty')
end if
if ( index(trim(filename_spec)," ") /= 0 )then
call shr_sys_abort ('restart_filename: filename specifier can not contain a space:'//trim(filename_spec))
end if
year = yr_spec
month = mon_spec
day = day_spec
ncsec = sec_spec
! Go through each character in the filename specifyer and interpret if special string
i = 1
restart_filename = ''
do while ( i <= len_trim(filename_spec) )
if ( filename_spec(i:i) == "%" )then
i = i + 1
select case( filename_spec(i:i) )
case( 'c' ) ! runid
string = trim(runid)
case( 'y' ) ! year
if ( year > 99999 ) then
format = '(i6.6)'
else if ( year > 9999 ) then
format = '(i5.5)'
else
format = '(i4.4)'
end if
write(string,format) year
case( 'm' ) ! month
write(string,'(i2.2)') month
case( 'd' ) ! day
write(string,'(i2.2)') day
case( 's' ) ! second
write(string,'(i5.5)') ncsec
case( '%' ) ! percent character
string = "%"
case default
call shr_sys_abort ('restart_filename: Invalid expansion character: '//filename_spec(i:i))
end select
else
n = index( filename_spec(i:), "%" )
if ( n == 0 ) n = len_trim( filename_spec(i:) ) + 1
if ( n == 0 ) exit
string = filename_spec(i:n+i-2)
i = n + i - 2
end if
if ( len_trim(restart_filename) == 0 )then
restart_filename = trim(string)
else
if ( (len_trim(restart_filename)+len_trim(string)) >= char_len_long )then
call shr_sys_abort ('restart_filename Resultant filename too long')
end if
restart_filename = trim(restart_filename) // trim(string)
end if
i = i + 1
end do
if ( len_trim(restart_filename) == 0 )then
call shr_sys_abort ('restart_filename: Resulting filename is empty')
end if
end function restart_filename
end module ice_comp_mct