! SVN:$Id: ice_history.F90 1304 2015-11-12 19:28:14Z afrobert@nps.edu $ !======================================================================= ! Driver for core history output ! ! The following variables are currently hard-wired as snapshots ! (instantaneous rather than time-averages): ! divu, shear, sig1, sig2, trsig, mlt_onset, frz_onset, hisnap, aisnap ! ! Options for histfreq: '1','h','d','m','y','x', where x means that ! output stream will not be used (recommended for efficiency). ! histfreq_n can be any nonnegative integer, where 0 means that the ! corresponding histfreq frequency will not be used. ! The flags (f_) can be set to '1','h','d','m','y' or 'x', where ! n means the field will not be written. To output the same field at ! more than one frequency, for instance monthy and daily, set ! f_ = 'md'. ! ! 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 CCSM 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 use ice_kinds_mod implicit none private public :: init_hist, accum_hist save !======================================================================= contains !======================================================================= ! Initialize history files ! ! authors Tony Craig, NCAR ! Elizabeth C. Hunke, LANL ! C.M. Bitz, UW ! Bruce P. Briegleb, NCAR ! William H. Lipscomb, LANL subroutine init_hist (dt) use ice_atmo, only: formdrag use ice_blocks, only: nx_block, ny_block use ice_broadcast, only: broadcast_scalar, broadcast_array use ice_communicate, only: my_task, master_task use ice_constants, only: c0, c1, c2, c100, mps_to_cmpdy, rhofresh, & Tffresh, kg_to_g, secday use ice_calendar, only: yday, days_per_year, histfreq, & histfreq_n, nstreams use ice_domain_size, only: max_blocks, max_nstrm use ice_dyn_shared, only: kdyn use ice_exit, only: abort_ice use ice_fileunits, only: nu_nml, nml_filename, nu_diag, & get_fileunit, release_fileunit use ice_flux, only: mlt_onset, frz_onset, albcnt use ice_history_shared ! everything use ice_history_mechred, only: init_hist_mechred_2D, init_hist_mechred_3Dc use ice_history_pond, only: init_hist_pond_2D, init_hist_pond_3Dc use ice_history_bgc, only: init_hist_bgc_2D, init_hist_bgc_3Dc, & init_hist_bgc_4Db use ice_history_drag, only: init_hist_drag_2D use ice_restart_shared, only: restart use ice_state, only: tr_iage, tr_FY, tr_lvl, tr_pond, tr_aero, tr_brine use ice_zbgc_shared, only: skl_bgc real (kind=dbl_kind), intent(in) :: & dt ! time step ! local variables integer (kind=int_kind) :: n, ns, ns1, ns2 integer (kind=int_kind), dimension(max_nstrm) :: & ntmp integer (kind=int_kind) :: nml_error ! namelist i/o error flag !----------------------------------------------------------------- ! read namelist !----------------------------------------------------------------- call get_fileunit(nu_nml) if (my_task == master_task) then open (nu_nml, file=nml_filename, status='old',iostat=nml_error) if (nml_error /= 0) then nml_error = -1 else nml_error = 1 endif do while (nml_error > 0) read(nu_nml, nml=icefields_nml,iostat=nml_error) end do if (nml_error == 0) close(nu_nml) endif call release_fileunit(nu_nml) call broadcast_scalar(nml_error, master_task) if (nml_error /= 0) then close (nu_nml) call abort_ice('ice: error reading icefields_nml') endif ! histfreq options ('1','h','d','m','y') nstreams = 0 do ns = 1, max_nstrm if (histfreq(ns) == '1' .or. histfreq(ns) == 'h' .or. & histfreq(ns) == 'd' .or. histfreq(ns) == 'm' .or. & histfreq(ns) == 'y') then nstreams = nstreams + 1 if (ns >= 2) then if (histfreq(ns-1) == 'x') then call abort_ice('ice: histfreq all non x must be at start of array') endif endif else if (histfreq(ns) /= 'x') then call abort_ice('ice: histfreq contains illegal element') endif enddo if (nstreams == 0) write (nu_diag,*) 'WARNING: No history output' do ns1 = 1, nstreams do ns2 = 1, nstreams if (histfreq(ns1) == histfreq(ns2) .and. ns1/=ns2 & .and. my_task == master_task) then call abort_ice('ice: histfreq elements must be unique') endif enddo enddo if (.not. tr_iage) then f_iage = 'x' f_dagedtt = 'x' f_dagedtd = 'x' endif if (.not. tr_FY) f_FY = 'x' if (kdyn /= 2) then f_a11 = 'x' f_a12 = 'x' f_e11 = 'x' f_e12 = 'x' f_e22 = 'x' f_s11 = 'x' f_s12 = 'x' f_s22 = 'x' f_yieldstress11 = 'x' f_yieldstress12 = 'x' f_yieldstress22 = 'x' endif ! these must be output at the same frequency because of ! cos(zenith angle) averaging if (f_albice(1:1) /= 'x' .and. f_albsni(1:1) /= 'x') f_albice = f_albsni if (f_albsno(1:1) /= 'x') f_albsno = f_albice if (f_albpnd(1:1) /= 'x') f_albpnd = f_albice if (f_coszen(1:1) /= 'x' .and. f_albice(1:1) /= 'x') f_coszen = f_albice if (f_coszen(1:1) /= 'x' .and. f_albsni(1:1) /= 'x') f_coszen = f_albsni ! to prevent array-out-of-bounds when aggregating if (f_fmeltt_ai(1:1) /= 'x') f_fmelttn_ai = f_fmeltt_ai #ifndef ncdf f_bounds = .false. #endif ! write dimensions for 3D or 4D history variables ! note: list of variables checked here is incomplete if (f_aicen(1:1) /= 'x' .or. f_vicen(1:1) /= 'x' .or. & f_Tinz (1:1) /= 'x' .or. f_Sinz (1:1) /= 'x') f_NCAT = .true. if (f_Tinz (1:1) /= 'x' .or. f_Sinz (1:1) /= 'x') f_VGRDi = .true. if (f_Tsnz (1:1) /= 'x') f_VGRDs = .true. call broadcast_scalar (f_tmask, master_task) call broadcast_scalar (f_blkmask, master_task) call broadcast_scalar (f_tarea, master_task) call broadcast_scalar (f_uarea, master_task) call broadcast_scalar (f_dxt, master_task) call broadcast_scalar (f_dyt, master_task) call broadcast_scalar (f_dxu, master_task) call broadcast_scalar (f_dyu, master_task) call broadcast_scalar (f_HTN, master_task) call broadcast_scalar (f_HTE, master_task) call broadcast_scalar (f_ANGLE, master_task) call broadcast_scalar (f_ANGLET, master_task) call broadcast_scalar (f_bounds, master_task) call broadcast_scalar (f_NCAT, master_task) call broadcast_scalar (f_VGRDi, master_task) call broadcast_scalar (f_VGRDs, master_task) call broadcast_scalar (f_VGRDb, master_task) ! call broadcast_scalar (f_example, master_task) call broadcast_scalar (f_hi, master_task) call broadcast_scalar (f_hs, master_task) call broadcast_scalar (f_Tsfc, master_task) call broadcast_scalar (f_aice, master_task) call broadcast_scalar (f_uvel, master_task) call broadcast_scalar (f_vvel, master_task) call broadcast_scalar (f_uatm, master_task) call broadcast_scalar (f_vatm, master_task) call broadcast_scalar (f_sice, master_task) call broadcast_scalar (f_fswdn, master_task) call broadcast_scalar (f_flwdn, master_task) call broadcast_scalar (f_snow, master_task) call broadcast_scalar (f_snow_ai, master_task) call broadcast_scalar (f_rain, master_task) call broadcast_scalar (f_rain_ai, master_task) call broadcast_scalar (f_sst, master_task) call broadcast_scalar (f_sss, master_task) call broadcast_scalar (f_uocn, master_task) call broadcast_scalar (f_vocn, master_task) call broadcast_scalar (f_frzmlt, master_task) call broadcast_scalar (f_fswfac, master_task) call broadcast_scalar (f_fswint_ai, master_task) call broadcast_scalar (f_fswabs, master_task) call broadcast_scalar (f_fswabs_ai, master_task) call broadcast_scalar (f_albsni, master_task) call broadcast_scalar (f_alvdr, master_task) call broadcast_scalar (f_alidr, master_task) call broadcast_scalar (f_alvdf, master_task) call broadcast_scalar (f_alidf, master_task) call broadcast_scalar (f_albice, master_task) call broadcast_scalar (f_albsno, master_task) call broadcast_scalar (f_albpnd, master_task) call broadcast_scalar (f_coszen, master_task) call broadcast_scalar (f_flat, master_task) call broadcast_scalar (f_flat_ai, master_task) call broadcast_scalar (f_fsens, master_task) call broadcast_scalar (f_fsens_ai, master_task) call broadcast_scalar (f_flwup, master_task) call broadcast_scalar (f_flwup_ai, master_task) call broadcast_scalar (f_evap, master_task) call broadcast_scalar (f_evap_ai, master_task) call broadcast_scalar (f_Tair, master_task) call broadcast_scalar (f_Tref, master_task) call broadcast_scalar (f_Qref, master_task) call broadcast_scalar (f_congel, master_task) call broadcast_scalar (f_frazil, master_task) call broadcast_scalar (f_snoice, master_task) call broadcast_scalar (f_dsnow, master_task) call broadcast_scalar (f_meltt, master_task) call broadcast_scalar (f_melts, master_task) call broadcast_scalar (f_meltb, master_task) call broadcast_scalar (f_meltl, master_task) call broadcast_scalar (f_fresh, master_task) call broadcast_scalar (f_fresh_ai, master_task) call broadcast_scalar (f_fsalt, master_task) call broadcast_scalar (f_fsalt_ai, master_task) call broadcast_scalar (f_fhocn, master_task) call broadcast_scalar (f_fhocn_ai, master_task) call broadcast_scalar (f_fswthru, master_task) call broadcast_scalar (f_fswthru_ai, master_task) call broadcast_scalar (f_strairx, master_task) call broadcast_scalar (f_strairy, master_task) call broadcast_scalar (f_strtltx, master_task) call broadcast_scalar (f_strtlty, master_task) call broadcast_scalar (f_strcorx, master_task) call broadcast_scalar (f_strcory, master_task) call broadcast_scalar (f_strocnx, master_task) call broadcast_scalar (f_strocny, master_task) call broadcast_scalar (f_strintx, master_task) call broadcast_scalar (f_strinty, master_task) call broadcast_scalar (f_strength, master_task) call broadcast_scalar (f_divu, master_task) call broadcast_scalar (f_shear, master_task) call broadcast_scalar (f_sig1, master_task) call broadcast_scalar (f_sig2, master_task) call broadcast_scalar (f_dvidtt, master_task) call broadcast_scalar (f_dvidtd, master_task) call broadcast_scalar (f_daidtt, master_task) call broadcast_scalar (f_daidtd, master_task) call broadcast_scalar (f_dagedtt, master_task) call broadcast_scalar (f_dagedtd, master_task) call broadcast_scalar (f_mlt_onset, master_task) call broadcast_scalar (f_frz_onset, master_task) call broadcast_scalar (f_aisnap, master_task) call broadcast_scalar (f_hisnap, master_task) call broadcast_scalar (f_aicen, master_task) call broadcast_scalar (f_vicen, master_task) call broadcast_scalar (f_vsnon, master_task) call broadcast_scalar (f_trsig, master_task) call broadcast_scalar (f_icepresent, master_task) call broadcast_scalar (f_fsurf_ai, master_task) call broadcast_scalar (f_fcondtop_ai, master_task) call broadcast_scalar (f_fmeltt_ai, master_task) call broadcast_scalar (f_fsurfn_ai, master_task) call broadcast_scalar (f_fcondtopn_ai, master_task) call broadcast_scalar (f_fmelttn_ai, master_task) call broadcast_scalar (f_flatn_ai, master_task) call broadcast_scalar (f_fsensn_ai, master_task) ! call broadcast_scalar (f_field3dz, master_task) call broadcast_scalar (f_keffn_top, master_task) call broadcast_scalar (f_Tinz, master_task) call broadcast_scalar (f_Sinz, master_task) call broadcast_scalar (f_Tsnz, master_task) call broadcast_scalar (f_iage, master_task) call broadcast_scalar (f_FY, master_task) call broadcast_scalar (f_a11, master_task) call broadcast_scalar (f_a12, master_task) call broadcast_scalar (f_e11, master_task) call broadcast_scalar (f_e12, master_task) call broadcast_scalar (f_e22, master_task) call broadcast_scalar (f_s11, master_task) call broadcast_scalar (f_s12, master_task) call broadcast_scalar (f_s22, master_task) call broadcast_scalar (f_yieldstress11, master_task) call broadcast_scalar (f_yieldstress12, master_task) call broadcast_scalar (f_yieldstress22, master_task) ! 2D variables do ns1 = 1, nstreams if (histfreq(ns1) /= 'x') then !!!!! begin example ! call define_hist_field(n_example,"example","m",tstr2D, tcstr, & ! "example: mean ice thickness", & ! "ice volume per unit grid cell area", c1, c0, & ! ns1, f_example) !!!!! end example call define_hist_field(n_hi,"hi","m",tstr2D, tcstr, & "grid cell mean ice thickness", & "ice volume per unit grid cell area", c1, c0, & ns1, f_hi) call define_hist_field(n_hs,"hs","m",tstr2D, tcstr, & "grid cell mean snow thickness", & "snow volume per unit grid cell area", c1, c0, & ns1, f_hs) call define_hist_field(n_Tsfc,"Tsfc","C",tstr2D, tcstr, & "snow/ice surface temperature", & "averaged with Tf if no ice is present", c1, c0, & ns1, f_Tsfc) call define_hist_field(n_aice,"aice","1",tstr2D, tcstr, & "ice area (aggregate)", & "none", c1, c0, & ns1, f_aice) call define_hist_field(n_uvel,"uvel","m/s",ustr2D, ucstr, & "ice velocity (x)", & "positive is x direction on U grid", c1, c0, & ns1, f_uvel) call define_hist_field(n_vvel,"vvel","m/s",ustr2D, ucstr, & "ice velocity (y)", & "positive is y direction on U grid", c1, c0, & ns1, f_vvel) call define_hist_field(n_uatm,"uatm","m/s",ustr2D, ucstr, & "atm velocity (x)", & "positive is x direction on U grid", c1, c0, & ns1, f_uatm) call define_hist_field(n_vatm,"vatm","m/s",ustr2D, ucstr, & "atm velocity (y)", & "positive is y direction on U grid", c1, c0, & ns1, f_vatm) call define_hist_field(n_sice,"sice","ppt",tstr2D, tcstr, & "bulk ice salinity", & "none", c1, c0, & ns1, f_sice) call define_hist_field(n_fswdn,"fswdn","W/m^2",tstr2D, tcstr, & "down solar flux", & "positive downward", c1, c0, & ns1, f_fswdn) call define_hist_field(n_flwdn,"flwdn","W/m^2",tstr2D, tcstr, & "down longwave flux", & "positive downward", c1, c0, & ns1, f_flwdn) call define_hist_field(n_snow,"snow","cm/day",tstr2D, tcstr, & "snowfall rate (cpl)", & "none", mps_to_cmpdy/rhofresh, c0, & ns1, f_snow) call define_hist_field(n_snow_ai,"snow_ai","cm/day",tstr2D, tcstr, & "snowfall rate", & "weighted by ice area", mps_to_cmpdy/rhofresh, c0, & ns1, f_snow_ai) call define_hist_field(n_rain,"rain","cm/day",tstr2D, tcstr, & "rainfall rate (cpl)", & "none", mps_to_cmpdy/rhofresh, c0, & ns1, f_rain) call define_hist_field(n_rain_ai,"rain_ai","cm/day",tstr2D, tcstr, & "rainfall rate", & "weighted by ice area", mps_to_cmpdy/rhofresh, c0, & ns1, f_rain_ai) call define_hist_field(n_sst,"sst","C",tstr2D, tcstr, & "sea surface temperature", & "none", c1, c0, & ns1, f_sst) call define_hist_field(n_sss,"sss","ppt",tstr2D, tcstr, & "sea surface salinity", & "none", c1, c0, & ns1, f_sss) call define_hist_field(n_uocn,"uocn","m/s",ustr2D, ucstr, & "ocean current (x)", & "positive is x direction on U grid", c1, c0, & ns1, f_uocn) call define_hist_field(n_vocn,"vocn","m/s",ustr2D, ucstr, & "ocean current (y)", & "positive is y direction on U grid", c1, c0, & ns1, f_vocn) call define_hist_field(n_frzmlt,"frzmlt","W/m^2",tstr2D, tcstr, & "freeze/melt potential", & "if >0, new ice forms; if <0, ice melts", c1, c0, & ns1, f_frzmlt) call define_hist_field(n_fswfac,"fswfac","1",tstr2D, tcstr, & "shortwave scaling factor", & "ratio of netsw new:old", c1, c0, & ns1, f_fswfac) call define_hist_field(n_fswint_ai,"fswint_ai","W/m^2",tstr2D, tcstr, & "shortwave absorbed in ice interior", & "does not include surface", c1, c0, & ns1, f_fswint_ai) call define_hist_field(n_fswabs,"fswabs","W/m^2",tstr2D, tcstr, & "snow/ice/ocn absorbed solar flux (cpl)", & "positive downward", c1, c0, & ns1, f_fswabs) call define_hist_field(n_fswabs_ai,"fswabs_ai","W/m^2",tstr2D, tcstr, & "snow/ice/ocn absorbed solar flux", & "weighted by ice area", c1, c0, & ns1, f_fswabs_ai) ! call define_hist_field(n_albsni,"albsni","%",tstr2D, tcstr, & ! "snow/ice broad band albedo", & ! "scaled (divided) by aice", c100, c0, & ! ns1, f_albsni) call define_hist_field(n_albsni,"albsni","%",tstr2D, tcstr, & "snow/ice broad band albedo", & "averaged for coszen>0, weighted by aice", c100, c0, & ns1, f_albsni) call define_hist_field(n_alvdr,"alvdr","%",tstr2D, tcstr, & "visible direct albedo", & "scaled (divided) by aice", c100, c0, & ns1, f_alvdr) call define_hist_field(n_alidr,"alidr","%",tstr2D, tcstr, & "near IR direct albedo", & "scaled (divided) by aice", c100, c0, & ns1, f_alidr) call define_hist_field(n_alvdf,"alvdf","%",tstr2D, tcstr, & "visible diffuse albedo", & "scaled (divided) by aice", c100, c0, & ns1, f_alvdf) call define_hist_field(n_alidf,"alidf","%",tstr2D, tcstr, & "near IR diffuse albedo", & "scaled (divided) by aice", c100, c0, & ns1, f_alidf) call define_hist_field(n_albice,"albice","%",tstr2D, tcstr, & "bare ice albedo", & "averaged for coszen>0, weighted by aice", c100, c0, & ns1, f_albice) call define_hist_field(n_albsno,"albsno","%",tstr2D, tcstr, & "snow albedo", & "averaged for coszen>0, weighted by aice", c100, c0, & ns1, f_albsno) call define_hist_field(n_albpnd,"albpnd","%",tstr2D, tcstr, & "melt pond albedo", & "averaged for coszen>0, weighted by aice", c100, c0, & ns1, f_albpnd) call define_hist_field(n_coszen,"coszen","radian",tstr2D, tcstr, & "cosine of the zenith angle", & "negative below horizon", c1, c0, & ns1, f_coszen) call define_hist_field(n_flat,"flat","W/m^2",tstr2D, tcstr, & "latent heat flux (cpl)", & "positive downward", c1, c0, & ns1, f_flat) call define_hist_field(n_flat_ai,"flat_ai","W/m^2",tstr2D, tcstr, & "latent heat flux", & "weighted by ice area", c1, c0, & ns1, f_flat_ai) call define_hist_field(n_fsens,"fsens","W/m^2",tstr2D, tcstr, & "sensible heat flux (cpl)", & "positive downward", c1, c0, & ns1, f_fsens) call define_hist_field(n_fsens_ai,"fsens_ai","W/m^2",tstr2D, tcstr, & "sensible heat flux", & "weighted by ice area", c1, c0, & ns1, f_fsens_ai) call define_hist_field(n_flwup,"flwup","W/m^2",tstr2D, tcstr, & "upward longwave flux (cpl)", & "positive downward", c1, c0, & ns1, f_flwup) call define_hist_field(n_flwup_ai,"flwup_ai","W/m^2",tstr2D, tcstr, & "upward longwave flux", & "weighted by ice area", c1, c0, & ns1, f_flwup_ai) call define_hist_field(n_evap,"evap","cm/day",tstr2D, tcstr, & "evaporative water flux (cpl)", & "none", mps_to_cmpdy/rhofresh, c0, & ns1, f_evap) call define_hist_field(n_evap_ai,"evap_ai","cm/day",tstr2D, tcstr, & "evaporative water flux", & "weighted by ice area", mps_to_cmpdy/rhofresh, c0, & ns1, f_evap_ai) call define_hist_field(n_Tair,"Tair","C",tstr2D, tcstr, & "air temperature", & "none", c1, -Tffresh, & ns1, f_Tair) call define_hist_field(n_Tref,"Tref","C",tstr2D, tcstr, & "2m reference temperature", & "none", c1, -Tffresh, & ns1, f_Tref) call define_hist_field(n_Qref,"Qref","g/kg",tstr2D, tcstr, & "2m reference specific humidity", & "none", kg_to_g, c0, & ns1, f_Qref) call define_hist_field(n_congel,"congel","cm/day",tstr2D, tcstr, & "congelation ice growth", & "none", mps_to_cmpdy/dt, c0, & ns1, f_congel) call define_hist_field(n_frazil,"frazil","cm/day",tstr2D, tcstr, & "frazil ice growth", & "none", mps_to_cmpdy/dt, c0, & ns1, f_frazil) call define_hist_field(n_snoice,"snoice","cm/day",tstr2D, tcstr, & "snow-ice formation", & "none", mps_to_cmpdy/dt, c0, & ns1, f_snoice) call define_hist_field(n_dsnow,"dsnow","cm/day",tstr2D, tcstr, & "snow formation", & "none", mps_to_cmpdy/dt, c0, & ns1, f_dsnow) call define_hist_field(n_meltt,"meltt","cm/day",tstr2D, tcstr, & "top ice melt", & "none", mps_to_cmpdy/dt, c0, & ns1, f_meltt) call define_hist_field(n_melts,"melts","cm/day",tstr2D, tcstr, & "top snow melt", & "none", mps_to_cmpdy/dt, c0, & ns1, f_melts) call define_hist_field(n_meltb,"meltb","cm/day",tstr2D, tcstr, & "basal ice melt", & "none", mps_to_cmpdy/dt, c0, & ns1, f_meltb) call define_hist_field(n_meltl,"meltl","cm/day",tstr2D, tcstr, & "lateral ice melt", & "none", mps_to_cmpdy/dt, c0, & ns1, f_meltl) call define_hist_field(n_fresh,"fresh","cm/day",tstr2D, tcstr, & "freshwtr flx ice to ocn (cpl)", & "if positive, ocean gains fresh water", & mps_to_cmpdy/rhofresh, c0, & ns1, f_fresh) call define_hist_field(n_fresh_ai,"fresh_ai","cm/day",tstr2D, tcstr, & "freshwtr flx ice to ocn", & "weighted by ice area", mps_to_cmpdy/rhofresh, c0, & ns1, f_fresh_ai) call define_hist_field(n_fsalt,"fsalt","kg/m^2/s",tstr2D, tcstr, & "salt flux ice to ocn (cpl)", & "if positive, ocean gains salt", c1, c0, & ns1, f_fsalt) call define_hist_field(n_fsalt_ai,"fsalt_ai","kg/m^2/s",tstr2D, tcstr, & "salt flux ice to ocean", & "weighted by ice area", c1, c0, & ns1, f_fsalt_ai) call define_hist_field(n_fhocn,"fhocn","W/m^2",tstr2D, tcstr, & "heat flux ice to ocn (cpl)", & "if positive, ocean gains heat", c1, c0, & ns1, f_fhocn) call define_hist_field(n_fhocn_ai,"fhocn_ai","W/m^2",tstr2D, tcstr, & "heat flux ice to ocean", & "weighted by ice area", c1, c0, & ns1, f_fhocn_ai) call define_hist_field(n_fswthru,"fswthru","W/m^2",tstr2D, tcstr, & "SW thru ice to ocean (cpl)", & "if positive, ocean gains heat", c1, c0, & ns1, f_fswthru) call define_hist_field(n_fswthru_ai,"fswthru_ai","W/m^2",tstr2D, tcstr,& "SW flux thru ice to ocean", & "weighted by ice area", c1, c0, & ns1, f_fswthru_ai) call define_hist_field(n_strairx,"strairx","N/m^2",ustr2D, ucstr, & "atm/ice stress (x)", & "positive is x direction on U grid", c1, c0, & ns1, f_strairx) call define_hist_field(n_strairy,"strairy","N/m^2",ustr2D, ucstr, & "atm/ice stress (y)", & "positive is y direction on U grid", c1, c0, & ns1, f_strairy) call define_hist_field(n_strtltx,"strtltx","N/m^2",ustr2D, ucstr, & "sea sfc tilt stress (x)", & "none", c1, c0, & ns1, f_strtltx) call define_hist_field(n_strtlty,"strtlty","N/m^2",ustr2D, ucstr, & "sea sfc tilt stress (y)", & "none", c1, c0, & ns1, f_strtlty) call define_hist_field(n_strcorx,"strcorx","N/m^2",ustr2D, ucstr, & "coriolis stress (x)", & "positive is x direction on U grid", c1, c0, & ns1, f_strcorx) call define_hist_field(n_strcory,"strcory","N/m^2",ustr2D, ucstr, & "coriolis stress (y)", & "positive is y direction on U grid", c1, c0, & ns1, f_strcory) call define_hist_field(n_strocnx,"strocnx","N/m^2",ustr2D, ucstr, & "ocean/ice stress (x)", & "positive is x direction on U grid", c1, c0, & ns1, f_strocnx) call define_hist_field(n_strocny,"strocny","N/m^2",ustr2D, ucstr, & "ocean/ice stress (y)", & "positive is y direction on U grid", c1, c0, & ns1, f_strocny) call define_hist_field(n_strintx,"strintx","N/m^2",ustr2D, ucstr, & "internal ice stress (x)", & "positive is x direction on U grid", c1, c0, & ns1, f_strintx) call define_hist_field(n_strinty,"strinty","N/m^2",ustr2D, ucstr, & "internal ice stress (y)", & "positive is y direction on U grid", c1, c0, & ns1, f_strinty) call define_hist_field(n_strength,"strength","N/m",tstr2D, tcstr, & "compressive ice strength", & "none", c1, c0, & ns1, f_strength) call define_hist_field(n_divu,"divu","%/day",tstr2D, tcstr, & "strain rate (divergence)", & "none", secday*c100, c0, & ns1, f_divu) call define_hist_field(n_shear,"shear","%/day",tstr2D, tcstr, & "strain rate (shear)", & "none", secday*c100, c0, & ns1, f_shear) call define_hist_field(n_sig1,"sig1","1",ustr2D, ucstr, & "norm. principal stress 1", & "sig1 is instantaneous", c1, c0, & ns1, f_sig1) call define_hist_field(n_sig2,"sig2","1",ustr2D, ucstr, & "norm. principal stress 2", & "sig2 is instantaneous", c1, c0, & ns1, f_sig2) call define_hist_field(n_dvidtt,"dvidtt","cm/day",tstr2D, tcstr, & "volume tendency thermo", & "none", mps_to_cmpdy, c0, & ns1, f_dvidtt) call define_hist_field(n_dvidtd,"dvidtd","cm/day",tstr2D, tcstr, & "volume tendency dynamics", & "none", mps_to_cmpdy, c0, & ns1, f_dvidtd) call define_hist_field(n_daidtt,"daidtt","%/day",tstr2D, tcstr, & "area tendency thermo", & "none", secday*c100, c0, & ns1, f_daidtt) call define_hist_field(n_daidtd,"daidtd","%/day",tstr2D, tcstr, & "area tendency dynamics", & "none", secday*c100, c0, & ns1, f_daidtd) call define_hist_field(n_dagedtt,"dagedtt","day/day",tstr2D, tcstr, & "age tendency thermo", & "excludes time step increment", c1, c0, & ns1, f_dagedtt) call define_hist_field(n_dagedtd,"dagedtd","day/day",tstr2D, tcstr, & "age tendency dynamics", & "excludes time step increment", c1, c0, & ns1, f_dagedtd) call define_hist_field(n_mlt_onset,"mlt_onset","day of year", & tstr2D, tcstr,"melt onset date", & "midyear restart gives erroneous dates", c1, c0, & ns1, f_mlt_onset) call define_hist_field(n_frz_onset,"frz_onset","day of year", & tstr2D, tcstr,"freeze onset date", & "midyear restart gives erroneous dates", c1, c0, & ns1, f_frz_onset) call define_hist_field(n_hisnap,"hisnap","m",tstr2D, tcstr, & "ice volume snapshot", & "none", c1, c0, & ns1, f_hisnap) call define_hist_field(n_aisnap,"aisnap","1",tstr2D, tcstr, & "ice area snapshot", & "none", c1, c0, & ns1, f_aisnap) call define_hist_field(n_trsig,"trsig","N/m^2",tstr2D, tcstr, & "internal stress tensor trace", & "ice strength approximation", c1, c0, & ns1, f_trsig) call define_hist_field(n_icepresent,"ice_present","1",tstr2D, tcstr, & "fraction of time-avg interval that ice is present", & "ice extent flag", c1, c0, & ns1, f_icepresent) call define_hist_field(n_fsurf_ai,"fsurf_ai","W/m^2",tstr2D, tcstr, & "net surface heat flux", & "positive downward, excludes conductive flux, weighted by ice area", & c1, c0, & ns1, f_fsurf_ai) call define_hist_field(n_fcondtop_ai,"fcondtop_ai","W/m^2", & tstr2D, tcstr,"top surface conductive heat flux", & "positive downward, weighted by ice area", c1, c0, & ns1, f_fcondtop_ai) call define_hist_field(n_fmeltt_ai,"fmeltt_ai","W/m^2",tstr2D, tcstr, & "net surface heat flux causing melt", & "always >= 0, weighted by ice area", c1, c0, & ns1, f_fmeltt_ai) call define_hist_field(n_a11,"a11"," ",tstr2D, tcstr, & "a11: component a11 of the structure tensor", & "none", c1, c0, & ns1, f_a11) call define_hist_field(n_a12,"a12"," ",tstr2D, tcstr, & "a12: component a12 of the structure tensor", & "none", c1, c0, & ns1, f_a12) call define_hist_field(n_e11,"e11","1/s",tstr2D, tcstr, & "e11: component e11 of the strain rate tensor", & "none", c1, c0, & ns1, f_e11) call define_hist_field(n_e12,"e12","1/s",tstr2D, tcstr, & "e12: component e12 of the strain rate tensor", & "none", c1, c0, & ns1, f_e12) call define_hist_field(n_e22,"e22","1/s",tstr2D, tcstr, & "e22: component e22 of the strain rate tensor", & "none", c1, c0, & ns1, f_e22) call define_hist_field(n_s11,"s11","kg/s^2",tstr2D, tcstr, & "s11: component s11 of the stress tensor", & "none", c1, c0, & ns1, f_s11) call define_hist_field(n_s12,"s12","kg/s^2",tstr2D, tcstr, & "s12: component s12 of the stress tensor", & "none", c1, c0, & ns1, f_s12) call define_hist_field(n_s22,"s22","kg/s^2",tstr2D, tcstr, & "s22: component s12 of the stress tensor", & "none", c1, c0, & ns1, f_s22) call define_hist_field(n_yieldstress11,"yieldstress11","kg/s^2",tstr2D, tcstr, & "yieldstress11: component 11 of the yieldstress tensor", & "none", c1, c0, & ns1, f_yieldstress11) call define_hist_field(n_yieldstress12,"yieldstress12","kg/s^2",tstr2D, tcstr, & "yieldstress12: component 12 of the yieldstress tensor", & "none", c1, c0, & ns1, f_yieldstress12) call define_hist_field(n_yieldstress22,"yieldstress22","kg/s^2",tstr2D, tcstr, & "yieldstress22: component 12 of the yieldstress tensor", & "none", c1, c0, & ns1, f_yieldstress22) ! Tracers ! Ice Age call define_hist_field(n_iage,"iage","years",tstr2D, tcstr, & "sea ice age", & "none", c1/(secday*days_per_year), c0, & ns1, f_iage) ! First Year Ice Area call define_hist_field(n_FY,"FYarea"," ",tstr2D, tcstr, & "first-year ice area", & "weighted by ice area", c1, c0, & ns1, f_FY) endif ! if (histfreq(ns1) /= 'x') then enddo ! ns1 ! other 2D history variables ! mechanical redistribution call init_hist_mechred_2D ! melt ponds if (tr_pond) call init_hist_pond_2D ! biogeochemistry if (tr_aero .or. tr_brine .or. skl_bgc) call init_hist_bgc_2D if (formdrag) call init_hist_drag_2D !----------------------------------------------------------------- ! 3D (category) variables looped separately for ordering !----------------------------------------------------------------- do ns1 = 1, nstreams if (histfreq(ns1) /= 'x') then call define_hist_field(n_aicen,"aicen","1",tstr3Dc, tcstr, & "ice area, categories","none", c1, c0, & ns1, f_aicen) call define_hist_field(n_vicen,"vicen","m",tstr3Dc, tcstr, & "ice volume, categories","none", c1, c0, & ns1, f_vicen) call define_hist_field(n_vsnon,"vsnon","m",tstr3Dc, tcstr, & "snow depth on ice, categories","volume per unit area of snow", c1, c0, & ns1, f_vsnon) call define_hist_field(n_fsurfn_ai,"fsurfn_ai","W/m^2",tstr3Dc, tcstr, & "net surface heat flux, categories","weighted by ice area", c1, c0, & ns1, f_fsurfn_ai) call define_hist_field(n_fcondtopn_ai,"fcondtopn_ai","W/m^2",tstr3Dc, tcstr, & "top sfc conductive heat flux, cat","weighted by ice area", c1, c0, & ns1, f_fcondtopn_ai) call define_hist_field(n_fmelttn_ai,"fmelttn_ai","W/m^2",tstr3Dc, tcstr, & "net sfc heat flux causing melt, cat","weighted by ice area", c1, c0, & ns1, f_fmelttn_ai) call define_hist_field(n_flatn_ai,"flatn_ai","W/m^2",tstr3Dc, tcstr, & "latent heat flux, category","weighted by ice area", c1, c0, & ns1, f_flatn_ai) call define_hist_field(n_fsensn_ai,"fsensn_ai","W/m^2",tstr3Dc, tcstr, & "sensible heat flux, category","weighted by ice area", c1, c0, & ns1, f_fsensn_ai) call define_hist_field(n_keffn_top,"keffn_top","W/m^2/K",tstr3Dc, tcstr, & "effective thermal conductivity of the top ice layer, categories", & "multilayer scheme", c1, c0, & ns1, f_keffn_top) endif ! if (histfreq(ns1) /= 'x') then enddo ! ns1 ! other 3D (category) history variables ! mechanical redistribution call init_hist_mechred_3Dc ! melt ponds if (tr_pond) call init_hist_pond_3Dc ! biogeochemistry if (tr_brine) call init_hist_bgc_3Dc !----------------------------------------------------------------- ! 3D (vertical) variables must be looped separately !----------------------------------------------------------------- ! do ns1 = 1, nstreams ! if (histfreq(ns1) /= 'x') then ! call define_hist_field(n_field3dz,"field3dz","1",tstr3Dz, tcstr, & ! "example 3dz field", & ! "vertical profile", c1, c0, & ! ns1, f_field3dz) ! endif ! if (histfreq(ns1) /= 'x') then ! enddo ! ns1 !----------------------------------------------------------------- ! 4D (categories, vertical) variables must be looped separately !----------------------------------------------------------------- do ns1 = 1, nstreams if (histfreq(ns1) /= 'x') then call define_hist_field(n_Tinz,"Tinz","C",tstr4Di, tcstr, & "ice internal temperatures on CICE grid", & "vertical profile", c1, c0, & ns1, f_Tinz) call define_hist_field(n_Sinz,"Sinz","ppt",tstr4Di, tcstr, & "ice internal bulk salinity", & "vertical profile", c1, c0, & ns1, f_Sinz) endif ! if (histfreq(ns1) /= 'x') then enddo ! ns1 do ns1 = 1, nstreams if (histfreq(ns1) /= 'x') then call define_hist_field(n_Tsnz,"Tsnz","C",tstr4Ds, tcstr, & "snow internal temperatures", & "vertical profile", c1, c0, & ns1, f_Tsnz) endif ! if (histfreq(ns1) /= 'x') then enddo if (f_Tinz (1:1) /= 'x') then if (allocated(Tinz4d)) deallocate(Tinz4d) allocate(Tinz4d(nx_block,ny_block,nzilyr,ncat_hist)) endif if (f_Sinz (1:1) /= 'x') then if (allocated(Sinz4d)) deallocate(Sinz4d) allocate(Sinz4d(nx_block,ny_block,nzilyr,ncat_hist)) endif if (f_Tsnz (1:1) /= 'x') then if (allocated(Tsnz4d)) deallocate(Tsnz4d) allocate(Tsnz4d(nx_block,ny_block,nzslyr,ncat_hist)) endif if (f_Sinz (1:1) /= 'x') then if (allocated(Sinz4d)) deallocate(Sinz4d) allocate(Sinz4d(nx_block,ny_block,nzilyr,ncat_hist)) endif ! other 4D history variables ! biogeochemistry if (tr_brine) call init_hist_bgc_4Db !----------------------------------------------------------------- ! fill igrd array with namelist values !----------------------------------------------------------------- igrd=.true. igrd(n_tmask ) = f_tmask igrd(n_blkmask ) = f_blkmask igrd(n_tarea ) = f_tarea igrd(n_uarea ) = f_uarea igrd(n_dxt ) = f_dxt igrd(n_dyt ) = f_dyt igrd(n_dxu ) = f_dxu igrd(n_dyu ) = f_dyu igrd(n_HTN ) = f_HTN igrd(n_HTE ) = f_HTE igrd(n_ANGLE ) = f_ANGLE igrd(n_ANGLET ) = f_ANGLET igrdz=.true. igrdz(n_NCAT ) = f_NCAT igrdz(n_VGRDi ) = f_VGRDi igrdz(n_VGRDs ) = f_VGRDs igrdz(n_VGRDb ) = f_VGRDb !----------------------------------------------------------------- ! diagnostic output !----------------------------------------------------------------- ntmp(:) = 0 if (my_task == master_task) then write(nu_diag,*) ' ' write(nu_diag,*) 'The following variables will be ', & 'written to the history tape: ' write(nu_diag,101) 'description','units','variable','frequency','x' do n=1,num_avail_hist_fields_tot if (avail_hist_fields(n)%vhistfreq_n /= 0) & write(nu_diag,100) avail_hist_fields(n)%vdesc, & avail_hist_fields(n)%vunit, avail_hist_fields(n)%vname, & avail_hist_fields(n)%vhistfreq,avail_hist_fields(n)%vhistfreq_n do ns = 1, nstreams if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) & ntmp(ns)=ntmp(ns)+1 enddo enddo ! num_avail_hist_fields_tot write(nu_diag,*) ' ' endif 100 format (1x,a40,2x,a16,2x,a12,1x,a1,2x,i6) 101 format (2x,a19,10x,a16,9x,a12,2x,a,3x,a1) call broadcast_array(ntmp, master_task) do ns = 1, nstreams if (ntmp(ns)==0) histfreq_n(ns) = 0 enddo !----------------------------------------------------------------- ! initialize the history arrays !----------------------------------------------------------------- if (allocated(a2D)) deallocate(a2D) if (num_avail_hist_fields_2D > 0) & allocate(a2D(nx_block,ny_block,num_avail_hist_fields_2D,max_blocks)) if (allocated(a3Dc)) deallocate(a3Dc) if (num_avail_hist_fields_3Dc > 0) & allocate(a3Dc(nx_block,ny_block,ncat_hist,num_avail_hist_fields_3Dc,max_blocks)) nzlyr = max(nzilyr, nzslyr) if (allocated(a3Dz)) deallocate(a3Dz) if (num_avail_hist_fields_3Dz > 0) & allocate(a3Dz(nx_block,ny_block,nzlyr,num_avail_hist_fields_3Dz,max_blocks)) nzlyrb = nzblyr if (allocated(a3Db)) deallocate(a3Db) if (num_avail_hist_fields_3Db > 0) & allocate(a3Db(nx_block,ny_block,nzlyrb,num_avail_hist_fields_3Db,max_blocks)) if (allocated(a4Di)) deallocate(a4Di) if (num_avail_hist_fields_4Di > 0) & allocate(a4Di(nx_block,ny_block,nzilyr,ncat_hist,num_avail_hist_fields_4Di,max_blocks)) if (allocated(a4Ds)) deallocate(a4Ds) if (num_avail_hist_fields_4Ds > 0) & allocate(a4Ds(nx_block,ny_block,nzslyr,ncat_hist,num_avail_hist_fields_4Ds,max_blocks)) if (allocated(a4Db)) deallocate(a4Db) if (num_avail_hist_fields_4Db > 0) & allocate(a4Db(nx_block,ny_block,nzblyr,ncat_hist,num_avail_hist_fields_4Db,max_blocks)) if (allocated(a2D)) a2D (:,:,:,:) = c0 if (allocated(a3Dc)) a3Dc(:,:,:,:,:) = c0 if (allocated(a3Dz)) a3Dz(:,:,:,:,:) = c0 if (allocated(a3Db)) a3Db(:,:,:,:,:) = c0 if (allocated(a4Di)) a4Di(:,:,:,:,:,:) = c0 if (allocated(a4Ds)) a4Ds(:,:,:,:,:,:) = c0 if (allocated(a4Db)) a4Db(:,:,:,:,:,:) = c0 avgct(:) = c0 albcnt(:,:,:,:) = c0 if (restart .and. yday >= c2) then ! restarting midyear gives erroneous onset dates mlt_onset = 999._dbl_kind frz_onset = 999._dbl_kind else mlt_onset = c0 frz_onset = c0 endif end subroutine init_hist !======================================================================= ! accumulate average ice quantities or snapshots ! ! author: Elizabeth C. Hunke, LANL subroutine accum_hist (dt) use ice_blocks, only: block, get_block, nx_block, ny_block use ice_constants, only: c0, c1, p25, puny, secday, depressT, & awtvdr, awtidr, awtvdf, awtidf, Lfresh, rhos, cp_ice, spval use ice_domain, only: blocks_ice, nblocks use ice_grid, only: tmask, lmask_n, lmask_s use ice_calendar, only: new_year, write_history, & write_ic, time, histfreq, nstreams, month, & new_month use ice_dyn_eap, only: a11, a12, e11, e12, e22, s11, s12, s22, & yieldstress11, yieldstress12, yieldstress22 use ice_dyn_shared, only: kdyn, principal_stress use ice_flux, only: fsw, flw, fsnow, frain, sst, sss, uocn, vocn, & frzmlt_init, fswfac, fswabs, fswthru, alvdr, alvdf, alidr, alidf, & albice, albsno, albpnd, coszen, flat, fsens, flwout, evap, & Tair, Tref, Qref, congel, frazil, snoice, dsnow, & melts, meltb, meltt, meltl, fresh, fsalt, fresh_ai, fsalt_ai, & fhocn, fhocn_ai, uatm, vatm, & fswthru_ai, strairx, strairy, strtltx, strtlty, strintx, strinty, & strocnx, strocny, fm, daidtt, dvidtt, daidtd, dvidtd, fsurf, & fcondtop, fsurfn, fcondtopn, flatn, fsensn, albcnt, prs_sig, & stressp_1, stressm_1, stress12_1, & stressp_2, stressm_2, stress12_2, & stressp_3, stressm_3, stress12_3, & stressp_4, stressm_4, stress12_4, sig1, sig2, & mlt_onset, frz_onset, dagedtt, dagedtd, fswint_ai, keffn_top use ice_atmo, only: formdrag use ice_history_shared ! almost everything use ice_history_write, only: ice_write_hist use ice_history_bgc, only: accum_hist_bgc use ice_history_mechred, only: accum_hist_mechred use ice_history_pond, only: accum_hist_pond use ice_history_drag, only: accum_hist_drag use ice_state ! almost everything use ice_therm_shared, only: calculate_Tin_from_qin, Tmlt, ktherm use ice_therm_mushy, only: temperature_mush, temperature_snow use ice_timers, only: ice_timer_start, ice_timer_stop, timer_readwrite use ice_zbgc_shared, only: skl_bgc real (kind=dbl_kind), intent(in) :: & dt ! time step ! local variables integer (kind=int_kind) :: & i,j,k,ic,n,ns,nn, & iblk , & ! block index ilo,ihi,jlo,jhi , & ! beginning and end of physical domain nstrm ! nstreams (1 if writing initial condition) real (kind=dbl_kind) :: & ravgct , & ! 1/avgct ravgctz ! 1/avgct real (kind=dbl_kind) :: & qn , & ! temporary variable for enthalpy Tmlts ! temporary variable for melting temperature real (kind=dbl_kind), dimension (nx_block,ny_block) :: & worka, workb type (block) :: & this_block ! block information for current block !--------------------------------------------------------------- ! increment step counter !--------------------------------------------------------------- n2D = num_avail_hist_fields_2D n3Dccum = n2D + num_avail_hist_fields_3Dc n3Dzcum = n3Dccum + num_avail_hist_fields_3Dz n3Dbcum = n3Dzcum + num_avail_hist_fields_3Db n4Dicum = n3Dbcum + num_avail_hist_fields_4Di n4Dscum = n4Dicum + num_avail_hist_fields_4Ds n4Dbcum = n4Dscum + num_avail_hist_fields_4Db ! should equal num_avail_hist_fields_tot do ns = 1,nstreams if (.not. hist_avg .or. histfreq(ns) == '1') then ! write snapshots do n = 1,n2D if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) & a2D(:,:,n,:) = c0 enddo do n = n2D + 1, n3Dccum nn = n - n2D if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) & a3Dc(:,:,:,nn,:) = c0 enddo do n = n3Dccum + 1, n3Dzcum nn = n - n3Dccum if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) & a3Dz(:,:,:,nn,:) = c0 enddo do n = n3Dzcum + 1, n3Dbcum nn = n - n3Dzcum if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) & a3Db(:,:,:,nn,:) = c0 enddo do n = n3Dbcum + 1, n4Dicum nn = n - n3Dbcum if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) & a4Di(:,:,:,:,nn,:) = c0 enddo do n = n4Dicum + 1, n4Dscum nn = n - n4Dicum if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) & a4Ds(:,:,:,:,nn,:) = c0 enddo do n = n4Dscum + 1, n4Dbcum nn = n - n4Dscum if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) & a4Db(:,:,:,:,nn,:) = c0 enddo avgct(ns) = c1 else ! write averages over time histfreq avgct(ns) = avgct(ns) + c1 ! if (avgct(ns) == c1) time_beg(ns) = (time-dt)/int(secday) if (avgct(ns) == c1) then time_beg(ns) = (time-dt)/int(secday) !sgl time_beg(ns) = real(time_beg(ns),kind=real_kind) endif endif enddo !--------------------------------------------------------------- ! increment field !--------------------------------------------------------------- !$OMP PARALLEL DO PRIVATE(iblk,i,j,ilo,ihi,jlo,jhi,this_block, & !$OMP k,n,qn,ns,worka,workb) 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 workb(:,:) = aice_init(:,:,iblk) ! if (f_example(1:1) /= 'x') & ! call accum_hist_field(n_example,iblk, vice(:,:,iblk), a2D) if (f_hi (1:1) /= 'x') & call accum_hist_field(n_hi, iblk, vice(:,:,iblk), a2D) if (f_hs (1:1) /= 'x') & call accum_hist_field(n_hs, iblk, vsno(:,:,iblk), a2D) if (f_Tsfc (1:1) /= 'x') & call accum_hist_field(n_Tsfc, iblk, trcr(:,:,nt_Tsfc,iblk), a2D) if (f_aice (1:1) /= 'x') & call accum_hist_field(n_aice, iblk, aice(:,:,iblk), a2D) if (f_uvel (1:1) /= 'x') & call accum_hist_field(n_uvel, iblk, uvel(:,:,iblk), a2D) if (f_vvel (1:1) /= 'x') & call accum_hist_field(n_vvel, iblk, vvel(:,:,iblk), a2D) if (f_uatm (1:1) /= 'x') & call accum_hist_field(n_uatm, iblk, uatm(:,:,iblk), a2D) if (f_vatm (1:1) /= 'x') & call accum_hist_field(n_vatm, iblk, vatm(:,:,iblk), a2D) if (f_sice (1:1) /= 'x') then do j = jlo, jhi do i = ilo, ihi worka(i,j) = c0 do k = 1, nzilyr worka(i,j) = worka(i,j) + trcr(i,j,nt_sice+k-1,iblk) enddo worka(i,j) = worka(i,j) / nzilyr enddo enddo call accum_hist_field(n_sice, iblk, worka(:,:), a2D) endif if (f_fswdn (1:1) /= 'x') & call accum_hist_field(n_fswdn, iblk, fsw(:,:,iblk), a2D) if (f_flwdn (1:1) /= 'x') & call accum_hist_field(n_flwdn, iblk, flw(:,:,iblk), a2D) if (f_snow (1:1) /= 'x') & call accum_hist_field(n_snow, iblk, fsnow(:,:,iblk), a2D) if (f_snow_ai(1:1) /= 'x') & call accum_hist_field(n_snow_ai,iblk, fsnow(:,:,iblk)*workb(:,:), a2D) if (f_rain (1:1) /= 'x') & call accum_hist_field(n_rain, iblk, frain(:,:,iblk), a2D) if (f_rain_ai(1:1) /= 'x') & call accum_hist_field(n_rain_ai,iblk, frain(:,:,iblk)*workb(:,:), a2D) if (f_sst (1:1) /= 'x') & call accum_hist_field(n_sst, iblk, sst(:,:,iblk), a2D) if (f_sss (1:1) /= 'x') & call accum_hist_field(n_sss, iblk, sss(:,:,iblk), a2D) if (f_uocn (1:1) /= 'x') & call accum_hist_field(n_uocn, iblk, uocn(:,:,iblk), a2D) if (f_vocn (1:1) /= 'x') & call accum_hist_field(n_vocn, iblk, vocn(:,:,iblk), a2D) if (f_frzmlt (1:1) /= 'x') & call accum_hist_field(n_frzmlt, iblk, frzmlt_init(:,:,iblk), a2D) if (f_fswfac (1:1) /= 'x') & call accum_hist_field(n_fswfac, iblk, fswfac(:,:,iblk), a2D) if (f_fswabs (1:1) /= 'x') & call accum_hist_field(n_fswabs, iblk, fswabs(:,:,iblk), a2D) if (f_fswint_ai (1:1) /= 'x') & call accum_hist_field(n_fswint_ai, iblk, fswint_ai(:,:,iblk), a2D) if (f_fswabs_ai(1:1)/= 'x') & call accum_hist_field(n_fswabs_ai, iblk, fswabs(:,:,iblk)*workb(:,:), a2D) if (f_albsni (1:1) /= 'x') & call accum_hist_field(n_albsni, iblk, & (awtvdr*alvdr(:,:,iblk) & + awtidr*alidr(:,:,iblk) & + awtvdf*alvdf(:,:,iblk) & + awtidf*alidf(:,:,iblk))*aice(:,:,iblk), a2D) ! awtvdr*alvdr(:,:,iblk) & ! + awtidr*alidr(:,:,iblk) & ! + awtvdf*alvdf(:,:,iblk) & ! + awtidf*alidf(:,:,iblk), a2D) if (f_alvdr (1:1) /= 'x') & call accum_hist_field(n_alvdr, iblk, alvdr(:,:,iblk), a2D) if (f_alidr (1:1) /= 'x') & call accum_hist_field(n_alidr, iblk, alidr(:,:,iblk), a2D) if (f_alvdf (1:1) /= 'x') & call accum_hist_field(n_alvdf, iblk, alvdf(:,:,iblk), a2D) if (f_alidf (1:1) /= 'x') & call accum_hist_field(n_alidf, iblk, alidf(:,:,iblk), a2D) if (f_albice (1:1) /= 'x') & call accum_hist_field(n_albice, iblk, albice(:,:,iblk), a2D) if (f_albsno (1:1) /= 'x') & call accum_hist_field(n_albsno, iblk, albsno(:,:,iblk), a2D) if (f_albpnd (1:1) /= 'x') & call accum_hist_field(n_albpnd, iblk, albpnd(:,:,iblk), a2D) if (f_coszen (1:1) /= 'x') & call accum_hist_field(n_coszen, iblk, coszen(:,:,iblk), a2D) if (f_flat (1:1) /= 'x') & call accum_hist_field(n_flat, iblk, flat(:,:,iblk), a2D) if (f_flat_ai(1:1) /= 'x') & call accum_hist_field(n_flat_ai,iblk, flat(:,:,iblk)*workb(:,:), a2D) if (f_fsens (1:1) /= 'x') & call accum_hist_field(n_fsens, iblk, fsens(:,:,iblk), a2D) if (f_fsens_ai(1:1)/= 'x') & call accum_hist_field(n_fsens_ai,iblk, fsens(:,:,iblk)*workb(:,:), a2D) if (f_flwup (1:1) /= 'x') & call accum_hist_field(n_flwup, iblk, flwout(:,:,iblk), a2D) if (f_flwup_ai(1:1)/= 'x') & call accum_hist_field(n_flwup_ai,iblk, flwout(:,:,iblk)*workb(:,:), a2D) if (f_evap (1:1) /= 'x') & call accum_hist_field(n_evap, iblk, evap(:,:,iblk), a2D) if (f_evap_ai(1:1) /= 'x') & call accum_hist_field(n_evap_ai,iblk, evap(:,:,iblk)*workb(:,:), a2D) if (f_Tair (1:1) /= 'x') & call accum_hist_field(n_Tair, iblk, Tair(:,:,iblk), a2D) if (f_Tref (1:1) /= 'x') & call accum_hist_field(n_Tref, iblk, Tref(:,:,iblk), a2D) if (f_Qref (1:1) /= 'x') & call accum_hist_field(n_Qref, iblk, Qref(:,:,iblk), a2D) if (f_congel (1:1) /= 'x') & call accum_hist_field(n_congel, iblk, congel(:,:,iblk), a2D) if (f_frazil (1:1) /= 'x') & call accum_hist_field(n_frazil, iblk, frazil(:,:,iblk), a2D) if (f_snoice (1:1) /= 'x') & call accum_hist_field(n_snoice, iblk, snoice(:,:,iblk), a2D) if (f_dsnow (1:1) /= 'x') & call accum_hist_field(n_dsnow, iblk, dsnow(:,:,iblk), a2D) if (f_meltt (1:1) /= 'x') & call accum_hist_field(n_meltt, iblk, meltt(:,:,iblk), a2D) if (f_melts (1:1) /= 'x') & call accum_hist_field(n_melts, iblk, melts(:,:,iblk), a2D) if (f_meltb (1:1) /= 'x') & call accum_hist_field(n_meltb, iblk, meltb(:,:,iblk), a2D) if (f_meltl (1:1) /= 'x') & call accum_hist_field(n_meltl, iblk, meltl(:,:,iblk), a2D) if (f_fresh (1:1) /= 'x') & call accum_hist_field(n_fresh, iblk, fresh(:,:,iblk), a2D) if (f_fresh_ai(1:1)/= 'x') & call accum_hist_field(n_fresh_ai,iblk, fresh_ai(:,:,iblk), a2D) if (f_fsalt (1:1) /= 'x') & call accum_hist_field(n_fsalt, iblk, fsalt(:,:,iblk), a2D) if (f_fsalt_ai(1:1)/= 'x') & call accum_hist_field(n_fsalt_ai,iblk, fsalt_ai(:,:,iblk), a2D) if (f_fhocn (1:1) /= 'x') & call accum_hist_field(n_fhocn, iblk, fhocn(:,:,iblk), a2D) if (f_fhocn_ai(1:1)/= 'x') & call accum_hist_field(n_fhocn_ai,iblk, fhocn_ai(:,:,iblk), a2D) if (f_fswthru(1:1) /= 'x') & call accum_hist_field(n_fswthru, iblk, fswthru(:,:,iblk), a2D) if (f_fswthru_ai(1:1)/= 'x') & call accum_hist_field(n_fswthru_ai,iblk, fswthru_ai(:,:,iblk), a2D) if (f_strairx(1:1) /= 'x') & call accum_hist_field(n_strairx, iblk, strairx(:,:,iblk), a2D) if (f_strairy(1:1) /= 'x') & call accum_hist_field(n_strairy, iblk, strairy(:,:,iblk), a2D) if (f_strtltx(1:1) /= 'x') & call accum_hist_field(n_strtltx, iblk, strtltx(:,:,iblk), a2D) if (f_strtlty(1:1) /= 'x') & call accum_hist_field(n_strtlty, iblk, strtlty(:,:,iblk), a2D) if (f_strcorx(1:1) /= 'x') & call accum_hist_field(n_strcorx, iblk, fm(:,:,iblk)*vvel(:,:,iblk), a2D) if (f_strcory(1:1) /= 'x') & call accum_hist_field(n_strcory, iblk,-fm(:,:,iblk)*uvel(:,:,iblk), a2D) if (f_strocnx(1:1) /= 'x') & call accum_hist_field(n_strocnx, iblk, strocnx(:,:,iblk), a2D) if (f_strocny(1:1) /= 'x') & call accum_hist_field(n_strocny, iblk, strocny(:,:,iblk), a2D) if (f_strintx(1:1) /= 'x') & call accum_hist_field(n_strintx, iblk, strintx(:,:,iblk), a2D) if (f_strinty(1:1) /= 'x') & call accum_hist_field(n_strinty, iblk, strinty(:,:,iblk), a2D) if (f_strength(1:1)/= 'x') & call accum_hist_field(n_strength,iblk, strength(:,:,iblk), a2D) ! The following fields (divu, shear, sig1, and sig2) will be smeared ! if averaged over more than a few days. ! Snapshots may be more useful (see below). ! if (f_divu (1:1) /= 'x') & ! call accum_hist_field(n_divu, iblk, divu(:,:,iblk), a2D) ! if (f_shear (1:1) /= 'x') & ! call accum_hist_field(n_shear, iblk, shear(:,:,iblk), a2D) ! if (f_sig1 (1:1) /= 'x') & ! call accum_hist_field(n_sig1, iblk, sig1(:,:,iblk), a2D) ! if (f_sig2 (1:1) /= 'x') & ! call accum_hist_field(n_sig2, iblk, sig2(:,:,iblk), a2D) ! if (f_trsig (1:1) /= 'x') & ! call accum_hist_field(n_trsig, iblk, trsig(:,:,iblk), a2D) if (f_dvidtt (1:1) /= 'x') & call accum_hist_field(n_dvidtt, iblk, dvidtt(:,:,iblk), a2D) if (f_dvidtd (1:1) /= 'x') & call accum_hist_field(n_dvidtd, iblk, dvidtd(:,:,iblk), a2D) if (f_daidtt (1:1) /= 'x') & call accum_hist_field(n_daidtt, iblk, daidtt(:,:,iblk), a2D) if (f_daidtd (1:1) /= 'x') & call accum_hist_field(n_daidtd, iblk, daidtd(:,:,iblk), a2D) if (f_dagedtt (1:1) /= 'x') & call accum_hist_field(n_dagedtt, iblk, dagedtt(:,:,iblk), a2D) if (f_dagedtd (1:1) /= 'x') & call accum_hist_field(n_dagedtd, iblk, dagedtd(:,:,iblk), a2D) if (f_fsurf_ai(1:1)/= 'x') & call accum_hist_field(n_fsurf_ai,iblk, fsurf(:,:,iblk)*workb(:,:), a2D) if (f_fcondtop_ai(1:1)/= 'x') & call accum_hist_field(n_fcondtop_ai, iblk, & fcondtop(:,:,iblk)*workb(:,:), a2D) if (f_icepresent(1:1) /= 'x') then worka(:,:) = c0 do j = jlo, jhi do i = ilo, ihi if (aice(i,j,iblk) > puny) worka(i,j) = c1 enddo enddo call accum_hist_field(n_icepresent, iblk, worka(:,:), a2D) endif ! 3D category fields if (f_aicen (1:1) /= 'x') & call accum_hist_field(n_aicen-n2D, iblk, ncat_hist, & aicen(:,:,1:ncat_hist,iblk), a3Dc) if (f_vicen (1:1) /= 'x') & call accum_hist_field(n_vicen-n2D, iblk, ncat_hist, & vicen(:,:,1:ncat_hist,iblk), a3Dc) if (f_vsnon (1:1) /= 'x') & call accum_hist_field(n_vsnon-n2D, iblk, ncat_hist, & vsnon(:,:,1:ncat_hist,iblk), a3Dc) if (f_keffn_top (1:1) /= 'x') & call accum_hist_field(n_keffn_top-n2D, iblk, ncat_hist, & keffn_top(:,:,1:ncat_hist,iblk), a3Dc) if (f_fsurfn_ai (1:1) /= 'x') & call accum_hist_field(n_fsurfn_ai-n2D, iblk, ncat_hist, & fsurfn(:,:,1:ncat_hist,iblk)*aicen_init(:,:,1:ncat_hist,iblk), a3Dc) if (f_fcondtopn_ai (1:1) /= 'x') & call accum_hist_field(n_fcondtopn_ai-n2D, iblk, ncat_hist, & fcondtopn(:,:,1:ncat_hist,iblk)*aicen_init(:,:,1:ncat_hist,iblk), a3Dc) if (f_flatn_ai (1:1) /= 'x') & call accum_hist_field(n_flatn_ai-n2D, iblk, ncat_hist, & flatn(:,:,1:ncat_hist,iblk)*aicen_init(:,:,1:ncat_hist,iblk), a3Dc) if (f_fsensn_ai (1:1) /= 'x') & call accum_hist_field(n_fsensn_ai-n2D, iblk, ncat_hist, & fsensn(:,:,1:ncat_hist,iblk)*aicen_init(:,:,1:ncat_hist,iblk), a3Dc) ! Calculate surface heat flux that causes melt (calculated by the ! atmos in HadGEM3 so needed for checking purposes) if (f_fmelttn_ai (1:1) /= 'x') & call accum_hist_field(n_fmelttn_ai-n2D, iblk, ncat_hist, & max(fsurfn(:,:,1:ncat_hist,iblk) - fcondtopn(:,:,1:ncat_hist,iblk),c0) & *aicen_init(:,:,1:ncat_hist,iblk), a3Dc) ! example for 3D field (x,y,z) ! if (f_field3dz (1:1) /= 'x') & ! call accum_hist_field(n_field3dz-n3Dccum, iblk, nzilyr, & ! field3dz(:,:,1:nzilyr,iblk), a3Dz) ! 4D category fields if (f_Tinz (1:1) /= 'x') then Tinz4d(:,:,:,:) = c0 if (ktherm == 2) then do n = 1, ncat_hist do j = jlo, jhi do i = ilo, ihi do k = 1, nzilyr Tinz4d(i,j,k,n) = temperature_mush( & trcrn(i,j,nt_qice+k-1,n,iblk), trcrn(i,j,nt_sice+k-1,n,iblk)) enddo enddo enddo enddo else do n = 1, ncat_hist do j = jlo, jhi do i = ilo, ihi do k = 1, nzilyr qn = trcrn(i,j,nt_qice+k-1,n,iblk) ! Tinz4d(i,j,k,n) = calculate_Tin_from_qin(qn,Tmlt(k)) Tmlts = -trcrn(i,j,nt_sice+k-1,n,iblk)*depressT Tinz4d(i,j,k,n) = calculate_Tin_from_qin(qn,Tmlts) enddo enddo enddo enddo endif call accum_hist_field(n_Tinz-n3Dbcum, iblk, nzilyr, ncat_hist, & Tinz4d(:,:,1:nzilyr,1:ncat_hist), a4Di) endif if (f_Sinz (1:1) /= 'x') then Sinz4d(:,:,:,:) = c0 do n = 1, ncat_hist do j = jlo, jhi do i = ilo, ihi if (vicen(i,j,n,iblk) > puny) then Sinz4d(i,j,1:nzilyr,n) = trcrn(i,j,nt_sice:nt_sice+nzilyr-1,n,iblk) endif enddo enddo enddo call accum_hist_field(n_Sinz-n3Dbcum, iblk, nzilyr, ncat_hist, & Sinz4d(:,:,1:nzilyr,1:ncat_hist), a4Di) endif if (f_Tsnz (1:1) /= 'x') then Tsnz4d(:,:,:,:) = c0 if (ktherm == 2) then do n = 1, ncat_hist do j = jlo, jhi do i = ilo, ihi do k = 1, nzslyr qn = trcrn(i,j,nt_qsno+k-1,n,iblk) Tsnz4d(i,j,k,n) = temperature_snow(trcrn(i,j,nt_qsno+k-1,n,iblk)) enddo enddo enddo enddo else do n = 1, ncat_hist do j = jlo, jhi do i = ilo, ihi do k = 1, nzslyr qn = trcrn(i,j,nt_qsno+k-1,n,iblk) Tsnz4d(i,j,k,n) = (Lfresh + qn/rhos)/cp_ice enddo enddo enddo enddo endif call accum_hist_field(n_Tsnz-n4Dicum, iblk, nzslyr, ncat_hist, & Tsnz4d(:,:,1:nzslyr,1:ncat_hist), a4Ds) endif ! Calculate aggregate surface melt flux by summing category values if (f_fmeltt_ai(1:1) /= 'x') then do ns = 1, nstreams if (n_fmeltt_ai(ns) /= 0) then worka(:,:) = c0 do j = jlo, jhi do i = ilo, ihi if (tmask(i,j,iblk)) then do n=1,ncat_hist worka(i,j) = worka(i,j) + a3Dc(i,j,n,n_fmelttn_ai(ns)-n2D,iblk) enddo ! n endif ! tmask enddo ! i enddo ! j a2D(:,:,n_fmeltt_ai(ns),iblk) = worka(:,:) endif enddo endif !--------------------------------------------------------------- ! accumulate other history output !--------------------------------------------------------------- ! mechanical redistribution call accum_hist_mechred (iblk) ! melt ponds if (tr_pond) call accum_hist_pond (iblk) ! biogeochemistry if (tr_aero .or. tr_brine .or. skl_bgc) call accum_hist_bgc (iblk) ! form drag if (formdrag) call accum_hist_drag (iblk) enddo ! iblk !$OMP END PARALLEL DO !--------------------------------------------------------------- ! Write output files at prescribed intervals !--------------------------------------------------------------- nstrm = nstreams if (write_ic) nstrm = 1 do ns = 1, nstrm if (write_history(ns) .or. write_ic) then !--------------------------------------------------------------- ! Mask out land points and convert units !--------------------------------------------------------------- ravgct = c1/avgct(ns) !$OMP PARALLEL DO PRIVATE(iblk,i,j,ilo,ihi,jlo,jhi,this_block, & !$OMP n,nn,ravgctz) 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 n = 1, num_avail_hist_fields_2D if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) then do j = jlo, jhi do i = ilo, ihi if (.not. tmask(i,j,iblk)) then ! mask out land points a2D(i,j,n,iblk) = spval else ! convert units a2D(i,j,n,iblk) = avail_hist_fields(n)%cona*a2D(i,j,n,iblk) & * ravgct + avail_hist_fields(n)%conb endif enddo ! i enddo ! j ! back out albedo/zenith angle dependence if (avail_hist_fields(n)%vname(1:6) == 'albice') then do j = jlo, jhi do i = ilo, ihi if (tmask(i,j,iblk)) then ravgctz = c0 if (albcnt(i,j,iblk,ns) > puny) & ravgctz = c1/albcnt(i,j,iblk,ns) if (f_albice (1:1) /= 'x' .and. n_albice(ns) /= 0) & a2D(i,j,n_albice(ns),iblk) = & a2D(i,j,n_albice(ns),iblk)*avgct(ns)*ravgctz if (f_albsno (1:1) /= 'x' .and. n_albsno(ns) /= 0) & a2D(i,j,n_albsno(ns),iblk) = & a2D(i,j,n_albsno(ns),iblk)*avgct(ns)*ravgctz if (f_albpnd (1:1) /= 'x' .and. n_albpnd(ns) /= 0) & a2D(i,j,n_albpnd(ns),iblk) = & a2D(i,j,n_albpnd(ns),iblk)*avgct(ns)*ravgctz endif enddo ! i enddo ! j endif if (avail_hist_fields(n)%vname(1:6) == 'albsni') then do j = jlo, jhi do i = ilo, ihi if (tmask(i,j,iblk)) then ravgctz = c0 if (albcnt(i,j,iblk,ns) > puny) & ravgctz = c1/albcnt(i,j,iblk,ns) if (f_albsni (1:1) /= 'x' .and. n_albsni(ns) /= 0) & a2D(i,j,n_albsni(ns),iblk) = & a2D(i,j,n_albsni(ns),iblk)*avgct(ns)*ravgctz endif enddo ! i enddo ! j endif endif enddo ! n do n = 1, num_avail_hist_fields_3Dc nn = n2D + n if (avail_hist_fields(nn)%vhistfreq == histfreq(ns)) then do k = 1, ncat_hist do j = jlo, jhi do i = ilo, ihi if (.not. tmask(i,j,iblk)) then ! mask out land points a3Dc(i,j,k,n,iblk) = spval else ! convert units a3Dc(i,j,k,n,iblk) = avail_hist_fields(nn)%cona*a3Dc(i,j,k,n,iblk) & * ravgct + avail_hist_fields(nn)%conb endif enddo ! i enddo ! j enddo ! k endif enddo ! n do n = 1, num_avail_hist_fields_3Dz nn = n3Dccum + n if (avail_hist_fields(nn)%vhistfreq == histfreq(ns)) then do k = 1, nzlyr do j = jlo, jhi do i = ilo, ihi if (.not. tmask(i,j,iblk)) then ! mask out land points a3Dz(i,j,k,n,iblk) = spval else ! convert units a3Dz(i,j,k,n,iblk) = avail_hist_fields(nn)%cona*a3Dz(i,j,k,n,iblk) & * ravgct + avail_hist_fields(nn)%conb endif enddo ! i enddo ! j enddo ! k endif enddo ! n do n = 1, num_avail_hist_fields_3Db nn = n3Dzcum + n if (avail_hist_fields(nn)%vhistfreq == histfreq(ns)) then do k = 1, nzblyr do j = jlo, jhi do i = ilo, ihi if (.not. tmask(i,j,iblk)) then ! mask out land points a3Db(i,j,k,n,iblk) = spval else ! convert units a3Db(i,j,k,n,iblk) = avail_hist_fields(nn)%cona*a3Db(i,j,k,n,iblk) & * ravgct + avail_hist_fields(nn)%conb endif enddo ! i enddo ! j enddo ! k endif enddo ! n do n = 1, num_avail_hist_fields_4Di nn = n3Dbcum + n if (avail_hist_fields(nn)%vhistfreq == histfreq(ns)) then do k = 1, nzilyr do ic = 1, ncat_hist do j = jlo, jhi do i = ilo, ihi if (.not. tmask(i,j,iblk)) then ! mask out land points a4Di(i,j,k,ic,n,iblk) = spval else ! convert units a4Di(i,j,k,ic,n,iblk) = avail_hist_fields(nn)%cona*a4Di(i,j,k,ic,n,iblk) & * ravgct + avail_hist_fields(nn)%conb endif enddo ! i enddo ! j enddo ! ic enddo ! k endif enddo ! n do n = 1, num_avail_hist_fields_4Ds nn = n4Dicum + n if (avail_hist_fields(nn)%vhistfreq == histfreq(ns)) then do k = 1, nzslyr do ic = 1, ncat_hist do j = jlo, jhi do i = ilo, ihi if (.not. tmask(i,j,iblk)) then ! mask out land points a4Ds(i,j,k,ic,n,iblk) = spval else ! convert units a4Ds(i,j,k,ic,n,iblk) = avail_hist_fields(nn)%cona*a4Ds(i,j,k,ic,n,iblk) & * ravgct + avail_hist_fields(nn)%conb endif enddo ! i enddo ! j enddo ! ic enddo ! k endif enddo ! n do n = 1, num_avail_hist_fields_4Db nn = n4Dscum + n if (avail_hist_fields(nn)%vhistfreq == histfreq(ns)) then do k = 1, nzblyr do ic = 1, ncat_hist do j = jlo, jhi do i = ilo, ihi if (.not. tmask(i,j,iblk)) then ! mask out land points a4Db(i,j,k,ic,n,iblk) = spval else ! convert units a4Db(i,j,k,ic,n,iblk) = avail_hist_fields(nn)%cona*a4Db(i,j,k,ic,n,iblk) & * ravgct + avail_hist_fields(nn)%conb endif enddo ! i enddo ! j enddo ! ic enddo ! k endif enddo ! n !--------------------------------------------------------------- ! snapshots !--------------------------------------------------------------- ! compute sig1 and sig2 call principal_stress (nx_block, ny_block, & stressp_1 (:,:,iblk), & stressm_1 (:,:,iblk), & stress12_1(:,:,iblk), & prs_sig (:,:,iblk), & sig1 (:,:,iblk), & sig2 (:,:,iblk)) do j = jlo, jhi do i = ilo, ihi if (.not. tmask(i,j,iblk)) then ! mask out land points if (n_divu (ns) /= 0) a2D(i,j,n_divu(ns), iblk) = spval if (n_shear (ns) /= 0) a2D(i,j,n_shear(ns), iblk) = spval if (n_sig1 (ns) /= 0) a2D(i,j,n_sig1(ns), iblk) = spval if (n_sig2 (ns) /= 0) a2D(i,j,n_sig2(ns), iblk) = spval if (n_mlt_onset(ns) /= 0) a2D(i,j,n_mlt_onset(ns),iblk) = spval if (n_frz_onset(ns) /= 0) a2D(i,j,n_frz_onset(ns),iblk) = spval if (n_hisnap (ns) /= 0) a2D(i,j,n_hisnap(ns), iblk) = spval if (n_aisnap (ns) /= 0) a2D(i,j,n_aisnap(ns), iblk) = spval if (n_trsig (ns) /= 0) a2D(i,j,n_trsig(ns), iblk) = spval if (n_iage (ns) /= 0) a2D(i,j,n_iage(ns), iblk) = spval if (n_FY (ns) /= 0) a2D(i,j,n_FY(ns), iblk) = spval if (n_a11 (ns) /= 0) a2D(i,j,n_a11(ns), iblk) = spval if (n_a12 (ns) /= 0) a2D(i,j,n_a12(ns), iblk) = spval if (n_e11 (ns) /= 0) a2D(i,j,n_e11(ns), iblk) = spval if (n_e12 (ns) /= 0) a2D(i,j,n_e12(ns), iblk) = spval if (n_e22 (ns) /= 0) a2D(i,j,n_e22(ns), iblk) = spval if (n_s11 (ns) /= 0) a2D(i,j,n_s11(ns), iblk) = spval if (n_s12 (ns) /= 0) a2D(i,j,n_s12(ns), iblk) = spval if (n_s22 (ns) /= 0) a2D(i,j,n_s22(ns), iblk) = spval if (n_yieldstress11 (ns) /= 0) a2D(i,j,n_yieldstress11(ns),iblk) = spval if (n_yieldstress12 (ns) /= 0) a2D(i,j,n_yieldstress12(ns),iblk) = spval if (n_yieldstress22 (ns) /= 0) a2D(i,j,n_yieldstress22(ns),iblk) = spval else if (n_divu (ns) /= 0) a2D(i,j,n_divu(ns),iblk) = & divu (i,j,iblk)*avail_hist_fields(n_divu(ns))%cona if (n_shear (ns) /= 0) a2D(i,j,n_shear(ns),iblk) = & shear(i,j,iblk)*avail_hist_fields(n_shear(ns))%cona if (n_sig1 (ns) /= 0) a2D(i,j,n_sig1(ns),iblk) = & sig1 (i,j,iblk)*avail_hist_fields(n_sig1(ns))%cona if (n_sig2 (ns) /= 0) a2D(i,j,n_sig2(ns),iblk) = & sig2 (i,j,iblk)*avail_hist_fields(n_sig2(ns))%cona if (n_mlt_onset(ns) /= 0) a2D(i,j,n_mlt_onset(ns),iblk) = & mlt_onset(i,j,iblk) if (n_frz_onset(ns) /= 0) a2D(i,j,n_frz_onset(ns),iblk) = & frz_onset(i,j,iblk) if (n_hisnap (ns) /= 0) a2D(i,j,n_hisnap(ns),iblk) = & vice(i,j,iblk) if (n_aisnap (ns) /= 0) a2D(i,j,n_aisnap(ns),iblk) = & aice(i,j,iblk) if (kdyn == 2) then ! for EAP dynamics different time of output if (n_trsig (ns) /= 0) a2D(i,j,n_trsig(ns),iblk ) = & prs_sig(i,j,iblk) else if (n_trsig (ns) /= 0) a2D(i,j,n_trsig(ns),iblk ) = & p25*(stressp_1(i,j,iblk) & + stressp_2(i,j,iblk) & + stressp_3(i,j,iblk) & + stressp_4(i,j,iblk)) endif if (n_iage (ns) /= 0) a2D(i,j,n_iage(ns),iblk) = & trcr(i,j,nt_iage,iblk)*avail_hist_fields(n_iage(ns))%cona if (n_FY (ns) /= 0) a2D(i,j,n_FY(ns),iblk) = & trcr(i,j,nt_FY,iblk)*avail_hist_fields(n_FY(ns))%cona if (n_a11 (ns) /= 0) a2D(i,j,n_a11(ns),iblk) = & a11 (i,j,iblk)*avail_hist_fields(n_a11(ns))%cona if (n_a12 (ns) /= 0) a2D(i,j,n_a12(ns),iblk) = & a12 (i,j,iblk)*avail_hist_fields(n_a12(ns))%cona if (n_e11 (ns) /= 0) a2D(i,j,n_e11(ns),iblk) = & e11 (i,j,iblk)*avail_hist_fields(n_e11(ns))%cona if (n_e12 (ns) /= 0) a2D(i,j,n_e12(ns),iblk) = & e12 (i,j,iblk)*avail_hist_fields(n_e12(ns))%cona if (n_e22 (ns) /= 0) a2D(i,j,n_e22(ns),iblk) = & e22 (i,j,iblk)*avail_hist_fields(n_e22(ns))%cona if (n_s11 (ns) /= 0) a2D(i,j,n_s11(ns),iblk) = & s11 (i,j,iblk)*avail_hist_fields(n_s11(ns))%cona if (n_s12 (ns) /= 0) a2D(i,j,n_s12(ns),iblk) = & s12 (i,j,iblk)*avail_hist_fields(n_s12(ns))%cona if (n_s22 (ns) /= 0) a2D(i,j,n_s22(ns),iblk) = & s22 (i,j,iblk)*avail_hist_fields(n_s22(ns))%cona if (n_yieldstress11 (ns) /= 0) a2D(i,j,n_yieldstress11(ns),iblk) = & yieldstress11 (i,j,iblk)*avail_hist_fields(n_yieldstress11(ns))%cona if (n_yieldstress12 (ns) /= 0) a2D(i,j,n_yieldstress12(ns),iblk) = & yieldstress12 (i,j,iblk)*avail_hist_fields(n_yieldstress12(ns))%cona if (n_yieldstress22 (ns) /= 0) a2D(i,j,n_yieldstress22(ns),iblk) = & yieldstress22 (i,j,iblk)*avail_hist_fields(n_yieldstress22(ns))%cona endif enddo ! i enddo ! j enddo ! iblk !$OMP END PARALLEL DO time_end(ns) = time/int(secday) !sgl time_end(ns) = real(time_end(ns),kind=real_kind) !--------------------------------------------------------------- ! write file !--------------------------------------------------------------- call ice_timer_start(timer_readwrite) ! reading/writing call ice_write_hist (ns) call ice_timer_stop(timer_readwrite) ! reading/writing !--------------------------------------------------------------- ! reset to zero !------------------------------------------------------------ if (write_ic) then if (allocated(a2D)) a2D (:,:,:,:) = c0 if (allocated(a3Dc)) a3Dc(:,:,:,:,:) = c0 if (allocated(a3Dz)) a3Dz(:,:,:,:,:) = c0 if (allocated(a3Db)) a3Db(:,:,:,:,:) = c0 if (allocated(a4Di)) a4Di(:,:,:,:,:,:) = c0 if (allocated(a4Ds)) a4Ds(:,:,:,:,:,:) = c0 if (allocated(a4Db)) a4Db(:,:,:,:,:,:) = c0 avgct(:) = c0 albcnt(:,:,:,:) = c0 write_ic = .false. ! write initial condition once at most else avgct(ns) = c0 albcnt(:,:,:,ns) = c0 endif ! if (write_history(ns)) albcnt(:,:,:,ns) = c0 do n = 1,n2D if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) a2D(:,:,n,:) = c0 enddo do n = n2D + 1, n3Dccum nn = n - n2D if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) a3Dc(:,:,:,nn,:) = c0 enddo do n = n3Dccum + 1, n3Dzcum nn = n - n3Dccum if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) a3Dz(:,:,:,nn,:) = c0 enddo do n = n3Dzcum + 1, n3Dbcum nn = n - n3Dzcum if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) a3Db(:,:,:,nn,:) = c0 enddo do n = n3Dbcum + 1, n4Dicum nn = n - n3Dbcum if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) a4Di(:,:,:,:,nn,:) = c0 enddo do n = n4Dicum + 1, n4Dscum nn = n - n4Dicum if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) a4Ds(:,:,:,:,nn,:) = c0 enddo do n = n4Dscum + 1, n4Dbcum nn = n - n4Dscum if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) a4Db(:,:,:,:,nn,:) = c0 enddo endif ! write_history or write_ic enddo ! nstreams !$OMP PARALLEL DO PRIVATE(iblk,i,j,ilo,ihi,jlo,jhi,this_block) 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 if (new_year) then do j=jlo,jhi do i=ilo,ihi ! reset NH Jan 1 if (lmask_n(i,j,iblk)) mlt_onset(i,j,iblk) = c0 ! reset SH Jan 1 if (lmask_s(i,j,iblk)) frz_onset(i,j,iblk) = c0 enddo enddo endif ! new_year if ( (month .eq. 7) .and. new_month ) then do j=jlo,jhi do i=ilo,ihi ! reset SH Jul 1 if (lmask_s(i,j,iblk)) mlt_onset(i,j,iblk) = c0 ! reset NH Jul 1 if (lmask_n(i,j,iblk)) frz_onset(i,j,iblk) = c0 enddo enddo endif ! 1st of July enddo ! iblk !$OMP END PARALLEL DO end subroutine accum_hist !======================================================================= end module ice_history !=======================================================================