From 52d2b10b54e2bd4be3f17a8f43834d1020df5919 Mon Sep 17 00:00:00 2001 From: apcraig Date: Wed, 23 Oct 2019 16:43:17 -0600 Subject: [PATCH 1/4] initial CICE6 port to CESM+NUOPC --- cicecore/cicedynB/general/ice_init.F90 | 14 +- .../io/io_pio2/ice_history_write.F90 | 1064 ++++++++++++++++ .../infrastructure/io/io_pio2/ice_pio.F90 | 365 ++++++ .../infrastructure/io/io_pio2/ice_restart.F90 | 876 +++++++++++++ cicecore/drivers/nuopc/CICE_FinalMod.F90 | 17 +- cicecore/drivers/nuopc/CICE_InitMod.F90 | 462 ++++--- cicecore/drivers/nuopc/CICE_RunMod.F90 | 1112 +++++++++-------- cicecore/drivers/nuopc/CICE_copyright.txt | 17 + cicecore/drivers/nuopc/ice_comp_nuopc.F90 | 99 +- cicecore/drivers/nuopc/ice_constants.F90 | 229 ---- cicecore/drivers/nuopc/ice_import_export.F90 | 47 +- cicecore/drivers/nuopc/ice_prescribed_mod.F90 | 121 +- cicecore/drivers/nuopc/ice_scam.F90 | 1 - cicecore/drivers/nuopc/ice_shr_methods.F90 | 8 +- cicecore/shared/ice_fileunits.F90 | 6 +- 15 files changed, 3464 insertions(+), 974 deletions(-) create mode 100644 cicecore/cicedynB/infrastructure/io/io_pio2/ice_history_write.F90 create mode 100644 cicecore/cicedynB/infrastructure/io/io_pio2/ice_pio.F90 create mode 100644 cicecore/cicedynB/infrastructure/io/io_pio2/ice_restart.F90 create mode 100644 cicecore/drivers/nuopc/CICE_copyright.txt delete mode 100644 cicecore/drivers/nuopc/ice_constants.F90 diff --git a/cicecore/cicedynB/general/ice_init.F90 b/cicecore/cicedynB/general/ice_init.F90 index 1d3aedf6e..4beb2c3e1 100644 --- a/cicecore/cicedynB/general/ice_init.F90 +++ b/cicecore/cicedynB/general/ice_init.F90 @@ -466,11 +466,15 @@ subroutine input_data history_file = trim(runid) // ".cice" // trim(inst_suffix) //".h" restart_file = trim(runid) // ".cice" // trim(inst_suffix) //".r" incond_file = trim(runid) // ".cice" // trim(inst_suffix) //".i" - inquire(file='ice_modelio.nml'//trim(inst_suffix),exist=exists) - if (exists) then - call get_fileUnit(nu_diag) - call shr_file_setIO('ice_modelio.nml'//trim(inst_suffix),nu_diag) - end if + ! Note the nuopc cap will set nu_diag before this point - so just + ! need to check that it is non-zero first + if (nu_diag == ice_stdout) then + inquire(file='ice_modelio.nml'//trim(inst_suffix),exist=exists) + if (exists) then + call get_fileUnit(nu_diag) + call shr_file_setIO('ice_modelio.nml'//trim(inst_suffix),nu_diag) + end if + endif else ! each task gets unique ice log filename when if test is true, for debugging if (1 == 0) then diff --git a/cicecore/cicedynB/infrastructure/io/io_pio2/ice_history_write.F90 b/cicecore/cicedynB/infrastructure/io/io_pio2/ice_history_write.F90 new file mode 100644 index 000000000..771d0e313 --- /dev/null +++ b/cicecore/cicedynB/infrastructure/io/io_pio2/ice_history_write.F90 @@ -0,0 +1,1064 @@ +!======================================================================= +! +! Writes history in netCDF format +! +! authors Tony Craig and Bruce Briegleb, NCAR +! Elizabeth C. Hunke and William H. Lipscomb, LANL +! C. M. Bitz, UW +! +! 2004 WHL: Block structure added +! 2006 ECH: Accepted some CESM code into mainstream CICE +! Added ice_present, aicen, vicen; removed aice1...10, vice1...1. +! Added histfreq_n and histfreq='h' options, removed histfreq='w' +! Converted to free source form (F90) +! Added option for binary output instead of netCDF +! 2009 D Bailey and ECH: Generalized for multiple frequency output +! 2010 Alison McLaren and ECH: Added 3D capability +! + module ice_history_write + + use ice_kinds_mod + use ice_fileunits, only: nu_diag + use ice_exit, only: abort_ice + use icepack_intfc, only: icepack_warnings_flush, icepack_warnings_aborted + use icepack_intfc, only: icepack_query_parameters + + implicit none + private + public :: ice_write_hist + +!======================================================================= + + contains + +!======================================================================= +! +! write average ice quantities or snapshots +! +! author: Elizabeth C. Hunke, LANL + + subroutine ice_write_hist (ns) + +#ifdef ncdf + use ice_blocks, only: nx_block, ny_block + use ice_broadcast, only: broadcast_scalar + use ice_calendar, only: time, sec, idate, idate0, write_ic, & + histfreq, dayyr, days_per_year, use_leap_years + use ice_communicate, only: my_task, master_task + use ice_constants, only: c0, c360, spval, spval_dbl + use ice_domain, only: distrb_info, nblocks + use ice_domain_size, only: nx_global, ny_global, max_blocks, max_nstrm + use ice_gather_scatter, only: gather_global + use ice_grid, only: TLON, TLAT, ULON, ULAT, hm, bm, tarea, uarea, & + dxu, dxt, dyu, dyt, HTN, HTE, ANGLE, ANGLET, tmask, & + lont_bounds, latt_bounds, lonu_bounds, latu_bounds + use ice_history_shared + use ice_arrays_column, only: hin_max + use ice_restart_shared, only: runid, lcdf64 + use netcdf +#endif + use ice_pio + use pio + + integer (kind=int_kind), intent(in) :: ns + + ! local variables + +#ifdef ncdf + integer (kind=int_kind) :: i,j,k,ic,n,nn, & + ncid,status,imtid,jmtid,kmtidi,kmtids,kmtidb, cmtid,timid, & + length,nvertexid,ivertex,kmtida + integer (kind=int_kind), dimension(2) :: dimid2 + integer (kind=int_kind), dimension(3) :: dimid3 + integer (kind=int_kind), dimension(4) :: dimidz + integer (kind=int_kind), dimension(5) :: dimidcz + integer (kind=int_kind), dimension(3) :: dimid_nverts + integer (kind=int_kind), dimension(5) :: dimidex + real (kind=real_kind) :: ltime + real (kind= dbl_kind) :: ltime2 + character (char_len) :: title + character (char_len_long) :: ncfile(max_nstrm) + + integer (kind=int_kind) :: iyear, imonth, iday + integer (kind=int_kind) :: icategory,ind,i_aice,boundid + + character (char_len) :: start_time,current_date,current_time + character (len=16) :: c_aice + character (len=8) :: cdate + + type(file_desc_t) :: File + type(io_desc_t) :: iodesc2d, & + iodesc3dc, iodesc3dv, iodesc3di, iodesc3db, iodesc3da, & + iodesc4di, iodesc4ds + type(var_desc_t) :: varid + + ! 4 coordinate variables: TLON, TLAT, ULON, ULAT + INTEGER (kind=int_kind), PARAMETER :: ncoord = 4 + + ! 4 vertices in each grid cell + INTEGER (kind=int_kind), PARAMETER :: nverts = 4 + + ! 4 variables describe T, U grid boundaries: + ! lont_bounds, latt_bounds, lonu_bounds, latu_bounds + INTEGER (kind=int_kind), PARAMETER :: nvar_verts = 4 + + TYPE coord_attributes ! netcdf coordinate attributes + character (len=11) :: short_name + character (len=45) :: long_name + character (len=20) :: units + END TYPE coord_attributes + + TYPE req_attributes ! req'd netcdf attributes + type (coord_attributes) :: req + character (len=20) :: coordinates + END TYPE req_attributes + + TYPE(req_attributes), dimension(nvar) :: var + TYPE(coord_attributes), dimension(ncoord) :: coord_var + TYPE(coord_attributes), dimension(nvar_verts) :: var_nverts + TYPE(coord_attributes), dimension(nvarz) :: var_nz + CHARACTER (char_len), dimension(ncoord) :: coord_bounds + + real (kind=dbl_kind), allocatable :: workr2(:,:,:) + real (kind=dbl_kind), allocatable :: workr3(:,:,:,:) + real (kind=dbl_kind), allocatable :: workr4(:,:,:,:,:) + real (kind=dbl_kind), allocatable :: workr3v(:,:,:,:) + + character(len=char_len_long) :: & + filename + + integer (kind=int_kind), dimension(1) :: & + tim_start,tim_length ! dimension quantities for netCDF + + integer (kind=int_kind), dimension(2) :: & + bnd_start,bnd_length ! dimension quantities for netCDF + + real (kind=dbl_kind) :: secday + real (kind=dbl_kind) :: rad_to_deg + + character(len=*), parameter :: subname = '(ice_write_hist)' + + call icepack_query_parameters(secday_out=secday) + call icepack_query_parameters(rad_to_deg_out=rad_to_deg) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + + if (my_task == master_task) then + call construct_filename(ncfile(ns),'nc',ns) + + ! add local directory path name to ncfile + if (write_ic) then + ncfile(ns) = trim(incond_dir)//ncfile(ns) + else + ncfile(ns) = trim(history_dir)//ncfile(ns) + endif + filename = ncfile(ns) + end if + call broadcast_scalar(filename, master_task) + + ! create file + + File%fh=-1 + call ice_pio_init(mode='write', filename=trim(filename), File=File, & + clobber=.true., cdf64=lcdf64) + + call ice_pio_initdecomp(iodesc=iodesc2d) + call ice_pio_initdecomp(ndim3=ncat_hist, iodesc=iodesc3dc) + call ice_pio_initdecomp(ndim3=nzilyr, iodesc=iodesc3di) + call ice_pio_initdecomp(ndim3=nzblyr, iodesc=iodesc3db) + call ice_pio_initdecomp(ndim3=nzalyr, iodesc=iodesc3da) + call ice_pio_initdecomp(ndim3=nverts, inner_dim=.true., iodesc=iodesc3dv) + call ice_pio_initdecomp(ndim3=nzilyr, ndim4=ncat_hist, iodesc=iodesc4di) + call ice_pio_initdecomp(ndim3=nzslyr, ndim4=ncat_hist, iodesc=iodesc4ds) + + ltime2 = time/int(secday) + ltime = real(time/int(secday),kind=real_kind) + + !----------------------------------------------------------------- + ! define dimensions + !----------------------------------------------------------------- + + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_def_dim(File,'d2',2,boundid) + endif + + status = pio_def_dim(File,'ni',nx_global,imtid) + status = pio_def_dim(File,'nj',ny_global,jmtid) + status = pio_def_dim(File,'nc',ncat_hist,cmtid) + status = pio_def_dim(File,'nkice',nzilyr,kmtidi) + status = pio_def_dim(File,'nksnow',nzslyr,kmtids) + status = pio_def_dim(File,'nkbio',nzblyr,kmtidb) + status = pio_def_dim(File,'nkaer',nzalyr,kmtida) + status = pio_def_dim(File,'time',PIO_UNLIMITED,timid) + status = pio_def_dim(File,'nvertices',nverts,nvertexid) + + !----------------------------------------------------------------- + ! define coordinate variables: time, time_bounds + !----------------------------------------------------------------- + +!sgl status = pio_def_var(File,'time',pio_real,(/timid/),varid) + status = pio_def_var(File,'time',pio_double,(/timid/),varid) + status = pio_put_att(File,varid,'long_name','model time') + + write(cdate,'(i8.8)') idate0 + write(title,'(a,a,a,a,a,a,a)') 'days since ', & + cdate(1:4),'-',cdate(5:6),'-',cdate(7:8),' 00:00:00' + status = pio_put_att(File,varid,'units',trim(title)) + + if (days_per_year == 360) then + status = pio_put_att(File,varid,'calendar','360_day') + elseif (days_per_year == 365 .and. .not.use_leap_years ) then + status = pio_put_att(File,varid,'calendar','NoLeap') + elseif (use_leap_years) then + status = pio_put_att(File,varid,'calendar','Gregorian') + else + call abort_ice(subname//'ERROR: invalid calendar settings') + endif + + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_put_att(File,varid,'bounds','time_bounds') + endif + + ! Define attributes for time_bounds if hist_avg is true + if (hist_avg .and. histfreq(ns) /= '1') then + dimid2(1) = boundid + dimid2(2) = timid +!sgl status = pio_def_var(File,'time_bounds',pio_real,dimid2,varid) + status = pio_def_var(File,'time_bounds',pio_double,dimid2,varid) + status = pio_put_att(File,varid,'long_name', & + 'boundaries for time-averaging interval') + write(cdate,'(i8.8)') idate0 + write(title,'(a,a,a,a,a,a,a,a)') 'days since ', & + cdate(1:4),'-',cdate(5:6),'-',cdate(7:8),' 00:00:00' + status = pio_put_att(File,varid,'units',trim(title)) + endif + + !----------------------------------------------------------------- + ! define information for required time-invariant variables + !----------------------------------------------------------------- + + ind = 0 + ind = ind + 1 + coord_var(ind) = coord_attributes('TLON', & + 'T grid center longitude', 'degrees_east') + coord_bounds(ind) = 'lont_bounds' + ind = ind + 1 + coord_var(ind) = coord_attributes('TLAT', & + 'T grid center latitude', 'degrees_north') + coord_bounds(ind) = 'latt_bounds' + ind = ind + 1 + coord_var(ind) = coord_attributes('ULON', & + 'U grid center longitude', 'degrees_east') + coord_bounds(ind) = 'lonu_bounds' + ind = ind + 1 + coord_var(ind) = coord_attributes('ULAT', & + 'U grid center latitude', 'degrees_north') + coord_bounds(ind) = 'latu_bounds' + + var_nz(1) = coord_attributes('NCAT', 'category maximum thickness', 'm') + var_nz(2) = coord_attributes('VGRDi', 'vertical ice levels', '1') + var_nz(3) = coord_attributes('VGRDs', 'vertical snow levels', '1') + var_nz(4) = coord_attributes('VGRDb', 'vertical ice-bio levels', '1') + var_nz(5) = coord_attributes('VGRDa', 'vertical snow-ice-bio levels', '1') + + !----------------------------------------------------------------- + ! define information for optional time-invariant variables + !----------------------------------------------------------------- + + var(n_tmask)%req = coord_attributes('tmask', & + 'ocean grid mask', ' ') + var(n_tmask)%coordinates = 'TLON TLAT' + + var(n_blkmask)%req = coord_attributes('blkmask', & + 'ice grid block mask', ' ') + var(n_blkmask)%coordinates = 'TLON TLAT' + + var(n_tarea)%req = coord_attributes('tarea', & + 'area of T grid cells', 'm^2') + var(n_tarea)%coordinates = 'TLON TLAT' + + var(n_uarea)%req = coord_attributes('uarea', & + 'area of U grid cells', 'm^2') + var(n_uarea)%coordinates = 'ULON ULAT' + var(n_dxt)%req = coord_attributes('dxt', & + 'T cell width through middle', 'm') + var(n_dxt)%coordinates = 'TLON TLAT' + var(n_dyt)%req = coord_attributes('dyt', & + 'T cell height through middle', 'm') + var(n_dyt)%coordinates = 'TLON TLAT' + var(n_dxu)%req = coord_attributes('dxu', & + 'U cell width through middle', 'm') + var(n_dxu)%coordinates = 'ULON ULAT' + var(n_dyu)%req = coord_attributes('dyu', & + 'U cell height through middle', 'm') + var(n_dyu)%coordinates = 'ULON ULAT' + var(n_HTN)%req = coord_attributes('HTN', & + 'T cell width on North side','m') + var(n_HTN)%coordinates = 'TLON TLAT' + var(n_HTE)%req = coord_attributes('HTE', & + 'T cell width on East side', 'm') + var(n_HTE)%coordinates = 'TLON TLAT' + var(n_ANGLE)%req = coord_attributes('ANGLE', & + 'angle grid makes with latitude line on U grid', & + 'radians') + var(n_ANGLE)%coordinates = 'ULON ULAT' + var(n_ANGLET)%req = coord_attributes('ANGLET', & + 'angle grid makes with latitude line on T grid', & + 'radians') + var(n_ANGLET)%coordinates = 'TLON TLAT' + + ! These fields are required for CF compliance + ! dimensions (nx,ny,nverts) + var_nverts(n_lont_bnds) = coord_attributes('lont_bounds', & + 'longitude boundaries of T cells', 'degrees_east') + var_nverts(n_latt_bnds) = coord_attributes('latt_bounds', & + 'latitude boundaries of T cells', 'degrees_north') + var_nverts(n_lonu_bnds) = coord_attributes('lonu_bounds', & + 'longitude boundaries of U cells', 'degrees_east') + var_nverts(n_latu_bnds) = coord_attributes('latu_bounds', & + 'latitude boundaries of U cells', 'degrees_north') + + !----------------------------------------------------------------- + ! define attributes for time-invariant variables + !----------------------------------------------------------------- + + dimid2(1) = imtid + dimid2(2) = jmtid + + do i = 1, ncoord + status = pio_def_var(File, trim(coord_var(i)%short_name), pio_real, & + dimid2, varid) + status = pio_put_att(File,varid,'long_name',trim(coord_var(i)%long_name)) + status = pio_put_att(File, varid, 'units', trim(coord_var(i)%units)) + status = pio_put_att(File, varid, 'missing_value', spval) + status = pio_put_att(File, varid,'_FillValue',spval) + if (coord_var(i)%short_name == 'ULAT') then + status = pio_put_att(File,varid,'comment', & + trim('Latitude of NE corner of T grid cell')) + endif + if (f_bounds) then + status = pio_put_att(File, varid, 'bounds', trim(coord_bounds(i))) + endif + enddo + + ! Extra dimensions (NCAT, NZILYR, NZSLYR, NZBLYR, NZALYR) + dimidex(1)=cmtid + dimidex(2)=kmtidi + dimidex(3)=kmtids + dimidex(4)=kmtidb + dimidex(5)=kmtida + + do i = 1, nvarz + if (igrdz(i)) then + status = pio_def_var(File, trim(var_nz(i)%short_name), pio_real, & + (/dimidex(i)/), varid) + status = pio_put_att(File, varid, 'long_name', var_nz(i)%long_name) + status = pio_put_att(File, varid, 'units' , var_nz(i)%units) + endif + enddo + + ! Attributes for tmask defined separately, since it has no units + if (igrd(n_tmask)) then + status = pio_def_var(File, 'tmask', pio_real, dimid2, varid) + status = pio_put_att(File,varid, 'long_name', 'ocean grid mask') + status = pio_put_att(File, varid, 'coordinates', 'TLON TLAT') + status = pio_put_att(File, varid, 'missing_value', spval) + status = pio_put_att(File, varid,'_FillValue',spval) + status = pio_put_att(File,varid,'comment', '0 = land, 1 = ocean') + endif + if (igrd(n_blkmask)) then + status = pio_def_var(File, 'blkmask', pio_real, dimid2, varid) + status = pio_put_att(File,varid, 'long_name', 'ice grid block mask') + status = pio_put_att(File, varid, 'coordinates', 'TLON TLAT') + status = pio_put_att(File,varid,'comment', 'mytask + iblk/100') + status = pio_put_att(File, varid, 'missing_value', spval) + status = pio_put_att(File, varid,'_FillValue',spval) + endif + + do i = 3, nvar ! note: n_tmask=1, n_blkmask=2 + if (igrd(i)) then + status = pio_def_var(File, trim(var(i)%req%short_name), & + pio_real, dimid2, varid) + status = pio_put_att(File,varid, 'long_name', trim(var(i)%req%long_name)) + status = pio_put_att(File, varid, 'units', trim(var(i)%req%units)) + status = pio_put_att(File, varid, 'coordinates', trim(var(i)%coordinates)) + status = pio_put_att(File, varid, 'missing_value', spval) + status = pio_put_att(File, varid,'_FillValue',spval) + endif + enddo + + ! Fields with dimensions (nverts,nx,ny) + dimid_nverts(1) = nvertexid + dimid_nverts(2) = imtid + dimid_nverts(3) = jmtid + do i = 1, nvar_verts + if (f_bounds) then + status = pio_def_var(File, trim(var_nverts(i)%short_name), & + pio_real,dimid_nverts, varid) + status = & + pio_put_att(File,varid, 'long_name', trim(var_nverts(i)%long_name)) + status = & + pio_put_att(File, varid, 'units', trim(var_nverts(i)%units)) + status = pio_put_att(File, varid, 'missing_value', spval) + status = pio_put_att(File, varid,'_FillValue',spval) + endif + enddo + + !----------------------------------------------------------------- + ! define attributes for time-variant variables + !----------------------------------------------------------------- + + !----------------------------------------------------------------- + ! 2D + !----------------------------------------------------------------- + + dimid3(1) = imtid + dimid3(2) = jmtid + dimid3(3) = timid + + do n=1,num_avail_hist_fields_2D + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & + pio_real, dimid3, varid) + status = pio_put_att(File,varid,'units', & + trim(avail_hist_fields(n)%vunit)) + status = pio_put_att(File,varid, 'long_name', & + trim(avail_hist_fields(n)%vdesc)) + status = pio_put_att(File,varid,'coordinates', & + trim(avail_hist_fields(n)%vcoord)) + status = pio_put_att(File,varid,'cell_measures', & + trim(avail_hist_fields(n)%vcellmeas)) + status = pio_put_att(File,varid,'missing_value',spval) + status = pio_put_att(File,varid,'_FillValue',spval) + + ! Add cell_methods attribute to variables if averaged + if (hist_avg .and. histfreq(ns) /= '1') then + if (TRIM(avail_hist_fields(n)%vname)/='sig1' & + .or.TRIM(avail_hist_fields(n)%vname)/='sig2' & + .or.TRIM(avail_hist_fields(n)%vname)/='sistreave' & + .or.TRIM(avail_hist_fields(n)%vname)/='sistremax' & + .or.TRIM(avail_hist_fields(n)%vname)/='sigP') then + status = pio_put_att(File,varid,'cell_methods','time: mean') + endif + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg & + .or. n==n_divu(ns) .or. n==n_shear(ns) & ! snapshots + .or. n==n_sig1(ns) .or. n==n_sig2(ns) & + .or. n==n_sigP(ns) .or. n==n_trsig(ns) & + .or. n==n_sistreave(ns) .or. n==n_sistremax(ns) & + .or. n==n_mlt_onset(ns) .or. n==n_frz_onset(ns) & + .or. n==n_hisnap(ns) .or. n==n_aisnap(ns)) then + status = pio_put_att(File,varid,'time_rep','instantaneous') + else + status = pio_put_att(File,varid,'time_rep','averaged') + endif + endif + enddo ! num_avail_hist_fields_2D + + !----------------------------------------------------------------- + ! 3D (category) + !----------------------------------------------------------------- + + dimidz(1) = imtid + dimidz(2) = jmtid + dimidz(3) = cmtid + dimidz(4) = timid + + do n = n2D + 1, n3Dccum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & + pio_real, dimidz, varid) + status = pio_put_att(File,varid,'units', & + trim(avail_hist_fields(n)%vunit)) + status = pio_put_att(File,varid, 'long_name', & + trim(avail_hist_fields(n)%vdesc)) + status = pio_put_att(File,varid,'coordinates', & + trim(avail_hist_fields(n)%vcoord)) + status = pio_put_att(File,varid,'cell_measures', & + trim(avail_hist_fields(n)%vcellmeas)) + status = pio_put_att(File,varid,'missing_value',spval) + status = pio_put_att(File,varid,'_FillValue',spval) + + ! Add cell_methods attribute to variables if averaged + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_put_att(File,varid,'cell_methods','time: mean') + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg) then + status = pio_put_att(File,varid,'time_rep','instantaneous') + else + status = pio_put_att(File,varid,'time_rep','averaged') + endif + endif + enddo ! num_avail_hist_fields_3Dc + + !----------------------------------------------------------------- + ! 3D (ice layers) + !----------------------------------------------------------------- + + dimidz(1) = imtid + dimidz(2) = jmtid + dimidz(3) = kmtidi + dimidz(4) = timid + + do n = n3Dccum + 1, n3Dzcum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & + pio_real, dimidz, varid) + status = pio_put_att(File,varid,'units', & + trim(avail_hist_fields(n)%vunit)) + status = pio_put_att(File,varid, 'long_name', & + trim(avail_hist_fields(n)%vdesc)) + status = pio_put_att(File,varid,'coordinates', & + trim(avail_hist_fields(n)%vcoord)) + status = pio_put_att(File,varid,'cell_measures', & + trim(avail_hist_fields(n)%vcellmeas)) + status = pio_put_att(File,varid,'missing_value',spval) + status = pio_put_att(File,varid,'_FillValue',spval) + + ! Add cell_methods attribute to variables if averaged + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_put_att(File,varid,'cell_methods','time: mean') + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg) then + status = pio_put_att(File,varid,'time_rep','instantaneous') + else + status = pio_put_att(File,varid,'time_rep','averaged') + endif + endif + enddo ! num_avail_hist_fields_3Dz + + !----------------------------------------------------------------- + ! 3D (biology ice layers) + !----------------------------------------------------------------- + + dimidz(1) = imtid + dimidz(2) = jmtid + dimidz(3) = kmtidb + dimidz(4) = timid + + do n = n3Dzcum + 1, n3Dbcum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & + pio_real, dimidz, varid) + status = pio_put_att(File,varid,'units', & + trim(avail_hist_fields(n)%vunit)) + status = pio_put_att(File,varid, 'long_name', & + trim(avail_hist_fields(n)%vdesc)) + status = pio_put_att(File,varid,'coordinates', & + trim(avail_hist_fields(n)%vcoord)) + status = pio_put_att(File,varid,'cell_measures', & + trim(avail_hist_fields(n)%vcellmeas)) + status = pio_put_att(File,varid,'missing_value',spval) + status = pio_put_att(File,varid,'_FillValue',spval) + + ! Add cell_methods attribute to variables if averaged + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_put_att(File,varid,'cell_methods','time: mean') + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg) then + status = pio_put_att(File,varid,'time_rep','instantaneous') + else + status = pio_put_att(File,varid,'time_rep','averaged') + endif + endif + enddo ! num_avail_hist_fields_3Db + + !----------------------------------------------------------------- + ! 3D (biology snow layers) + !----------------------------------------------------------------- + + dimidz(1) = imtid + dimidz(2) = jmtid + dimidz(3) = kmtida + dimidz(4) = timid + + do n = n3Dbcum + 1, n3Dacum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & + pio_real, dimidz, varid) + status = pio_put_att(File,varid,'units', & + trim(avail_hist_fields(n)%vunit)) + status = pio_put_att(File,varid, 'long_name', & + trim(avail_hist_fields(n)%vdesc)) + status = pio_put_att(File,varid,'coordinates', & + trim(avail_hist_fields(n)%vcoord)) + status = pio_put_att(File,varid,'cell_measures', & + trim(avail_hist_fields(n)%vcellmeas)) + status = pio_put_att(File,varid,'missing_value',spval) + status = pio_put_att(File,varid,'_FillValue',spval) + + ! Add cell_methods attribute to variables if averaged + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_put_att(File,varid,'cell_methods','time: mean') + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg) then + status = pio_put_att(File,varid,'time_rep','instantaneous') + else + status = pio_put_att(File,varid,'time_rep','averaged') + endif + endif + enddo ! num_avail_hist_fields_3Da + + !----------------------------------------------------------------- + ! define attributes for 4D variables + ! time coordinate is dropped + !----------------------------------------------------------------- + + !----------------------------------------------------------------- + ! 4D (ice categories) + !----------------------------------------------------------------- + + dimidcz(1) = imtid + dimidcz(2) = jmtid + dimidcz(3) = kmtidi + dimidcz(4) = cmtid + dimidcz(5) = timid + + do n = n3Dacum + 1, n4Dicum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & + pio_real, dimidcz, varid) + status = pio_put_att(File,varid,'units', & + trim(avail_hist_fields(n)%vunit)) + status = pio_put_att(File,varid, 'long_name', & + trim(avail_hist_fields(n)%vdesc)) + status = pio_put_att(File,varid,'coordinates', & + trim(avail_hist_fields(n)%vcoord)) + status = pio_put_att(File,varid,'cell_measures', & + trim(avail_hist_fields(n)%vcellmeas)) + status = pio_put_att(File,varid,'missing_value',spval) + status = pio_put_att(File,varid,'_FillValue',spval) + + ! Add cell_methods attribute to variables if averaged + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_put_att(File,varid,'cell_methods','time: mean') + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg) then + status = pio_put_att(File,varid,'time_rep','instantaneous') + else + status = pio_put_att(File,varid,'time_rep','averaged') + endif + endif + enddo ! num_avail_hist_fields_4Di + + !----------------------------------------------------------------- + ! 4D (snow layers) + !----------------------------------------------------------------- + + dimidcz(1) = imtid + dimidcz(2) = jmtid + dimidcz(3) = kmtids + dimidcz(4) = cmtid + dimidcz(5) = timid + + do n = n4Dicum + 1, n4Dscum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_def_var(File, trim(avail_hist_fields(n)%vname), & + pio_real, dimidcz, varid) + status = pio_put_att(File,varid,'units', & + trim(avail_hist_fields(n)%vunit)) + status = pio_put_att(File,varid, 'long_name', & + trim(avail_hist_fields(n)%vdesc)) + status = pio_put_att(File,varid,'coordinates', & + trim(avail_hist_fields(n)%vcoord)) + status = pio_put_att(File,varid,'cell_measures', & + trim(avail_hist_fields(n)%vcellmeas)) + status = pio_put_att(File,varid,'missing_value',spval) + status = pio_put_att(File,varid,'_FillValue',spval) + + ! Add cell_methods attribute to variables if averaged + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_put_att(File,varid,'cell_methods','time: mean') + endif + + if (histfreq(ns) == '1' .or. .not. hist_avg) then + status = pio_put_att(File,varid,'time_rep','instantaneous') + else + status = pio_put_att(File,varid,'time_rep','averaged') + endif + endif + enddo ! num_avail_hist_fields_4Ds + + !----------------------------------------------------------------- + ! global attributes + !----------------------------------------------------------------- + ! ... the user should change these to something useful ... + !----------------------------------------------------------------- +#ifdef CESMCOUPLED + status = pio_put_att(File,pio_global,'title',runid) +#else + title = 'sea ice model output for CICE' + status = pio_put_att(File,pio_global,'title',trim(title)) +#endif + title = 'Diagnostic and Prognostic Variables' + status = pio_put_att(File,pio_global,'contents',trim(title)) + + write(title,'(2a)') 'Los Alamos Sea Ice Model, ', trim(version_name) + status = pio_put_att(File,pio_global,'source',trim(title)) + + if (use_leap_years) then + write(title,'(a,i3,a)') 'This year has ',int(dayyr),' days' + else + write(title,'(a,i3,a)') 'All years have exactly ',int(dayyr),' days' + endif + status = pio_put_att(File,pio_global,'comment',trim(title)) + + write(title,'(a,i8.8)') 'File written on model date ',idate + status = pio_put_att(File,pio_global,'comment2',trim(title)) + + write(title,'(a,i6)') 'seconds elapsed into model date: ',sec + status = pio_put_att(File,pio_global,'comment3',trim(title)) + + title = 'CF-1.0' + status = & + pio_put_att(File,pio_global,'conventions',trim(title)) + + call date_and_time(date=current_date, time=current_time) + write(start_time,1000) current_date(1:4), current_date(5:6), & + current_date(7:8), current_time(1:2), & + current_time(3:4) +1000 format('This dataset was created on ', & + a,'-',a,'-',a,' at ',a,':',a) + status = pio_put_att(File,pio_global,'history',trim(start_time)) + + status = pio_put_att(File,pio_global,'io_flavor','io_pio') + + !----------------------------------------------------------------- + ! end define mode + !----------------------------------------------------------------- + + status = pio_enddef(File) + + !----------------------------------------------------------------- + ! write time variable + !----------------------------------------------------------------- + + status = pio_inq_varid(File,'time',varid) +!sgl status = pio_put_var(File,varid,(/1/),ltime) + status = pio_put_var(File,varid,(/1/),ltime2) + + !----------------------------------------------------------------- + ! write time_bounds info + !----------------------------------------------------------------- + + if (hist_avg .and. histfreq(ns) /= '1') then + status = pio_inq_varid(File,'time_bounds',varid) + time_bounds=(/time_beg(ns),time_end(ns)/) + bnd_start = (/1,1/) + bnd_length = (/2,1/) + status = pio_put_var(File,varid,ival=time_bounds, & + start=bnd_start(:),count=bnd_length(:)) + endif + + !----------------------------------------------------------------- + ! write coordinate variables + !----------------------------------------------------------------- + + allocate(workr2(nx_block,ny_block,nblocks)) + + do i = 1,ncoord + status = pio_inq_varid(File, coord_var(i)%short_name, varid) + SELECT CASE (coord_var(i)%short_name) + CASE ('TLON') + ! Convert T grid longitude from -180 -> 180 to 0 to 360 + workr2(:,:,:) = mod(tlon(:,:,1:nblocks)*rad_to_deg + c360, c360) + CASE ('TLAT') + workr2(:,:,:) = tlat(:,:,1:nblocks)*rad_to_deg + CASE ('ULON') + workr2(:,:,:) = ulon(:,:,1:nblocks)*rad_to_deg + CASE ('ULAT') + workr2(:,:,:) = ulat(:,:,1:nblocks)*rad_to_deg + END SELECT + call pio_write_darray(File, varid, iodesc2d, & + workr2, status, fillval=spval_dbl) + enddo + + ! Extra dimensions (NCAT, VGRD*) + + do i = 1, nvarz + if (igrdz(i)) then + status = pio_inq_varid(File, var_nz(i)%short_name, varid) + SELECT CASE (var_nz(i)%short_name) + CASE ('NCAT') + status = pio_put_var(File, varid, hin_max(1:ncat_hist)) + CASE ('VGRDi') + status = pio_put_var(File, varid, (/(k, k=1,nzilyr)/)) + CASE ('VGRDs') + status = pio_put_var(File, varid, (/(k, k=1,nzslyr)/)) + CASE ('VGRDb') + status = pio_put_var(File, varid, (/(k, k=1,nzblyr)/)) + CASE ('VGRDa') + status = pio_put_var(File, varid, (/(k, k=1,nzalyr)/)) + END SELECT + endif + enddo + + !----------------------------------------------------------------- + ! write grid masks, area and rotation angle + !----------------------------------------------------------------- + +! if (igrd(n_tmask)) then +! status = pio_inq_varid(File, 'tmask', varid) +! call pio_write_darray(File, varid, iodesc2d, & +! hm(:,:,1:nblocks), status, fillval=spval_dbl) +! endif +! if (igrd(n_blkmask)) then +! status = pio_inq_varid(File, 'blkmask', varid) +! call pio_write_darray(File, varid, iodesc2d, & +! bm(:,:,1:nblocks), status, fillval=spval_dbl) +! endif + + do i = 1, nvar ! note: n_tmask=1, n_blkmask=2 + if (igrd(i)) then + SELECT CASE (var(i)%req%short_name) + CASE ('tmask') + workr2 = hm(:,:,1:nblocks) + CASE ('blkmask') + workr2 = bm(:,:,1:nblocks) + CASE ('tarea') + workr2 = tarea(:,:,1:nblocks) + CASE ('uarea') + workr2 = uarea(:,:,1:nblocks) + CASE ('dxu') + workr2 = dxu(:,:,1:nblocks) + CASE ('dyu') + workr2 = dyu(:,:,1:nblocks) + CASE ('dxt') + workr2 = dxt(:,:,1:nblocks) + CASE ('dyt') + workr2 = dyt(:,:,1:nblocks) + CASE ('HTN') + workr2 = HTN(:,:,1:nblocks) + CASE ('HTE') + workr2 = HTE(:,:,1:nblocks) + CASE ('ANGLE') + workr2 = ANGLE(:,:,1:nblocks) + CASE ('ANGLET') + workr2 = ANGLET(:,:,1:nblocks) + END SELECT + status = pio_inq_varid(File, var(i)%req%short_name, varid) + call pio_write_darray(File, varid, iodesc2d, & + workr2, status, fillval=spval_dbl) + endif + enddo + + !---------------------------------------------------------------- + ! Write coordinates of grid box vertices + !---------------------------------------------------------------- + + if (f_bounds) then + allocate(workr3v(nverts,nx_block,ny_block,nblocks)) + workr3v (:,:,:,:) = c0 + do i = 1, nvar_verts + SELECT CASE (var_nverts(i)%short_name) + CASE ('lont_bounds') + do ivertex = 1, nverts + workr3v(ivertex,:,:,:) = lont_bounds(ivertex,:,:,1:nblocks) + enddo + CASE ('latt_bounds') + do ivertex = 1, nverts + workr3v(ivertex,:,:,:) = latt_bounds(ivertex,:,:,1:nblocks) + enddo + CASE ('lonu_bounds') + do ivertex = 1, nverts + workr3v(ivertex,:,:,:) = lonu_bounds(ivertex,:,:,1:nblocks) + enddo + CASE ('latu_bounds') + do ivertex = 1, nverts + workr3v(ivertex,:,:,:) = latu_bounds(ivertex,:,:,1:nblocks) + enddo + END SELECT + + status = pio_inq_varid(File, var_nverts(i)%short_name, varid) + call pio_write_darray(File, varid, iodesc3dv, & + workr3v, status, fillval=spval_dbl) + enddo + deallocate(workr3v) + endif ! f_bounds + + + !----------------------------------------------------------------- + ! write variable data + !----------------------------------------------------------------- + + ! 2D + do n=1,num_avail_hist_fields_2D + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_inq_varid(File,avail_hist_fields(n)%vname,varid) + if (status /= pio_noerr) call abort_ice(subname// & + 'ERROR getting varid for '//avail_hist_fields(n)%vname) + workr2(:,:,:) = a2D(:,:,n,1:nblocks) + call pio_setframe(File, varid, int(1,kind=PIO_OFFSET_KIND)) + call pio_write_darray(File, varid, iodesc2d,& + workr2, status, fillval=spval_dbl) + endif + enddo ! num_avail_hist_fields_2D + + deallocate(workr2) + + ! 3D (category) + allocate(workr3(nx_block,ny_block,nblocks,ncat_hist)) + do n = n2D + 1, n3Dccum + nn = n - n2D + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_inq_varid(File,avail_hist_fields(n)%vname,varid) + if (status /= pio_noerr) call abort_ice(subname// & + 'ERROR: getting varid for '//avail_hist_fields(n)%vname) + do j = 1, nblocks + do i = 1, ncat_hist + workr3(:,:,j,i) = a3Dc(:,:,i,nn,j) + enddo + enddo + call pio_setframe(File, varid, int(1,kind=PIO_OFFSET_KIND)) + call pio_write_darray(File, varid, iodesc3dc,& + workr3, status, fillval=spval_dbl) + endif + enddo ! num_avail_hist_fields_3Dc + deallocate(workr3) + + ! 3D (vertical ice) + allocate(workr3(nx_block,ny_block,nblocks,nzilyr)) + do n = n3Dccum+1, n3Dzcum + nn = n - n3Dccum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_inq_varid(File,avail_hist_fields(n)%vname,varid) + if (status /= pio_noerr) call abort_ice(subname// & + 'ERROR: getting varid for '//avail_hist_fields(n)%vname) + do j = 1, nblocks + do i = 1, nzilyr + workr3(:,:,j,i) = a3Dz(:,:,i,nn,j) + enddo + enddo + call pio_setframe(File, varid, int(1,kind=PIO_OFFSET_KIND)) + call pio_write_darray(File, varid, iodesc3di,& + workr3, status, fillval=spval_dbl) + endif + enddo ! num_avail_hist_fields_3Dz + deallocate(workr3) + + ! 3D (vertical ice biology) + allocate(workr3(nx_block,ny_block,nblocks,nzblyr)) + do n = n3Dzcum+1, n3Dbcum + nn = n - n3Dzcum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_inq_varid(File,avail_hist_fields(n)%vname,varid) + if (status /= pio_noerr) call abort_ice(subname// & + 'ERROR: getting varid for '//avail_hist_fields(n)%vname) + do j = 1, nblocks + do i = 1, nzblyr + workr3(:,:,j,i) = a3Db(:,:,i,nn,j) + enddo + enddo + call pio_setframe(File, varid, int(1,kind=PIO_OFFSET_KIND)) + call pio_write_darray(File, varid, iodesc3db,& + workr3, status, fillval=spval_dbl) + endif + enddo ! num_avail_hist_fields_3Db + deallocate(workr3) + + ! 3D (vertical snow biology) + allocate(workr3(nx_block,ny_block,nblocks,nzalyr)) + do n = n3Dbcum+1, n3Dacum + nn = n - n3Dbcum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_inq_varid(File,avail_hist_fields(n)%vname,varid) + if (status /= pio_noerr) call abort_ice(subname// & + 'ERROR: getting varid for '//avail_hist_fields(n)%vname) + do j = 1, nblocks + do i = 1, nzalyr + workr3(:,:,j,i) = a3Da(:,:,i,nn,j) + enddo + enddo + call pio_setframe(File, varid, int(1,kind=PIO_OFFSET_KIND)) + call pio_write_darray(File, varid, iodesc3da,& + workr3, status, fillval=spval_dbl) + endif + enddo ! num_avail_hist_fields_3Db + deallocate(workr3) + + allocate(workr4(nx_block,ny_block,nblocks,ncat_hist,nzilyr)) + ! 4D (categories, vertical ice) + do n = n3Dacum+1, n4Dicum + nn = n - n3Dacum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_inq_varid(File,avail_hist_fields(n)%vname,varid) + if (status /= pio_noerr) call abort_ice(subname// & + 'ERROR: getting varid for '//avail_hist_fields(n)%vname) + do j = 1, nblocks + do i = 1, ncat_hist + do k = 1, nzilyr + workr4(:,:,j,i,k) = a4Di(:,:,k,i,nn,j) + enddo ! k + enddo ! i + enddo ! j + call pio_setframe(File, varid, int(1,kind=PIO_OFFSET_KIND)) + call pio_write_darray(File, varid, iodesc4di,& + workr4, status, fillval=spval_dbl) + endif + enddo ! num_avail_hist_fields_4Di + deallocate(workr4) + + allocate(workr4(nx_block,ny_block,nblocks,ncat_hist,nzslyr)) + ! 4D (categories, vertical ice) + do n = n4Dicum+1, n4Dscum + nn = n - n4Dicum + if (avail_hist_fields(n)%vhistfreq == histfreq(ns) .or. write_ic) then + status = pio_inq_varid(File,avail_hist_fields(n)%vname,varid) + if (status /= pio_noerr) call abort_ice(subname// & + 'ERROR: getting varid for '//avail_hist_fields(n)%vname) + do j = 1, nblocks + do i = 1, ncat_hist + do k = 1, nzslyr + workr4(:,:,j,i,k) = a4Ds(:,:,k,i,nn,j) + enddo ! k + enddo ! i + enddo ! j + call pio_setframe(File, varid, int(1,kind=PIO_OFFSET_KIND)) + call pio_write_darray(File, varid, iodesc4ds,& + workr4, status, fillval=spval_dbl) + endif + enddo ! num_avail_hist_fields_4Di + deallocate(workr4) + +! similarly for num_avail_hist_fields_4Db (define workr4b, iodesc4db) + + + !----------------------------------------------------------------- + ! clean-up PIO descriptors + !----------------------------------------------------------------- + + call pio_freedecomp(File,iodesc2d) + call pio_freedecomp(File,iodesc3dv) + call pio_freedecomp(File,iodesc3dc) + call pio_freedecomp(File,iodesc3di) + call pio_freedecomp(File,iodesc3db) + call pio_freedecomp(File,iodesc3da) + call pio_freedecomp(File,iodesc4di) + call pio_freedecomp(File,iodesc4ds) + + !----------------------------------------------------------------- + ! close output dataset + !----------------------------------------------------------------- + + call pio_closefile(File) + if (my_task == master_task) then + write(nu_diag,*) ' ' + write(nu_diag,*) 'Finished writing ',trim(ncfile(ns)) + endif + +#endif + + end subroutine ice_write_hist + +!======================================================================= + + end module ice_history_write + +!======================================================================= diff --git a/cicecore/cicedynB/infrastructure/io/io_pio2/ice_pio.F90 b/cicecore/cicedynB/infrastructure/io/io_pio2/ice_pio.F90 new file mode 100644 index 000000000..5fff64944 --- /dev/null +++ b/cicecore/cicedynB/infrastructure/io/io_pio2/ice_pio.F90 @@ -0,0 +1,365 @@ +!============================================================================ +! Writes netcdf files +! Created by Mariana Vertenstein, June 2009 + + module ice_pio + + use shr_kind_mod, only: r8 => shr_kind_r8, in=>shr_kind_in + use shr_kind_mod, only: cl => shr_kind_cl + use shr_sys_mod , only: shr_sys_flush + use ice_kinds_mod + use ice_blocks + use ice_broadcast + use ice_communicate + use ice_domain, only : nblocks, blocks_ice + use ice_domain_size + use ice_fileunits + use ice_exit, only: abort_ice + use icepack_intfc, only: icepack_warnings_flush, icepack_warnings_aborted + use pio + + implicit none + private + + interface ice_pio_initdecomp + module procedure ice_pio_initdecomp_2d + module procedure ice_pio_initdecomp_3d + module procedure ice_pio_initdecomp_4d + module procedure ice_pio_initdecomp_3d_inner + end interface + + public ice_pio_init + public ice_pio_initdecomp + + type(iosystem_desc_t), pointer, public :: ice_pio_subsystem + +!=============================================================================== + + contains + +!=============================================================================== + +! Initialize the io subsystem +! 2009-Feb-17 - J. Edwards - initial version + + subroutine ice_pio_init(mode, filename, File, clobber, cdf64) + + use shr_pio_mod, only: shr_pio_getiosys, shr_pio_getiotype + + implicit none + character(len=*) , intent(in), optional :: mode + character(len=*) , intent(in), optional :: filename + type(file_desc_t) , intent(inout), optional :: File + logical , intent(in), optional :: clobber + logical , intent(in), optional :: cdf64 + + ! local variables + + integer (int_kind) :: & + nml_error ! namelist read error flag + + integer :: pio_iotype + logical :: exists + logical :: lclobber + logical :: lcdf64 + integer :: status + integer :: nmode + character(len=*), parameter :: subname = '(ice_pio_init)' + logical, save :: first_call = .true. + + ice_pio_subsystem => shr_pio_getiosys(inst_name) + pio_iotype = shr_pio_getiotype(inst_name) + + if (present(mode) .and. present(filename) .and. present(File)) then + + if (trim(mode) == 'write') then + lclobber = .false. + if (present(clobber)) lclobber=clobber + + lcdf64 = .false. + if (present(cdf64)) lcdf64=cdf64 + + if (File%fh<0) then + ! filename not open + inquire(file=trim(filename),exist=exists) + if (exists) then + if (lclobber) then + nmode = pio_clobber + if (lcdf64) nmode = ior(nmode,PIO_64BIT_OFFSET) + status = pio_createfile(ice_pio_subsystem, File, pio_iotype, trim(filename), nmode) + if (my_task == master_task) then + write(nu_diag,*) subname,' create file ',trim(filename) + end if + else + nmode = pio_write + status = pio_openfile(ice_pio_subsystem, File, pio_iotype, trim(filename), nmode) + if (my_task == master_task) then + write(nu_diag,*) subname,' open file ',trim(filename) + end if + endif + else + nmode = pio_noclobber + if (lcdf64) nmode = ior(nmode,PIO_64BIT_OFFSET) + status = pio_createfile(ice_pio_subsystem, File, pio_iotype, trim(filename), nmode) + if (my_task == master_task) then + write(nu_diag,*) subname,' create file ',trim(filename) + end if + endif + else + ! filename is already open, just return + endif + end if + + if (trim(mode) == 'read') then + inquire(file=trim(filename),exist=exists) + if (exists) then + status = pio_openfile(ice_pio_subsystem, File, pio_iotype, trim(filename), pio_nowrite) + else + if(my_task==master_task) then + write(nu_diag,*) 'ice_pio_ropen ERROR: file invalid ',trim(filename) + end if + call abort_ice(subname//'ERROR: aborting with invalid file') + endif + end if + + end if + + end subroutine ice_pio_init + +!================================================================================ + + subroutine ice_pio_initdecomp_2d(iodesc) + + type(io_desc_t), intent(out) :: iodesc + + integer (kind=int_kind) :: & + iblk,ilo,ihi,jlo,jhi,lon,lat,i,j,n,k + + type(block) :: this_block + + integer(kind=int_kind), pointer :: dof2d(:) + character(len=*), parameter :: subname = '(ice_pio_initdecomp_2d)' + + allocate(dof2d(nx_block*ny_block*nblocks)) + + 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=1,ny_block + do i=1,nx_block + n = n+1 + if (j < jlo .or. j>jhi) then + dof2d(n) = 0 + else if (i < ilo .or. i > ihi) then + dof2d(n) = 0 + else + lon = this_block%i_glob(i) + lat = this_block%j_glob(j) + dof2d(n) = (lat-1)*nx_global + lon + endif + enddo !i + enddo !j + end do + + call pio_initdecomp(ice_pio_subsystem, pio_double, (/nx_global,ny_global/), & + dof2d, iodesc) + + deallocate(dof2d) + + end subroutine ice_pio_initdecomp_2d + +!================================================================================ + + subroutine ice_pio_initdecomp_3d (ndim3, iodesc, remap) + + integer(kind=int_kind), intent(in) :: ndim3 + type(io_desc_t), intent(out) :: iodesc + logical, optional :: remap + integer (kind=int_kind) :: & + iblk,ilo,ihi,jlo,jhi,lon,lat,i,j,n,k + + type(block) :: this_block + logical :: lremap + integer(kind=int_kind), pointer :: dof3d(:) + character(len=*), parameter :: subname = '(ice_pio_initdecomp_2d)' + + allocate(dof3d(nx_block*ny_block*nblocks*ndim3)) + lremap=.false. + if (present(remap)) lremap=remap + if (lremap) then + ! Reorder the ndim3 and nblocks loops to avoid a temporary array in restart read/write + 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 k=1,ndim3 + do j=1,ny_block + do i=1,nx_block + n = n+1 + if (j < jlo .or. j>jhi) then + dof3d(n)=0 + else if (i < ilo .or. i > ihi) then + dof3d(n) = 0 + else + lon = this_block%i_glob(i) + lat = this_block%j_glob(j) + dof3d(n) = ((lat-1)*nx_global + lon) + (k-1)*nx_global*ny_global + endif + enddo !i + enddo !j + enddo !ndim3 + enddo ! iblk + else + n=0 + do k=1,ndim3 + 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=1,ny_block + do i=1,nx_block + n = n+1 + if (j < jlo .or. j>jhi) then + dof3d(n)=0 + else if (i < ilo .or. i > ihi) then + dof3d(n) = 0 + else + lon = this_block%i_glob(i) + lat = this_block%j_glob(j) + dof3d(n) = ((lat-1)*nx_global + lon) + (k-1)*nx_global*ny_global + endif + enddo !i + enddo !j + enddo ! iblk + enddo !ndim3 + endif + + call pio_initdecomp(ice_pio_subsystem, pio_double, (/nx_global,ny_global,ndim3/), & + dof3d, iodesc) + + deallocate(dof3d) + + end subroutine ice_pio_initdecomp_3d + +!================================================================================ + + subroutine ice_pio_initdecomp_3d_inner(ndim3, inner_dim, iodesc) + + integer(kind=int_kind), intent(in) :: ndim3 + logical, intent(in) :: inner_dim + type(io_desc_t), intent(out) :: iodesc + + integer (kind=int_kind) :: & + iblk,ilo,ihi,jlo,jhi,lon,lat,i,j,n,k + + type(block) :: this_block + + integer(kind=int_kind), pointer :: dof3d(:) + + character(len=*), parameter :: subname = '(ice_pio_initdecomp_3d_inner)' + + allocate(dof3d(nx_block*ny_block*nblocks*ndim3)) + + 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=1,ny_block + do i=1,nx_block + do k=1,ndim3 + n = n+1 + if (j < jlo .or. j>jhi) then + dof3d(n) = 0 + else if (i < ilo .or. i > ihi) then + dof3d(n) = 0 + else + lon = this_block%i_glob(i) + lat = this_block%j_glob(j) + dof3d(n) = k + ((lon-1) + (lat-1)*nx_global)*ndim3 + endif + end do !ndim3 + enddo !i + enddo !j + end do !iblk + + call pio_initdecomp(ice_pio_subsystem, pio_double, (/ndim3,nx_global,ny_global/), & + dof3d, iodesc) + + deallocate(dof3d) + + end subroutine ice_pio_initdecomp_3d_inner + +!================================================================================ + + subroutine ice_pio_initdecomp_4d (ndim3, ndim4, iodesc) + + integer(kind=int_kind), intent(in) :: ndim3, ndim4 + type(io_desc_t), intent(out) :: iodesc + + integer (kind=int_kind) :: & + iblk,ilo,ihi,jlo,jhi,lon,lat,i,j,n,k,l + + type(block) :: this_block + + integer(kind=int_kind), pointer :: dof4d(:) + + character(len=*), parameter :: subname = '(ice_pio_initdecomp_4d)' + + allocate(dof4d(nx_block*ny_block*nblocks*ndim3*ndim4)) + + n=0 + do l=1,ndim4 + do k=1,ndim3 + 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=1,ny_block + do i=1,nx_block + n = n+1 + if (j < jlo .or. j>jhi) then + dof4d(n)=0 + else if (i < ilo .or. i > ihi) then + dof4d(n) = 0 + else + lon = this_block%i_glob(i) + lat = this_block%j_glob(j) + dof4d(n) = ((lat-1)*nx_global + lon) & + + (k-1)*nx_global*ny_global & + + (l-1)*nx_global*ny_global*ndim3 + endif + enddo !i + enddo !j + enddo ! iblk + enddo !ndim3 + enddo !ndim4 + + call pio_initdecomp(ice_pio_subsystem, pio_double, & + (/nx_global,ny_global,ndim3,ndim4/), dof4d, iodesc) + + deallocate(dof4d) + + end subroutine ice_pio_initdecomp_4d + +!================================================================================ + + end module ice_pio + +!================================================================================ diff --git a/cicecore/cicedynB/infrastructure/io/io_pio2/ice_restart.F90 b/cicecore/cicedynB/infrastructure/io/io_pio2/ice_restart.F90 new file mode 100644 index 000000000..f00501bfc --- /dev/null +++ b/cicecore/cicedynB/infrastructure/io/io_pio2/ice_restart.F90 @@ -0,0 +1,876 @@ +!======================================================================= +! +! Read and write ice model restart files using pio interfaces. +! authors David A Bailey, NCAR + + module ice_restart + + use ice_broadcast + use ice_exit, only: abort_ice + use ice_fileunits, only: nu_diag, nu_restart, nu_rst_pointer + use ice_kinds_mod + use ice_restart_shared, only: & + restart, restart_ext, restart_dir, restart_file, pointer_file, & + runid, runtype, use_restart_time, restart_format, lcdf64, lenstr + use ice_pio + use pio + use icepack_intfc, only: icepack_warnings_flush, icepack_warnings_aborted + use icepack_intfc, only: icepack_query_tracer_flags, icepack_query_tracer_indices + use icepack_intfc, only: icepack_query_parameters, icepack_query_tracer_numbers + + implicit none + private + public :: init_restart_write, init_restart_read, & + read_restart_field, write_restart_field, final_restart + + type(file_desc_t) :: File + type(var_desc_t) :: vardesc + + type(io_desc_t) :: iodesc2d + type(io_desc_t) :: iodesc3d_ncat + +!======================================================================= + + contains + +!======================================================================= + +! Sets up restart file for reading. +! author David A Bailey, NCAR + + subroutine init_restart_read(ice_ic) + + use ice_calendar, only: istep0, istep1, time, time_forc, nyr, month, & + mday, sec, npt + use ice_communicate, only: my_task, master_task + use ice_domain_size, only: ncat + use ice_read_write, only: ice_open + + character(len=char_len_long), intent(in), optional :: ice_ic + + ! local variables + + character(len=char_len_long) :: & + filename, filename0 + + integer (kind=int_kind) :: status + + character(len=*), parameter :: subname = '(init_restart_read)' + + if (present(ice_ic)) then + filename = trim(ice_ic) + else + if (my_task == master_task) then + open(nu_rst_pointer,file=pointer_file) + read(nu_rst_pointer,'(a)') filename0 + filename = trim(filename0) + close(nu_rst_pointer) + write(nu_diag,*) 'Read ',pointer_file(1:lenstr(pointer_file)) + endif + call broadcast_scalar(filename, master_task) + endif + + if (my_task == master_task) then + write(nu_diag,*) 'Using restart dump=', trim(filename) + end if + + if (restart_format == 'pio') then + File%fh=-1 + call ice_pio_init(mode='read', filename=trim(filename), File=File) + + call ice_pio_initdecomp(iodesc=iodesc2d) + call ice_pio_initdecomp(ndim3=ncat , iodesc=iodesc3d_ncat,remap=.true.) + + if (use_restart_time) then + status = pio_get_att(File, pio_global, 'istep1', istep0) + status = pio_get_att(File, pio_global, 'time', time) + status = pio_get_att(File, pio_global, 'time_forc', time_forc) + call pio_seterrorhandling(File, PIO_BCAST_ERROR) + status = pio_get_att(File, pio_global, 'nyr', nyr) + call pio_seterrorhandling(File, PIO_INTERNAL_ERROR) + if (status == PIO_noerr) then + status = pio_get_att(File, pio_global, 'month', month) + status = pio_get_att(File, pio_global, 'mday', mday) + status = pio_get_att(File, pio_global, 'sec', sec) + endif + endif ! use namelist values if use_restart_time = F + endif + + if (my_task == master_task) then + write(nu_diag,*) 'Restart read at istep=',istep0,time,time_forc + endif + + call broadcast_scalar(istep0,master_task) + call broadcast_scalar(time,master_task) + call broadcast_scalar(time_forc,master_task) + call broadcast_scalar(nyr,master_task) + + istep1 = istep0 + + ! if runid is bering then need to correct npt for istep0 + if (trim(runid) == 'bering') then + npt = npt - istep0 + endif + + end subroutine init_restart_read + +!======================================================================= + +! Sets up restart file for writing. +! author David A Bailey, NCAR + + subroutine init_restart_write(filename_spec) + + use ice_calendar, only: sec, month, mday, nyr, istep1, & + time, time_forc, year_init + use ice_communicate, only: my_task, master_task + use ice_domain_size, only: nx_global, ny_global, ncat, nilyr, nslyr, & + n_aero, nblyr, n_zaero, n_algae, n_doc, & + n_dic, n_don, n_fed, n_fep + use ice_dyn_shared, only: kdyn + use ice_arrays_column, only: oceanmixed_ice + + logical (kind=log_kind) :: & + solve_zsal, skl_bgc, z_tracers + + logical (kind=log_kind) :: & + tr_iage, tr_FY, tr_lvl, tr_aero, tr_pond_cesm, & + tr_pond_topo, tr_pond_lvl, tr_brine, & + tr_bgc_N, tr_bgc_C, tr_bgc_Nit, & + tr_bgc_Sil, tr_bgc_DMS, & + tr_bgc_chl, tr_bgc_Am, & + tr_bgc_PON, tr_bgc_DON, & + tr_zaero, tr_bgc_Fe, & + tr_bgc_hum + + integer (kind=int_kind) :: & + nbtrcr + + character(len=char_len_long), intent(in), optional :: filename_spec + + ! local variables + + integer (kind=int_kind) :: & + iyear, imonth, iday ! year, month, day + + character(len=char_len_long) :: filename + + integer (kind=int_kind) :: dimid_ni, dimid_nj, dimid_ncat, & + dimid_nilyr, dimid_nslyr, dimid_naero + + integer (kind=int_kind), allocatable :: dims(:) + + integer (kind=int_kind) :: & + k, n, & ! loop index + status ! status variable from netCDF routine + + character (len=3) :: nchar, ncharb + + character(len=*), parameter :: subname = '(init_restart_write)' + + call icepack_query_tracer_numbers(nbtrcr_out=nbtrcr) + call icepack_query_tracer_flags( & + tr_iage_out=tr_iage, tr_FY_out=tr_FY, tr_lvl_out=tr_lvl, & + tr_aero_out=tr_aero, tr_pond_cesm_out=tr_pond_cesm, & + tr_pond_topo_out=tr_pond_topo, tr_pond_lvl_out=tr_pond_lvl, tr_brine_out=tr_brine, & + tr_bgc_N_out=tr_bgc_N, tr_bgc_C_out=tr_bgc_C, tr_bgc_Nit_out=tr_bgc_Nit, & + tr_bgc_Sil_out=tr_bgc_Sil, tr_bgc_DMS_out=tr_bgc_DMS, & + tr_bgc_chl_out=tr_bgc_chl, tr_bgc_Am_out=tr_bgc_Am, & + tr_bgc_PON_out=tr_bgc_PON, tr_bgc_DON_out=tr_bgc_DON, & + tr_zaero_out=tr_zaero, tr_bgc_Fe_out=tr_bgc_Fe, & + tr_bgc_hum_out=tr_bgc_hum) + call icepack_query_parameters(solve_zsal_out=solve_zsal, skl_bgc_out=skl_bgc, & + z_tracers_out=z_tracers) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + + ! construct path/file + if (present(filename_spec)) then + filename = trim(filename_spec) + else + iyear = nyr + year_init - 1 + imonth = month + iday = mday + + write(filename,'(a,a,a,i4.4,a,i2.2,a,i2.2,a,i5.5)') & + restart_dir(1:lenstr(restart_dir)), & + restart_file(1:lenstr(restart_file)),'.', & + iyear,'-',month,'-',mday,'-',sec + end if + + if (restart_format /= 'bin') filename = trim(filename) // '.nc' + + ! write pointer (path/file) + if (my_task == master_task) then + open(nu_rst_pointer,file=pointer_file) + write(nu_rst_pointer,'(a)') filename + close(nu_rst_pointer) + endif + + if (restart_format == 'pio') then + + File%fh=-1 + call ice_pio_init(mode='write',filename=trim(filename), File=File, & + clobber=.true., cdf64=lcdf64 ) + + status = pio_put_att(File,pio_global,'istep1',istep1) + status = pio_put_att(File,pio_global,'time',time) + status = pio_put_att(File,pio_global,'time_forc',time_forc) + status = pio_put_att(File,pio_global,'nyr',nyr) + status = pio_put_att(File,pio_global,'month',month) + status = pio_put_att(File,pio_global,'mday',mday) + status = pio_put_att(File,pio_global,'sec',sec) + + status = pio_def_dim(File,'ni',nx_global,dimid_ni) + status = pio_def_dim(File,'nj',ny_global,dimid_nj) + status = pio_def_dim(File,'ncat',ncat,dimid_ncat) + + !----------------------------------------------------------------- + ! 2D restart fields + !----------------------------------------------------------------- + + allocate(dims(2)) + + dims(1) = dimid_ni + dims(2) = dimid_nj + + call define_rest_field(File,'uvel',dims) + call define_rest_field(File,'vvel',dims) + +#ifdef CESMCOUPLED + call define_rest_field(File,'coszen',dims) +#endif + call define_rest_field(File,'scale_factor',dims) + call define_rest_field(File,'swvdr',dims) + call define_rest_field(File,'swvdf',dims) + call define_rest_field(File,'swidr',dims) + call define_rest_field(File,'swidf',dims) + + call define_rest_field(File,'strocnxT',dims) + call define_rest_field(File,'strocnyT',dims) + + call define_rest_field(File,'stressp_1',dims) + call define_rest_field(File,'stressp_2',dims) + call define_rest_field(File,'stressp_3',dims) + call define_rest_field(File,'stressp_4',dims) + + call define_rest_field(File,'stressm_1',dims) + call define_rest_field(File,'stressm_2',dims) + call define_rest_field(File,'stressm_3',dims) + call define_rest_field(File,'stressm_4',dims) + + call define_rest_field(File,'stress12_1',dims) + call define_rest_field(File,'stress12_2',dims) + call define_rest_field(File,'stress12_3',dims) + call define_rest_field(File,'stress12_4',dims) + + call define_rest_field(File,'iceumask',dims) + + if (oceanmixed_ice) then + call define_rest_field(File,'sst',dims) + call define_rest_field(File,'frzmlt',dims) + endif + + if (tr_FY) then + call define_rest_field(File,'frz_onset',dims) + end if + + if (kdyn == 2) then + call define_rest_field(File,'a11_1',dims) + call define_rest_field(File,'a11_2',dims) + call define_rest_field(File,'a11_3',dims) + call define_rest_field(File,'a11_4',dims) + call define_rest_field(File,'a12_1',dims) + call define_rest_field(File,'a12_2',dims) + call define_rest_field(File,'a12_3',dims) + call define_rest_field(File,'a12_4',dims) + endif + + if (tr_pond_lvl) then + call define_rest_field(File,'fsnow',dims) + endif + + if (nbtrcr > 0) then + if (tr_bgc_N) then + do k=1,n_algae + write(nchar,'(i3.3)') k + call define_rest_field(File,'algalN'//trim(nchar),dims) + enddo + endif + if (tr_bgc_C) then + do k=1,n_doc + write(nchar,'(i3.3)') k + call define_rest_field(File,'doc'//trim(nchar),dims) + enddo + do k=1,n_dic + write(nchar,'(i3.3)') k + call define_rest_field(File,'dic'//trim(nchar),dims) + enddo + endif + call define_rest_field(File,'nit' ,dims) + if (tr_bgc_Am) & + call define_rest_field(File,'amm' ,dims) + if (tr_bgc_Sil) & + call define_rest_field(File,'sil' ,dims) + if (tr_bgc_hum) & + call define_rest_field(File,'hum' ,dims) + if (tr_bgc_DMS) then + call define_rest_field(File,'dmsp' ,dims) + call define_rest_field(File,'dms' ,dims) + endif + if (tr_bgc_DON) then + do k=1,n_don + write(nchar,'(i3.3)') k + call define_rest_field(File,'don'//trim(nchar),dims) + enddo + endif + if (tr_bgc_Fe ) then + do k=1,n_fed + write(nchar,'(i3.3)') k + call define_rest_field(File,'fed'//trim(nchar),dims) + enddo + do k=1,n_fep + write(nchar,'(i3.3)') k + call define_rest_field(File,'fep'//trim(nchar),dims) + enddo + endif + if (tr_zaero) then + do k=1,n_zaero + write(nchar,'(i3.3)') k + call define_rest_field(File,'zaeros'//trim(nchar),dims) + enddo + endif + endif !nbtrcr + + if (solve_zsal) call define_rest_field(File,'sss',dims) + + deallocate(dims) + + !----------------------------------------------------------------- + ! 3D restart fields (ncat) + !----------------------------------------------------------------- + + allocate(dims(3)) + + dims(1) = dimid_ni + dims(2) = dimid_nj + dims(3) = dimid_ncat + + call define_rest_field(File,'aicen',dims) + call define_rest_field(File,'vicen',dims) + call define_rest_field(File,'vsnon',dims) + call define_rest_field(File,'Tsfcn',dims) + + if (tr_iage) then + call define_rest_field(File,'iage',dims) + end if + + if (tr_FY) then + call define_rest_field(File,'FY',dims) + end if + + if (tr_lvl) then + call define_rest_field(File,'alvl',dims) + call define_rest_field(File,'vlvl',dims) + end if + + if (tr_pond_cesm) then + call define_rest_field(File,'apnd',dims) + call define_rest_field(File,'hpnd',dims) + end if + + if (tr_pond_topo) then + call define_rest_field(File,'apnd',dims) + call define_rest_field(File,'hpnd',dims) + call define_rest_field(File,'ipnd',dims) + end if + + if (tr_pond_lvl) then + call define_rest_field(File,'apnd',dims) + call define_rest_field(File,'hpnd',dims) + call define_rest_field(File,'ipnd',dims) + call define_rest_field(File,'dhs',dims) + call define_rest_field(File,'ffrac',dims) + end if + + if (tr_brine) then + call define_rest_field(File,'fbrn',dims) + call define_rest_field(File,'first_ice',dims) + endif + + if (skl_bgc) then + do k = 1, n_algae + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_N'//trim(nchar) ,dims) + enddo + if (tr_bgc_C) then + ! do k = 1, n_algae + ! write(nchar,'(i3.3)') k + ! call define_rest_field(File,'bgc_C'//trim(nchar) ,dims) + ! enddo + do k = 1, n_doc + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_DOC'//trim(nchar) ,dims) + enddo + do k = 1, n_dic + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_DIC'//trim(nchar) ,dims) + enddo + endif + if (tr_bgc_chl) then + do k = 1, n_algae + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_chl'//trim(nchar) ,dims) + enddo + endif + call define_rest_field(File,'bgc_Nit' ,dims) + if (tr_bgc_Am) & + call define_rest_field(File,'bgc_Am' ,dims) + if (tr_bgc_Sil) & + call define_rest_field(File,'bgc_Sil' ,dims) + if (tr_bgc_hum) & + call define_rest_field(File,'bgc_hum' ,dims) + if (tr_bgc_DMS) then + call define_rest_field(File,'bgc_DMSPp',dims) + call define_rest_field(File,'bgc_DMSPd',dims) + call define_rest_field(File,'bgc_DMS' ,dims) + endif + if (tr_bgc_PON) & + call define_rest_field(File,'bgc_PON' ,dims) + if (tr_bgc_DON) then + do k = 1, n_don + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_DON'//trim(nchar) ,dims) + enddo + endif + if (tr_bgc_Fe ) then + do k = 1, n_fed + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_Fed'//trim(nchar) ,dims) + enddo + do k = 1, n_fep + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_Fep'//trim(nchar) ,dims) + enddo + endif + endif !skl_bgc + if (solve_zsal) & + call define_rest_field(File,'Rayleigh',dims) + + !----------------------------------------------------------------- + ! 4D restart fields, written as layers of 3D + !----------------------------------------------------------------- + + do k=1,nilyr + write(nchar,'(i3.3)') k + call define_rest_field(File,'sice'//trim(nchar),dims) + call define_rest_field(File,'qice'//trim(nchar),dims) + enddo + + do k=1,nslyr + write(nchar,'(i3.3)') k + call define_rest_field(File,'qsno'//trim(nchar),dims) + enddo + + if (tr_aero) then + do k=1,n_aero + write(nchar,'(i3.3)') k + call define_rest_field(File,'aerosnossl'//nchar, dims) + call define_rest_field(File,'aerosnoint'//nchar, dims) + call define_rest_field(File,'aeroicessl'//nchar, dims) + call define_rest_field(File,'aeroiceint'//nchar, dims) + enddo + endif + + if (solve_zsal) then + do k = 1, nblyr + write(nchar,'(i3.3)') k + call define_rest_field(File,'zSalinity'//trim(nchar),dims) + enddo + endif + if (z_tracers) then + if (tr_zaero) then + do n = 1, n_zaero + write(ncharb,'(i3.3)') n + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'zaero'//trim(ncharb)//trim(nchar),dims) + enddo !k + enddo !n + endif !tr_zaero + if (tr_bgc_Nit) then + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_Nit'//trim(nchar),dims) + enddo + endif + if (tr_bgc_N) then + do n = 1, n_algae + write(ncharb,'(i3.3)') n + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_N'//trim(ncharb)//trim(nchar),dims) + enddo + enddo + endif + if (tr_bgc_C) then + ! do n = 1, n_algae + ! write(ncharb,'(i3.3)') n + ! do k = 1, nblyr+3 + ! write(nchar,'(i3.3)') k + ! call + ! define_rest_field(File,'bgc_C'//trim(ncharb)//trim(nchar),dims) + ! enddo + ! enddo + do n = 1, n_doc + write(ncharb,'(i3.3)') n + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_DOC'//trim(ncharb)//trim(nchar),dims) + enddo + enddo + do n = 1, n_dic + write(ncharb,'(i3.3)') n + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_DIC'//trim(ncharb)//trim(nchar),dims) + enddo + enddo + endif + if (tr_bgc_chl) then + do n = 1, n_algae + write(ncharb,'(i3.3)') n + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_chl'//trim(ncharb)//trim(nchar),dims) + enddo + enddo + endif + if (tr_bgc_Am) then + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_Am'//trim(nchar),dims) + enddo + endif + if (tr_bgc_Sil) then + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_Sil'//trim(nchar),dims) + enddo + endif + if (tr_bgc_hum) then + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_hum'//trim(nchar),dims) + enddo + endif + if (tr_bgc_DMS) then + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_DMSPp'//trim(nchar),dims) + call define_rest_field(File,'bgc_DMSPd'//trim(nchar),dims) + call define_rest_field(File,'bgc_DMS'//trim(nchar),dims) + enddo + endif + if (tr_bgc_PON) then + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_PON'//trim(nchar),dims) + enddo + endif + if (tr_bgc_DON) then + do n = 1, n_don + write(ncharb,'(i3.3)') n + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_DON'//trim(ncharb)//trim(nchar),dims) + enddo + enddo + endif + if (tr_bgc_Fe ) then + do n = 1, n_fed + write(ncharb,'(i3.3)') n + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_Fed'//trim(ncharb)//trim(nchar),dims) + enddo + enddo + do n = 1, n_fep + write(ncharb,'(i3.3)') n + do k = 1, nblyr+3 + write(nchar,'(i3.3)') k + call define_rest_field(File,'bgc_Fep'//trim(ncharb)//trim(nchar),dims) + enddo + enddo + endif + do k = 1, nbtrcr + write(nchar,'(i3.3)') k + call define_rest_field(File,'zbgc_frac'//trim(nchar),dims) + enddo + endif !z_tracers + + deallocate(dims) + status = pio_enddef(File) + + call ice_pio_initdecomp(iodesc=iodesc2d) + call ice_pio_initdecomp(ndim3=ncat , iodesc=iodesc3d_ncat, remap=.true.) + + endif + + if (my_task == master_task) then + write(nu_diag,*) 'Writing ',filename(1:lenstr(filename)) + endif + + end subroutine init_restart_write + +!======================================================================= + +! Reads a single restart field +! author David A Bailey, NCAR + + subroutine read_restart_field(nu,nrec,work,atype,vname,ndim3,diag, & + field_loc, field_type) + + use ice_blocks, only: nx_block, ny_block + use ice_communicate, only: my_task, master_task + use ice_constants, only: c0, field_loc_center + use ice_boundary, only: ice_HaloUpdate + use ice_domain, only: halo_info, distrb_info, nblocks + use ice_domain_size, only: max_blocks, ncat + use ice_global_reductions, only: global_minval, global_maxval, global_sum + + integer (kind=int_kind), intent(in) :: & + nu , & ! unit number (not used for netcdf) + ndim3 , & ! third dimension + nrec ! record number (0 for sequential access) + + real (kind=dbl_kind), dimension(nx_block,ny_block,ndim3,max_blocks), intent(inout) :: & + work ! input array (real, 8-byte) + + character (len=4), intent(in) :: & + atype ! format for output array + ! (real/integer, 4-byte/8-byte) + + logical (kind=log_kind), intent(in) :: & + diag ! if true, write diagnostic output + + character (len=*), intent(in) :: vname + + integer (kind=int_kind), optional, intent(in) :: & + field_loc, & ! location of field on staggered grid + field_type ! type of field (scalar, vector, angle) + + ! local variables + + integer (kind=int_kind) :: & + j, & ! dimension counter + n, & ! number of dimensions for variable + ndims, & ! number of variable dimensions + status ! status variable from netCDF routine + + real (kind=dbl_kind) :: amin,amax,asum + + character(len=*), parameter :: subname = '(read_restart_field)' + + if (restart_format == "pio") then + if (my_task == master_task) & + write(nu_diag,*)'Parallel restart file read: ',vname + + call pio_seterrorhandling(File, PIO_BCAST_ERROR) + + status = pio_inq_varid(File,trim(vname),vardesc) + + if (status /= 0) then + call abort_ice(subname//"ERROR: CICE restart? Missing variable: "//trim(vname)) + endif + + status = pio_inq_varndims(File, vardesc, ndims) + + call pio_seterrorhandling(File, PIO_INTERNAL_ERROR) + +! if (ndim3 == ncat .and. ncat>1) then + if (ndim3 == ncat .and. ndims == 3) then + call pio_read_darray(File, vardesc, iodesc3d_ncat, work, status) + if (present(field_loc)) then + do n=1,ndim3 + call ice_HaloUpdate (work(:,:,n,:), halo_info, & + field_loc, field_type) + enddo + endif +! elseif (ndim3 == 1) then + elseif (ndim3 == 1 .and. ndims == 2) then + call pio_read_darray(File, vardesc, iodesc2d, work, status) + if (present(field_loc)) then + call ice_HaloUpdate (work(:,:,1,:), halo_info, & + field_loc, field_type) + endif + else + write(nu_diag,*) "ndim3 not supported ",ndim3 + endif + + if (diag) then + if (ndim3 > 1) then + do n=1,ndim3 + amin = global_minval(work(:,:,n,:),distrb_info) + amax = global_maxval(work(:,:,n,:),distrb_info) + asum = global_sum(work(:,:,n,:), distrb_info, field_loc_center) + if (my_task == master_task) then + write(nu_diag,*) ' min and max =', amin, amax + write(nu_diag,*) ' sum =',asum + endif + enddo + else + amin = global_minval(work(:,:,1,:),distrb_info) + amax = global_maxval(work(:,:,1,:),distrb_info) + asum = global_sum(work(:,:,1,:), distrb_info, field_loc_center) + if (my_task == master_task) then + write(nu_diag,*) ' min and max =', amin, amax + write(nu_diag,*) ' sum =',asum + write(nu_diag,*) '' + endif + endif + + endif + else + call abort_ice(subname//"ERROR: Invalid restart_format: "//trim(restart_format)) + endif + + end subroutine read_restart_field + +!======================================================================= + +! Writes a single restart field. +! author David A Bailey, NCAR + + subroutine write_restart_field(nu,nrec,work,atype,vname,ndim3,diag) + + use ice_blocks, only: nx_block, ny_block + use ice_communicate, only: my_task, master_task + use ice_constants, only: c0, field_loc_center + use ice_domain, only: distrb_info, nblocks + use ice_domain_size, only: max_blocks, ncat + use ice_global_reductions, only: global_minval, global_maxval, global_sum + + integer (kind=int_kind), intent(in) :: & + nu , & ! unit number + ndim3 , & ! third dimension + nrec ! record number (0 for sequential access) + + real (kind=dbl_kind), dimension(nx_block,ny_block,ndim3,max_blocks), intent(in) :: & + work ! input array (real, 8-byte) + + character (len=4), intent(in) :: & + atype ! format for output array + ! (real/integer, 4-byte/8-byte) + + logical (kind=log_kind), intent(in) :: & + diag ! if true, write diagnostic output + + character (len=*), intent(in) :: vname + + ! local variables + + integer (kind=int_kind) :: & + j, & ! dimension counter + n, & ! dimension counter + ndims, & ! number of variable dimensions + status ! status variable from netCDF routine + + real (kind=dbl_kind) :: amin,amax,asum + + character(len=*), parameter :: subname = '(write_restart_field)' + + if (restart_format == "pio") then + if (my_task == master_task) & + write(nu_diag,*)'Parallel restart file write: ',vname + + status = pio_inq_varid(File,trim(vname),vardesc) + + status = pio_inq_varndims(File, vardesc, ndims) + + if (ndims==3) then + call pio_write_darray(File, vardesc, iodesc3d_ncat,work(:,:,:,1:nblocks), & + status, fillval=c0) + elseif (ndims == 2) then + call pio_write_darray(File, vardesc, iodesc2d, work(:,:,1,1:nblocks), & + status, fillval=c0) + else + write(nu_diag,*) "ndims not supported",ndims,ndim3 + endif + + if (diag) then + if (ndim3 > 1) then + do n=1,ndim3 + amin = global_minval(work(:,:,n,:),distrb_info) + amax = global_maxval(work(:,:,n,:),distrb_info) + asum = global_sum(work(:,:,n,:), distrb_info, field_loc_center) + if (my_task == master_task) then + write(nu_diag,*) ' min and max =', amin, amax + write(nu_diag,*) ' sum =',asum + endif + enddo + else + amin = global_minval(work(:,:,1,:),distrb_info) + amax = global_maxval(work(:,:,1,:),distrb_info) + asum = global_sum(work(:,:,1,:), distrb_info, field_loc_center) + if (my_task == master_task) then + write(nu_diag,*) ' min and max =', amin, amax + write(nu_diag,*) ' sum =',asum + endif + endif + endif + else + call abort_ice(subname//"ERROR: Invalid restart_format: "//trim(restart_format)) + endif + + end subroutine write_restart_field + +!======================================================================= + +! Finalize the restart file. +! author David A Bailey, NCAR + + subroutine final_restart() + + use ice_calendar, only: istep1, time, time_forc + use ice_communicate, only: my_task, master_task + + character(len=*), parameter :: subname = '(final_restart)' + + if (restart_format == 'pio') then + call PIO_freeDecomp(File,iodesc2d) + call PIO_freeDecomp(File,iodesc3d_ncat) + call pio_closefile(File) + endif + + if (my_task == master_task) & + write(nu_diag,*) 'Restart read/written ',istep1,time,time_forc + + end subroutine final_restart + +!======================================================================= + +! Defines a restart field +! author David A Bailey, NCAR + + subroutine define_rest_field(File, vname, dims) + + type(file_desc_t) , intent(in) :: File + character (len=*) , intent(in) :: vname + integer (kind=int_kind), intent(in) :: dims(:) + + integer (kind=int_kind) :: & + status ! status variable from netCDF routine + + character(len=*), parameter :: subname = '(define_rest_field)' + + status = pio_def_var(File,trim(vname),pio_double,dims,vardesc) + + end subroutine define_rest_field + +!======================================================================= + + end module ice_restart + +!======================================================================= diff --git a/cicecore/drivers/nuopc/CICE_FinalMod.F90 b/cicecore/drivers/nuopc/CICE_FinalMod.F90 index a3b4e6084..c2331e4e5 100644 --- a/cicecore/drivers/nuopc/CICE_FinalMod.F90 +++ b/cicecore/drivers/nuopc/CICE_FinalMod.F90 @@ -1,4 +1,3 @@ -! SVN:$Id: CICE_FinalMod.F90 744 2013-09-27 22:53:24Z eclare $ !======================================================================= ! ! This module contains routines for the final exit of the CICE model, @@ -12,6 +11,10 @@ module CICE_FinalMod use ice_kinds_mod + use ice_communicate, only: my_task, master_task + use ice_exit, only: end_run, abort_ice + use ice_fileunits, only: nu_diag, release_all_fileunits + use icepack_intfc, only: icepack_warnings_flush, icepack_warnings_aborted implicit none private @@ -27,16 +30,16 @@ module CICE_FinalMod subroutine CICE_Finalize - use ice_exit, only: end_run - use ice_fileunits, only: nu_diag, release_all_fileunits use ice_restart_shared, only: runid use ice_timers, only: ice_timer_stop, ice_timer_print_all, timer_total + character(len=*), parameter :: subname = '(CICE_Finalize)' + !------------------------------------------------------------------- ! stop timers and print timer info !------------------------------------------------------------------- -! call ice_timer_stop(timer_total) ! stop timing entire run + call ice_timer_stop(timer_total) ! stop timing entire run call ice_timer_print_all(stats=.false.) ! print timing information !echmod if (nu_diag /= 6) close (nu_diag) ! diagnostic output @@ -52,10 +55,6 @@ subroutine CICE_Finalize ! quit MPI !------------------------------------------------------------------- -#ifdef CESMCOUPLED -#define coupled -#endif - #ifndef coupled call end_run ! quit MPI #endif @@ -72,9 +71,9 @@ end subroutine CICE_Finalize subroutine writeout_finished_file() use ice_restart_shared, only: restart_dir - use ice_communicate, only: my_task, master_task character(len=char_len_long) :: filename + character(len=*), parameter :: subname = '(writeout_finished_file)' if (my_task == master_task) then diff --git a/cicecore/drivers/nuopc/CICE_InitMod.F90 b/cicecore/drivers/nuopc/CICE_InitMod.F90 index 2e845f2c7..56d30cb72 100644 --- a/cicecore/drivers/nuopc/CICE_InitMod.F90 +++ b/cicecore/drivers/nuopc/CICE_InitMod.F90 @@ -1,4 +1,3 @@ -! SVN:$Id: CICE_InitMod.F90 746 2013-09-28 22:47:56Z eclare $ !======================================================================= ! ! This module contains the CICE initialization routine that sets model @@ -14,154 +13,260 @@ module CICE_InitMod use ice_kinds_mod + use ice_exit, only: abort_ice + use ice_fileunits, only: init_fileunits, nu_diag + use icepack_intfc, only: icepack_aggregate + use icepack_intfc, only: icepack_init_itd, icepack_init_itd_hist + use icepack_intfc, only: icepack_configure + use icepack_intfc, only: icepack_warnings_flush, icepack_warnings_aborted + use icepack_intfc, only: icepack_query_parameters, icepack_query_tracer_flags, & + icepack_query_tracer_indices, icepack_query_tracer_numbers implicit none private - public :: cice_init - -#define coupled + public :: CICE_Initialize, cice_init !======================================================================= contains +!======================================================================= + +! Initialize the basic state, grid and all necessary parameters for +! running the CICE model. Return the initial state in routine +! export state. +! Note: This initialization driver is designed for standalone and +! CESM-coupled applications. For other +! applications (e.g., standalone CAM), this driver would be +! replaced by a different driver that calls subroutine cice_init, +! where most of the work is done. + + subroutine CICE_Initialize + + character(len=*), parameter :: subname='(CICE_Initialize)' + + !-------------------------------------------------------------------- + ! model initialization + !-------------------------------------------------------------------- + + call cice_init + + end subroutine CICE_Initialize + !======================================================================= ! ! Initialize CICE model. subroutine cice_init(mpicom_ice) - ! Initialize the basic state, grid and all necessary parameters for - ! running the CICE model. Return the initial state in routine - ! export state. - - use ice_aerosol , only: faero_default - use ice_algae , only: get_forcing_bgc - use ice_calendar , only: dt, dt_dyn, write_ic, init_calendar, calendar, time - use ice_communicate , only: init_communicate - use ice_diagnostics , only: init_diags - use ice_domain , only: init_domain_blocks - use ice_dyn_eap , only: init_eap - use ice_dyn_shared , only: kdyn, init_evp - use ice_fileunits , only: init_fileunits - use ice_flux , only: init_coupler_flux, init_history_therm - use ice_flux , only: init_history_dyn, init_flux_atm, init_flux_ocn - use ice_forcing , only: init_forcing_ocn, init_forcing_atmo - use ice_forcing , only: get_forcing_atmo, get_forcing_ocn - use ice_grid , only: init_grid1, init_grid2 - use ice_history , only: init_hist, accum_hist - use ice_restart_shared , only: restart, runid, runtype - use ice_init , only: input_data, init_state - use ice_itd , only: init_itd + use ice_arrays_column, only: hin_max, c_hi_range, alloc_arrays_column + use ice_state, only: alloc_state + use ice_flux_bgc, only: alloc_flux_bgc + use ice_calendar, only: dt, dt_dyn, time, istep, istep1, write_ic, & + init_calendar, calendar + use ice_communicate, only: init_communicate, my_task, master_task + use ice_diagnostics, only: init_diags + use ice_domain, only: init_domain_blocks + use ice_domain_size, only: ncat + use ice_dyn_eap, only: init_eap, alloc_dyn_eap + use ice_dyn_shared, only: kdyn, init_evp, alloc_dyn_shared + use ice_flux, only: init_coupler_flux, init_history_therm, & + init_history_dyn, init_flux_atm, init_flux_ocn, alloc_flux + use ice_forcing, only: init_forcing_ocn, init_forcing_atmo, & + get_forcing_atmo, get_forcing_ocn, alloc_forcing + use ice_forcing_bgc, only: get_forcing_bgc, get_atm_bgc, & + faero_default, faero_optics, alloc_forcing_bgc + use ice_grid, only: init_grid1, init_grid2, alloc_grid + use ice_history, only: init_hist, accum_hist + use ice_restart_shared, only: restart, runid, runtype + use ice_init, only: input_data, init_state + use ice_init_column, only: init_thermo_vertical, init_shortwave, init_zbgc, input_zbgc, count_tracers use ice_kinds_mod - use ice_restoring , only: ice_HaloRestore_init - use ice_shortwave , only: init_shortwave - use ice_state , only: tr_aero - use ice_therm_vertical , only: init_thermo_vertical - use ice_timers , only: timer_total, init_ice_timers, ice_timer_start, ice_timer_stop - use ice_transport_driver , only: init_transport - use ice_zbgc , only: init_zbgc - use ice_zbgc_shared , only: skl_bgc - - ! !INPUT/OUTPUT PARAMETERS: - integer (kind=int_kind), optional, intent(in) :: mpicom_ice ! communicator for sequential ccsm - !---------------------------------------------------- - - call init_communicate(mpicom_ice) ! initial setup for message passing - call init_fileunits ! unit numbers - call input_data ! namelist variables - if (trim(runid) == 'bering') then - call check_finished_file - end if - call init_zbgc ! vertical biogeochemistry namelist - - call init_domain_blocks ! set up block decomposition - call init_grid1 ! domain distribution - call init_ice_timers ! initialize all timers - call ice_timer_start(timer_total) ! start timing entire run - call init_grid2 ! grid variables - - call init_calendar ! initialize some calendar stuff - call init_hist (dt) ! initialize output history file + use ice_restoring, only: ice_HaloRestore_init + use ice_timers, only: timer_total, init_ice_timers, ice_timer_start + use ice_transport_driver, only: init_transport +#ifdef popcice + use drv_forcing, only: sst_sss +#endif + + integer (kind=int_kind), optional, intent(in) :: & + mpicom_ice ! communicator for sequential ccsm + + logical(kind=log_kind) :: tr_aero, tr_zaero, skl_bgc, z_tracers + character(len=*), parameter :: subname = '(cice_init)' + + call init_communicate(mpicom_ice) ! initial setup for message passing + call init_fileunits ! unit numbers + + call icepack_configure() ! initialize icepack + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(trim(subname), & + file=__FILE__,line= __LINE__) + + call input_data ! namelist variables + call input_zbgc ! vertical biogeochemistry namelist + call count_tracers ! count tracers + + call init_domain_blocks ! set up block decomposition + call init_grid1 ! domain distribution + call alloc_grid ! allocate grid arrays + call alloc_arrays_column ! allocate column arrays + call alloc_state ! allocate state arrays + call alloc_dyn_shared ! allocate dyn shared arrays + call alloc_flux_bgc ! allocate flux_bgc arrays + call alloc_flux ! allocate flux arrays + call init_ice_timers ! initialize all timers + call ice_timer_start(timer_total) ! start timing entire run + call init_grid2 ! grid variables + call init_zbgc ! vertical biogeochemistry initialization + + call init_calendar ! initialize some calendar stuff + call init_hist (dt) ! initialize output history file if (kdyn == 2) then - call init_eap (dt_dyn) ! define eap dynamics parameters, variables - else ! for both kdyn = 0 or 1 - call init_evp (dt_dyn) ! define evp dynamics parameters, variables + call alloc_dyn_eap ! allocate dyn_eap arrays + call init_eap (dt_dyn) ! define eap dynamics parameters, variables + else ! for both kdyn = 0 or 1 + call init_evp (dt_dyn) ! define evp dynamics parameters, variables + endif + + call init_coupler_flux ! initialize fluxes exchanged with coupler +#ifdef popcice + call sst_sss ! POP data for CICE initialization +#endif + call init_thermo_vertical ! initialize vertical thermodynamics + + call icepack_init_itd(ncat, hin_max) ! ice thickness distribution + if (my_task == master_task) then + call icepack_init_itd_hist(ncat, hin_max, c_hi_range) ! output endif + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + + call calendar(time) ! determine the initial date + + call init_forcing_ocn(dt) ! initialize sss and sst from data + call init_state ! initialize the ice state + call init_transport ! initialize horizontal transport + call ice_HaloRestore_init ! restored boundary conditions - call init_coupler_flux ! initialize fluxes exchanged with coupler - call init_thermo_vertical ! initialize vertical thermodynamics - call init_itd ! initialize ice thickness distribution - call calendar(time) ! determine the initial date + call icepack_query_parameters(skl_bgc_out=skl_bgc, z_tracers_out=z_tracers) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(trim(subname), & + file=__FILE__,line= __LINE__) - call init_forcing_ocn(dt) ! initialize sss and sst from data - call init_state ! initialize the ice state - call init_transport ! initialize horizontal transport - call ice_HaloRestore_init ! restored boundary conditions + if (skl_bgc .or. z_tracers) call alloc_forcing_bgc ! allocate biogeochemistry arrays - call init_restart ! initialize restart variables + call init_restart ! initialize restart variables + call init_diags ! initialize diagnostic output points + call init_history_therm ! initialize thermo history variables + call init_history_dyn ! initialize dynamic history variables - call init_diags ! initialize diagnostic output points - call init_history_therm ! initialize thermo history variables - call init_history_dyn ! initialize dynamic history variables + call icepack_query_tracer_flags(tr_aero_out=tr_aero, tr_zaero_out=tr_zaero) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(trim(subname), & + file=__FILE__,line= __LINE__) + + if (tr_aero .or. tr_zaero) call faero_optics !initialize aerosol optical + !property tables ! Initialize shortwave components using swdn from previous timestep ! if restarting. These components will be scaled to current forcing ! in prep_radiation. - - if (trim(runtype) == 'continue' .or. restart) then + if (trim(runtype) == 'continue' .or. restart) & call init_shortwave ! initialize radiative transfer - end if - !-------------------------------------------------------------------- - ! coupler communication or forcing data initialization - !-------------------------------------------------------------------- + !-------------------------------------------------------------------- + ! coupler communication or forcing data initialization + !-------------------------------------------------------------------- call init_forcing_atmo ! initialize atmospheric forcing (standalone) - if (runtype == 'initial' .and. .not. restart) then +#ifndef coupled +#ifndef CESMCOUPLED + call get_forcing_atmo ! atmospheric forcing from data + call get_forcing_ocn(dt) ! ocean forcing from data + + ! aerosols + ! if (tr_aero) call faero_data ! data file + ! if (tr_zaero) call fzaero_data ! data file (gx1) + if (tr_aero .or. tr_zaero) call faero_default ! default values + if (skl_bgc .or. z_tracers) call get_forcing_bgc ! biogeochemistry +#endif +#endif + if (z_tracers) call get_atm_bgc ! biogeochemistry + + if (runtype == 'initial' .and. .not. restart) & call init_shortwave ! initialize radiative transfer using current swdn - end if call init_flux_atm ! initialize atmosphere fluxes sent to coupler call init_flux_ocn ! initialize ocean fluxes sent to coupler - ! if (write_ic) call accum_hist(dt) ! write initial conditions - - call ice_timer_stop(timer_total) ! stop timing entire run +! if (write_ic) call accum_hist(dt) ! write initial conditions end subroutine cice_init !======================================================================= - subroutine init_restart() - - use ice_aerosol , only: init_aerosol - use ice_age , only: init_age, restart_age, read_restart_age - use ice_blocks , only: nx_block, ny_block - use ice_brine , only: init_hbrine - use ice_calendar , only: time, calendar - use ice_domain , only: nblocks - use ice_domain_size , only: ncat, max_ntrcr - use ice_dyn_eap , only: read_restart_eap - use ice_dyn_shared , only: kdyn - use ice_firstyear , only: init_fy, restart_FY, read_restart_FY - use ice_flux , only: sss - use ice_grid , only: tmask - use ice_init , only: ice_ic - use ice_itd , only: aggregate - use ice_lvl , only: init_lvl, restart_lvl, read_restart_lvl - use ice_meltpond_cesm , only: init_meltponds_cesm, restart_pond_cesm, read_restart_pond_cesm - use ice_meltpond_lvl , only: init_meltponds_lvl, restart_pond_lvl, read_restart_pond_lvl, dhsn - use ice_meltpond_topo , only: init_meltponds_topo, restart_pond_topo, read_restart_pond_topo - use ice_restart_shared , only: runtype, restart - use ice_restart_driver , only: restartfile, restartfile_v4 - use ice_state ! almost everything - use ice_zbgc , only: init_bgc - use ice_zbgc_shared , only: skl_bgc - - integer(kind=int_kind) :: iblk, ltmp + subroutine init_restart + + use ice_arrays_column, only: dhsn + use ice_blocks, only: nx_block, ny_block + use ice_calendar, only: time, calendar + use ice_constants, only: c0 + use ice_domain, only: nblocks + use ice_domain_size, only: ncat, n_aero + use ice_dyn_eap, only: read_restart_eap + use ice_dyn_shared, only: kdyn + use ice_grid, only: tmask + use ice_init, only: ice_ic + use ice_init_column, only: init_age, init_FY, init_lvl, & + init_meltponds_cesm, init_meltponds_lvl, init_meltponds_topo, & + init_aerosol, init_hbrine, init_bgc + use ice_restart_column, only: restart_age, read_restart_age, & + restart_FY, read_restart_FY, restart_lvl, read_restart_lvl, & + restart_pond_cesm, read_restart_pond_cesm, & + restart_pond_lvl, read_restart_pond_lvl, & + restart_pond_topo, read_restart_pond_topo, & + restart_aero, read_restart_aero, & + restart_hbrine, read_restart_hbrine, & + restart_zsal, restart_bgc + use ice_restart_driver, only: restartfile + use ice_restart_shared, only: runtype, restart + use ice_state ! almost everything + + integer(kind=int_kind) :: & + i, j , & ! horizontal indices + iblk ! block index + logical(kind=log_kind) :: & + tr_iage, tr_FY, tr_lvl, tr_pond_cesm, tr_pond_lvl, & + tr_pond_topo, tr_aero, tr_brine, & + skl_bgc, z_tracers, solve_zsal + integer(kind=int_kind) :: & + ntrcr + integer(kind=int_kind) :: & + nt_alvl, nt_vlvl, nt_apnd, nt_hpnd, nt_ipnd, & + nt_iage, nt_FY, nt_aero + + character(len=*), parameter :: subname = '(init_restart)' + + call icepack_query_tracer_numbers(ntrcr_out=ntrcr) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + + call icepack_query_parameters(skl_bgc_out=skl_bgc, & + z_tracers_out=z_tracers, solve_zsal_out=solve_zsal) + call icepack_query_tracer_flags(tr_iage_out=tr_iage, tr_FY_out=tr_FY, & + tr_lvl_out=tr_lvl, tr_pond_cesm_out=tr_pond_cesm, tr_pond_lvl_out=tr_pond_lvl, & + tr_pond_topo_out=tr_pond_topo, tr_aero_out=tr_aero, tr_brine_out=tr_brine) + call icepack_query_tracer_indices(nt_alvl_out=nt_alvl, nt_vlvl_out=nt_vlvl, & + nt_apnd_out=nt_apnd, nt_hpnd_out=nt_hpnd, nt_ipnd_out=nt_ipnd, & + nt_iage_out=nt_iage, nt_FY_out=nt_FY, nt_aero_out=nt_aero) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) if (trim(runtype) == 'continue') then ! start from core restart file @@ -169,14 +274,11 @@ subroutine init_restart() call calendar(time) ! update time parameters if (kdyn == 2) call read_restart_eap ! EAP else if (restart) then ! ice_ic = core restart file - ltmp = len_trim(ice_ic) - if (ice_ic(ltmp-2:ltmp) == '.nc') then - call restartfile (ice_ic) ! or 'default' or 'none' - else - call restartfile_v4 (ice_ic) ! CICE v4.1 binary restart file - !!! uncomment if EAP restart data exists - ! if (kdyn == 2) call read_restart_eap - endif + call restartfile (ice_ic) ! or 'default' or 'none' + !!! uncomment to create netcdf + ! call restartfile_v4 (ice_ic) ! CICE v4.1 binary restart file + !!! uncomment if EAP restart data exists + ! if (kdyn == 2) call read_restart_eap endif ! tracers @@ -188,7 +290,7 @@ subroutine init_restart() call read_restart_age else do iblk = 1, nblocks - call init_age(nx_block, ny_block, ncat, trcrn(:,:,nt_iage,:,iblk)) + call init_age(trcrn(:,:,nt_iage,:,iblk)) enddo ! iblk endif endif @@ -199,7 +301,7 @@ subroutine init_restart() call read_restart_FY else do iblk = 1, nblocks - call init_FY(nx_block, ny_block, ncat, trcrn(:,:,nt_FY,:,iblk)) + call init_FY(trcrn(:,:,nt_FY,:,iblk)) enddo ! iblk endif endif @@ -210,8 +312,8 @@ subroutine init_restart() call read_restart_lvl else do iblk = 1, nblocks - call init_lvl(nx_block, ny_block, ncat, & - trcrn(:,:,nt_alvl,:,iblk), trcrn(:,:,nt_vlvl,:,iblk)) + call init_lvl(iblk,trcrn(:,:,nt_alvl,:,iblk), & + trcrn(:,:,nt_vlvl,:,iblk)) enddo ! iblk endif endif @@ -223,8 +325,8 @@ subroutine init_restart() call read_restart_pond_cesm else do iblk = 1, nblocks - call init_meltponds_cesm(nx_block, ny_block, ncat, & - trcrn(:,:,nt_apnd,:,iblk), trcrn(:,:,nt_hpnd,:,iblk)) + call init_meltponds_cesm(trcrn(:,:,nt_apnd,:,iblk), & + trcrn(:,:,nt_hpnd,:,iblk)) enddo ! iblk endif endif @@ -236,9 +338,10 @@ subroutine init_restart() call read_restart_pond_lvl else do iblk = 1, nblocks - call init_meltponds_lvl(nx_block, ny_block, ncat, & - trcrn(:,:,nt_apnd,:,iblk), trcrn(:,:,nt_hpnd,:,iblk), & - trcrn(:,:,nt_ipnd,:,iblk), dhsn(:,:,:,iblk)) + call init_meltponds_lvl(trcrn(:,:,nt_apnd,:,iblk), & + trcrn(:,:,nt_hpnd,:,iblk), & + trcrn(:,:,nt_ipnd,:,iblk), & + dhsn(:,:,:,iblk)) enddo ! iblk endif endif @@ -250,15 +353,38 @@ subroutine init_restart() call read_restart_pond_topo else do iblk = 1, nblocks - call init_meltponds_topo(nx_block, ny_block, ncat, & - trcrn(:,:,nt_apnd,:,iblk), trcrn(:,:,nt_hpnd,:,iblk), & - trcrn(:,:,nt_ipnd,:,iblk)) + call init_meltponds_topo(trcrn(:,:,nt_apnd,:,iblk), & + trcrn(:,:,nt_hpnd,:,iblk), & + trcrn(:,:,nt_ipnd,:,iblk)) + enddo ! iblk + endif ! .not. restart_pond + endif + if (tr_aero) then ! ice aerosol + if (trim(runtype) == 'continue') restart_aero = .true. + if (restart_aero) then + call read_restart_aero + else + do iblk = 1, nblocks + call init_aerosol(trcrn(:,:,nt_aero:nt_aero+4*n_aero-1,:,iblk)) enddo ! iblk - endif ! .not restart_pond + endif ! .not. restart_aero endif - if (tr_aero) call init_aerosol ! ice aerosol - if (tr_brine) call init_hbrine ! brine height tracer - if (skl_bgc) call init_bgc ! biogeochemistry + + if (trim(runtype) == 'continue') then + if (tr_brine) & + restart_hbrine = .true. + if (solve_zsal) & + restart_zsal = .true. + if (skl_bgc .or. z_tracers) & + restart_bgc = .true. + endif + + if (tr_brine .or. skl_bgc) then ! brine height tracer + call init_hbrine + if (tr_brine .and. restart_hbrine) call read_restart_hbrine + endif + + if (solve_zsal .or. skl_bgc .or. z_tracers) call init_bgc ! biogeochemistry !----------------------------------------------------------------- ! aggregate tracers @@ -266,54 +392,38 @@ subroutine init_restart() !$OMP PARALLEL DO PRIVATE(iblk) do iblk = 1, nblocks - - call aggregate (nx_block, ny_block, & - aicen(:,:,:,iblk), & - trcrn(:,:,:,:,iblk),& - vicen(:,:,:,iblk), & - vsnon(:,:,:,iblk), & - aice (:,:, iblk), & - trcr (:,:,:,iblk), & - vice (:,:, iblk), & - vsno (:,:, iblk), & - aice0(:,:, iblk), & - tmask(:,:, iblk), & - max_ntrcr, & - trcr_depend) - + do j = 1, ny_block + do i = 1, nx_block + if (tmask(i,j,iblk)) then + call icepack_aggregate (ncat, & + aicen(i,j,:,iblk), & + trcrn(i,j,:,:,iblk),& + vicen(i,j,:,iblk), & + vsnon(i,j,:,iblk), & + aice (i,j, iblk), & + trcr (i,j,:,iblk), & + vice (i,j, iblk), & + vsno (i,j, iblk), & + aice0(i,j, iblk), & + ntrcr, & + trcr_depend, & + trcr_base, & + n_trcr_strata, & + nt_strata) + else + ! tcraig, reset all tracer values on land to zero + trcrn(i,j,:,:,iblk) = c0 + endif + enddo + enddo enddo !$OMP END PARALLEL DO - end subroutine init_restart - -!======================================================================= -! -! Check whether a file indicating that the previous run finished cleanly -! If so, then do not continue the current restart. This is needed only -! for runs on machine 'bering' (set using runid = 'bering'). -! -! author: Adrian Turner, LANL - - subroutine check_finished_file() - - use ice_communicate , only: my_task, master_task - use ice_exit , only: abort_ice - use ice_restart_shared , only: restart_dir + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) - character(len=char_len_long) :: filename - logical :: lexist = .false. - - if (my_task == master_task) then - - filename = trim(restart_dir)//"finished" - inquire(file=filename, exist=lexist) - if (lexist) then - call abort_ice("Found already finished file - quitting") - end if - - endif - - end subroutine check_finished_file + end subroutine init_restart !======================================================================= diff --git a/cicecore/drivers/nuopc/CICE_RunMod.F90 b/cicecore/drivers/nuopc/CICE_RunMod.F90 index cc9fab18c..9ed65826b 100644 --- a/cicecore/drivers/nuopc/CICE_RunMod.F90 +++ b/cicecore/drivers/nuopc/CICE_RunMod.F90 @@ -7,529 +7,653 @@ ! William H. Lipscomb, LANL ! ! 2006 ECH: moved exit timeLoop to prevent execution of unnecessary timestep -! 2006 ECH: Streamlined for efficiency +! 2006 ECH: Streamlined for efficiency ! 2006 ECH: Converted to free source form (F90) ! 2007 BPB: Modified Delta-Eddington shortwave interface ! 2008 ECH: moved ESMF code to its own driver -module CICE_RunMod + module CICE_RunMod - use ice_kinds_mod - use perf_mod , only: t_startf, t_stopf, t_barrierf - use ice_fileunits , only: nu_diag + use ice_kinds_mod + use perf_mod, only : t_startf, t_stopf, t_barrierf + use ice_fileunits, only: nu_diag + use ice_arrays_column, only: oceanmixed_ice + use ice_constants, only: c0, c1 + use ice_constants, only: field_loc_center, field_type_scalar + use ice_exit, only: abort_ice + use icepack_intfc, only: icepack_warnings_flush, icepack_warnings_aborted + use icepack_intfc, only: icepack_max_aero + use icepack_intfc, only: icepack_query_parameters + use icepack_intfc, only: icepack_query_tracer_flags, icepack_query_tracer_numbers - implicit none - private + implicit none + private + public :: CICE_Run, ice_step - public :: CICE_Run - public :: ice_step +!======================================================================= - private :: coupling_prep - private :: sfcflux_to_ocn + contains !======================================================================= -contains -!======================================================================= +! +! This is the main driver routine for advancing CICE forward in time. +! +! author Elizabeth C. Hunke, LANL +! Philip W. Jones, LANL +! William H. Lipscomb, LANL + + subroutine CICE_Run + + use ice_calendar, only: istep, istep1, time, dt, stop_now, calendar + use ice_forcing, only: get_forcing_atmo, get_forcing_ocn, atm_data_type + use ice_forcing_bgc, only: get_forcing_bgc, get_atm_bgc, & + faero_default + use ice_flux, only: init_flux_atm, init_flux_ocn + use ice_timers, only: ice_timer_start, ice_timer_stop, & + timer_couple, timer_step + + logical (kind=log_kind) :: & + tr_aero, tr_zaero, skl_bgc, z_tracers + character(len=*), parameter :: subname = '(CICE_Run)' + + !-------------------------------------------------------------------- + ! initialize error code and step timer + !-------------------------------------------------------------------- + + call ice_timer_start(timer_step) ! start timing entire run + + call icepack_query_parameters(skl_bgc_out=skl_bgc, z_tracers_out=z_tracers) + call icepack_query_tracer_flags(tr_aero_out=tr_aero, tr_zaero_out=tr_zaero) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + + !-------------------------------------------------------------------- + ! timestep loop + !-------------------------------------------------------------------- + +! timeLoop: do + + istep = istep + 1 ! update time step counters + istep1 = istep1 + 1 + time = time + dt ! determine the time and date + + call ice_timer_start(timer_couple) ! atm/ocn coupling + +#ifndef coupled +#ifndef CESMCOUPLED + call get_forcing_atmo ! atmospheric forcing from data + call get_forcing_ocn(dt) ! ocean forcing from data + + ! aerosols + ! if (tr_aero) call faero_data ! data file + ! if (tr_zaero) call fzaero_data ! data file (gx1) + if (tr_aero .or. tr_zaero) call faero_default ! default values + + if (skl_bgc .or. z_tracers) call get_forcing_bgc ! biogeochemistry +#endif +#endif + if (z_tracers) call get_atm_bgc ! biogeochemistry + + call init_flux_atm ! Initialize atmosphere fluxes sent to coupler + call init_flux_ocn ! initialize ocean fluxes sent to coupler + + call calendar(time) ! at the end of the timestep + + call ice_timer_stop(timer_couple) ! atm/ocn coupling + + call ice_step + +! if (stop_now >= 1) exit timeLoop +! enddo timeLoop + + !-------------------------------------------------------------------- + ! end of timestep loop + !-------------------------------------------------------------------- + + call ice_timer_stop(timer_step) ! end timestepping loop timer - subroutine CICE_Run(restart_filename) - - ! This is the main driver routine for advancing CICE forward in time. - ! - ! author Elizabeth C. Hunke, LANL - ! Philip W. Jones, LANL - ! William H. Lipscomb, LANL - - use ice_aerosol , only: faero_default - use ice_algae , only: get_forcing_bgc - use ice_calendar , only: istep, istep1, time, dt, stop_now, calendar - use ice_forcing , only: get_forcing_atmo, get_forcing_ocn - use ice_flux , only: init_flux_atm, init_flux_ocn - use ice_state , only: tr_aero - use ice_timers , only: ice_timer_start, ice_timer_stop - use ice_timers , only: timer_couple, timer_step - use ice_zbgc_shared , only: skl_bgc - - character(len=*), intent(in), optional :: restart_filename - !-------------------------------------------------------------------- - - !-------------------------------------------------------------------- - ! timestep loop - !-------------------------------------------------------------------- - - call ice_timer_start(timer_step) ! start timing entire run - - istep = istep + 1 ! update time step counters - istep1 = istep1 + 1 - time = time + dt ! determine the time and date - - call init_flux_atm ! initialize atmosphere fluxes sent to coupler - call init_flux_ocn ! initialize ocean fluxes sent to coupler - - call calendar(time) ! at the end of the timestep - - if (present(restart_filename)) then - call ice_step(restart_filename) - else - call ice_step() - end if - - call ice_timer_stop(timer_step) ! end timestepping loop timer - - end subroutine CICE_Run - - !======================================================================= - - - subroutine ice_step(restart_filename) - - ! Calls drivers for physics components, some initialization, and output - ! - ! author Elizabeth C. Hunke, LANL - ! William H. Lipscomb, LANL - - use ice_age , only: write_restart_age - use ice_aerosol , only: write_restart_aero - use ice_boundary , only: ice_HaloUpdate - use ice_brine , only: hbrine_diags, write_restart_hbrine - use ice_calendar , only: dt, dt_dyn, ndtd, diagfreq, write_restart, istep, idate, sec - use ice_constants , only: field_loc_center, field_type_scalar - use ice_diagnostics , only: init_mass_diags, runtime_diags, print_points_state - use ice_domain , only: halo_info, nblocks - use ice_domain_size , only: nslyr - use ice_dyn_eap , only: write_restart_eap - use ice_dyn_shared , only: kdyn - use ice_firstyear , only: write_restart_FY - use ice_flux , only: scale_factor, init_history_therm - use ice_history , only: accum_hist - use ice_lvl , only: write_restart_lvl - use ice_restart , only: final_restart - use ice_restart_driver , only: dumpfile - use ice_meltpond_cesm , only: write_restart_pond_cesm - use ice_meltpond_lvl , only: write_restart_pond_lvl - use ice_meltpond_topo , only: write_restart_pond_topo - use ice_restoring , only: restore_ice, ice_HaloRestore - use ice_state , only: nt_qsno, trcrn, tr_iage, tr_FY, tr_lvl - use ice_state , only: nt_qsno, trcrn, tr_iage, tr_FY, tr_lvl - use ice_state , only: tr_pond_cesm, tr_pond_lvl, tr_pond_topo, tr_brine, tr_aero - use ice_step_mod , only: prep_radiation, step_therm1, step_therm2 - use ice_step_mod , only: post_thermo, step_dynamics, step_radiation - use ice_therm_shared , only: calc_Tsfc - use ice_timers , only: ice_timer_start, ice_timer_stop - use ice_timers , only: timer_diags, timer_column, timer_thermo, timer_bound - use ice_timers , only: timer_hist, timer_readwrite - use ice_algae , only: bgc_diags, write_restart_bgc - use ice_zbgc , only: init_history_bgc, biogeochemistry - use ice_zbgc_shared , only: skl_bgc - use ice_communicate , only: MPI_COMM_ICE - use ice_prescribed_mod - - character(len=*), intent(in), optional :: restart_filename - - integer (kind=int_kind) :: & - iblk , & ! block index - k ! dynamics supercycling index - - !----------------------------------------------------------------- - ! restoring on grid boundaries - !----------------------------------------------------------------- - - if (restore_ice) call ice_HaloRestore - - !----------------------------------------------------------------- - ! initialize diagnostics - !----------------------------------------------------------------- - - call ice_timer_start(timer_diags) ! diagnostics/history - call init_mass_diags ! diagnostics per timestep - call init_history_therm - call init_history_bgc - call ice_timer_stop(timer_diags) ! diagnostics/history - - if(prescribed_ice) then ! read prescribed ice - call t_barrierf('cice_run_presc_BARRIER',MPI_COMM_ICE) - call t_startf ('cice_run_presc') - call ice_prescribed_run(idate, sec) - call t_stopf ('cice_run_presc') - endif - - call ice_timer_start(timer_column) ! column physics - call ice_timer_start(timer_thermo) ! thermodynamics - - !$OMP PARALLEL DO PRIVATE(iblk) - do iblk = 1, nblocks - - !----------------------------------------------------------------- - ! Scale radiation fields - !----------------------------------------------------------------- - - if (calc_Tsfc) call prep_radiation (dt, iblk) - - !----------------------------------------------------------------- - ! thermodynamics - !----------------------------------------------------------------- - - call step_therm1 (dt, iblk) ! vertical thermodynamics - call biogeochemistry (dt, iblk) ! biogeochemistry - if (.not.prescribed_ice) & - call step_therm2 (dt, iblk) ! ice thickness distribution thermo - - enddo ! iblk - !$OMP END PARALLEL DO - - call post_thermo (dt) ! finalize thermo update - - call ice_timer_stop(timer_thermo) ! thermodynamics - call ice_timer_stop(timer_column) ! column physics - - !----------------------------------------------------------------- - ! dynamics, transport, ridging - !----------------------------------------------------------------- - - if (.not.prescribed_ice .and. kdyn>0) then - do k = 1, ndtd - call step_dynamics (dt_dyn, ndtd) - enddo - endif - - !----------------------------------------------------------------- - ! albedo, shortwave radiation - !----------------------------------------------------------------- - - call ice_timer_start(timer_column) ! column physics - call ice_timer_start(timer_thermo) ! thermodynamics - - !$OMP PARALLEL DO PRIVATE(iblk) - do iblk = 1, nblocks - - call step_radiation (dt, iblk) - - !----------------------------------------------------------------- - ! get ready for coupling and the next time step - !----------------------------------------------------------------- - - call coupling_prep (iblk) - - enddo ! iblk - !$OMP END PARALLEL DO - - call ice_timer_start(timer_bound) - call ice_HaloUpdate (scale_factor, halo_info, & - field_loc_center, field_type_scalar) - call ice_timer_stop(timer_bound) - - call ice_timer_stop(timer_thermo) ! thermodynamics - call ice_timer_stop(timer_column) ! column physics - - !----------------------------------------------------------------- - ! write data - !----------------------------------------------------------------- - - call ice_timer_start(timer_diags) ! diagnostics - if (mod(istep,diagfreq) == 0) then - call runtime_diags(dt) ! log file - if (skl_bgc) call bgc_diags (dt) - if (tr_brine) call hbrine_diags (dt) - endif - call ice_timer_stop(timer_diags) ! diagnostics - - call ice_timer_start(timer_hist) ! history - call accum_hist (dt) ! history file - call ice_timer_stop(timer_hist) ! history - - call ice_timer_start(timer_readwrite) ! reading/writing - if (write_restart == 1) then - if (present(restart_filename)) then - call dumpfile(filename_spec=restart_filename) ! core variables for restarting - else - call dumpfile() - end if - if (tr_iage) call write_restart_age - if (tr_FY) call write_restart_FY - if (tr_lvl) call write_restart_lvl - if (tr_pond_cesm) call write_restart_pond_cesm - if (tr_pond_lvl) call write_restart_pond_lvl - if (tr_pond_topo) call write_restart_pond_topo - if (tr_aero) call write_restart_aero - if (skl_bgc) call write_restart_bgc - if (tr_brine) call write_restart_hbrine - if (kdyn == 2) call write_restart_eap - call final_restart - endif - - call ice_timer_stop(timer_readwrite) ! reading/writing - - end subroutine ice_step - - !======================================================================= - - subroutine coupling_prep (iblk) - - ! - ! Prepare for coupling - ! - ! authors: Elizabeth C. Hunke, LANL - - use ice_blocks , only: block, nx_block, ny_block - use ice_calendar , only: dt, nstreams - use ice_constants , only: c0, c1, puny, rhofresh - use ice_domain_size , only: ncat - use ice_flux , only: alvdf, alidf, alvdr, alidr, albice, albsno - use ice_flux , only: albpnd, albcnt, apeff_ai, coszen, fpond, fresh, l_mpond_fresh - use ice_flux , only: alvdf_ai, alidf_ai, alvdr_ai, alidr_ai, fhocn_ai - use ice_flux , only: fresh_ai, fsalt_ai, fsalt - use ice_flux , only: fswthru_ai, fhocn, fswthru, scale_factor - use ice_flux , only: fswthruvdr, fswthruvdf, fswthruidr, fswthruidf - use ice_flux , only: swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyt - use ice_flux , only: fsens, flat, fswabs, flwout, evap, Tref, Qref, Uref, faero_ocn - use ice_flux , only: fsurfn_f, flatn_f, scale_fluxes, frzmlt_init, frzmlt, wind - use ice_flux , only: snowfrac - use ice_grid , only: tmask - use ice_ocean , only: oceanmixed_ice, ocean_mixed_layer - use ice_shortwave , only: alvdfn, alidfn, alvdrn, alidrn - use ice_shortwave , only: albicen, albsnon, albpndn, apeffn, snowfracn - use ice_state , only: aicen, aice, aice_init, nbtrcr - use ice_therm_shared , only: calc_Tsfc - use ice_timers , only: timer_couple, ice_timer_start, ice_timer_stop - use ice_zbgc_shared , only: flux_bio, flux_bio_ai - - integer (kind=int_kind), intent(in) :: & - iblk ! block index - - ! local variables - - integer (kind=int_kind) :: & + end subroutine CICE_Run + +!======================================================================= +! +! Calls drivers for physics components, some initialization, and output +! +! author Elizabeth C. Hunke, LANL +! William H. Lipscomb, LANL + + subroutine ice_step + + use ice_boundary, only: ice_HaloUpdate + use ice_calendar, only: dt, dt_dyn, ndtd, diagfreq, write_restart, istep, idate, sec + use ice_diagnostics, only: init_mass_diags, runtime_diags + use ice_diagnostics_bgc, only: hbrine_diags, zsal_diags, bgc_diags + use ice_domain, only: halo_info, nblocks + use ice_dyn_eap, only: write_restart_eap + use ice_dyn_shared, only: kdyn, kridge + use ice_flux, only: scale_factor, init_history_therm, & + daidtt, daidtd, dvidtt, dvidtd, dagedtt, dagedtd + use ice_history, only: accum_hist + use ice_history_bgc, only: init_history_bgc + use ice_restart, only: final_restart + use ice_restart_column, only: write_restart_age, write_restart_FY, & + write_restart_lvl, write_restart_pond_cesm, write_restart_pond_lvl, & + write_restart_pond_topo, write_restart_aero, & + write_restart_bgc, write_restart_hbrine + use ice_restart_driver, only: dumpfile + use ice_restoring, only: restore_ice, ice_HaloRestore + use ice_step_mod, only: prep_radiation, step_therm1, step_therm2, & + update_state, step_dyn_horiz, step_dyn_ridge, step_radiation, & + biogeochemistry, save_init + use ice_timers, only: ice_timer_start, ice_timer_stop, & + timer_diags, timer_column, timer_thermo, timer_bound, & + timer_hist, timer_readwrite + use ice_communicate, only: MPI_COMM_ICE + use ice_prescribed_mod + + integer (kind=int_kind) :: & + iblk , & ! block index + k , & ! dynamics supercycling index + ktherm ! thermodynamics is off when ktherm = -1 + + real (kind=dbl_kind) :: & + offset ! d(age)/dt time offset + + logical (kind=log_kind) :: & + tr_iage, tr_FY, tr_lvl, & + tr_pond_cesm, tr_pond_lvl, tr_pond_topo, tr_brine, tr_aero, & + calc_Tsfc, skl_bgc, solve_zsal, z_tracers + + character(len=*), parameter :: subname = '(ice_step)' + + call icepack_query_parameters(calc_Tsfc_out=calc_Tsfc, skl_bgc_out=skl_bgc, & + solve_zsal_out=solve_zsal, z_tracers_out=z_tracers, ktherm_out=ktherm) + call icepack_query_tracer_flags(tr_iage_out=tr_iage, tr_FY_out=tr_FY, & + tr_lvl_out=tr_lvl, tr_pond_cesm_out=tr_pond_cesm, tr_pond_lvl_out=tr_pond_lvl, & + tr_pond_topo_out=tr_pond_topo, tr_brine_out=tr_brine, tr_aero_out=tr_aero) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + + !----------------------------------------------------------------- + ! restoring on grid boundaries + !----------------------------------------------------------------- + + if (restore_ice) call ice_HaloRestore + + !----------------------------------------------------------------- + ! initialize diagnostics and save initial state values + !----------------------------------------------------------------- + + call ice_timer_start(timer_diags) ! diagnostics/history + call init_mass_diags ! diagnostics per timestep + call init_history_therm + call init_history_bgc + call ice_timer_stop(timer_diags) ! diagnostics/history + + if (prescribed_ice) then ! read prescribed ice + call t_barrierf('cice_run_presc_BARRIER',MPI_COMM_ICE) + call t_startf ('cice_run_presc') + call ice_prescribed_run(idate, sec) + call t_stopf ('cice_run_presc') + endif + + call save_init + + call ice_timer_start(timer_column) ! column physics + call ice_timer_start(timer_thermo) ! thermodynamics + + !$OMP PARALLEL DO PRIVATE(iblk) + do iblk = 1, nblocks + + if (ktherm >= 0) then + + !----------------------------------------------------------------- + ! scale radiation fields + !----------------------------------------------------------------- + + if (calc_Tsfc) call prep_radiation (iblk) + + !----------------------------------------------------------------- + ! thermodynamics and biogeochemistry + !----------------------------------------------------------------- + + call step_therm1 (dt, iblk) ! vertical thermodynamics + call biogeochemistry (dt, iblk) ! biogeochemistry + if (.not.prescribed_ice) & + call step_therm2 (dt, iblk) ! ice thickness distribution thermo + + endif + + enddo ! iblk + !$OMP END PARALLEL DO + + ! clean up, update tendency diagnostics + offset = dt + call update_state (dt, daidtt, dvidtt, dagedtt, offset) + + call ice_timer_stop(timer_thermo) ! thermodynamics + call ice_timer_stop(timer_column) ! column physics + + !----------------------------------------------------------------- + ! dynamics, transport, ridging + !----------------------------------------------------------------- + + if (.not.prescribed_ice) then + do k = 1, ndtd + + ! momentum, stress, transport + call step_dyn_horiz (dt_dyn) + + ! ridging + !$OMP PARALLEL DO PRIVATE(iblk) + do iblk = 1, nblocks + if (kridge > 0) call step_dyn_ridge (dt_dyn, ndtd, iblk) + enddo + !$OMP END PARALLEL DO + + ! clean up, update tendency diagnostics + offset = c0 + call update_state (dt_dyn, daidtd, dvidtd, dagedtd, offset) + + enddo + endif + + !----------------------------------------------------------------- + ! albedo, shortwave radiation + !----------------------------------------------------------------- + + call ice_timer_start(timer_column) ! column physics + call ice_timer_start(timer_thermo) ! thermodynamics + + !$OMP PARALLEL DO PRIVATE(iblk) + do iblk = 1, nblocks + + if (ktherm >= 0) call step_radiation (dt, iblk) + + !----------------------------------------------------------------- + ! get ready for coupling and the next time step + !----------------------------------------------------------------- + + call coupling_prep (iblk) + + enddo ! iblk + !$OMP END PARALLEL DO + + call ice_timer_start(timer_bound) + call ice_HaloUpdate (scale_factor, halo_info, & + field_loc_center, field_type_scalar) + call ice_timer_stop(timer_bound) + + call ice_timer_stop(timer_thermo) ! thermodynamics + call ice_timer_stop(timer_column) ! column physics + + !----------------------------------------------------------------- + ! write data + !----------------------------------------------------------------- + + call ice_timer_start(timer_diags) ! diagnostics + if (mod(istep,diagfreq) == 0) then + call runtime_diags(dt) ! log file + if (solve_zsal) call zsal_diags + if (skl_bgc .or. z_tracers) call bgc_diags + if (tr_brine) call hbrine_diags + endif + call ice_timer_stop(timer_diags) ! diagnostics + + call ice_timer_start(timer_hist) ! history + call accum_hist (dt) ! history file + call ice_timer_stop(timer_hist) ! history + + call ice_timer_start(timer_readwrite) ! reading/writing + if (write_restart == 1) then + call dumpfile ! core variables for restarting + if (tr_iage) call write_restart_age + if (tr_FY) call write_restart_FY + if (tr_lvl) call write_restart_lvl + if (tr_pond_cesm) call write_restart_pond_cesm + if (tr_pond_lvl) call write_restart_pond_lvl + if (tr_pond_topo) call write_restart_pond_topo + if (tr_aero) call write_restart_aero + if (solve_zsal .or. skl_bgc .or. z_tracers) & + call write_restart_bgc + if (tr_brine) call write_restart_hbrine + if (kdyn == 2) call write_restart_eap + call final_restart + endif + + call ice_timer_stop(timer_readwrite) ! reading/writing + + end subroutine ice_step + +!======================================================================= +! +! Prepare for coupling +! +! authors: Elizabeth C. Hunke, LANL + + subroutine coupling_prep (iblk) + + use ice_arrays_column, only: alvdfn, alidfn, alvdrn, alidrn, & + albicen, albsnon, albpndn, apeffn, fzsal_g, fzsal, snowfracn + use ice_blocks, only: nx_block, ny_block, get_block, block + use ice_domain, only: blocks_ice + use ice_calendar, only: dt, nstreams + use ice_domain_size, only: ncat + use ice_flux, only: alvdf, alidf, alvdr, alidr, albice, albsno, & + albpnd, albcnt, apeff_ai, fpond, fresh, l_mpond_fresh, & + alvdf_ai, alidf_ai, alvdr_ai, alidr_ai, fhocn_ai, & + fresh_ai, fsalt_ai, fsalt, & + fswthru_ai, fhocn, fswthru, scale_factor, snowfrac, & + swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyT, & + fsens, flat, fswabs, flwout, evap, Tref, Qref, & + scale_fluxes, frzmlt_init, frzmlt, Uref, wind, fsurfn_f, flatn_f + use ice_flux_bgc, only: faero_ocn, fzsal_ai, fzsal_g_ai, flux_bio, flux_bio_ai, & + fnit, fsil, famm, fdmsp, fdms, fhum, fdust, falgalN, & + fdoc, fdic, fdon, ffep, ffed, bgcflux_ice_to_ocn + use ice_grid, only: tmask + use ice_state, only: aicen, aice, aice_init + use ice_step_mod, only: ocean_mixed_layer + use ice_timers, only: timer_couple, ice_timer_start, ice_timer_stop + + integer (kind=int_kind), intent(in) :: & + iblk ! block index + + ! local variables + + integer (kind=int_kind) :: & + ilo,ihi,jlo,jhi, & ! beginning and end of physical domain n , & ! thickness category index i,j , & ! horizontal indices - k ! tracer index + k , & ! tracer index + nbtrcr ! - real (kind=dbl_kind) :: cszn ! counter for history averaging + type (block) :: & + this_block ! block information for current block - real (kind=dbl_kind) :: netsw + logical (kind=log_kind) :: & + skl_bgc , & ! + calc_Tsfc ! - !----------------------------------------------------------------- - ! Save current value of frzmlt for diagnostics. - ! Update mixed layer with heat and radiation from ice. - !----------------------------------------------------------------- + real (kind=dbl_kind) :: & + cszn , & ! counter for history averaging + puny , & ! + rhofresh , & ! + netsw ! flag for shortwave radiation presence - do j = 1, ny_block - do i = 1, nx_block - frzmlt_init (i,j,iblk) = frzmlt(i,j,iblk) - enddo - enddo + character(len=*), parameter :: subname = '(coupling_prep)' - call ice_timer_start(timer_couple,iblk) ! atm/ocn coupling + !----------------------------------------------------------------- - if (oceanmixed_ice) & + call icepack_query_parameters(puny_out=puny, rhofresh_out=rhofresh) + call icepack_query_parameters(skl_bgc_out=skl_bgc) + call icepack_query_tracer_numbers(nbtrcr_out=nbtrcr) + call icepack_query_parameters(calc_Tsfc_out=calc_Tsfc) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + + !----------------------------------------------------------------- + ! Save current value of frzmlt for diagnostics. + ! Update mixed layer with heat and radiation from ice. + !----------------------------------------------------------------- + + do j = 1, ny_block + do i = 1, nx_block + frzmlt_init (i,j,iblk) = frzmlt(i,j,iblk) + enddo + enddo + + call ice_timer_start(timer_couple,iblk) ! atm/ocn coupling + + if (oceanmixed_ice) & call ocean_mixed_layer (dt,iblk) ! ocean surface fluxes and sst - !----------------------------------------------------------------- - ! Aggregate albedos - !----------------------------------------------------------------- - - do j = 1, ny_block - do i = 1, nx_block - alvdf(i,j,iblk) = c0 - alidf(i,j,iblk) = c0 - alvdr(i,j,iblk) = c0 - alidr(i,j,iblk) = c0 - - albice(i,j,iblk) = c0 - albsno(i,j,iblk) = c0 - albpnd(i,j,iblk) = c0 - apeff_ai(i,j,iblk) = c0 - snowfrac(i,j,iblk) = c0 - - ! for history averaging - cszn = c0 - netsw = swvdr(i,j,iblk)+swidr(i,j,iblk)+swvdf(i,j,iblk)+swidf(i,j,iblk) - if (netsw > puny) cszn = c1 - do n = 1, nstreams - albcnt(i,j,iblk,n) = albcnt(i,j,iblk,n) + cszn - enddo - enddo - enddo - do n = 1, ncat - do j = 1, ny_block - do i = 1, nx_block - alvdf(i,j,iblk) = alvdf(i,j,iblk) & - + alvdfn(i,j,n,iblk)*aicen(i,j,n,iblk) - alidf(i,j,iblk) = alidf(i,j,iblk) & - + alidfn(i,j,n,iblk)*aicen(i,j,n,iblk) - alvdr(i,j,iblk) = alvdr(i,j,iblk) & - + alvdrn(i,j,n,iblk)*aicen(i,j,n,iblk) - alidr(i,j,iblk) = alidr(i,j,iblk) & - + alidrn(i,j,n,iblk)*aicen(i,j,n,iblk) - - if (coszen(i,j,iblk) > puny) then ! sun above horizon - albice(i,j,iblk) = albice(i,j,iblk) & - + albicen(i,j,n,iblk)*aicen(i,j,n,iblk) - albsno(i,j,iblk) = albsno(i,j,iblk) & - + albsnon(i,j,n,iblk)*aicen(i,j,n,iblk) - albpnd(i,j,iblk) = albpnd(i,j,iblk) & - + albpndn(i,j,n,iblk)*aicen(i,j,n,iblk) - endif - - apeff_ai(i,j,iblk) = apeff_ai(i,j,iblk) & ! for history - + apeffn(i,j,n,iblk)*aicen(i,j,n,iblk) - snowfrac(i,j,iblk) = snowfrac(i,j,iblk) & ! for history - + snowfracn(i,j,n,iblk)*aicen(i,j,n,iblk) - enddo - enddo - enddo - - do j = 1, ny_block - do i = 1, nx_block - - !----------------------------------------------------------------- - ! reduce fresh by fpond for coupling - !----------------------------------------------------------------- - - if (l_mpond_fresh) then - fpond(i,j,iblk) = fpond(i,j,iblk) * rhofresh/dt - fresh(i,j,iblk) = fresh(i,j,iblk) - fpond(i,j,iblk) - endif - - !---------------------------------------------------------------- - ! Store grid box mean albedos and fluxes before scaling by aice - !---------------------------------------------------------------- - - alvdf_ai (i,j,iblk) = alvdf (i,j,iblk) - alidf_ai (i,j,iblk) = alidf (i,j,iblk) - alvdr_ai (i,j,iblk) = alvdr (i,j,iblk) - alidr_ai (i,j,iblk) = alidr (i,j,iblk) - fresh_ai (i,j,iblk) = fresh (i,j,iblk) - fsalt_ai (i,j,iblk) = fsalt (i,j,iblk) - fhocn_ai (i,j,iblk) = fhocn (i,j,iblk) - fswthru_ai(i,j,iblk) = fswthru(i,j,iblk) - - if (nbtrcr > 0) then - do k = 1, nbtrcr - flux_bio_ai (i,j,k,iblk) = flux_bio (i,j,k,iblk) - enddo - endif - - !----------------------------------------------------------------- - ! Save net shortwave for scaling factor in scale_factor - !----------------------------------------------------------------- - scale_factor(i,j,iblk) = & - swvdr(i,j,iblk)*(c1 - alvdr_ai(i,j,iblk)) & - + swvdf(i,j,iblk)*(c1 - alvdf_ai(i,j,iblk)) & - + swidr(i,j,iblk)*(c1 - alidr_ai(i,j,iblk)) & - + swidf(i,j,iblk)*(c1 - alidf_ai(i,j,iblk)) - - enddo - enddo - - !----------------------------------------------------------------- - ! Divide fluxes by ice area - ! - the CCSM coupler assumes fluxes are per unit ice area - ! - also needed for global budget in diagnostics - !----------------------------------------------------------------- - - call scale_fluxes (nx_block, ny_block, & - tmask (:,:,iblk) , nbtrcr, & - aice (:,:,iblk) , Tf (:,:,iblk), & - Tair (:,:,iblk) , Qa (:,:,iblk), & - strairxT (:,:,iblk) , strairyT(:,:,iblk), & - fsens (:,:,iblk) , flat (:,:,iblk), & - fswabs (:,:,iblk) , flwout (:,:,iblk), & - evap (:,:,iblk) , & - Tref (:,:,iblk) , Qref (:,:,iblk), & - fresh (:,:,iblk) , fsalt (:,:,iblk), & - fhocn (:,:,iblk) , fswthru (:,:,iblk), & - fswthruvdr(:,:,iblk), fswthruvdf(:,:,iblk), & - fswthruidr(:,:,iblk) , fswthruidf(:,:,iblk), & - faero_ocn (:,:,:,iblk), & - alvdr (:,:,iblk) , alidr (:,:,iblk), & - alvdf (:,:,iblk) , alidf (:,:,iblk), & - flux_bio (:,:,1:nbtrcr,iblk), & - Uref=Uref (:,:,iblk), wind=wind(:,:,iblk) ) - - !echmod - comment this out for efficiency, if .not. calc_Tsfc - if (.not. calc_Tsfc) then + !----------------------------------------------------------------- + ! Aggregate albedos + !----------------------------------------------------------------- + + do j = 1, ny_block + do i = 1, nx_block + alvdf(i,j,iblk) = c0 + alidf(i,j,iblk) = c0 + alvdr(i,j,iblk) = c0 + alidr(i,j,iblk) = c0 + + albice(i,j,iblk) = c0 + albsno(i,j,iblk) = c0 + albpnd(i,j,iblk) = c0 + apeff_ai(i,j,iblk) = c0 + snowfrac(i,j,iblk) = c0 + + ! for history averaging + cszn = c0 + netsw = swvdr(i,j,iblk)+swidr(i,j,iblk)+swvdf(i,j,iblk)+swidf(i,j,iblk) + if (netsw > puny) cszn = c1 + do n = 1, nstreams + albcnt(i,j,iblk,n) = albcnt(i,j,iblk,n) + cszn + enddo + enddo + enddo + + 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 n = 1, ncat + do j = jlo, jhi + do i = ilo, ihi + if (aicen(i,j,n,iblk) > puny) then + + alvdf(i,j,iblk) = alvdf(i,j,iblk) & + + alvdfn(i,j,n,iblk)*aicen(i,j,n,iblk) + alidf(i,j,iblk) = alidf(i,j,iblk) & + + alidfn(i,j,n,iblk)*aicen(i,j,n,iblk) + alvdr(i,j,iblk) = alvdr(i,j,iblk) & + + alvdrn(i,j,n,iblk)*aicen(i,j,n,iblk) + alidr(i,j,iblk) = alidr(i,j,iblk) & + + alidrn(i,j,n,iblk)*aicen(i,j,n,iblk) + + netsw = swvdr(i,j,iblk) + swidr(i,j,iblk) & + + swvdf(i,j,iblk) + swidf(i,j,iblk) + if (netsw > puny) then ! sun above horizon + albice(i,j,iblk) = albice(i,j,iblk) & + + albicen(i,j,n,iblk)*aicen(i,j,n,iblk) + albsno(i,j,iblk) = albsno(i,j,iblk) & + + albsnon(i,j,n,iblk)*aicen(i,j,n,iblk) + albpnd(i,j,iblk) = albpnd(i,j,iblk) & + + albpndn(i,j,n,iblk)*aicen(i,j,n,iblk) + endif + + apeff_ai(i,j,iblk) = apeff_ai(i,j,iblk) & ! for history + + apeffn(i,j,n,iblk)*aicen(i,j,n,iblk) + snowfrac(i,j,iblk) = snowfrac(i,j,iblk) & ! for history + + snowfracn(i,j,n,iblk)*aicen(i,j,n,iblk) + + endif ! aicen > puny + enddo + enddo + enddo + + do j = 1, ny_block + do i = 1, nx_block + + !----------------------------------------------------------------- + ! reduce fresh by fpond for coupling + !----------------------------------------------------------------- + + if (l_mpond_fresh) then + fpond(i,j,iblk) = fpond(i,j,iblk) * rhofresh/dt + fresh(i,j,iblk) = fresh(i,j,iblk) - fpond(i,j,iblk) + endif + + !---------------------------------------------------------------- + ! Store grid box mean albedos and fluxes before scaling by aice + !---------------------------------------------------------------- + + alvdf_ai (i,j,iblk) = alvdf (i,j,iblk) + alidf_ai (i,j,iblk) = alidf (i,j,iblk) + alvdr_ai (i,j,iblk) = alvdr (i,j,iblk) + alidr_ai (i,j,iblk) = alidr (i,j,iblk) + fresh_ai (i,j,iblk) = fresh (i,j,iblk) + fsalt_ai (i,j,iblk) = fsalt (i,j,iblk) + fhocn_ai (i,j,iblk) = fhocn (i,j,iblk) + fswthru_ai(i,j,iblk) = fswthru(i,j,iblk) + fzsal_ai (i,j,iblk) = fzsal (i,j,iblk) + fzsal_g_ai(i,j,iblk) = fzsal_g(i,j,iblk) + + if (nbtrcr > 0) then + do k = 1, nbtrcr + flux_bio_ai (i,j,k,iblk) = flux_bio (i,j,k,iblk) + enddo + endif + + !----------------------------------------------------------------- + ! Save net shortwave for scaling factor in scale_factor + !----------------------------------------------------------------- + scale_factor(i,j,iblk) = & + swvdr(i,j,iblk)*(c1 - alvdr_ai(i,j,iblk)) & + + swvdf(i,j,iblk)*(c1 - alvdf_ai(i,j,iblk)) & + + swidr(i,j,iblk)*(c1 - alidr_ai(i,j,iblk)) & + + swidf(i,j,iblk)*(c1 - alidf_ai(i,j,iblk)) + + enddo + enddo + + !----------------------------------------------------------------- + ! Divide fluxes by ice area + ! - the CESM coupler assumes fluxes are per unit ice area + ! - also needed for global budget in diagnostics + !----------------------------------------------------------------- + + call scale_fluxes (nx_block, ny_block, & + tmask (:,:,iblk), nbtrcr, icepack_max_aero, & + aice (:,:,iblk), Tf (:,:,iblk), & + Tair (:,:,iblk), Qa (:,:,iblk), & + strairxT (:,:,iblk), strairyT(:,:,iblk), & + fsens (:,:,iblk), flat (:,:,iblk), & + fswabs (:,:,iblk), flwout (:,:,iblk), & + evap (:,:,iblk), & + Tref (:,:,iblk), Qref (:,:,iblk), & + fresh (:,:,iblk), fsalt (:,:,iblk), & + fhocn (:,:,iblk), fswthru (:,:,iblk), & + faero_ocn(:,:,:,iblk), & + alvdr (:,:,iblk), alidr (:,:,iblk), & + alvdf (:,:,iblk), alidf (:,:,iblk), & + fzsal (:,:,iblk), fzsal_g (:,:,iblk), & + flux_bio(:,:,1:nbtrcr,iblk), & + Uref=Uref(:,:,iblk), wind=wind(:,:,iblk) ) + + !----------------------------------------------------------------- + ! Define ice-ocean bgc fluxes + !----------------------------------------------------------------- + + if (nbtrcr > 0 .or. skl_bgc) then + call bgcflux_ice_to_ocn (nx_block, ny_block, & + flux_bio(:,:,1:nbtrcr,iblk), & + fnit(:,:,iblk), fsil(:,:,iblk), & + famm(:,:,iblk), fdmsp(:,:,iblk), & + fdms(:,:,iblk), fhum(:,:,iblk), & + fdust(:,:,iblk), falgalN(:,:,:,iblk), & + fdoc(:,:,:,iblk), fdic(:,:,:,iblk), & + fdon(:,:,:,iblk), ffep(:,:,:,iblk), & + ffed(:,:,:,iblk)) + endif + +!echmod - comment this out for efficiency, if .not. calc_Tsfc + if (.not. calc_Tsfc) then !--------------------------------------------------------------- - ! If surface fluxes were provided, conserve these fluxes at ice - ! free points by passing to ocean. + ! If surface fluxes were provided, conserve these fluxes at ice + ! free points by passing to ocean. !--------------------------------------------------------------- - call sfcflux_to_ocn & - (nx_block, ny_block, & - tmask (:,:,iblk), aice_init(:,:,iblk), & - fsurfn_f (:,:,:,iblk), flatn_f(:,:,:,iblk), & - fresh (:,:,iblk), fhocn (:,:,iblk)) - endif - !echmod - - call ice_timer_stop(timer_couple,iblk) ! atm/ocn coupling - - end subroutine coupling_prep - - !======================================================================= - ! - ! If surface heat fluxes are provided to CICE instead of CICE calculating - ! them internally (i.e. .not. calc_Tsfc), then these heat fluxes can - ! be provided at points which do not have ice. (This is could be due to - ! the heat fluxes being calculated on a lower resolution grid or the - ! heat fluxes not recalculated at every CICE timestep.) At ice free points, - ! conserve energy and water by passing these fluxes to the ocean. - ! - ! author: A. McLaren, Met Office - - subroutine sfcflux_to_ocn(nx_block, ny_block, & - tmask, aice, & - fsurfn_f, flatn_f, & - fresh, fhocn) - - use ice_domain_size, only: ncat - - integer (kind=int_kind), intent(in) :: & - nx_block, ny_block ! block dimensions - - logical (kind=log_kind), dimension (nx_block,ny_block), & - intent(in) :: & - tmask ! land/boundary mask, thickness (T-cell) - - real (kind=dbl_kind), dimension(nx_block,ny_block), & - intent(in):: & - aice ! initial ice concentration - - real (kind=dbl_kind), dimension(nx_block,ny_block,ncat), & - intent(in) :: & - fsurfn_f, & ! net surface heat flux (provided as forcing) - flatn_f ! latent heat flux (provided as forcing) - - real (kind=dbl_kind), dimension(nx_block,ny_block), & - intent(inout):: & - fresh , & ! fresh water flux to ocean (kg/m2/s) - fhocn ! actual ocn/ice heat flx (W/m**2) - - ! local variables - integer (kind=int_kind) :: & - i, j, n ! horizontal indices - - real (kind=dbl_kind) :: & - rLsub ! 1/Lsub - !-------------------------------------------------- + call sfcflux_to_ocn & + (nx_block, ny_block, & + tmask (:,:,iblk), aice_init(:,:,iblk), & + fsurfn_f (:,:,:,iblk), flatn_f(:,:,:,iblk), & + fresh (:,:,iblk), fhocn (:,:,iblk)) + endif +!echmod + + call ice_timer_stop(timer_couple,iblk) ! atm/ocn coupling + + end subroutine coupling_prep + +!======================================================================= +! +! If surface heat fluxes are provided to CICE instead of CICE calculating +! them internally (i.e. .not. calc_Tsfc), then these heat fluxes can +! be provided at points which do not have ice. (This is could be due to +! the heat fluxes being calculated on a lower resolution grid or the +! heat fluxes not recalculated at every CICE timestep.) At ice free points, +! conserve energy and water by passing these fluxes to the ocean. +! +! author: A. McLaren, Met Office + + subroutine sfcflux_to_ocn(nx_block, ny_block, & + tmask, aice, & + fsurfn_f, flatn_f, & + fresh, fhocn) + + use ice_domain_size, only: ncat + + integer (kind=int_kind), intent(in) :: & + nx_block, ny_block ! block dimensions + + logical (kind=log_kind), dimension (nx_block,ny_block), intent(in) :: & + tmask ! land/boundary mask, thickness (T-cell) + + real (kind=dbl_kind), dimension(nx_block,ny_block), intent(in):: & + aice ! initial ice concentration + + real (kind=dbl_kind), dimension(nx_block,ny_block,ncat), intent(in) :: & + fsurfn_f, & ! net surface heat flux (provided as forcing) + flatn_f ! latent heat flux (provided as forcing) + + real (kind=dbl_kind), dimension(nx_block,ny_block), intent(inout):: & + fresh , & ! fresh water flux to ocean (kg/m2/s) + fhocn ! actual ocn/ice heat flx (W/m**2) #ifdef CICE_IN_NEMO - rLsub = c1 / Lsub - - do n = 1, ncat - do j = 1, ny_block - do i = 1, nx_block - if (tmask(i,j) .and. aice(i,j) <= puny) then - fhocn(i,j) = fhocn(i,j) + fsurfn_f(i,j,n) + flatn_f(i,j,n) - fresh(i,j) = fresh(i,j) + flatn_f(i,j,n) * rLsub - endif - enddo ! i - enddo ! j - enddo ! n -#endif + ! local variables + integer (kind=int_kind) :: & + i, j, n ! horizontal indices + + real (kind=dbl_kind) :: & + puny, & ! + rLsub ! 1/Lsub + + character(len=*), parameter :: subname = '(sfcflux_to_ocn)' + + call icepack_query_parameters(puny_out=puny) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + rLsub = c1 / Lsub + + do n = 1, ncat + do j = 1, ny_block + do i = 1, nx_block + if (tmask(i,j) .and. aice(i,j) <= puny) then + fhocn(i,j) = fhocn(i,j) & + + fsurfn_f(i,j,n) + flatn_f(i,j,n) + fresh(i,j) = fresh(i,j) & + + flatn_f(i,j,n) * rLsub + endif + enddo ! i + enddo ! j + enddo ! n + +#endif + + end subroutine sfcflux_to_ocn - end subroutine sfcflux_to_ocn +!======================================================================= -end module CICE_RunMod + end module CICE_RunMod !======================================================================= diff --git a/cicecore/drivers/nuopc/CICE_copyright.txt b/cicecore/drivers/nuopc/CICE_copyright.txt new file mode 100644 index 000000000..959a3ce32 --- /dev/null +++ b/cicecore/drivers/nuopc/CICE_copyright.txt @@ -0,0 +1,17 @@ +! Copyright (c) 2019, Triad National Security, LLC +! All rights reserved. +! +! Copyright 2019. Triad National Security, LLC. This software was +! produced under U.S. Government contract DE-AC52-06NA25396 for Los +! Alamos National Laboratory (LANL), which is operated by Triad +! National Security, LLC for the U.S. Department of Energy. The U.S. +! Government has rights to use, reproduce, and distribute this software. +! NEITHER THE GOVERNMENT NOR TRIAD NATIONAL SECURITY, LLC MAKES ANY +! WARRANTY, EXPRESS OR IMPLIED, OR ASSUMES ANY LIABILITY FOR THE USE OF +! THIS SOFTWARE. If software is modified to produce derivative works, +! such modified software should be clearly marked, so as not to confuse +! it with the version available from LANL. +! +! The full license and distribution policy are available from +! https://github.com/CICE-Consortium +! diff --git a/cicecore/drivers/nuopc/ice_comp_nuopc.F90 b/cicecore/drivers/nuopc/ice_comp_nuopc.F90 index 8ccbed03e..52dd67c55 100644 --- a/cicecore/drivers/nuopc/ice_comp_nuopc.F90 +++ b/cicecore/drivers/nuopc/ice_comp_nuopc.F90 @@ -20,8 +20,9 @@ module ice_comp_nuopc use shr_file_mod , only : shr_file_getlogunit, shr_file_setlogunit use shr_string_mod , only : shr_string_listGetNum use shr_orb_mod , only : shr_orb_decl - use shr_const_mod , only : shr_const_pi + use shr_const_mod use shr_cal_mod , only : shr_cal_noleap, shr_cal_gregorian, shr_cal_ymd2date + use ice_constants , only : ice_init_constants use ice_shr_methods , only : chkerr, state_setscalar, state_getscalar, state_diagnose, alarmInit use ice_shr_methods , only : set_component_logging, get_component_instance use ice_shr_methods , only : state_flddebug @@ -33,25 +34,29 @@ module ice_comp_nuopc use ice_blocks , only : nblocks_tot, get_block_parameter use ice_distribution , only : ice_distributiongetblockloc use ice_grid , only : tlon, tlat, hm, tarea, ULON, ULAT - use ice_constants , only : rad_to_deg use ice_communicate , only : my_task, master_task, mpi_comm_ice use ice_calendar , only : force_restart_now, write_ic use ice_calendar , only : idate, mday, time, month, daycal, time2sec, year_init use ice_calendar , only : sec, dt, calendar, calendar_type, nextsw_cday, istep - use ice_orbital , only : eccen, obliqr, lambm0, mvelpp - use ice_kinds_mod , only : dbl_kind, int_kind + use ice_kinds_mod , only : dbl_kind, int_kind, char_len use ice_scam , only : scmlat, scmlon, single_column use ice_fileunits , only : nu_diag, inst_index, inst_name, inst_suffix, release_all_fileunits - use ice_ocean , only : tfrz_option - use ice_therm_shared , only : ktherm use ice_restart_shared , only : runid, runtype, restart_dir, restart_file use ice_history , only : accum_hist +#if (defined NEWCODE) use ice_history_shared , only : model_doi_url ! TODO: add this functionality +#endif use ice_prescribed_mod , only : ice_prescribed_init +#if (defined NEWCODE) use ice_atmo , only : flux_convergence_tolerance, flux_convergence_max_iteration use ice_atmo , only : use_coldair_outbreak_mod +#endif use CICE_InitMod , only : CICE_Init use CICE_RunMod , only : CICE_Run + use ice_exit , only : abort_ice + use icepack_intfc , only : icepack_warnings_flush, icepack_warnings_aborted + use icepack_intfc , only : icepack_init_orbit, icepack_init_parameters + use icepack_intfc , only : icepack_query_tracer_flags, icepack_query_parameters use perf_mod , only : t_startf, t_stopf, t_barrierf use ice_timers @@ -253,7 +258,9 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) integer , allocatable :: gindex_elim(:) integer , allocatable :: gindex(:) integer :: globalID - character(CL) :: cvalue + character(ESMF_MAXSTR) :: cvalue + real(r8) :: eccen, obliqr, lambm0, mvelpp + character(len=char_len) :: tfrz_option character(ESMF_MAXSTR) :: convCIM, purpComp type(ESMF_VM) :: vm type(ESMF_Time) :: currTime ! Current time @@ -299,6 +306,8 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) integer :: num_elim_blocks ! local number of eliminated blocks integer :: num_total_blocks integer :: my_elim_start, my_elim_end + real(dbl_kind) :: rad_to_deg + integer(int_kind) :: ktherm character(*), parameter :: F00 = "('(ice_comp_nuopc) ',2a,1x,d21.14)" character(len=*), parameter :: subname=trim(modName)//':(InitializeRealize) ' !-------------------------------- @@ -331,6 +340,41 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) call t_startf ('cice_init_total') + !---------------------------------------------------------------------------- + ! Initialize constants + !---------------------------------------------------------------------------- + + call ice_init_constants(omega_in=SHR_CONST_OMEGA, radius_in=SHR_CONST_REARTH, & + spval_dbl_in=SHR_CONST_SPVAL) + call icepack_init_parameters( & + secday_in = SHR_CONST_CDAY, & + rhoi_in = SHR_CONST_RHOICE, & + rhow_in = SHR_CONST_RHOSW, & + cp_air_in = SHR_CONST_CPDAIR, & + cp_ice_in = SHR_CONST_CPICE, & + cp_ocn_in = SHR_CONST_CPSW, & + gravit_in = SHR_CONST_G, & + rhofresh_in = SHR_CONST_RHOFW, & + zvir_in = SHR_CONST_ZVIR, & + vonkar_in = SHR_CONST_KARMAN, & + cp_wv_in = SHR_CONST_CPWV, & + stefan_boltzmann_in = SHR_CONST_STEBOL, & + Tffresh_in= SHR_CONST_TKFRZ, & + Lsub_in = SHR_CONST_LATSUB, & + Lvap_in = SHR_CONST_LATVAP, & +! Lfresh_in = SHR_CONST_LATICE, & ! computed in init_parameters as Lsub-Lvap + Timelt_in = SHR_CONST_TKFRZ-SHR_CONST_TKFRZ, & + Tsmelt_in = SHR_CONST_TKFRZ-SHR_CONST_TKFRZ, & + ice_ref_salinity_in = SHR_CONST_ICE_REF_SAL, & + depressT_in = 0.054_dbl_kind, & + Tocnfrz_in= -34.0_dbl_kind*0.054_dbl_kind, & + pi_in = SHR_CONST_PI, & + snowpatch_in = 0.005_dbl_kind, & + dragio_in = 0.00962_dbl_kind) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + !---------------------------------------------------------------------------- ! Determine attributes - also needed in realize phase to get grid information !---------------------------------------------------------------------------- @@ -359,6 +403,12 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) read(cvalue,*) mvelpp end if + call icepack_init_orbit(eccen_in=eccen, mvelpp_in=mvelpp, & + lambm0_in=lambm0, obliqr_in=obliqr) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + ! Determine runtype and possibly nextsw_cday call NUOPC_CompAttributeGet(gcomp, name='start_type', value=cvalue, isPresent=isPresent, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return @@ -430,7 +480,12 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) if (.not. isPresent) then tfrz_option = 'linear_salt' ! TODO: is this right? This must be the same as mom is using for the calculation. end if + call icepack_init_parameters(tfrz_option_in=tfrz_option) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) +#if (defined NEWCODE) call NUOPC_CompAttributeGet(gcomp, name="flux_convergence", value=cvalue, isPresent=isPresent, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return if (isPresent) then @@ -454,6 +509,7 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) else use_coldair_outbreak_mod = .false. end if +#endif ! Get clock information before call to cice_init @@ -520,6 +576,11 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) ! reset shr logging to my log file !---------------------------------------------------------------------------- + call icepack_query_parameters(ktherm_out=ktherm) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + ! Now write output to nu_diag - this must happen AFTER call to cice_init if (localPet == 0) then write(nu_diag,F00) trim(subname),' cice init nextsw_cday = ',nextsw_cday @@ -530,8 +591,10 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) write(nu_diag,*) trim(subname),' inst_name = ',trim(inst_name) write(nu_diag,*) trim(subname),' inst_index = ',inst_index write(nu_diag,*) trim(subname),' inst_suffix = ',trim(inst_suffix) +#if (defined NEWCODE) write(nu_diag,*) trim(subname),' flux_convergence = ', flux_convergence_tolerance write(nu_diag,*) trim(subname),' flux_convergence_max_iteration = ', flux_convergence_max_iteration +#endif endif !--------------------------------------------------------------------------- @@ -760,6 +823,11 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) latMesh(n) = ownedElemCoords(2*n) end do + call icepack_query_parameters(rad_to_deg_out=rad_to_deg) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + ! obtain internally generated cice lats and lons for error checks allocate(lon(lsize)) allocate(lat(lsize)) @@ -895,6 +963,7 @@ subroutine ModelAdvance(gcomp, rc) type(ESMF_Time) :: nextTime type(ESMF_State) :: importState, exportState character(ESMF_MAXSTR) :: cvalue + real(r8) :: eccen, obliqr, lambm0, mvelpp integer :: shrlogunit ! original log unit integer :: k,n ! index logical :: stop_now ! .true. ==> stop at the end of this run phase @@ -912,7 +981,7 @@ subroutine ModelAdvance(gcomp, rc) character(CL) :: restart_filename logical :: isPresent character(*) , parameter :: F00 = "('(ice_comp_nuopc) ',2a,i8,d21.14)" - character(len=*),parameter :: subname=trim(modName)//':(ModelAdvance) ' + character(len=*),parameter :: subname=trim(modName)//':(ModelAdvance) ' !-------------------------------- rc = ESMF_SUCCESS @@ -976,6 +1045,12 @@ subroutine ModelAdvance(gcomp, rc) read(cvalue,*) mvelpp end if + call icepack_init_orbit(eccen_in=eccen, mvelpp_in=mvelpp, & + lambm0_in=lambm0, obliqr_in=obliqr) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + !-------------------------------- ! check that cice internal time is in sync with master clock before timestep update !-------------------------------- @@ -1054,11 +1129,11 @@ subroutine ModelAdvance(gcomp, rc) ! Advance cice and timestep update !-------------------------------- - if (force_restart_now) then - call CICE_Run(restart_filename=restart_filename) - else +!tcraig if (force_restart_now) then +! call CICE_Run(restart_filename=restart_filename) +! else call CICE_Run() - end if +! end if !-------------------------------- ! Create export state diff --git a/cicecore/drivers/nuopc/ice_constants.F90 b/cicecore/drivers/nuopc/ice_constants.F90 deleted file mode 100644 index fce26c1d6..000000000 --- a/cicecore/drivers/nuopc/ice_constants.F90 +++ /dev/null @@ -1,229 +0,0 @@ -!======================================================================= -!BOP -! -! !MODULE: ice_constants - sets physical constants -! -! !DESCRIPTION: -! -! This module defines a variety of physical and numerical constants -! used throughout the ice model \\ -! -! Code originally based on constants.F in POP -! -! !REVISION HISTORY: -! SVN:$Id: ice_constants.F90 37 2006-11-29 18:06:44Z eclare $ -! -! author Elizabeth C. Hunke, LANL -! -! !INTERFACE: - - module ice_constants -! -! !USES: -! - use shr_const_mod - use ice_kinds_mod -! -!EOP -! - implicit none - save - - public - - !----------------------------------------------------------------- - ! physical constants - !----------------------------------------------------------------- - - real (kind=dbl_kind), parameter :: & - rhos = 330.0_dbl_kind ,&! density of snow (kg/m^3) - rhoi = SHR_CONST_RHOICE ,&! density of ice (kg/m^3) - rhow = SHR_CONST_RHOSW ,&! density of seawater (kg/m^3) - cp_air = SHR_CONST_CPDAIR ,&! specific heat of air (J/kg/K) - ! (Briegleb JGR 97 11475-11485 July 1992) - emissivity = 0.95_dbl_kind ,&! emissivity of snow and ice - cp_ice = SHR_CONST_CPICE ,&! specific heat of fresh ice (J/kg/K) - cp_ocn = SHR_CONST_CPSW ,&! specific heat of ocn (J/kg/K) - depressT = 0.054_dbl_kind ,&! Tf:brine salinity ratio (C/ppt) - dragio = 0.00536_dbl_kind ,&! ice-ocn drag coefficient - albocn = 0.06_dbl_kind ! ocean albedo - - real (kind=dbl_kind), parameter :: & - gravit = SHR_CONST_G ,&! gravitational acceleration (m/s^2) - omega = SHR_CONST_OMEGA ,&! angular velocity of earth (rad/sec) - radius = SHR_CONST_REARTH ! earth radius (m) - - real (kind=dbl_kind), parameter :: & - secday = SHR_CONST_CDAY ,&! seconds in calendar day - viscosity_dyn = 1.79e-3_dbl_kind,&! dynamic viscosity of brine (kg/m/s) - Tocnfrz= -34.0_dbl_kind*depressT,&! freezing temp of seawater (C), - ! used as Tsfcn for open water - rhofresh = SHR_CONST_RHOFW ,&! density of fresh water (kg/m^3) - zvir = SHR_CONST_ZVIR ,&! rh2o/rair - 1.0 - vonkar = SHR_CONST_KARMAN,&! von Karman constant - cp_wv = SHR_CONST_CPWV ,&! specific heat of water vapor (J/kg/K) - stefan_boltzmann = SHR_CONST_STEBOL,&! W/m^2/K^4 - Tffresh = SHR_CONST_TKFRZ ,&! freezing temp of fresh ice (K) - Lsub = SHR_CONST_LATSUB,&! latent heat, sublimation freshwater (J/kg) - Lvap = SHR_CONST_LATVAP,&! latent heat, vaporization freshwater (J/kg) - Lfresh = SHR_CONST_LATICE,&! latent heat of melting of fresh ice (J/kg) - Timelt = SHR_CONST_TKFRZ-SHR_CONST_TKFRZ,&! melting temp. ice top surface (C) - Tsmelt = SHR_CONST_TKFRZ-SHR_CONST_TKFRZ,&! melting temp. snow top surface (C) - ice_ref_salinity = SHR_CONST_ICE_REF_SAL ,&! (psu) -! ocn_ref_salinity = SHR_CONST_OCN_REF_SAL ,&! (psu) -! rho_air = SHR_CONST_RHODAIR,&! ambient air density (kg/m^3) - spval_dbl = SHR_CONST_SPVAL ! special value - - real (kind=real_kind), parameter :: & - spval = 1.0e30_real_kind ! special value for netCDF output - - real (kind=dbl_kind), parameter :: & - iceruf = 0.0005_dbl_kind ,&! default ice surface roughness (m) - - ! (Ebert, Schramm and Curry JGR 100 15965-15975 Aug 1995) - kappav = 1.4_dbl_kind ,&! vis extnctn coef in ice, wvlngth<700nm (1/m) - !kappan = 17.6_dbl_kind,&! vis extnctn coef in ice, wvlngth<700nm (1/m) - - kice = 2.03_dbl_kind ,&! thermal conductivity of fresh ice(W/m/deg) - kseaice= 2.00_dbl_kind ,&! thermal conductivity of sea ice (W/m/deg) - ! (used in zero layer thermodynamics option) - ksno = 0.30_dbl_kind ,&! thermal conductivity of snow (W/m/deg) - zref = 10._dbl_kind ,&! reference height for stability (m) - ! hs0 = 0.03_dbl_kind, &! parameter for delta-Eddington snow frac - ! hsmin = 0.0001_dbl_kind, &! minimum snow thickness for dEdd - hs_min = 1.e-4_dbl_kind ,&! min snow thickness for computing zTsn (m) - snowpatch = 0.005_dbl_kind ! parameter for fractional snow area (m) - ! tcx note cice snowpatch = 0.02 - - ! weights for albedos - ! 4 Jan 2007 BPB Following are appropriate for complete cloud - ! in a summer polar atmosphere with 1.5m bare sea ice surface: - ! .636/.364 vis/nir with only 0.5% direct for each band. - real (kind=dbl_kind), parameter :: & ! currently used only - awtvdr = 0.00318_dbl_kind, &! visible, direct ! for history and - awtidr = 0.00182_dbl_kind, &! near IR, direct ! diagnostics - awtvdf = 0.63282_dbl_kind, &! visible, diffuse - awtidf = 0.36218_dbl_kind ! near IR, diffuse - - real (kind=dbl_kind), parameter :: & - qqqice = 11637800._dbl_kind ,&! for qsat over ice - TTTice = 5897.8_dbl_kind ,&! for qsat over ice - qqqocn = 627572.4_dbl_kind ,&! for qsat over ocn - TTTocn = 5107.4_dbl_kind ! for qsat over ocn - - ! these are currently set so as to have no effect on the decomposition - real (kind=dbl_kind), parameter :: & - shlat = 30.0_dbl_kind ,&! artificial masking edge (deg) - nhlat = -30.0_dbl_kind ! artificial masking edge (deg) - - !----------------------------------------------------------------- - ! numbers - !----------------------------------------------------------------- - - real (kind=dbl_kind), parameter :: & - c0 = 0.0_dbl_kind, & - c1 = 1.0_dbl_kind, & - c1p5 = 1.5_dbl_kind, & - c2 = 2.0_dbl_kind, & - c3 = 3.0_dbl_kind, & - c4 = 4.0_dbl_kind, & - c5 = 5.0_dbl_kind, & - c6 = 6.0_dbl_kind, & - c7 = 7.0_dbl_kind, & - c8 = 8.0_dbl_kind, & - c9 = 9.0_dbl_kind, & - c10 = 10.0_dbl_kind, & - c12 = 12.0_dbl_kind, & - c15 = 15.0_dbl_kind, & - c16 = 16.0_dbl_kind, & - c20 = 20.0_dbl_kind, & - c25 = 25.0_dbl_kind, & - c30 = 30.0_dbl_kind, & - c90 = 90.0_dbl_kind, & - c100 = 100.0_dbl_kind, & - c180 = 180.0_dbl_kind, & - c360 = 360.0_dbl_kind, & - c365 = 365.0_dbl_kind, & - c400 = 400.0_dbl_kind, & - c3600= 3600.0_dbl_kind, & - c1000= 1000.0_dbl_kind, & - p001 = 0.001_dbl_kind, & - p01 = 0.01_dbl_kind, & - p025 = 0.025_dbl_kind, & - p1 = 0.1_dbl_kind, & - p2 = 0.2_dbl_kind, & - p4 = 0.4_dbl_kind, & - p5 = 0.5_dbl_kind, & - p6 = 0.6_dbl_kind, & - p05 = 0.05_dbl_kind, & - p15 = 0.15_dbl_kind, & - p25 = 0.25_dbl_kind, & - p75 = 0.75_dbl_kind, & - p166 = c1/c6, & - p333 = c1/c3, & - p666 = c2/c3, & - p111 = c1/c9, & - p055 = p111*p5, & - p027 = p055*p5, & - p222 = c2/c9, & - eps04 = 1.0e-4_dbl_kind, & - eps11 = 1.0e-11_dbl_kind, & - eps12 = 1.0e-12_dbl_kind, & - eps13 = 1.0e-13_dbl_kind, & - eps15 = 1.0e-15_dbl_kind, & - eps16 = 1.0e-16_dbl_kind, & - puny = eps11, & - bignum = 1.0e+30_dbl_kind, & - pi = SHR_CONST_PI ,&! pi - pih = p5*pi, & - piq = p5*pih, & - pi2 = c2*pi - - !----------------------------------------------------------------- - ! location of fields for staggered grids - !----------------------------------------------------------------- - - integer (int_kind), parameter :: & - field_loc_unknown = 0, & - field_loc_noupdate = -1, & - field_loc_center = 1, & - field_loc_NEcorner = 2, & - field_loc_Nface = 3, & - field_loc_Eface = 4, & - field_loc_Wface = 5 - - - !----------------------------------------------------------------- - ! field type attribute - necessary for handling - ! changes of direction across tripole boundary - !----------------------------------------------------------------- - - integer (int_kind), parameter :: & - field_type_unknown = 0, & - field_type_noupdate = -1, & - field_type_scalar = 1, & - field_type_vector = 2, & - field_type_angle = 3 - - !----------------------------------------------------------------- - ! conversion factors - !----------------------------------------------------------------- - - real (kind=dbl_kind), parameter :: & - cm_to_m = 0.01_dbl_kind ,&! cm to meters - m_to_cm = 100._dbl_kind ,&! meters to cm - m2_to_km2 = 1.e-6_dbl_kind ,&! m^2 to km^2 - kg_to_g = 1000._dbl_kind ,&! kilograms to grams - mps_to_cmpdy = 8.64e6_dbl_kind ,&! m per s to cm per day - rad_to_deg = 180._dbl_kind/pi ! degree-radian conversion - - ! useful for debugging - integer (kind=int_kind), parameter :: & - mtest = -999, itest = 1, jtest = 1, ntest = 1, btest = 1 -! mtest = 2, itest = 50, jtest = 53, ntest = 1, btest = 1 - -!======================================================================= - - end module ice_constants - -!======================================================================= diff --git a/cicecore/drivers/nuopc/ice_import_export.F90 b/cicecore/drivers/nuopc/ice_import_export.F90 index e6ec87e97..0fe2510aa 100644 --- a/cicecore/drivers/nuopc/ice_import_export.F90 +++ b/cicecore/drivers/nuopc/ice_import_export.F90 @@ -7,23 +7,28 @@ module ice_import_export use shr_frz_mod , only : shr_frz_freezetemp use shr_kind_mod , only : r8 => shr_kind_r8, cl=>shr_kind_cl, cs=>shr_kind_cs use ice_kinds_mod , only : int_kind, dbl_kind, char_len_long, log_kind - use ice_constants , only : c0, c1, tffresh, spval_dbl + use ice_constants , only : c0, c1, spval_dbl use ice_constants , only : field_loc_center, field_type_scalar, field_type_vector use ice_blocks , only : block, get_block, nx_block, ny_block use ice_domain , only : nblocks, blocks_ice, halo_info, distrb_info use ice_domain_size , only : nx_global, ny_global, block_size_x, block_size_y, max_blocks, ncat + use ice_exit , only : abort_ice use ice_flux , only : strairxt, strairyt, strocnxt, strocnyt use ice_flux , only : alvdr, alidr, alvdf, alidf, Tref, Qref, Uref use ice_flux , only : flat, fsens, flwout, evap, fswabs, fhocn, fswthru +#if (defined NEWCODE) use ice_flux , only : fswthruvdr, fswthruvdf, fswthruidr, fswthruidf use ice_flux , only : send_i2x_per_cat, fswthrun_ai +#endif use ice_flux , only : fresh, fsalt, zlvl, uatm, vatm, potT, Tair, Qa use ice_flux , only : rhoa, swvdr, swvdf, swidr, swidf, flw, frain use ice_flux , only : fsnow, uocn, vocn, sst, ss_tltx, ss_tlty, frzmlt use ice_flux , only : sss, tf, wind, fsw +#if (defined NEWCODE) use ice_flux , only : faero_atm, faero_ocn use ice_flux , only : fiso_atm, fiso_ocn, fiso_rain, fiso_evap use ice_flux , only : Qa_iso, Qref_iso, HDO_ocn, H2_18O_ocn, H2_16O_ocn +#endif use ice_state , only : vice, vsno, aice, aicen_init, trcr use ice_grid , only : tlon, tlat, tarea, tmask, anglet, hm, ocn_gridcell_frac use ice_grid , only : grid_type, t2ugrid_vector @@ -32,6 +37,8 @@ module ice_import_export use ice_communicate , only : my_task, master_task, MPI_COMM_ICE use ice_prescribed_mod , only : prescribed_ice use ice_shr_methods , only : chkerr, state_reset + use icepack_intfc , only : icepack_warnings_flush, icepack_warnings_aborted + use icepack_intfc , only : icepack_query_parameters, icepack_query_tracer_flags use perf_mod , only : t_startf, t_stopf, t_barrierf implicit none @@ -115,10 +122,12 @@ subroutine ice_advertise_fields(gcomp, importState, exportState, flds_scalar_nam read(cvalue,*) flds_wiso call ESMF_LogWrite('flds_wiso = '// trim(cvalue), ESMF_LOGMSG_INFO) +#if (defined NEWCODE) call NUOPC_CompAttributeGet(gcomp, name='flds_i2o_per_cat', value=cvalue, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return read(cvalue,*) send_i2x_per_cat call ESMF_LogWrite('flds_i2o_per_cat = '// trim(cvalue), ESMF_LOGMSG_INFO) +#endif !----------------- ! advertise import fields @@ -189,10 +198,12 @@ subroutine ice_advertise_fields(gcomp, importState, exportState, flds_scalar_nam call fldlist_add(fldsFrIce_num, fldsFrIce, 'inst_ice_ir_dir_albedo' ) call fldlist_add(fldsFrIce_num, fldsFrIce, 'inst_ice_vis_dif_albedo' ) call fldlist_add(fldsFrIce_num, fldsFrIce, 'inst_ice_ir_dif_albedo' ) +#if (defined NEWCODE) if (send_i2x_per_cat) then call fldlist_add(fldsFrIce_num, fldsFrIce, 'ice_fraction_n', & ungridded_lbound=1, ungridded_ubound=ncat) end if +#endif ! ice/atm fluxes computed by ice call fldlist_add(fldsFrIce_num, fldsFrIce, 'stress_on_air_ice_zonal' ) @@ -210,10 +221,12 @@ subroutine ice_advertise_fields(gcomp, importState, exportState, flds_scalar_nam call fldlist_add(fldsFrIce_num, fldsFrIce, 'mean_sw_pen_to_ocn_vis_dif_flx' ) call fldlist_add(fldsFrIce_num, fldsFrIce, 'mean_sw_pen_to_ocn_ir_dir_flx' ) call fldlist_add(fldsFrIce_num, fldsFrIce, 'mean_sw_pen_to_ocn_ir_dif_flx' ) +#if (defined NEWCODE) if (send_i2x_per_cat) then call fldlist_add(fldsFrIce_num, fldsFrIce, 'mean_sw_pen_to_ocn_ifrac_n', & ungridded_lbound=1, ungridded_ubound=ncat) end if +#endif call fldlist_add(fldsFrIce_num , fldsFrIce, 'mean_fresh_water_to_ocean_rate' ) call fldlist_add(fldsFrIce_num , fldsFrIce, 'mean_salt_rate' ) call fldlist_add(fldsFrIce_num , fldsFrIce, 'stress_on_ocn_ice_zonal' ) @@ -332,9 +345,21 @@ subroutine ice_import( importState, rc ) real (kind=dbl_kind),allocatable :: aflds(:,:,:,:) real (kind=dbl_kind) :: workx, worky real (kind=dbl_kind) :: MIN_RAIN_TEMP, MAX_SNOW_TEMP + real (kind=dbl_kind) :: tffresh character(len=*), parameter :: subname = 'ice_import' !----------------------------------------------------- + call icepack_query_parameters(Tffresh_out=Tffresh) +! call icepack_query_parameters(tfrz_option_out=tfrz_option, & +! modal_aero_out=modal_aero, z_tracers_out=z_tracers, skl_bgc_out=skl_bgc, & +! Tffresh_out=Tffresh) +! call icepack_query_tracer_flags(tr_aero_out=tr_aero, tr_iage_out=tr_iage, & +! tr_FY_out=tr_FY, tr_pond_out=tr_pond, tr_lvl_out=tr_lvl, & +! tr_zaero_out=tr_zaero, tr_bgc_Nit_out=tr_bgc_Nit) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=u_FILE_u, line=__LINE__) + ! Note that the precipitation fluxes received from the mediator ! are in units of kg/s/m^2 which is what CICE requires. ! Note also that the read in below includes only values needed @@ -486,6 +511,7 @@ subroutine ice_import( importState, rc ) ! Get aerosols from mediator !------------------------------------------------------- +#if (defined NEWCODE) if (State_FldChk(importState, 'Faxa_bcph')) then ! the following indices are based on what the atmosphere is sending ! bcphidry ungridded_index=1 @@ -521,6 +547,7 @@ subroutine ice_import( importState, rc ) call state_getimport(importState, 'Faxa_dstdry', output=faero_atm, index=3, do_sum=.true., ungridded_index=4, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return end if +#endif !------------------------------------------------------- ! Water isotopes from the mediator @@ -530,6 +557,7 @@ subroutine ice_import( importState, rc ) ! 18O => ungridded_index=2 ! HDO => ungridded_index=3 +#if (defined NEWCODE) if (State_FldChk(importState, 'shum_wiso')) then call state_getimport(importState, 'inst_spec_humid_height_lowest_wiso', output=Qa_iso, index=1, ungridded_index=3, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return @@ -559,6 +587,7 @@ subroutine ice_import( importState, rc ) call state_getimport(importState, 'So_roce_wiso', output=H2_18O_ocn, ungridded_index=2, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return end if +#endif !----------------------------------------------------------------- ! rotate zonal/meridional vectors to local coordinates @@ -671,12 +700,24 @@ subroutine ice_export( exportState, rc ) real (kind=dbl_kind) :: tauyo (nx_block,ny_block,max_blocks) ! ice/ocean stress real (kind=dbl_kind) :: ailohi(nx_block,ny_block,max_blocks) ! fractional ice area real (kind=dbl_kind), allocatable :: tempfld(:,:,:) + real (kind=dbl_kind) :: tffresh character(len=*),parameter :: subname = 'ice_export' !----------------------------------------------------- rc = ESMF_SUCCESS if (dbug > 5) call ESMF_LogWrite(subname//' called', ESMF_LOGMSG_INFO) + call icepack_query_parameters(Tffresh_out=Tffresh) +! call icepack_query_parameters(tfrz_option_out=tfrz_option, & +! modal_aero_out=modal_aero, z_tracers_out=z_tracers, skl_bgc_out=skl_bgc, & +! Tffresh_out=Tffresh) +! call icepack_query_tracer_flags(tr_aero_out=tr_aero, tr_iage_out=tr_iage, & +! tr_FY_out=tr_FY, tr_pond_out=tr_pond, tr_lvl_out=tr_lvl, & +! tr_zaero_out=tr_zaero, tr_bgc_Nit_out=tr_bgc_Nit) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=u_FILE_u, line=__LINE__) + !calculate ice thickness from aice and vice. Also !create Tsrf from the first tracer (trcr) in ice_state.F @@ -890,6 +931,7 @@ subroutine ice_export( exportState, rc ) call state_setexport(exportState, 'mean_sw_pen_to_ocn' , input=fswthru, lmask=tmask, ifrac=ailohi, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return +#if (defined NEWCODE) ! flux of vis dir shortwave through ice to ocean call state_setexport(exportState, 'mean_sw_pen_to_ocn_vis_dir_flx' , input=fswthruvdr, lmask=tmask, ifrac=ailohi, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return @@ -905,6 +947,7 @@ subroutine ice_export( exportState, rc ) ! flux of ir dif shortwave through ice to ocean call state_setexport(exportState, 'mean_sw_pen_to_ocn_ir_dif_flx' , input=fswthruidf, lmask=tmask, ifrac=ailohi, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return +#endif ! heat exchange with ocean call state_setexport(exportState, 'net_heat_flx_to_ocn' , input=fhocn, lmask=tmask, ifrac=ailohi, rc=rc) @@ -926,6 +969,7 @@ subroutine ice_export( exportState, rc ) call state_setexport(exportState, 'stress_on_ocn_ice_merid' , input=tauyo, lmask=tmask, ifrac=ailohi, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return +#if (defined NEWCODE) ! ------ ! optional aerosol fluxes to ocean ! ------ @@ -995,6 +1039,7 @@ subroutine ice_export( exportState, rc ) lmask=tmask, ifrac=ailohi, ungridded_index=2, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return endif +#endif ! ------ ! optional short wave penetration to ocean ice category diff --git a/cicecore/drivers/nuopc/ice_prescribed_mod.F90 b/cicecore/drivers/nuopc/ice_prescribed_mod.F90 index 8183eb524..8d994e690 100644 --- a/cicecore/drivers/nuopc/ice_prescribed_mod.F90 +++ b/cicecore/drivers/nuopc/ice_prescribed_mod.F90 @@ -21,18 +21,22 @@ module ice_prescribed_mod use pio use ice_broadcast - use ice_communicate , only : my_task, master_task, MPI_COMM_ICE + use ice_communicate , only : my_task, master_task, MPI_COMM_ICE use ice_kinds_mod use ice_fileunits - use ice_exit , only : abort_ice - use ice_domain_size , only : nx_global, ny_global, ncat, nilyr, nslyr, max_blocks + use ice_exit , only : abort_ice + use ice_domain_size , only : nx_global, ny_global, ncat, nilyr, nslyr, max_blocks use ice_constants - use ice_blocks , only : nx_block, ny_block, block, get_block - use ice_domain , only : nblocks, distrb_info, blocks_ice - use ice_grid , only : TLAT, TLON, hm, tmask, tarea, grid_type, ocn_gridcell_frac - use ice_calendar , only : idate, sec, calendar_type - use ice_itd , only : hin_max + use ice_blocks , only : nx_block, ny_block, block, get_block + use ice_domain , only : nblocks, distrb_info, blocks_ice + use ice_grid , only : TLAT, TLON, hm, tmask, tarea, grid_type, ocn_gridcell_frac + use ice_calendar , only : idate, sec, calendar_type + use ice_arrays_column, only : hin_max use ice_read_write + use ice_exit, only: abort_ice + use icepack_intfc, only: icepack_warnings_flush, icepack_warnings_aborted + use icepack_intfc, only: icepack_query_tracer_indices, icepack_query_tracer_numbers + use icepack_intfc, only: icepack_query_parameters implicit none private ! except @@ -41,7 +45,6 @@ module ice_prescribed_mod public :: ice_prescribed_init ! initialize input data stream public :: ice_prescribed_run ! get time slices and time interp public :: ice_prescribed_phys ! set prescribed ice state and fluxes - private :: ice_set_domain ! !PUBLIC DATA MEMBERS: logical(kind=log_kind), public :: prescribed_ice ! true if prescribed ice @@ -65,19 +68,19 @@ module ice_prescribed_mod type(shr_strdata_type) :: sdat ! prescribed data stream character(len=char_len_long) :: fldList ! list of fields in data stream - real(kind=dbl_kind) :: ice_cov(nx_block,ny_block,max_blocks) ! ice cover - - real (kind=dbl_kind), parameter :: & - cp_sno = 0.0_dbl_kind & ! specific heat of snow (J/kg/K) - , rLfi = Lfresh*rhoi & ! latent heat of fusion ice (J/m^3) - , rLfs = Lfresh*rhos & ! latent heat of fusion snow (J/m^3) - , rLvi = Lvap*rhoi & ! latent heat of vapor*rhoice (J/m^3) - , rLvs = Lvap*rhos & ! latent heat of vapor*rhosno (J/m^3) - , rcpi = cp_ice*rhoi & ! heat capacity of fresh ice (J/m^3) - , rcps = cp_sno*rhos & ! heat capacity of snow (J/m^3) - , rcpidepressT = rcpi*depressT & ! param for finding T(z) from q (J/m^3) - , rLfidepressT = rLfi*depressT ! param for heat capacity (J deg/m^3) - ! heat capacity of sea ice, rhoi*C=rcpi+rLfidepressT*salinity/T^2 + real(kind=dbl_kind),allocatable :: ice_cov(:,:,:) ! ice cover + +! real (kind=dbl_kind), parameter :: & +! cp_sno = 0.0_dbl_kind & ! specific heat of snow (J/kg/K) +! , rLfi = Lfresh*rhoi & ! latent heat of fusion ice (J/m^3) +! , rLfs = Lfresh*rhos & ! latent heat of fusion snow (J/m^3) +! , rLvi = Lvap*rhoi & ! latent heat of vapor*rhoice (J/m^3) +! , rLvs = Lvap*rhos & ! latent heat of vapor*rhosno (J/m^3) +! , rcpi = cp_ice*rhoi & ! heat capacity of fresh ice (J/m^3) +! , rcps = cp_sno*rhos & ! heat capacity of snow (J/m^3) +! , rcpidepressT = rcpi*depressT & ! param for finding T(z) from q (J/m^3) +! , rLfidepressT = rLfi*depressT ! param for heat capacity (J deg/m^3) +! ! heat capacity of sea ice, rhoi*C=rcpi+rLfidepressT*salinity/T^2 !======================================================================= contains @@ -106,8 +109,7 @@ subroutine ice_prescribed_init(mpicom, compid, gindex) integer(kind=int_kind) :: nml_error ! namelist i/o error flag integer(kind=int_kind) :: n, nFile, ierr character(len=8) :: fillalgo - character(*),parameter :: subName = "('ice_prescribed_init2')" - character(*),parameter :: F00 = "('(ice_prescribed_init2) ',4a)" + character(*),parameter :: subName = '(ice_prescribed_init)' namelist /ice_prescribed_nml/ & prescribed_ice, & @@ -208,7 +210,7 @@ subroutine ice_prescribed_init(mpicom, compid, gindex) gsize = nx_global*ny_global lsize = size(gindex) call mct_gsMap_init( gsmap_ice, gindex, MPI_COMM_ICE, compid, lsize, gsize) - call ice_set_domain( lsize, MPI_COMM_ICE, gsmap_ice, dom_ice ) + call ice_prescribed_set_domain( lsize, MPI_COMM_ICE, gsmap_ice, dom_ice ) call shr_strdata_create(sdat,name="prescribed_ice", & mpicom=MPI_COMM_ICE, compid=compid, & @@ -264,8 +266,8 @@ subroutine ice_prescribed_run(mDateIn, secIn) type (block) :: this_block real(kind=dbl_kind) :: aice_max ! maximun ice concentration logical, save :: first_time = .true. - character(*),parameter :: subName = "('ice_prescribed_run')" - character(*),parameter :: F00 = "('(ice_prescribed_run) ',a,2g20.13)" + character(*),parameter :: subName = '(ice_prescribed_run)' + character(*),parameter :: F00 = "(a,2g20.13)" !------------------------------------------------------------------------ ! Interpolate to new ice coverage @@ -273,6 +275,10 @@ subroutine ice_prescribed_run(mDateIn, secIn) call shr_strdata_advance(sdat,mDateIn,SecIn,MPI_COMM_ICE,'cice_pice') + if (first_time) then + allocate(ice_cov(nx_block,ny_block,max_blocks)) + endif + ice_cov(:,:,:) = c0 ! This initializes ghost cells as well n=0 @@ -298,7 +304,7 @@ subroutine ice_prescribed_run(mDateIn, secIn) aice_max = maxval(ice_cov) if (aice_max > c10) then - write(nu_diag,F00) "ERROR: Ice conc data must be in fraction, aice_max= ",& + write(nu_diag,F00) subname//" ERROR: Ice conc data must be in fraction, aice_max= ",& aice_max call abort_ice(subName) end if @@ -336,7 +342,7 @@ subroutine ice_prescribed_phys ! !USES: use ice_flux use ice_state - use ice_itd, only : aggregate + use icepack_intfc, only : icepack_aggregate use ice_dyn_evp implicit none @@ -345,6 +351,7 @@ subroutine ice_prescribed_phys integer(kind=int_kind) :: nc ! ice category index integer(kind=int_kind) :: i,j,k ! longitude, latitude and level indices integer(kind=int_kind) :: iblk + integer(kind=int_kind) :: nt_Tsfc, nt_sice, nt_qice, nt_qsno, ntrcr real(kind=dbl_kind) :: slope ! diff in underlying ocean tmp and ice surface tmp real(kind=dbl_kind) :: Ti ! ice level temperature @@ -355,10 +362,23 @@ subroutine ice_prescribed_phys real(kind=dbl_kind) :: hs ! snow thickness real(kind=dbl_kind) :: zn ! normalized ice thickness real(kind=dbl_kind) :: salin(nilyr) ! salinity (ppt) + real(kind=dbl_kind) :: rad_to_deg, pi, puny + real(kind=dbl_kind) :: rhoi, rhos, cp_ice, cp_ocn, lfresh, depressT real(kind=dbl_kind), parameter :: nsal = 0.407_dbl_kind real(kind=dbl_kind), parameter :: msal = 0.573_dbl_kind real(kind=dbl_kind), parameter :: saltmax = 3.2_dbl_kind ! max salinity at ice base (ppm) + character(*),parameter :: subName = '(ice_prescribed_phys)' + + call icepack_query_tracer_indices(nt_Tsfc_out=nt_Tsfc, nt_sice_out=nt_sice, & + nt_qice_out=nt_qice, nt_qsno_out=nt_qsno) + call icepack_query_tracer_numbers(ntrcr_out=ntrcr) + call icepack_query_parameters(rad_to_deg_out=rad_to_deg, pi_out=pi, & + puny_out=puny, rhoi_out=rhoi, rhos_out=rhos, cp_ice_out=cp_ice, cp_ocn_out=cp_ocn, & + lfresh_out=lfresh, depressT_out=depressT) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) !----------------------------------------------------------------- ! Initialize ice state @@ -468,23 +488,27 @@ subroutine ice_prescribed_phys trcrn(i,j,nt_qice:nt_qice+nilyr-1,:,iblk) = c0 trcrn(i,j,nt_qsno:nt_qsno+nslyr-1,:,iblk) = c0 end if ! ice_cov >= eps04 + + !-------------------------------------------------------------------- + ! compute aggregate ice state and open water area + !-------------------------------------------------------------------- + call icepack_aggregate (ncat, & + aicen(i,j,:,iblk), & + trcrn(i,j,1:ntrcr,:,iblk), & + vicen(i,j,:,iblk), vsnon(i,j, :,iblk), & + aice (i,j, iblk), & + trcr (i,j,1:ntrcr, iblk), & + vice (i,j, iblk), vsno (i,j, iblk), & + aice0(i,j, iblk), & + ntrcr, & + trcr_depend(1:ntrcr), & + trcr_base(1:ntrcr,:), & + n_trcr_strata(1:ntrcr), & + nt_strata(1:ntrcr,:)) + end if ! tmask enddo ! i enddo ! j - - !-------------------------------------------------------------------- - ! compute aggregate ice state and open water area - !-------------------------------------------------------------------- - call aggregate (nx_block, ny_block, & - aicen(:,:,:,iblk), & - trcrn(:,:,1:ntrcr,:,iblk), & - vicen(:,:,:,iblk), vsnon(:,:, :,iblk), & - aice (:,:, iblk), & - trcr (:,:,1:ntrcr, iblk), & - vice (:,:, iblk), vsno (:,:, iblk), & - aice0(:,:, iblk), tmask(:,:, iblk), & - ntrcr, trcr_depend(1:ntrcr)) - enddo ! iblk do iblk = 1, nblocks @@ -515,7 +539,7 @@ end subroutine ice_prescribed_phys !=============================================================================== - subroutine ice_set_domain( lsize, mpicom, gsmap_i, dom_i ) + subroutine ice_prescribed_set_domain( lsize, mpicom, gsmap_i, dom_i ) ! Arguments integer , intent(in) :: lsize @@ -533,9 +557,16 @@ subroutine ice_set_domain( lsize, mpicom, gsmap_i, dom_i ) real(dbl_kind), pointer :: data5(:) ! temporary real(dbl_kind), pointer :: data6(:) ! temporary integer , pointer :: idata(:) ! temporary + real(kind=dbl_kind) :: rad_to_deg type(block) :: this_block ! block information for current block + character(*),parameter :: subName = '(ice_prescribed_set_domain)' !-------------------------------- + call icepack_query_parameters(rad_to_deg_out=rad_to_deg) + call icepack_warnings_flush(nu_diag) + if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & + file=__FILE__, line=__LINE__) + ! Initialize mct domain type call mct_gGrid_init(GGrid=dom_i, & CoordChars='lat:lon:hgt', OtherChars='area:aream:mask:frac', lsize=lsize ) @@ -604,7 +635,7 @@ subroutine ice_set_domain( lsize, mpicom, gsmap_i, dom_i ) deallocate(data1, data2, data3, data4, data5, data6) - end subroutine ice_set_domain + end subroutine ice_prescribed_set_domain #endif diff --git a/cicecore/drivers/nuopc/ice_scam.F90 b/cicecore/drivers/nuopc/ice_scam.F90 index 86a56d19c..f5280b259 100644 --- a/cicecore/drivers/nuopc/ice_scam.F90 +++ b/cicecore/drivers/nuopc/ice_scam.F90 @@ -12,4 +12,3 @@ module ice_scam end module ice_scam - diff --git a/cicecore/drivers/nuopc/ice_shr_methods.F90 b/cicecore/drivers/nuopc/ice_shr_methods.F90 index 8604e1dba..24a4226e5 100644 --- a/cicecore/drivers/nuopc/ice_shr_methods.F90 +++ b/cicecore/drivers/nuopc/ice_shr_methods.F90 @@ -86,6 +86,7 @@ subroutine memcheck(string, level, mastertask) ! local variables integer :: ierr integer, external :: GPTLprint_memusage + character(len=*), parameter :: subname='(memcheck)' !----------------------------------------------------------------------- if ((mastertask .and. memdebug_level > level) .or. memdebug_level > level+1) then @@ -107,6 +108,7 @@ subroutine get_component_instance(gcomp, inst_suffix, inst_index, rc) ! local variables logical :: isPresent character(len=4) :: cvalue + character(len=*), parameter :: subname='(get_component_instance)' !----------------------------------------------------------------------- rc = ESMF_SUCCESS @@ -140,6 +142,7 @@ subroutine set_component_logging(gcomp, mastertask, logunit, shrlogunit, rc) ! local variables character(len=CL) :: diro character(len=CL) :: logfile + character(len=*), parameter :: subname='(set_component_logging)' !----------------------------------------------------------------------- rc = ESMF_SUCCESS @@ -173,6 +176,7 @@ subroutine log_clock_advance(clock, component, logunit, rc) ! local variables character(len=CL) :: cvalue, prestring + character(len=*), parameter :: subname='(log_clock_advance)' !----------------------------------------------------------------------- rc = ESMF_SUCCESS @@ -371,6 +375,7 @@ subroutine state_flddebug(state, flds_scalar_name, prefix, ymd, tod, logunit, rc type(ESMF_Field) , allocatable :: lfields(:) integer , allocatable :: dimCounts(:) character(len=ESMF_MAXSTR) , allocatable :: fieldNameList(:) + character(len=*), parameter :: subname='(state_flddebug)' !----------------------------------------------------- ! Determine the list of fields and the dimension count for each field @@ -658,7 +663,7 @@ subroutine alarmInit( clock, alarm, option, & type(ESMF_Time) :: NextAlarm ! Next restart alarm time type(ESMF_TimeInterval) :: AlarmInterval ! Alarm interval integer :: sec - character(len=*), parameter :: subname = '(set_alarmInit): ' + character(len=*), parameter :: subname = '(alarmInit): ' !------------------------------------------------------------------------------- rc = ESMF_SUCCESS @@ -982,6 +987,7 @@ logical function chkerr(rc, line, file) character(len=*), intent(in) :: file integer :: lrc + character(len=*), parameter :: subname='(chkerr)' chkerr = .false. lrc = rc diff --git a/cicecore/shared/ice_fileunits.F90 b/cicecore/shared/ice_fileunits.F90 index 2620aa499..554f4c61a 100644 --- a/cicecore/shared/ice_fileunits.F90 +++ b/cicecore/shared/ice_fileunits.F90 @@ -99,7 +99,11 @@ subroutine init_fileunits character(len=*),parameter :: subname='(init_fileunits)' - nu_diag = ice_stdout ! default + ! Note - the nuopc cap sets nu_diag there - and so the following + ! should not be set for the nuopc cap + if (nu_diag == 0) then + nu_diag = ice_stdout ! default + end if allocate(ice_IOUnitsInUse(ice_IOUnitsMaxUnit)) ice_IOUnitsInUse = .false. From b54bca96dcc7fbc676ddb9981e3d30c7047c086e Mon Sep 17 00:00:00 2001 From: apcraig Date: Fri, 25 Oct 2019 16:38:04 -0600 Subject: [PATCH 2/4] clean up nu_diag initialization due to complexity added by nuopc cap --- cicecore/cicedynB/general/ice_init.F90 | 6 ++++++ cicecore/shared/ice_fileunits.F90 | 15 +++++++-------- 2 files changed, 13 insertions(+), 8 deletions(-) diff --git a/cicecore/cicedynB/general/ice_init.F90 b/cicecore/cicedynB/general/ice_init.F90 index 4beb2c3e1..053056847 100644 --- a/cicecore/cicedynB/general/ice_init.F90 +++ b/cicecore/cicedynB/general/ice_init.F90 @@ -468,6 +468,12 @@ subroutine input_data incond_file = trim(runid) // ".cice" // trim(inst_suffix) //".i" ! Note the nuopc cap will set nu_diag before this point - so just ! need to check that it is non-zero first + ! Note by tcraig - this is a terrible kludge and exists because with the nuopc + ! implementation this information is now defined in the cap, but the cice + ! model here doesn't care/know about the cap. In addition, the implementation + ! doesn't have the needed namelist and so shr_file_setIO aborts rather + ! than exit gracefully and return. An issue has been created in cime to + ! bring up this issue. In the long term, this if check should change. if (nu_diag == ice_stdout) then inquire(file='ice_modelio.nml'//trim(inst_suffix),exist=exists) if (exists) then diff --git a/cicecore/shared/ice_fileunits.F90 b/cicecore/shared/ice_fileunits.F90 index 554f4c61a..1174e3830 100644 --- a/cicecore/shared/ice_fileunits.F90 +++ b/cicecore/shared/ice_fileunits.F90 @@ -62,8 +62,7 @@ module ice_fileunits nu_restart_eap, & ! restart input file for eap dynamics nu_rst_pointer, & ! pointer to latest restart file nu_history , & ! binary history output file - nu_hdr , & ! header file for binary history output - nu_diag ! diagnostics output file + nu_hdr ! header file for binary history output character (32), public :: & nml_filename = 'ice_in' ! namelist input file name @@ -73,6 +72,9 @@ module ice_fileunits ice_stdout = 6, & ! reserved unit for standard output ice_stderr = 6 ! reserved unit for standard error + integer (kind=int_kind), public :: & + nu_diag = ice_stdout ! diagnostics output file, unit number may be overwritten + integer (kind=int_kind), public :: & ice_IOUnitsMinUnit = 11, & ! do not use unit numbers below ice_IOUnitsMaxUnit = 99 ! or above, set by setup_nml @@ -99,17 +101,14 @@ subroutine init_fileunits character(len=*),parameter :: subname='(init_fileunits)' - ! Note - the nuopc cap sets nu_diag there - and so the following - ! should not be set for the nuopc cap - if (nu_diag == 0) then - nu_diag = ice_stdout ! default - end if - allocate(ice_IOUnitsInUse(ice_IOUnitsMaxUnit)) ice_IOUnitsInUse = .false. + ice_IOUnitsInUse(ice_stdin) = .true. ! reserve unit 5 ice_IOUnitsInUse(ice_stdout) = .true. ! reserve unit 6 ice_IOUnitsInUse(ice_stderr) = .true. + if (nu_diag >= 1 .and. nu_diag <= ice_IOUnitsMaxUnit) & + ice_IOUnitsInUse(nu_diag) = .true. ! reserve unit nu_diag call get_fileunit(nu_grid) call get_fileunit(nu_kmt) From 2c4ffa0cd5db7ca4584f8f2df4ba4eb092e1ab3e Mon Sep 17 00:00:00 2001 From: apcraig Date: Fri, 25 Oct 2019 17:00:55 -0600 Subject: [PATCH 3/4] add nu_diag_set flag to work around the nu_diag setting in the nuopc cap --- cicecore/cicedynB/general/ice_init.F90 | 16 ++++++---------- cicecore/drivers/nuopc/ice_comp_nuopc.F90 | 6 +++--- cicecore/shared/ice_fileunits.F90 | 3 +++ 3 files changed, 12 insertions(+), 13 deletions(-) diff --git a/cicecore/cicedynB/general/ice_init.F90 b/cicecore/cicedynB/general/ice_init.F90 index 053056847..239a8c867 100644 --- a/cicecore/cicedynB/general/ice_init.F90 +++ b/cicecore/cicedynB/general/ice_init.F90 @@ -16,7 +16,7 @@ module ice_init use ice_communicate, only: my_task, master_task, ice_barrier use ice_constants, only: c0, c1, c2, c3, p2, p5 use ice_exit, only: abort_ice - use ice_fileunits, only: nu_nml, nu_diag, nml_filename, diag_type, & + use ice_fileunits, only: nu_nml, nu_diag, nu_diag_set, nml_filename, diag_type, & ice_stdout, get_fileunit, release_fileunit, bfbflag, flush_fileunit, & ice_IOUnitsMinUnit, ice_IOUnitsMaxUnit #ifdef CESMCOUPLED @@ -466,15 +466,11 @@ subroutine input_data history_file = trim(runid) // ".cice" // trim(inst_suffix) //".h" restart_file = trim(runid) // ".cice" // trim(inst_suffix) //".r" incond_file = trim(runid) // ".cice" // trim(inst_suffix) //".i" - ! Note the nuopc cap will set nu_diag before this point - so just - ! need to check that it is non-zero first - ! Note by tcraig - this is a terrible kludge and exists because with the nuopc - ! implementation this information is now defined in the cap, but the cice - ! model here doesn't care/know about the cap. In addition, the implementation - ! doesn't have the needed namelist and so shr_file_setIO aborts rather - ! than exit gracefully and return. An issue has been created in cime to - ! bring up this issue. In the long term, this if check should change. - if (nu_diag == ice_stdout) then + ! Note by tcraig - this if test is needed because the nuopc cap sets + ! nu_diag before this routine is called. This creates a conflict. + ! In addition, in the nuopc cap, shr_file_setIO will fail if the + ! needed namelist is missing (which it is in the CIME nuopc implementation) + if (.not. nu_diag_set) then inquire(file='ice_modelio.nml'//trim(inst_suffix),exist=exists) if (exists) then call get_fileUnit(nu_diag) diff --git a/cicecore/drivers/nuopc/ice_comp_nuopc.F90 b/cicecore/drivers/nuopc/ice_comp_nuopc.F90 index 52dd67c55..29cd34320 100644 --- a/cicecore/drivers/nuopc/ice_comp_nuopc.F90 +++ b/cicecore/drivers/nuopc/ice_comp_nuopc.F90 @@ -40,7 +40,7 @@ module ice_comp_nuopc use ice_calendar , only : sec, dt, calendar, calendar_type, nextsw_cday, istep use ice_kinds_mod , only : dbl_kind, int_kind, char_len use ice_scam , only : scmlat, scmlon, single_column - use ice_fileunits , only : nu_diag, inst_index, inst_name, inst_suffix, release_all_fileunits + use ice_fileunits , only : nu_diag, nu_diag_set, inst_index, inst_name, inst_suffix, release_all_fileunits use ice_restart_shared , only : runid, runtype, restart_dir, restart_file use ice_history , only : accum_hist #if (defined NEWCODE) @@ -553,11 +553,11 @@ subroutine InitializeRealize(gcomp, importState, exportState, clock, rc) ! Set cice logging !---------------------------------------------------------------------------- ! Note that sets the nu_diag module variable in ice_fileunits - ! nu_diag in this module is initialized to 0 in the module, and if this reset does not - ! happen here - then ice_init.F90 will obtain it from the input file ice_modelio.nml + ! Set the nu_diag_set flag so it's not reset later call set_component_logging(gcomp, my_task==master_task, nu_diag, shrlogunit, rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return + nu_diag_set = .true. call shr_file_setLogUnit (shrlogunit) diff --git a/cicecore/shared/ice_fileunits.F90 b/cicecore/shared/ice_fileunits.F90 index 1174e3830..6914c6b1d 100644 --- a/cicecore/shared/ice_fileunits.F90 +++ b/cicecore/shared/ice_fileunits.F90 @@ -75,6 +75,9 @@ module ice_fileunits integer (kind=int_kind), public :: & nu_diag = ice_stdout ! diagnostics output file, unit number may be overwritten + logical (kind=log_kind), public :: & + nu_diag_set = .false. ! flag to indicate whether nu_diag is already set + integer (kind=int_kind), public :: & ice_IOUnitsMinUnit = 11, & ! do not use unit numbers below ice_IOUnitsMaxUnit = 99 ! or above, set by setup_nml From 659aac8915d203c0917e537ed6cd92c920aa4e1f Mon Sep 17 00:00:00 2001 From: apcraig Date: Mon, 28 Oct 2019 17:59:40 -0600 Subject: [PATCH 4/4] update cice driver directories --- cicecore/drivers/{cesm => mct/cesm1}/CICE_FinalMod.F90 | 0 cicecore/drivers/{cesm => mct/cesm1}/CICE_InitMod.F90 | 0 cicecore/drivers/{cesm => mct/cesm1}/CICE_RunMod.F90 | 0 cicecore/drivers/{cesm => mct/cesm1}/CICE_RunMod.F90_debug | 0 cicecore/drivers/{cesm => mct/cesm1}/CICE_copyright.txt | 0 cicecore/drivers/{cesm => mct/cesm1}/ice_comp_esmf.F90 | 0 cicecore/drivers/{cesm => mct/cesm1}/ice_comp_mct.F90 | 0 cicecore/drivers/{cesm => mct/cesm1}/ice_cpl_indices.F90 | 0 cicecore/drivers/{cesm => mct/cesm1}/ice_import_export.F90 | 0 cicecore/drivers/{cesm => mct/cesm1}/ice_prescribed_mod.F90 | 0 cicecore/drivers/{cesm => mct/cesm1}/ice_scam.F90 | 0 cicecore/drivers/nuopc/{ => cmeps}/CICE_FinalMod.F90 | 0 cicecore/drivers/nuopc/{ => cmeps}/CICE_InitMod.F90 | 0 cicecore/drivers/nuopc/{ => cmeps}/CICE_RunMod.F90 | 0 cicecore/drivers/nuopc/{ => cmeps}/CICE_copyright.txt | 0 cicecore/drivers/nuopc/{ => cmeps}/ice_comp_nuopc.F90 | 0 cicecore/drivers/nuopc/{ => cmeps}/ice_import_export.F90 | 0 cicecore/drivers/nuopc/{ => cmeps}/ice_prescribed_mod.F90 | 0 cicecore/drivers/nuopc/{ => cmeps}/ice_scam.F90 | 0 cicecore/drivers/nuopc/{ => cmeps}/ice_shr_methods.F90 | 0 cicecore/drivers/{ => standalone}/cice/CICE.F90 | 0 cicecore/drivers/{ => standalone}/cice/CICE_FinalMod.F90 | 0 cicecore/drivers/{ => standalone}/cice/CICE_InitMod.F90 | 0 cicecore/drivers/{ => standalone}/cice/CICE_RunMod.F90 | 0 cicecore/drivers/{ => standalone}/cice/CICE_RunMod.F90_debug | 0 cicecore/drivers/{ => subroutine}/hadgem3/CICE.F90 | 0 cicecore/drivers/{ => subroutine}/hadgem3/CICE_FinalMod.F90 | 0 cicecore/drivers/{ => subroutine}/hadgem3/CICE_InitMod.F90 | 0 cicecore/drivers/{ => subroutine}/hadgem3/CICE_RunMod.F90 | 0 configuration/scripts/cice.settings | 2 +- 30 files changed, 1 insertion(+), 1 deletion(-) rename cicecore/drivers/{cesm => mct/cesm1}/CICE_FinalMod.F90 (100%) rename cicecore/drivers/{cesm => mct/cesm1}/CICE_InitMod.F90 (100%) rename cicecore/drivers/{cesm => mct/cesm1}/CICE_RunMod.F90 (100%) rename cicecore/drivers/{cesm => mct/cesm1}/CICE_RunMod.F90_debug (100%) rename cicecore/drivers/{cesm => mct/cesm1}/CICE_copyright.txt (100%) rename cicecore/drivers/{cesm => mct/cesm1}/ice_comp_esmf.F90 (100%) rename cicecore/drivers/{cesm => mct/cesm1}/ice_comp_mct.F90 (100%) rename cicecore/drivers/{cesm => mct/cesm1}/ice_cpl_indices.F90 (100%) rename cicecore/drivers/{cesm => mct/cesm1}/ice_import_export.F90 (100%) rename cicecore/drivers/{cesm => mct/cesm1}/ice_prescribed_mod.F90 (100%) rename cicecore/drivers/{cesm => mct/cesm1}/ice_scam.F90 (100%) rename cicecore/drivers/nuopc/{ => cmeps}/CICE_FinalMod.F90 (100%) rename cicecore/drivers/nuopc/{ => cmeps}/CICE_InitMod.F90 (100%) rename cicecore/drivers/nuopc/{ => cmeps}/CICE_RunMod.F90 (100%) rename cicecore/drivers/nuopc/{ => cmeps}/CICE_copyright.txt (100%) rename cicecore/drivers/nuopc/{ => cmeps}/ice_comp_nuopc.F90 (100%) rename cicecore/drivers/nuopc/{ => cmeps}/ice_import_export.F90 (100%) rename cicecore/drivers/nuopc/{ => cmeps}/ice_prescribed_mod.F90 (100%) rename cicecore/drivers/nuopc/{ => cmeps}/ice_scam.F90 (100%) rename cicecore/drivers/nuopc/{ => cmeps}/ice_shr_methods.F90 (100%) rename cicecore/drivers/{ => standalone}/cice/CICE.F90 (100%) rename cicecore/drivers/{ => standalone}/cice/CICE_FinalMod.F90 (100%) rename cicecore/drivers/{ => standalone}/cice/CICE_InitMod.F90 (100%) rename cicecore/drivers/{ => standalone}/cice/CICE_RunMod.F90 (100%) rename cicecore/drivers/{ => standalone}/cice/CICE_RunMod.F90_debug (100%) rename cicecore/drivers/{ => subroutine}/hadgem3/CICE.F90 (100%) rename cicecore/drivers/{ => subroutine}/hadgem3/CICE_FinalMod.F90 (100%) rename cicecore/drivers/{ => subroutine}/hadgem3/CICE_InitMod.F90 (100%) rename cicecore/drivers/{ => subroutine}/hadgem3/CICE_RunMod.F90 (100%) diff --git a/cicecore/drivers/cesm/CICE_FinalMod.F90 b/cicecore/drivers/mct/cesm1/CICE_FinalMod.F90 similarity index 100% rename from cicecore/drivers/cesm/CICE_FinalMod.F90 rename to cicecore/drivers/mct/cesm1/CICE_FinalMod.F90 diff --git a/cicecore/drivers/cesm/CICE_InitMod.F90 b/cicecore/drivers/mct/cesm1/CICE_InitMod.F90 similarity index 100% rename from cicecore/drivers/cesm/CICE_InitMod.F90 rename to cicecore/drivers/mct/cesm1/CICE_InitMod.F90 diff --git a/cicecore/drivers/cesm/CICE_RunMod.F90 b/cicecore/drivers/mct/cesm1/CICE_RunMod.F90 similarity index 100% rename from cicecore/drivers/cesm/CICE_RunMod.F90 rename to cicecore/drivers/mct/cesm1/CICE_RunMod.F90 diff --git a/cicecore/drivers/cesm/CICE_RunMod.F90_debug b/cicecore/drivers/mct/cesm1/CICE_RunMod.F90_debug similarity index 100% rename from cicecore/drivers/cesm/CICE_RunMod.F90_debug rename to cicecore/drivers/mct/cesm1/CICE_RunMod.F90_debug diff --git a/cicecore/drivers/cesm/CICE_copyright.txt b/cicecore/drivers/mct/cesm1/CICE_copyright.txt similarity index 100% rename from cicecore/drivers/cesm/CICE_copyright.txt rename to cicecore/drivers/mct/cesm1/CICE_copyright.txt diff --git a/cicecore/drivers/cesm/ice_comp_esmf.F90 b/cicecore/drivers/mct/cesm1/ice_comp_esmf.F90 similarity index 100% rename from cicecore/drivers/cesm/ice_comp_esmf.F90 rename to cicecore/drivers/mct/cesm1/ice_comp_esmf.F90 diff --git a/cicecore/drivers/cesm/ice_comp_mct.F90 b/cicecore/drivers/mct/cesm1/ice_comp_mct.F90 similarity index 100% rename from cicecore/drivers/cesm/ice_comp_mct.F90 rename to cicecore/drivers/mct/cesm1/ice_comp_mct.F90 diff --git a/cicecore/drivers/cesm/ice_cpl_indices.F90 b/cicecore/drivers/mct/cesm1/ice_cpl_indices.F90 similarity index 100% rename from cicecore/drivers/cesm/ice_cpl_indices.F90 rename to cicecore/drivers/mct/cesm1/ice_cpl_indices.F90 diff --git a/cicecore/drivers/cesm/ice_import_export.F90 b/cicecore/drivers/mct/cesm1/ice_import_export.F90 similarity index 100% rename from cicecore/drivers/cesm/ice_import_export.F90 rename to cicecore/drivers/mct/cesm1/ice_import_export.F90 diff --git a/cicecore/drivers/cesm/ice_prescribed_mod.F90 b/cicecore/drivers/mct/cesm1/ice_prescribed_mod.F90 similarity index 100% rename from cicecore/drivers/cesm/ice_prescribed_mod.F90 rename to cicecore/drivers/mct/cesm1/ice_prescribed_mod.F90 diff --git a/cicecore/drivers/cesm/ice_scam.F90 b/cicecore/drivers/mct/cesm1/ice_scam.F90 similarity index 100% rename from cicecore/drivers/cesm/ice_scam.F90 rename to cicecore/drivers/mct/cesm1/ice_scam.F90 diff --git a/cicecore/drivers/nuopc/CICE_FinalMod.F90 b/cicecore/drivers/nuopc/cmeps/CICE_FinalMod.F90 similarity index 100% rename from cicecore/drivers/nuopc/CICE_FinalMod.F90 rename to cicecore/drivers/nuopc/cmeps/CICE_FinalMod.F90 diff --git a/cicecore/drivers/nuopc/CICE_InitMod.F90 b/cicecore/drivers/nuopc/cmeps/CICE_InitMod.F90 similarity index 100% rename from cicecore/drivers/nuopc/CICE_InitMod.F90 rename to cicecore/drivers/nuopc/cmeps/CICE_InitMod.F90 diff --git a/cicecore/drivers/nuopc/CICE_RunMod.F90 b/cicecore/drivers/nuopc/cmeps/CICE_RunMod.F90 similarity index 100% rename from cicecore/drivers/nuopc/CICE_RunMod.F90 rename to cicecore/drivers/nuopc/cmeps/CICE_RunMod.F90 diff --git a/cicecore/drivers/nuopc/CICE_copyright.txt b/cicecore/drivers/nuopc/cmeps/CICE_copyright.txt similarity index 100% rename from cicecore/drivers/nuopc/CICE_copyright.txt rename to cicecore/drivers/nuopc/cmeps/CICE_copyright.txt diff --git a/cicecore/drivers/nuopc/ice_comp_nuopc.F90 b/cicecore/drivers/nuopc/cmeps/ice_comp_nuopc.F90 similarity index 100% rename from cicecore/drivers/nuopc/ice_comp_nuopc.F90 rename to cicecore/drivers/nuopc/cmeps/ice_comp_nuopc.F90 diff --git a/cicecore/drivers/nuopc/ice_import_export.F90 b/cicecore/drivers/nuopc/cmeps/ice_import_export.F90 similarity index 100% rename from cicecore/drivers/nuopc/ice_import_export.F90 rename to cicecore/drivers/nuopc/cmeps/ice_import_export.F90 diff --git a/cicecore/drivers/nuopc/ice_prescribed_mod.F90 b/cicecore/drivers/nuopc/cmeps/ice_prescribed_mod.F90 similarity index 100% rename from cicecore/drivers/nuopc/ice_prescribed_mod.F90 rename to cicecore/drivers/nuopc/cmeps/ice_prescribed_mod.F90 diff --git a/cicecore/drivers/nuopc/ice_scam.F90 b/cicecore/drivers/nuopc/cmeps/ice_scam.F90 similarity index 100% rename from cicecore/drivers/nuopc/ice_scam.F90 rename to cicecore/drivers/nuopc/cmeps/ice_scam.F90 diff --git a/cicecore/drivers/nuopc/ice_shr_methods.F90 b/cicecore/drivers/nuopc/cmeps/ice_shr_methods.F90 similarity index 100% rename from cicecore/drivers/nuopc/ice_shr_methods.F90 rename to cicecore/drivers/nuopc/cmeps/ice_shr_methods.F90 diff --git a/cicecore/drivers/cice/CICE.F90 b/cicecore/drivers/standalone/cice/CICE.F90 similarity index 100% rename from cicecore/drivers/cice/CICE.F90 rename to cicecore/drivers/standalone/cice/CICE.F90 diff --git a/cicecore/drivers/cice/CICE_FinalMod.F90 b/cicecore/drivers/standalone/cice/CICE_FinalMod.F90 similarity index 100% rename from cicecore/drivers/cice/CICE_FinalMod.F90 rename to cicecore/drivers/standalone/cice/CICE_FinalMod.F90 diff --git a/cicecore/drivers/cice/CICE_InitMod.F90 b/cicecore/drivers/standalone/cice/CICE_InitMod.F90 similarity index 100% rename from cicecore/drivers/cice/CICE_InitMod.F90 rename to cicecore/drivers/standalone/cice/CICE_InitMod.F90 diff --git a/cicecore/drivers/cice/CICE_RunMod.F90 b/cicecore/drivers/standalone/cice/CICE_RunMod.F90 similarity index 100% rename from cicecore/drivers/cice/CICE_RunMod.F90 rename to cicecore/drivers/standalone/cice/CICE_RunMod.F90 diff --git a/cicecore/drivers/cice/CICE_RunMod.F90_debug b/cicecore/drivers/standalone/cice/CICE_RunMod.F90_debug similarity index 100% rename from cicecore/drivers/cice/CICE_RunMod.F90_debug rename to cicecore/drivers/standalone/cice/CICE_RunMod.F90_debug diff --git a/cicecore/drivers/hadgem3/CICE.F90 b/cicecore/drivers/subroutine/hadgem3/CICE.F90 similarity index 100% rename from cicecore/drivers/hadgem3/CICE.F90 rename to cicecore/drivers/subroutine/hadgem3/CICE.F90 diff --git a/cicecore/drivers/hadgem3/CICE_FinalMod.F90 b/cicecore/drivers/subroutine/hadgem3/CICE_FinalMod.F90 similarity index 100% rename from cicecore/drivers/hadgem3/CICE_FinalMod.F90 rename to cicecore/drivers/subroutine/hadgem3/CICE_FinalMod.F90 diff --git a/cicecore/drivers/hadgem3/CICE_InitMod.F90 b/cicecore/drivers/subroutine/hadgem3/CICE_InitMod.F90 similarity index 100% rename from cicecore/drivers/hadgem3/CICE_InitMod.F90 rename to cicecore/drivers/subroutine/hadgem3/CICE_InitMod.F90 diff --git a/cicecore/drivers/hadgem3/CICE_RunMod.F90 b/cicecore/drivers/subroutine/hadgem3/CICE_RunMod.F90 similarity index 100% rename from cicecore/drivers/hadgem3/CICE_RunMod.F90 rename to cicecore/drivers/subroutine/hadgem3/CICE_RunMod.F90 diff --git a/configuration/scripts/cice.settings b/configuration/scripts/cice.settings index 0513f1dd1..8bb860916 100755 --- a/configuration/scripts/cice.settings +++ b/configuration/scripts/cice.settings @@ -12,7 +12,7 @@ setenv ICE_OBJDIR ${ICE_RUNDIR}/compile setenv ICE_RSTDIR ${ICE_RUNDIR}/restart setenv ICE_HSTDIR ${ICE_RUNDIR}/history setenv ICE_LOGDIR ${ICE_CASEDIR}/logs -setenv ICE_DRVOPT cice +setenv ICE_DRVOPT standalone/cice setenv ICE_IOTYPE netcdf # set to none if netcdf library is unavailable setenv ICE_CLEANBUILD true setenv ICE_QUIETMODE false