00001 module ice_comp_mct
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 use shr_kind_mod, only : r8 => shr_kind_r8
00014 use shr_sys_mod, only : shr_sys_abort, shr_sys_flush
00015
00016 use shr_file_mod, only : shr_file_getlogunit, shr_file_getloglevel, &
00017 shr_file_setloglevel, shr_file_setlogunit
00018 use mct_mod
00019 use esmf_mod, only : ESMF_Clock
00020
00021 use seq_flds_mod
00022 use seq_flds_indices
00023 use seq_cdata_mod, only : seq_cdata, seq_cdata_setptrs
00024 use seq_infodata_mod,only : seq_infodata_type, seq_infodata_getdata, &
00025 seq_infodata_putdata, seq_infodata_start_type_cont, &
00026 seq_infodata_start_type_brnch, seq_infodata_start_type_start
00027 use seq_timemgr_mod, only : seq_timemgr_eclockgetdata, &
00028 seq_timemgr_restartalarmison, &
00029 seq_timemgr_eclockdateinsync, &
00030 seq_timemgr_stopalarmison
00031 use perf_mod, only : t_startf, t_stopf
00032
00033 use ice_flux, only : strairxt, strairyt, strocnxt, strocnyt, &
00034 alvdr, alidr, alvdf, alidf, tref, qref, flat, &
00035 fsens, flwout, evap, fswabs, fhocn, fswthru, &
00036 fresh, fsalt, zlvl, uatm, vatm, potT, Tair, Qa, &
00037 rhoa, swvdr, swvdf, swidr, swidf, flw, frain, &
00038 fsnow, uocn, vocn, sst, ss_tltx, ss_tlty, frzmlt,&
00039 sss, tf, wind, fsw, init_flux_atm, init_flux_ocn,&
00040 faero
00041 use ice_state, only : vice, aice, trcr, filename_aero, filename_iage, &
00042 filename_volpn, filename_FY, &
00043 tr_aero, tr_iage, tr_FY, tr_pond
00044 use ice_domain_size, only : nx_global, ny_global, block_size_x, block_size_y, max_blocks
00045 use ice_domain, only : nblocks, blocks_ice, halo_info, distrb_info, profile_barrier
00046 use ice_blocks, only : block, get_block, nx_block, ny_block
00047 use ice_grid, only : tlon, tlat, tarea, tmask, anglet, hm, ocn_gridcell_frac, &
00048 grid_type, t2ugrid_vector
00049 use ice_constants, only : c0, c1, puny, tffresh, spval_dbl, rad_to_deg, radius, &
00050 field_loc_center, field_type_scalar, field_type_vector, c100
00051 use ice_communicate, only : my_task, master_task, lprint_stats, MPI_COMM_ICE
00052 use ice_calendar, only : idate, mday, time, month, daycal, secday, &
00053 sec, dt, dt_dyn, xndt_dyn, calendar, &
00054 calendar_type, nextsw_cday, days_per_year,&
00055 get_daycal, leap_year_count
00056 use ice_timers
00057 use ice_probability, only : init_numIceCells, print_numIceCells, &
00058 write_numIceCells, accum_numIceCells2
00059
00060 use ice_kinds_mod, only : int_kind, dbl_kind, char_len_long, log_kind
00061
00062 use ice_boundary, only : ice_HaloUpdate
00063 use ice_scam, only : scmlat, scmlon, single_column
00064 use ice_fileunits, only : nu_diag
00065 use ice_dyn_evp, only: kdyn
00066 use ice_prescribed_mod, only : prescribed_ice, ice_prescribed_run
00067 use ice_prescaero_mod
00068 use ice_step_mod
00069 use CICE_RunMod
00070 use ice_global_reductions
00071 use ice_broadcast
00072
00073
00074 implicit none
00075 public :: ice_init_mct
00076 public :: ice_run_mct
00077 public :: ice_final_mct
00078 SAVE
00079 private
00080
00081
00082
00083
00084
00085
00086
00087
00088 private :: ice_export_mct
00089 private :: ice_import_mct
00090 private :: ice_SetGSMap_mct
00091 private :: ice_domain_mct
00092
00093
00094
00095
00096
00097
00098 contains
00099
00100
00101
00102
00103
00104
00105
00106 subroutine ice_init_mct( EClock, cdata_i, x2i_i, i2x_i, NLFilename )
00107
00108
00109
00110
00111
00112
00113
00114 use CICE_InitMod
00115 use ice_restart, only: runid, runtype, restart_dir, restart_format
00116 use ice_history, only: history_dir, history_file
00117
00118
00119 type(ESMF_Clock) , intent(in) :: EClock
00120 type(seq_cdata) , intent(inout) :: cdata_i
00121 type(mct_aVect) , intent(inout) :: x2i_i, i2x_i
00122 character(len=*), optional , intent(in) :: NLFilename
00123
00124
00125
00126 integer :: ICEID
00127 integer :: mpicom_ice
00128 type(mct_gsMap) , pointer :: gsMap_ice
00129 type(mct_gGrid) , pointer :: dom_i
00130 type(seq_infodata_type) , pointer :: infodata
00131 integer :: lsize
00132
00133 character(len=256) :: drvarchdir
00134 character(len=32) :: starttype
00135 integer :: start_ymd
00136 integer :: start_tod
00137 integer :: ref_ymd
00138 integer :: ref_tod
00139 integer :: iyear
00140 integer :: dtime
00141 integer :: shrlogunit,shrloglev
00142 integer :: iam,ierr
00143 integer :: lbnum
00144 integer :: daycal(13)
00145 integer :: nleaps
00146
00147 real(r8) :: mrss, mrss0,msize,msize0
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158 call seq_cdata_setptrs(cdata_i, ID=ICEID, mpicom=mpicom_ice, &
00159 gsMap=gsMap_ice, dom=dom_i, infodata=infodata)
00160
00161
00162 call seq_infodata_GetData(infodata, nextsw_cday=nextsw_cday )
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172 single_column = .false.
00173 scmlat = -999.
00174 scmlon = -999.
00175
00176 call seq_infodata_GetData( infodata, case_name=runid , &
00177 single_column=single_column ,scmlat=scmlat,scmlon=scmlon)
00178 call seq_infodata_GetData( infodata, start_type=starttype)
00179
00180 if ( trim(starttype) == trim(seq_infodata_start_type_start)) then
00181 runtype = "initial"
00182 else if (trim(starttype) == trim(seq_infodata_start_type_cont) ) then
00183 runtype = "continue"
00184 else if (trim(starttype) == trim(seq_infodata_start_type_brnch)) then
00185 runtype = "branch"
00186 else
00187 write(nu_diag,*) 'ice_comp_mct ERROR: unknown starttype'
00188 call shr_sys_abort()
00189 end if
00190
00191
00192
00193 if (trim(runtype) /= 'initial') nextsw_cday = -1
00194
00195
00196
00197
00198
00199 call seq_timemgr_EClockGetData(EClock, dtime=dtime, calendar=calendar_type)
00200 dt = real(dtime)
00201
00202
00203
00204
00205
00206
00207 call t_startf ('cice_init')
00208 call cice_init( mpicom_ice )
00209 call t_stopf ('cice_init')
00210 call init_numIceCells
00211
00212
00213
00214
00215
00216 call shr_file_getLogUnit (shrlogunit)
00217 call shr_file_getLogLevel(shrloglev)
00218 call shr_file_setLogUnit (nu_diag)
00219
00220
00221
00222
00223
00224
00225
00226
00227
00228
00229
00230
00231
00232
00233
00234
00235
00236
00237
00238
00239 if (runtype == 'initial') then
00240 call seq_timemgr_EClockGetData(EClock, &
00241 start_ymd=start_ymd, start_tod=start_tod, &
00242 ref_ymd=ref_ymd, ref_tod=ref_tod)
00243
00244 if (ref_ymd /= start_ymd .or. ref_tod /= start_tod) then
00245 if (my_task == master_task) then
00246 write(nu_diag,*) 'ice_comp_mct: ref_ymd ',ref_ymd, &
00247 ' must equal start_ymd ',start_ymd
00248 write(nu_diag,*) 'ice_comp_mct: ref_ymd ',ref_tod, &
00249 ' must equal start_ymd ',start_tod
00250 end if
00251 call shr_sys_abort()
00252 end if
00253
00254 if (my_task == master_task) then
00255 write(nu_diag,*) '(ice_init_mct) idate from sync clock = ', &
00256 start_ymd
00257 write(nu_diag,*) '(ice_init_mct) tod from sync clock = ', &
00258 start_tod
00259 write(nu_diag,*) &
00260 '(ice_init_mct) resetting idate to match sync clock'
00261 end if
00262
00263 idate = start_ymd
00264 iyear = (idate/10000)
00265 month = (idate-iyear*10000)/100
00266 mday = idate-iyear*10000-month*100-1
00267
00268
00269 call get_daycal(year=iyear,days_per_year_in=days_per_year, &
00270 daycal_out=daycal)
00271
00272 nleaps = leap_year_count(iyear)
00273 time = (((iyear)*days_per_year + nleaps + daycal(month)+mday)*secday) &
00274 + start_tod
00275
00276 call shr_sys_flush(nu_diag)
00277 end if
00278
00279 call calendar(time)
00280
00281
00282
00283
00284
00285 call t_startf ('cice_mct_init')
00286
00287
00288
00289 call ice_SetGSMap_mct( mpicom_ice, ICEID, GSMap_ice )
00290 lsize = mct_gsMap_lsize(gsMap_ice, mpicom_ice)
00291
00292
00293
00294 call ice_domain_mct( lsize, gsMap_ice, dom_i )
00295
00296
00297
00298 call mct_aVect_init(x2i_i, rList=seq_flds_x2i_fields, lsize=lsize)
00299 call mct_aVect_zero(x2i_i)
00300
00301 call mct_aVect_init(i2x_i, rList=seq_flds_i2x_fields, lsize=lsize)
00302 call mct_aVect_zero(i2x_i)
00303
00304
00305
00306
00307
00308 call coupling_prep
00309
00310
00311
00312
00313
00314 call ice_export_mct (i2x_i)
00315 call seq_infodata_PutData( infodata, ice_prognostic=.true., &
00316 ice_nx = nx_global, ice_ny = ny_global )
00317 call t_stopf ('cice_mct_init')
00318
00319
00320
00321
00322
00323 call shr_file_setLogUnit (shrlogunit)
00324 call shr_file_setLogLevel(shrloglev)
00325
00326 call ice_timer_stop(timer_total)
00327
00328
00329
00330
00331
00332
00333
00334
00335 105 format( A, 2i8, A, f10.2, A, f10.2, A)
00336
00337 end subroutine ice_init_mct
00338
00339
00340
00341
00342
00343
00344
00345 subroutine ice_run_mct( EClock, cdata_i, x2i_i, i2x_i )
00346
00347
00348
00349
00350
00351 use ice_history
00352 use ice_restart
00353 use ice_diagnostics
00354 use ice_aerosol, only: write_restart_aero
00355 use ice_age, only: write_restart_age
00356 use ice_meltpond, only: write_restart_pond
00357 use ice_FY, only: write_restart_FY
00358 use ice_restoring, only: restore_ice, ice_HaloRestore
00359 use ice_shortwave, only: init_shortwave
00360
00361
00362 type(ESMF_Clock),intent(in) :: EClock
00363 type(seq_cdata), intent(inout) :: cdata_i
00364 type(mct_aVect), intent(inout) :: x2i_i
00365 type(mct_aVect), intent(inout) :: i2x_i
00366
00367
00368 integer :: k
00369 logical :: rstwr
00370 logical :: stop_now
00371 integer :: ymd
00372 integer :: tod
00373 integer :: yr_sync
00374 integer :: mon_sync
00375 integer :: day_sync
00376 integer :: tod_sync
00377 integer :: ymd_sync
00378 integer :: shrlogunit,shrloglev
00379 integer :: lbnum
00380 integer :: n
00381 type(mct_gGrid) , pointer :: dom_i
00382 type(seq_infodata_type) , pointer :: infodata
00383 character(len=char_len_long) :: fname
00384 character(len=char_len_long) :: string1, string2
00385 character(len=*), parameter :: SubName = "ice_run_mct"
00386
00387 real(r8) :: mrss, mrss0,msize,msize0
00388
00389
00390
00391
00392
00393
00394
00395
00396 call ice_timer_start(timer_total)
00397
00398
00399
00400
00401
00402 call shr_file_getLogUnit (shrlogunit)
00403 call shr_file_getLogLevel(shrloglev)
00404 call shr_file_setLogUnit (nu_diag)
00405
00406
00407
00408 call seq_cdata_setptrs(cdata_i, infodata=infodata, dom=dom_i)
00409 call seq_infodata_GetData(infodata, nextsw_cday=nextsw_cday )
00410
00411
00412
00413
00414
00415 call t_startf ('cice_import')
00416 call ice_timer_start(timer_cplrecv)
00417 call ice_import_mct( x2i_i )
00418 call ice_timer_stop(timer_cplrecv)
00419 call t_stopf ('cice_import')
00420
00421
00422
00423
00424
00425 call ice_timer_start(timer_step)
00426
00427 istep = istep + 1
00428 istep1 = istep1 + 1
00429 time = time + dt
00430 call calendar(time)
00431
00432 if (resttype == 'old' .and. istep == 1) &
00433 call init_shortwave
00434
00435
00436
00437
00438
00439 if (restore_ice) call ice_HaloRestore
00440
00441 call t_startf ('cice_initmd')
00442 call init_mass_diags
00443 call t_stopf ('cice_initmd')
00444
00445 if(prescribed_ice) then
00446 call t_startf ('cice_presc')
00447 call ice_prescribed_run(idate, sec)
00448 call t_stopf ('cice_presc')
00449 endif
00450
00451 if(prescribed_aero) then
00452 call ice_prescaero_run(idate, sec)
00453 endif
00454
00455 call init_flux_atm
00456 call init_flux_ocn
00457
00458
00459
00460
00461
00462 call t_startf ('cice_prep_radiation')
00463 call ice_timer_start(timer_sw)
00464 call prep_radiation(dt)
00465 call ice_timer_stop(timer_sw)
00466 call t_stopf ('cice_prep_radiation')
00467
00468
00469
00470
00471
00472 call t_startf ('cice_therm1')
00473 call step_therm1(dt)
00474 call t_stopf ('cice_therm1')
00475
00476
00477
00478
00479
00480 if (.not.prescribed_ice) then
00481 call t_startf ('cice_therm2')
00482 call step_therm2 (dt)
00483 call t_stopf ('cice_therm2')
00484 end if
00485
00486
00487
00488
00489
00490 if (.not.prescribed_ice .and. kdyn>0) then
00491 if (xndt_dyn > c1) then
00492 call t_startf ('cice_dyn')
00493 do k = 1, nint(xndt_dyn)
00494 call step_dynamics(dt_dyn,dt)
00495 enddo
00496 call t_stopf ('cice_dyn')
00497 else
00498 if (mod(time, dt_dyn) == c0) then
00499 call t_startf ('cice_dyn')
00500 call step_dynamics(dt_dyn,dt)
00501 call t_stopf ('cice_dyn')
00502 endif
00503 endif
00504 endif
00505
00506
00507
00508
00509
00510 call t_startf ('cice_radiation')
00511 call ice_timer_start(timer_sw)
00512 call step_radiation(dt)
00513 call ice_timer_stop(timer_sw)
00514 call t_stopf ('cice_radiation')
00515
00516
00517
00518
00519
00520 call coupling_prep
00521
00522 call ice_timer_stop(timer_step)
00523
00524
00525
00526
00527
00528 call t_startf ('cice_diag')
00529 call ice_timer_start(timer_diags)
00530 if (mod(istep,diagfreq) == 0) call runtime_diags(dt)
00531 call ice_timer_stop(timer_diags)
00532 call t_stopf ('cice_diag')
00533
00534 call t_startf ('cice_hist')
00535 call ice_timer_start(timer_hist)
00536 #if (defined _NOIO)
00537
00538
00539 #else
00540 call ice_write_hist (dt)
00541 #endif
00542 call ice_timer_stop(timer_hist)
00543 call t_stopf ('cice_hist')
00544
00545
00546
00547
00548 call accum_numIceCells2(aice)
00549
00550 rstwr = seq_timemgr_RestartAlarmIsOn(EClock)
00551 if (rstwr) then
00552 call seq_timemgr_EClockGetData(EClock, curr_ymd=ymd_sync, curr_tod=tod_sync, &
00553 curr_yr=yr_sync,curr_mon=mon_sync,curr_day=day_sync)
00554 fname = restart_filename(yr_sync, mon_sync, day_sync, tod_sync)
00555
00556 if (my_task == master_task) then
00557 write(nu_diag,*) &
00558 'ice_comp_mct: calling dumpfile for restart filename= ', fname
00559 endif
00560
00561 if (restart_format /= 'nc') then
00562
00563 if (tr_pond) then
00564 n = index(fname,'cice.r') + 6
00565 string1 = trim(fname(1:n-1))
00566 string2 = trim(fname(n:lenstr(fname)))
00567 write(filename_volpn,'(a,a,a,a)') &
00568 string1(1:lenstr(string1)),'.volpn', &
00569 string2(1:lenstr(string2))
00570 endif
00571 if (tr_aero) then
00572 n = index(fname,'cice.r') + 6
00573 string1 = trim(fname(1:n-1))
00574 string2 = trim(fname(n:lenstr(fname)))
00575 write(filename_aero,'(a,a,a,a)') &
00576 string1(1:lenstr(string1)),'.aero', &
00577 string2(1:lenstr(string2))
00578 endif
00579 if (tr_iage) then
00580 n = index(fname,'cice.r') + 6
00581 string1 = trim(fname(1:n-1))
00582 string2 = trim(fname(n:lenstr(fname)))
00583 write(filename_iage,'(a,a,a,a)') &
00584 string1(1:lenstr(string1)),'.age', &
00585 string2(1:lenstr(string2))
00586 endif
00587 if (tr_FY) then
00588 n = index(fname,'cice.r') + 6
00589 string1 = trim(fname(1:n-1))
00590 string2 = trim(fname(n:lenstr(fname)))
00591 write(filename_FY,'(a,a,a,a)') &
00592 string1(1:lenstr(string1)),'.FY', &
00593 string2(1:lenstr(string2))
00594 endif
00595
00596 endif
00597
00598 #if (defined _NOIO)
00599
00600
00601 #else
00602 call ice_timer_start(timer_readwrite)
00603 call dumpfile(fname)
00604 if (restart_format /= 'nc') then
00605 if (tr_aero) call write_restart_aero(filename_aero)
00606 if (tr_iage) call write_restart_age(filename_iage)
00607 if (tr_pond) call write_restart_pond(filename_volpn)
00608 if (tr_FY) call write_restart_FY(filename_FY)
00609 endif
00610 call ice_timer_stop(timer_readwrite)
00611 #endif
00612 end if
00613
00614
00615
00616
00617
00618 call t_startf ('cice_export')
00619 call ice_timer_start(timer_cplsend)
00620 call ice_export_mct ( i2x_i )
00621 call ice_timer_stop(timer_cplsend)
00622 call t_stopf ('cice_export')
00623
00624
00625
00626
00627
00628 tod = sec
00629 ymd = idate
00630 if (.not. seq_timemgr_EClockDateInSync( EClock, ymd, tod )) then
00631 call seq_timemgr_EClockGetData( EClock, curr_ymd=ymd_sync, &
00632 curr_tod=tod_sync )
00633 write(nu_diag,*)' cice ymd=',ymd ,' cice tod= ',tod
00634 write(nu_diag,*)' sync ymd=',ymd_sync,' sync tod= ',tod_sync
00635 call shr_sys_abort( SubName// &
00636 ":: Internal sea-ice clock not in sync with Sync Clock")
00637 end if
00638
00639
00640
00641 call shr_file_setLogUnit (shrlogunit)
00642 call shr_file_setLogLevel(shrloglev)
00643
00644
00645
00646
00647
00648
00649
00650
00651
00652 call ice_timer_stop(timer_total)
00653
00654 stop_now = seq_timemgr_StopAlarmIsOn( EClock )
00655 if (stop_now) then
00656 call ice_timer_print_all(stats=.true.)
00657 if(lprint_stats) then
00658 call write_numIceCells()
00659 endif
00660 call release_all_fileunits
00661 end if
00662
00663
00664
00665
00666
00667
00668
00669
00670
00671
00672
00673 105 format( A, 2i8, A, f10.2, A, f10.2, A)
00674
00675 end subroutine ice_run_mct
00676
00677
00678
00679
00680
00681
00682
00683 subroutine ice_final_mct( )
00684
00685
00686
00687
00688
00689
00690
00691
00692
00693
00694
00695
00696
00697
00698
00699
00700 end subroutine ice_final_mct
00701
00702
00703
00704 subroutine ice_SetGSMap_mct( mpicom_ice, ICEID, gsMap_ice )
00705
00706
00707
00708
00709
00710 implicit none
00711 integer , intent(in) :: mpicom_ice
00712 integer , intent(in) :: ICEID
00713 type(mct_gsMap), intent(inout) :: gsMap_ice
00714
00715
00716
00717 integer,allocatable :: gindex(:)
00718 integer :: lat
00719 integer :: lon
00720 integer :: i, j, iblk, n, gi
00721 integer :: lsize,gsize
00722 integer :: ier
00723 integer :: ilo, ihi, jlo, jhi
00724 type(block) :: this_block
00725
00726
00727
00728
00729
00730
00731
00732
00733
00734 n=0
00735 do iblk = 1, nblocks
00736 this_block = get_block(blocks_ice(iblk),iblk)
00737 ilo = this_block%ilo
00738 ihi = this_block%ihi
00739 jlo = this_block%jlo
00740 jhi = this_block%jhi
00741
00742 do j = jlo, jhi
00743 do i = ilo, ihi
00744 n = n+1
00745 enddo
00746 enddo
00747 enddo
00748 lsize = n
00749
00750
00751
00752 gsize = nx_global*ny_global
00753
00754 allocate(gindex(lsize),stat=ier)
00755 n=0
00756 do iblk = 1, nblocks
00757 this_block = get_block(blocks_ice(iblk),iblk)
00758 ilo = this_block%ilo
00759 ihi = this_block%ihi
00760 jlo = this_block%jlo
00761 jhi = this_block%jhi
00762
00763 do j = jlo, jhi
00764 do i = ilo, ihi
00765 n = n+1
00766 lon = this_block%i_glob(i)
00767 lat = this_block%j_glob(j)
00768 gi = (lat-1)*nx_global + lon
00769 gindex(n) = gi
00770 enddo
00771 enddo
00772 enddo
00773
00774 call mct_gsMap_init( gsMap_ice, gindex, mpicom_ice, ICEID, lsize, gsize )
00775
00776 deallocate(gindex)
00777
00778 end subroutine ice_SetGSMap_mct
00779
00780
00781
00782
00783 subroutine ice_export_mct( i2x_i )
00784
00785
00786 type(mct_aVect) , intent(inout) :: i2x_i
00787
00788 integer :: i, j, iblk, n, ij
00789 integer :: ilo, ihi, jlo, jhi
00790 integer (kind=int_kind) :: icells
00791 integer (kind=int_kind), dimension (nx_block*ny_block) :: indxi
00792 integer (kind=int_kind), dimension (nx_block*ny_block) :: indxj
00793
00794 real (kind=dbl_kind), dimension(nx_block,ny_block,max_blocks) ::
00795 Tsrf
00796 , tauxa
00797 , tauya
00798 , tauxo
00799 , tauyo
00800 , sicthk
00801 , ailohi
00802
00803 real (kind=dbl_kind) ::
00804 gsum, workx, worky
00805
00806 type(block) :: this_block
00807 logical :: flag
00808
00809
00810 flag=.false.
00811
00812
00813
00814
00815
00816 do iblk = 1, nblocks
00817 do j = 1, ny_block
00818 do i = 1, nx_block
00819
00820
00821 sicthk(i,j,iblk) = vice(i,j,iblk)/(aice(i,j,iblk)+puny)
00822
00823
00824 ailohi(i,j,iblk) = min(aice(i,j,iblk), c1)
00825
00826
00827 Tsrf(i,j,iblk) = Tffresh + trcr(i,j,1,iblk)
00828
00829
00830 workx = strairxT(i,j,iblk)
00831 worky = strairyT(i,j,iblk)
00832 tauxa(i,j,iblk) = workx*cos(ANGLET(i,j,iblk)) &
00833 - worky*sin(ANGLET(i,j,iblk))
00834 tauya(i,j,iblk) = worky*cos(ANGLET(i,j,iblk)) &
00835 + workx*sin(ANGLET(i,j,iblk))
00836
00837
00838 workx = -strocnxT(i,j,iblk)
00839 worky = -strocnyT(i,j,iblk)
00840 tauxo(i,j,iblk) = workx*cos(ANGLET(i,j,iblk)) &
00841 - worky*sin(ANGLET(i,j,iblk))
00842 tauyo(i,j,iblk) = worky*cos(ANGLET(i,j,iblk)) &
00843 + workx*sin(ANGLET(i,j,iblk))
00844
00845 enddo
00846 enddo
00847 enddo
00848
00849
00850 do iblk = 1, nblocks
00851 do j = 1, ny_block
00852 do i = 1, nx_block
00853 if (tmask(i,j,iblk) .and. ailohi(i,j,iblk) < c0 ) then
00854 flag = .true.
00855 endif
00856 end do
00857 end do
00858 end do
00859 if (flag) then
00860 do iblk = 1, nblocks
00861 do j = 1, ny_block
00862 do i = 1, nx_block
00863 if (tmask(i,j,iblk) .and. ailohi(i,j,iblk) < c0 ) then
00864 write(nu_diag,*) &
00865 ' (ice) send: ERROR ailohi < 0.0 ',i,j,ailohi(i,j,iblk)
00866 call shr_sys_flush(nu_diag)
00867 endif
00868 end do
00869 end do
00870 end do
00871 endif
00872
00873
00874
00875 i2x_i%rAttr(:,:) = spval_dbl
00876
00877 n=0
00878 do iblk = 1, nblocks
00879 this_block = get_block(blocks_ice(iblk),iblk)
00880 ilo = this_block%ilo
00881 ihi = this_block%ihi
00882 jlo = this_block%jlo
00883 jhi = this_block%jhi
00884
00885 do j = jlo, jhi
00886 do i = ilo, ihi
00887
00888 n = n+1
00889
00890
00891 i2x_i%rAttr(index_i2x_Si_sicthk,n) = sicthk(i,j,iblk)
00892 i2x_i%rAttr(index_i2x_Si_ifrac ,n) = ailohi(i,j,iblk)
00893
00894 if ( tmask(i,j,iblk) .and. ailohi(i,j,iblk) > c0 ) then
00895
00896 i2x_i%rAttr(index_i2x_Si_t ,n) = Tsrf(i,j,iblk)
00897 i2x_i%rAttr(index_i2x_Si_avsdr ,n) = alvdr(i,j,iblk)
00898 i2x_i%rAttr(index_i2x_Si_anidr ,n) = alidr(i,j,iblk)
00899 i2x_i%rAttr(index_i2x_Si_avsdf ,n) = alvdf(i,j,iblk)
00900 i2x_i%rAttr(index_i2x_Si_anidf ,n) = alidf(i,j,iblk)
00901 i2x_i%rAttr(index_i2x_Si_tref ,n) = Tref(i,j,iblk)
00902 i2x_i%rAttr(index_i2x_Si_qref ,n) = Qref(i,j,iblk)
00903
00904
00905
00906
00907 i2x_i%rAttr(index_i2x_Faii_taux ,n) = tauxa(i,j,iblk)
00908 i2x_i%rAttr(index_i2x_Faii_tauy ,n) = tauya(i,j,iblk)
00909 i2x_i%rAttr(index_i2x_Faii_lat ,n) = flat(i,j,iblk)
00910 i2x_i%rAttr(index_i2x_Faii_sen ,n) = fsens(i,j,iblk)
00911 i2x_i%rAttr(index_i2x_Faii_lwup ,n) = flwout(i,j,iblk)
00912 i2x_i%rAttr(index_i2x_Faii_evap ,n) = evap(i,j,iblk)
00913 i2x_i%rAttr(index_i2x_Faii_swnet,n) = fswabs(i,j,iblk)
00914
00915
00916 i2x_i%rAttr(index_i2x_Fioi_melth,n) = fhocn(i,j,iblk)
00917 i2x_i%rAttr(index_i2x_Fioi_swpen,n) = fswthru(i,j,iblk)
00918 i2x_i%rAttr(index_i2x_Fioi_meltw,n) = fresh(i,j,iblk)
00919 i2x_i%rAttr(index_i2x_Fioi_salt ,n) = fsalt(i,j,iblk)
00920 i2x_i%rAttr(index_i2x_Fioi_taux ,n) = tauxo(i,j,iblk)
00921 i2x_i%rAttr(index_i2x_Fioi_tauy ,n) = tauyo(i,j,iblk)
00922 end if
00923 enddo
00924 enddo
00925 enddo
00926
00927 end subroutine ice_export_mct
00928
00929
00930
00931 subroutine ice_import_mct( x2i_i )
00932
00933
00934
00935 implicit none
00936 type(mct_aVect) , intent(inout) :: x2i_i
00937
00938 integer :: i, j, iblk, n
00939 integer :: ilo, ihi, jlo, jhi
00940 type(block) :: this_block
00941
00942 integer,parameter :: nflds=15,nfldv=6
00943 real (kind=dbl_kind),allocatable :: aflds(:,:,:,:)
00944
00945 real (kind=dbl_kind) ::
00946 gsum, workx, worky, maxwork
00947 real (kind=dbl_kind), allocatable, dimension (:,:,:) ::
00948 work
00949 logical (kind=log_kind) ::
00950 first_call = .true.
00951 logical (kind=log_kind) ::
00952 prescribed_aero_tmp
00953
00954
00955
00956
00957
00958
00959
00960
00961
00962
00963
00964
00965
00966
00967
00968 allocate(aflds(nx_block,ny_block,nflds,nblocks))
00969 aflds = c0
00970
00971 n=0
00972 do iblk = 1, nblocks
00973 this_block = get_block(blocks_ice(iblk),iblk)
00974 ilo = this_block%ilo
00975 ihi = this_block%ihi
00976 jlo = this_block%jlo
00977 jhi = this_block%jhi
00978
00979 do j = jlo, jhi
00980 do i = ilo, ihi
00981
00982 n = n+1
00983 aflds(i,j, 1,iblk) = x2i_i%rAttr(index_x2i_So_t,n)
00984 aflds(i,j, 2,iblk) = x2i_i%rAttr(index_x2i_So_s,n)
00985 aflds(i,j, 3,iblk) = x2i_i%rAttr(index_x2i_Sa_z,n)
00986 aflds(i,j, 4,iblk) = x2i_i%rAttr(index_x2i_Sa_ptem,n)
00987 aflds(i,j, 5,iblk) = x2i_i%rAttr(index_x2i_Sa_tbot,n)
00988 aflds(i,j, 6,iblk) = x2i_i%rAttr(index_x2i_Sa_shum,n)
00989 aflds(i,j, 7,iblk) = x2i_i%rAttr(index_x2i_Sa_dens,n)
00990 aflds(i,j, 8,iblk) = x2i_i%rAttr(index_x2i_Fioo_q,n)
00991 aflds(i,j, 9,iblk) = x2i_i%rAttr(index_x2i_Faxa_swvdr,n)
00992 aflds(i,j,10,iblk) = x2i_i%rAttr(index_x2i_Faxa_swndr,n)
00993 aflds(i,j,11,iblk) = x2i_i%rAttr(index_x2i_Faxa_swvdf,n)
00994 aflds(i,j,12,iblk) = x2i_i%rAttr(index_x2i_Faxa_swndf,n)
00995 aflds(i,j,13,iblk) = x2i_i%rAttr(index_x2i_Faxa_lwdn,n)
00996 aflds(i,j,14,iblk) = x2i_i%rAttr(index_x2i_Faxa_rain,n)
00997 aflds(i,j,15,iblk) = x2i_i%rAttr(index_x2i_Faxa_snow,n)
00998
00999 enddo
01000 enddo
01001
01002 enddo
01003
01004 if (.not.prescribed_ice) then
01005 call t_startf ('cice_imp_halo')
01006 call ice_HaloUpdate(aflds, halo_info, field_loc_center, &
01007 field_type_scalar)
01008 call t_stopf ('cice_imp_halo')
01009 endif
01010
01011
01012 do iblk = 1, nblocks
01013 do j = 1,ny_block
01014 do i = 1,nx_block
01015 sst (i,j,iblk) = aflds(i,j, 1,iblk)
01016 sss (i,j,iblk) = aflds(i,j, 2,iblk)
01017 zlvl (i,j,iblk) = aflds(i,j, 3,iblk)
01018 potT (i,j,iblk) = aflds(i,j, 4,iblk)
01019 Tair (i,j,iblk) = aflds(i,j, 5,iblk)
01020 Qa (i,j,iblk) = aflds(i,j, 6,iblk)
01021 rhoa (i,j,iblk) = aflds(i,j, 7,iblk)
01022 frzmlt (i,j,iblk) = aflds(i,j, 8,iblk)
01023 swvdr(i,j,iblk) = aflds(i,j, 9,iblk)
01024 swidr(i,j,iblk) = aflds(i,j,10,iblk)
01025 swvdf(i,j,iblk) = aflds(i,j,11,iblk)
01026 swidf(i,j,iblk) = aflds(i,j,12,iblk)
01027 flw (i,j,iblk) = aflds(i,j,13,iblk)
01028 frain(i,j,iblk) = aflds(i,j,14,iblk)
01029 fsnow(i,j,iblk) = aflds(i,j,15,iblk)
01030 enddo
01031 enddo
01032 enddo
01033
01034
01035 deallocate(aflds)
01036 allocate(aflds(nx_block,ny_block,nfldv,nblocks))
01037 aflds = c0
01038
01039 n=0
01040 do iblk = 1, nblocks
01041 this_block = get_block(blocks_ice(iblk),iblk)
01042 ilo = this_block%ilo
01043 ihi = this_block%ihi
01044 jlo = this_block%jlo
01045 jhi = this_block%jhi
01046
01047 do j = jlo, jhi
01048 do i = ilo, ihi
01049 n = n+1
01050 aflds(i,j, 1,iblk) = x2i_i%rAttr(index_x2i_So_u,n)
01051 aflds(i,j, 2,iblk) = x2i_i%rAttr(index_x2i_So_v,n)
01052 aflds(i,j, 3,iblk) = x2i_i%rAttr(index_x2i_Sa_u,n)
01053 aflds(i,j, 4,iblk) = x2i_i%rAttr(index_x2i_Sa_v,n)
01054 aflds(i,j, 5,iblk) = x2i_i%rAttr(index_x2i_So_dhdx,n)
01055 aflds(i,j, 6,iblk) = x2i_i%rAttr(index_x2i_So_dhdy,n)
01056 enddo
01057 enddo
01058 enddo
01059
01060 if (.not.prescribed_ice) then
01061 call t_startf ('cice_imp_halo')
01062 call ice_HaloUpdate(aflds, halo_info, field_loc_center, &
01063 field_type_vector)
01064 call t_stopf ('cice_imp_halo')
01065 endif
01066
01067
01068 do iblk = 1, nblocks
01069 do j = 1,ny_block
01070 do i = 1,nx_block
01071 uocn (i,j,iblk) = aflds(i,j, 1,iblk)
01072 vocn (i,j,iblk) = aflds(i,j, 2,iblk)
01073 uatm (i,j,iblk) = aflds(i,j, 3,iblk)
01074 vatm (i,j,iblk) = aflds(i,j, 4,iblk)
01075 ss_tltx(i,j,iblk) = aflds(i,j, 5,iblk)
01076 ss_tlty(i,j,iblk) = aflds(i,j, 6,iblk)
01077 enddo
01078 enddo
01079 enddo
01080
01081
01082 deallocate(aflds)
01083
01084
01085
01086
01087
01088
01089
01090 if (tr_aero .and. first_call) then
01091 allocate(work(nx_block,ny_block,max_blocks))
01092 n=0
01093 do iblk = 1, nblocks
01094 this_block = get_block(blocks_ice(iblk),iblk)
01095 ilo = this_block%ilo
01096 ihi = this_block%ihi
01097 jlo = this_block%jlo
01098 jhi = this_block%jhi
01099 do j = jlo, jhi
01100 do i = ilo, ihi
01101 n = n+1
01102 work(i,j,iblk) = x2i_i%rAttr(index_x2i_Faxa_bcphodry,n)
01103 end do
01104 end do
01105 end do
01106 maxwork = global_maxval(work,distrb_info)
01107 call broadcast_scalar(maxwork, master_task)
01108 if (abs(maxwork) < c100) then
01109 prescribed_aero = .false.
01110 else
01111 prescribed_aero_tmp = .true.
01112 call ice_prescaero_init(prescribed_aero_tmp)
01113 end if
01114 deallocate(work)
01115 first_call = .false.
01116 endif
01117
01118 if (tr_aero .and. .not. prescribed_aero) then
01119
01120 n=0
01121 do iblk = 1, nblocks
01122 this_block = get_block(blocks_ice(iblk),iblk)
01123 ilo = this_block%ilo
01124 ihi = this_block%ihi
01125 jlo = this_block%jlo
01126 jhi = this_block%jhi
01127
01128 do j = jlo, jhi
01129 do i = ilo, ihi
01130
01131 n = n+1
01132 faero(i,j,1,iblk) = x2i_i%rAttr(index_x2i_Faxa_bcphodry,n)
01133
01134 faero(i,j,2,iblk) = x2i_i%rAttr(index_x2i_Faxa_bcphidry,n) &
01135 + x2i_i%rAttr(index_x2i_Faxa_bcphiwet,n)
01136
01137 faero(i,j,3,iblk) = x2i_i%rAttr(index_x2i_Faxa_dstwet1,n) &
01138 + x2i_i%rAttr(index_x2i_Faxa_dstdry1,n) &
01139 + x2i_i%rAttr(index_x2i_Faxa_dstwet2,n) &
01140 + x2i_i%rAttr(index_x2i_Faxa_dstdry2,n) &
01141 + x2i_i%rAttr(index_x2i_Faxa_dstwet3,n) &
01142 + x2i_i%rAttr(index_x2i_Faxa_dstdry3,n) &
01143 + x2i_i%rAttr(index_x2i_Faxa_dstwet4,n) &
01144 + x2i_i%rAttr(index_x2i_Faxa_dstdry4,n)
01145
01146 enddo
01147 enddo
01148
01149 enddo
01150
01151 endif
01152
01153
01154
01155
01156
01157
01158
01159
01160
01161
01162
01163
01164 call t_startf ('cice_imp_ocn')
01165
01166 do iblk = 1, nblocks
01167
01168 do j = 1,ny_block
01169 do i = 1,nx_block
01170
01171
01172 workx = uocn (i,j,iblk)
01173 worky = vocn (i,j,iblk)
01174 uocn(i,j,iblk) = workx*cos(ANGLET(i,j,iblk)) &
01175 + worky*sin(ANGLET(i,j,iblk))
01176 vocn(i,j,iblk) = worky*cos(ANGLET(i,j,iblk)) &
01177 - workx*sin(ANGLET(i,j,iblk))
01178
01179 workx = ss_tltx (i,j,iblk)
01180 worky = ss_tlty (i,j,iblk)
01181 ss_tltx(i,j,iblk) = workx*cos(ANGLET(i,j,iblk)) &
01182 + worky*sin(ANGLET(i,j,iblk))
01183 ss_tlty(i,j,iblk) = worky*cos(ANGLET(i,j,iblk)) &
01184 - workx*sin(ANGLET(i,j,iblk))
01185
01186 sst(i,j,iblk) = sst(i,j,iblk) - Tffresh
01187 Tf (i,j,iblk) = -1.8_dbl_kind
01188
01189
01190
01191 enddo
01192 enddo
01193 enddo
01194
01195 call t_stopf ('cice_imp_ocn')
01196
01197
01198
01199
01200 if (.not.prescribed_ice) then
01201 call t_startf ('cice_imp_t2u')
01202 call t2ugrid_vector(uocn)
01203 call t2ugrid_vector(vocn)
01204 call t2ugrid_vector(ss_tltx)
01205 call t2ugrid_vector(ss_tlty)
01206 call t_stopf ('cice_imp_t2u')
01207 end if
01208
01209
01210
01211
01212
01213 call t_startf ('cice_imp_atm')
01214
01215 do iblk = 1, nblocks
01216 do j = 1, ny_block
01217 do i = 1, nx_block
01218
01219
01220 workx = uatm(i,j,iblk)
01221 worky = vatm(i,j,iblk)
01222 uatm (i,j,iblk) = workx*cos(ANGLET(i,j,iblk)) &
01223 + worky*sin(ANGLET(i,j,iblk))
01224 vatm (i,j,iblk) = worky*cos(ANGLET(i,j,iblk)) &
01225 - workx*sin(ANGLET(i,j,iblk))
01226
01227 wind (i,j,iblk) = sqrt(uatm(i,j,iblk)**2 + vatm(i,j,iblk)**2)
01228 fsw (i,j,iblk) = swvdr(i,j,iblk) + swvdf(i,j,iblk) &
01229 + swidr(i,j,iblk) + swidf(i,j,iblk)
01230 enddo
01231 enddo
01232 enddo
01233
01234 call t_stopf ('cice_imp_atm')
01235
01236 end subroutine ice_import_mct
01237
01238
01239
01240 subroutine ice_domain_mct( lsize, gsMap_i, dom_i )
01241
01242
01243
01244
01245
01246 integer , intent(in) :: lsize
01247 type(mct_gsMap), intent(in) :: gsMap_i
01248 type(mct_ggrid), intent(inout) :: dom_i
01249
01250
01251
01252 integer :: i, j, iblk, n, gi
01253 integer :: ilo, ihi, jlo, jhi
01254 real(dbl_kind), pointer :: work_dom(:)
01255 real(dbl_kind), pointer :: data(:)
01256 integer , pointer :: idata(:)
01257 type(block) :: this_block
01258
01259
01260
01261
01262
01263 call mct_gGrid_init( GGrid=dom_i, CoordChars=trim(seq_flds_dom_coord), &
01264 OtherChars=trim(seq_flds_dom_other), lsize=lsize )
01265 call mct_aVect_zero(dom_i%data)
01266
01267 allocate(data(lsize))
01268
01269
01270
01271 call mct_gsMap_orderedPoints(gsMap_i, my_task, idata)
01272 call mct_gGrid_importIAttr(dom_i,'GlobGridNum',idata,lsize)
01273
01274
01275
01276
01277 data(:) = -9999.0_R8
01278 call mct_gGrid_importRAttr(dom_i,"lat" ,data,lsize)
01279 call mct_gGrid_importRAttr(dom_i,"lon" ,data,lsize)
01280 call mct_gGrid_importRAttr(dom_i,"area" ,data,lsize)
01281 call mct_gGrid_importRAttr(dom_i,"aream",data,lsize)
01282 data(:) = 0.0_R8
01283 call mct_gGrid_importRAttr(dom_i,"mask",data,lsize)
01284 call mct_gGrid_importRAttr(dom_i,"frac",data,lsize)
01285
01286
01287
01288 allocate(work_dom(lsize))
01289 work_dom(:) = 0.0_dbl_kind
01290
01291 n=0
01292 do iblk = 1, nblocks
01293 this_block = get_block(blocks_ice(iblk),iblk)
01294 ilo = this_block%ilo
01295 ihi = this_block%ihi
01296 jlo = this_block%jlo
01297 jhi = this_block%jhi
01298
01299 do j = jlo, jhi
01300 do i = ilo, ihi
01301 n = n+1
01302 data(n) = TLON(i,j,iblk)*rad_to_deg
01303 enddo
01304 enddo
01305 enddo
01306 call mct_gGrid_importRattr(dom_i,"lon",data,lsize)
01307
01308 n=0
01309 do iblk = 1, nblocks
01310 this_block = get_block(blocks_ice(iblk),iblk)
01311 ilo = this_block%ilo
01312 ihi = this_block%ihi
01313 jlo = this_block%jlo
01314 jhi = this_block%jhi
01315
01316 do j = jlo, jhi
01317 do i = ilo, ihi
01318 n = n+1
01319 data(n) = TLAT(i,j,iblk)*rad_to_deg
01320 enddo
01321 enddo
01322 enddo
01323 call mct_gGrid_importRattr(dom_i,"lat",data,lsize)
01324
01325 n=0
01326 do iblk = 1, nblocks
01327 this_block = get_block(blocks_ice(iblk),iblk)
01328 ilo = this_block%ilo
01329 ihi = this_block%ihi
01330 jlo = this_block%jlo
01331 jhi = this_block%jhi
01332
01333 do j = jlo, jhi
01334 do i = ilo, ihi
01335 n = n+1
01336 data(n) = tarea(i,j,iblk)/(radius*radius)
01337 enddo
01338 enddo
01339 enddo
01340 call mct_gGrid_importRattr(dom_i,"area",data,lsize)
01341
01342 n=0
01343 do iblk = 1, nblocks
01344 this_block = get_block(blocks_ice(iblk),iblk)
01345 ilo = this_block%ilo
01346 ihi = this_block%ihi
01347 jlo = this_block%jlo
01348 jhi = this_block%jhi
01349
01350 do j = jlo, jhi
01351 do i = ilo, ihi
01352 n = n+1
01353 data(n) = real(nint(hm(i,j,iblk)),kind=dbl_kind)
01354 enddo
01355 enddo
01356 enddo
01357 call mct_gGrid_importRattr(dom_i,"mask",data,lsize)
01358
01359 n=0
01360 do iblk = 1, nblocks
01361 this_block = get_block(blocks_ice(iblk),iblk)
01362 ilo = this_block%ilo
01363 ihi = this_block%ihi
01364 jlo = this_block%jlo
01365 jhi = this_block%jhi
01366
01367 do j = jlo, jhi
01368 do i = ilo, ihi
01369 n = n+1
01370 if (trim(grid_type) == 'latlon') then
01371 data(n) = ocn_gridcell_frac(i,j,iblk)
01372 else
01373 data(n) = real(nint(hm(i,j,iblk)),kind=dbl_kind)
01374 end if
01375 enddo
01376 enddo
01377 enddo
01378 call mct_gGrid_importRattr(dom_i,"frac",data,lsize)
01379
01380 deallocate(data)
01381 deallocate(idata)
01382 deallocate(work_dom)
01383
01384 end subroutine ice_domain_mct
01385
01386
01387
01388
01389
01390
01391
01392 character(len=char_len_long) function restart_filename( yr_spec, mon_spec, day_spec, sec_spec )
01393
01394
01395
01396
01397
01398
01399
01400
01401
01402
01403
01404
01405
01406
01407 use ice_restart, only: runid
01408
01409
01410 integer , intent(in) :: yr_spec
01411 integer , intent(in) :: mon_spec
01412 integer , intent(in) :: day_spec
01413 integer , intent(in) :: sec_spec
01414
01415
01416
01417 integer :: i, n
01418 integer :: year
01419 integer :: month
01420 integer :: day
01421 integer :: ncsec
01422 character(len=char_len_long) :: string
01423 character(len=char_len_long) :: format
01424 character(len=char_len_long) :: filename_spec = '%c.cice.r.%y-%m-%d-%s'
01425
01426
01427
01428
01429
01430 if ( len_trim(filename_spec) == 0 )then
01431 call shr_sys_abort ('restart_filename: filename specifier is empty')
01432 end if
01433 if ( index(trim(filename_spec)," ") /= 0 )then
01434 call shr_sys_abort ('restart_filename: filename specifier can not contain a space:'//trim(filename_spec))
01435 end if
01436
01437 year = yr_spec
01438 month = mon_spec
01439 day = day_spec
01440 ncsec = sec_spec
01441
01442
01443
01444 i = 1
01445 restart_filename = ''
01446 do while ( i <= len_trim(filename_spec) )
01447 if ( filename_spec(i:i) == "%" )then
01448 i = i + 1
01449 select case( filename_spec(i:i) )
01450 case( 'c' )
01451 string = trim(runid)
01452 case( 'y' )
01453 if ( year > 99999 ) then
01454 format = '(i6.6)'
01455 else if ( year > 9999 ) then
01456 format = '(i5.5)'
01457 else
01458 format = '(i4.4)'
01459 end if
01460 write(string,format) year
01461 case( 'm' )
01462 write(string,'(i2.2)') month
01463 case( 'd' )
01464 write(string,'(i2.2)') day
01465 case( 's' )
01466 write(string,'(i5.5)') ncsec
01467 case( '%' )
01468 string = "%"
01469 case default
01470 call shr_sys_abort ('restart_filename: Invalid expansion character: '//filename_spec(i:i))
01471 end select
01472 else
01473 n = index( filename_spec(i:), "%" )
01474 if ( n == 0 ) n = len_trim( filename_spec(i:) ) + 1
01475 if ( n == 0 ) exit
01476 string = filename_spec(i:n+i-2)
01477 i = n + i - 2
01478 end if
01479 if ( len_trim(restart_filename) == 0 )then
01480 restart_filename = trim(string)
01481 else
01482 if ( (len_trim(restart_filename)+len_trim(string)) >= char_len_long )then
01483 call shr_sys_abort ('restart_filename Resultant filename too long')
01484 end if
01485 restart_filename = trim(restart_filename) // trim(string)
01486 end if
01487 i = i + 1
01488 end do
01489 if ( len_trim(restart_filename) == 0 )then
01490 call shr_sys_abort ('restart_filename: Resulting filename is empty')
01491 end if
01492
01493 end function restart_filename
01494
01495 end module ice_comp_mct
01496