diff --git a/cicecore/cicedyn/analysis/ice_history.F90 b/cicecore/cicedyn/analysis/ice_history.F90 index 99b04afb9..72939f31a 100644 --- a/cicecore/cicedyn/analysis/ice_history.F90 +++ b/cicecore/cicedyn/analysis/ice_history.F90 @@ -273,7 +273,7 @@ subroutine init_hist (dt) 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 + histfreq(ns) == 'y' .or. histfreq(ns) == 'n') then nstreams = nstreams + 1 if (ns >= 2) then if (histfreq(ns-1) == 'x') then @@ -337,90 +337,6 @@ subroutine init_hist (dt) ! to prevent array-out-of-bounds when aggregating if (f_fmeltt_ai(1:1) /= 'x') f_fmelttn_ai = f_fmeltt_ai - ! Turn on all CMIP fields in one go. - - if (f_CMIP(1:1) /= 'x') then - f_sithick = 'mxxxx' - f_sisnthick = 'mxxxx' - f_siage = 'mxxxx' - f_sitemptop = 'mxxxx' - f_sitempsnic = 'mxxxx' - f_sitempbot = 'mxxxx' - f_sispeed = 'mxxxx' - f_siu = 'mxxxx' - f_siv = 'mxxxx' - f_sidmasstranx = 'mxxxx' - f_sidmasstrany = 'mxxxx' - f_sistrxdtop = 'mxxxx' - f_sistrydtop = 'mxxxx' - f_sistrxubot = 'mxxxx' - f_sistryubot = 'mxxxx' - f_sicompstren = 'mxxxx' - f_sialb = 'mxxxx' - f_sihc = 'mxxxx' - f_sisnhc = 'mxxxx' - f_sidconcth = 'mxxxx' - f_sidconcdyn = 'mxxxx' - f_sidmassth = 'mxxxx' - f_sidmassdyn = 'mxxxx' - f_sidmassgrowthwat = 'mxxxx' - f_sidmassgrowthbot = 'mxxxx' - f_sidmasssi = 'mxxxx' - f_sidmassevapsubl = 'mxxxx' - f_sndmasssubl = 'mxxxx' - f_sidmassmelttop = 'mxxxx' - f_sidmassmeltbot = 'mxxxx' - f_sidmasslat = 'mxxxx' - f_sndmasssnf = 'mxxxx' - f_sndmassmelt = 'mxxxx' - f_sndmassdyn = 'mxxxx' - f_siflswdtop = 'mxxxx' - f_siflswutop = 'mxxxx' - f_siflswdbot = 'mxxxx' - f_sifllwdtop = 'mxxxx' - f_sifllwutop = 'mxxxx' - f_siflsenstop = 'mxxxx' - f_siflsensupbot = 'mxxxx' - f_sifllatstop = 'mxxxx' - f_siflcondtop = 'mxxxx' - f_siflcondbot = 'mxxxx' - f_sipr = 'mxxxx' - f_sifb = 'mxxxx' - f_siflsaltbot = 'mxxxx' - f_siflfwbot = 'mxxxx' - f_siflfwdrain = 'mxxxx' - f_siforcetiltx = 'mxxxx' - f_siforcetilty = 'mxxxx' - f_siforcecoriolx = 'mxxxx' - f_siforcecorioly = 'mxxxx' - f_siforceintstrx = 'mxxxx' - f_siforceintstry = 'mxxxx' - f_sidragtop = 'mxxxx' - f_sistreave = 'mxxxx' - f_sistremax = 'mxxxx' - f_sirdgthick = 'mxxxx' - f_siitdconc = 'mxxxx' - f_siitdthick = 'mxxxx' - f_siitdsnthick = 'mxxxx' - f_aicen = 'mxxxx' - endif - - if (f_CMIP(2:2) == 'd') then - f_icepresent = f_CMIP - f_aice = f_CMIP - f_sithick = f_CMIP - f_sisnthick = f_CMIP - f_sitemptop = f_CMIP - f_siu = f_CMIP - f_siv = f_CMIP - f_sispeed = f_CMIP - f_sndmasssubl = f_CMIP - f_sndmasssnf = f_CMIP - f_sndmassmelt = f_CMIP - f_sndmassdyn = f_CMIP - f_sidmasssi = f_CMIP - endif - if (grid_ice == 'CD' .or. grid_ice == 'C') then f_uvelE = f_uvel f_vvelE = f_vvel @@ -494,12 +410,14 @@ subroutine init_hist (dt) call broadcast_scalar (f_NFSD, master_task) ! call broadcast_scalar (f_example, master_task) + call broadcast_scalar (f_CMIP, master_task) call broadcast_scalar (f_hi, master_task) call broadcast_scalar (f_hs, master_task) call broadcast_scalar (f_snowfrac, master_task) call broadcast_scalar (f_snowfracn, master_task) call broadcast_scalar (f_Tsfc, master_task) call broadcast_scalar (f_aice, master_task) + call broadcast_scalar (f_aice_init, master_task) call broadcast_scalar (f_uvel, master_task) call broadcast_scalar (f_vvel, master_task) call broadcast_scalar (f_icespd, master_task) @@ -629,8 +547,15 @@ subroutine init_hist (dt) 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_sitimefrac, master_task) call broadcast_scalar (f_sithick, master_task) call broadcast_scalar (f_siage, master_task) + call broadcast_scalar (f_siconc, master_task) + call broadcast_scalar (f_sisnconc, master_task) + call broadcast_scalar (f_sisnmass, master_task) + call broadcast_scalar (f_sivol, master_task) + call broadcast_scalar (f_simass, master_task) + call broadcast_scalar (f_sisaltmass, master_task) call broadcast_scalar (f_sisnthick, master_task) call broadcast_scalar (f_sitemptop, master_task) call broadcast_scalar (f_sitempsnic, master_task) @@ -645,8 +570,8 @@ subroutine init_hist (dt) call broadcast_scalar (f_sistryubot, master_task) call broadcast_scalar (f_sicompstren, master_task) call broadcast_scalar (f_sispeed, master_task) - call broadcast_scalar (f_sidir, master_task) - call broadcast_scalar (f_sialb, master_task) +! call broadcast_scalar (f_sidir, master_task) + call broadcast_scalar (f_sisali, master_task) call broadcast_scalar (f_sihc, master_task) call broadcast_scalar (f_sisnhc, master_task) call broadcast_scalar (f_sidconcth, master_task) @@ -655,23 +580,24 @@ subroutine init_hist (dt) call broadcast_scalar (f_sidmassdyn, master_task) call broadcast_scalar (f_sidmassgrowthwat, master_task) call broadcast_scalar (f_sidmassgrowthbot, master_task) - call broadcast_scalar (f_sidmasssi, master_task) + call broadcast_scalar (f_sidmassgrowthsi, master_task) call broadcast_scalar (f_sidmassevapsubl, master_task) - call broadcast_scalar (f_sndmasssubl, master_task) + call broadcast_scalar (f_sisndmasssubl, master_task) call broadcast_scalar (f_sidmassmelttop, master_task) call broadcast_scalar (f_sidmassmeltbot, master_task) - call broadcast_scalar (f_sidmasslat, master_task) - call broadcast_scalar (f_sndmasssnf, master_task) - call broadcast_scalar (f_sndmassmelt, master_task) - call broadcast_scalar (f_sndmassdyn, master_task) + call broadcast_scalar (f_sidmassmeltlat, master_task) + call broadcast_scalar (f_sisndmasssnf, master_task) + call broadcast_scalar (f_sisndmassmelt, master_task) + call broadcast_scalar (f_sisndmassdyn, master_task) + call broadcast_scalar (f_sisndmasssi, master_task) call broadcast_scalar (f_siflswdtop, master_task) call broadcast_scalar (f_siflswutop, master_task) call broadcast_scalar (f_siflswdbot, master_task) call broadcast_scalar (f_sifllwdtop, master_task) call broadcast_scalar (f_sifllwutop, master_task) call broadcast_scalar (f_siflsenstop, master_task) - call broadcast_scalar (f_siflsensupbot, master_task) - call broadcast_scalar (f_sifllatstop, master_task) + call broadcast_scalar (f_siflsensbot, master_task) + call broadcast_scalar (f_sifllattop, master_task) call broadcast_scalar (f_siflcondtop, master_task) call broadcast_scalar (f_siflcondbot, master_task) call broadcast_scalar (f_sipr, master_task) @@ -686,12 +612,15 @@ subroutine init_hist (dt) call broadcast_scalar (f_siforceintstrx, master_task) call broadcast_scalar (f_siforceintstry, master_task) call broadcast_scalar (f_siitdconc, master_task) + call broadcast_scalar (f_siitdsnconc, master_task) call broadcast_scalar (f_siitdthick, master_task) call broadcast_scalar (f_siitdsnthick, master_task) call broadcast_scalar (f_sidragtop, master_task) - call broadcast_scalar (f_sistreave, master_task) - call broadcast_scalar (f_sistremax, master_task) - call broadcast_scalar (f_sirdgthick, master_task) + call broadcast_scalar (f_sidragbot, master_task) + call broadcast_scalar (f_sistressave, master_task) + call broadcast_scalar (f_sistressmax, master_task) + call broadcast_scalar (f_sidivvel, master_task) + call broadcast_scalar (f_sishearvel, master_task) call broadcast_scalar (f_aicen, master_task) call broadcast_scalar (f_vicen, master_task) @@ -728,6 +657,11 @@ subroutine init_hist (dt) call broadcast_scalar (f_yieldstress12, master_task) call broadcast_scalar (f_yieldstress22, master_task) + if (f_CMIP(1:1) /= 'x') then + if (my_task == master_task) write(nu_diag,*) subname, & + 'WARNING: f_CMIP has been deprecated. Please use the set_nml.cmip namelist option' + endif + ! 2D variables do ns1 = 1, nstreams if (histfreq(ns1) /= 'x') then @@ -739,189 +673,194 @@ subroutine init_hist (dt) ! 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, & + 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, & + 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_snowfrac,"snowfrac","1",tstr2D, tcstr, & - "grid cell mean snow fraction", & - "snow fraction per unit grid cell area", c1, c0, & + "snow fraction of sea ice", & + "none", c1, c0, & ns1, f_snowfrac) - call define_hist_field(n_Tsfc,"Tsfc","C",tstr2D, tcstr, & - "snow/ice surface temperature", & - "averaged with Tf if no ice is present", c1, c0, & + 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, & + call define_hist_field(n_aice,"aice","1",tstr2D, tcstr, & + "ice area (aggregate)", & + "none", c1, c0, & ns1, f_aice) - call define_hist_field(n_uvelE,"uvelE","m/s",estr2D, ecstr, & - "ice velocity (x)", & - "positive is x direction on E grid", c1, c0, & + call define_hist_field(n_aice_init,"aice_init","1",tstr2D, tcstr, & + "ice area (aggregate) initial", & + "none", c1, c0, & + ns1, f_aice_init) + + call define_hist_field(n_uvelE,"uvelE","m/s",estr2D, ecstr, & + "ice velocity (x)", & + "positive is x direction on E grid", c1, c0, & ns1, f_uvelE) - call define_hist_field(n_vvelE,"vvelE","m/s",estr2D, ecstr, & - "ice velocity (y)", & - "positive is y direction on E grid", c1, c0, & + call define_hist_field(n_vvelE,"vvelE","m/s",estr2D, ecstr, & + "ice velocity (y)", & + "positive is y direction on E grid", c1, c0, & ns1, f_vvelE) call define_hist_field(n_icespdE,"icespdE","m/s",estr2D, ecstr, & - "sea ice speed", & + "sea ice speed", & "vector magnitude on E grid", c1, c0, & ns1, f_icespdE) call define_hist_field(n_icedirE,"icedirE","deg",estr2D, ecstr, & - "sea ice direction", & + "sea ice direction", & "vector direction - coming from on E grid", c1, c0, & ns1, f_icedirE) - call define_hist_field(n_uvelN,"uvelN","m/s",nstr2D, ncstr, & - "ice velocity (x)", & - "positive is x direction on N grid", c1, c0, & + call define_hist_field(n_uvelN,"uvelN","m/s",nstr2D, ncstr, & + "ice velocity (x)", & + "positive is x direction on N grid", c1, c0, & ns1, f_uvelN) - call define_hist_field(n_vvelN,"vvelN","m/s",nstr2D, ncstr, & - "ice velocity (y)", & - "positive is y direction on N grid", c1, c0, & + call define_hist_field(n_vvelN,"vvelN","m/s",nstr2D, ncstr, & + "ice velocity (y)", & + "positive is y direction on N grid", c1, c0, & ns1, f_vvelN) call define_hist_field(n_icespdN,"icespdN","m/s",nstr2D, ncstr, & - "sea ice speed", & + "sea ice speed", & "vector magnitude on N grid", c1, c0, & ns1, f_icespdN) call define_hist_field(n_icedirN,"icedirN","deg",nstr2D, ncstr, & - "sea ice direction", & + "sea ice direction", & "vector direction - coming from on N grid", c1, c0, & ns1, f_icedirN) - call define_hist_field(n_uvel,"uvel","m/s",ustr2D, ucstr, & - "ice velocity (x)", & - "positive is x direction on U grid", c1, c0, & + 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, & + 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_icespd,"icespd","m/s",ustr2D, ucstr, & - "sea ice speed", & + "sea ice speed", & "vector magnitude", c1, c0, & ns1, f_icespd) call define_hist_field(n_icedir,"icedir","deg",ustr2D, ucstr, & - "sea ice direction", & + "sea ice direction", & "vector direction - coming from", c1, c0, & ns1, f_icedir) - call define_hist_field(n_uatm,"uatm","m/s",str2D_gau, cstr_gau, & - "atm velocity (x)", & - "positive is x direction on U grid", c1, c0, & + call define_hist_field(n_uatm,"uatm","m/s",str2D_gau, cstr_gau, & + "atm velocity (x)", & + "positive is x direction on U grid", c1, c0, & ns1, f_uatm) - call define_hist_field(n_vatm,"vatm","m/s",str2D_gav, cstr_gav, & - "atm velocity (y)", & - "positive is y direction on U grid", c1, c0, & + call define_hist_field(n_vatm,"vatm","m/s",str2D_gav, cstr_gav, & + "atm velocity (y)", & + "positive is y direction on U grid", c1, c0, & ns1, f_vatm) call define_hist_field(n_atmspd,"atmspd","m/s",str2D_gau, cstr_gau, & - "atmosphere wind speed", & - "vector magnitude", c1, c0, & + "atmosphere wind speed", & + "vector magnitude", c1, c0, & ns1, f_atmspd) call define_hist_field(n_atmdir,"atmdir","deg",str2D_gau, cstr_gau, & - "atmosphere wind direction", & - "vector direction - coming from", c1, c0, & + "atmosphere wind direction", & + "vector direction - coming from", c1, c0, & ns1, f_atmdir) - call define_hist_field(n_sice,"sice","ppt",tstr2D, tcstr, & - "bulk ice salinity", & - "none", c1, c0, & + call define_hist_field(n_sice,"sice","ppt",tstr2D, tcstr, & + "bulk ice salinity", & + "none", c1, c0, & ns1, f_sice) call define_hist_field(n_fswup,"fswup","W/m^2",tstr2D, tcstr, & "upward solar flux", & - "positive upward", c1, c0, & + "positive upward", c1, c0, & ns1, f_fswup) call define_hist_field(n_fswdn,"fswdn","W/m^2",tstr2D, tcstr, & - "down solar flux", & - "positive downward", c1, c0, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "sea surface salinity", & + "none", c1, c0, & ns1, f_sss) call define_hist_field(n_uocn,"uocn","m/s",str2D_gou, cstr_gou, & - "ocean current (x)", & - "positive is x direction on U grid", c1, c0, & + "ocean current (x)", & + "positive is x direction on U grid", c1, c0, & ns1, f_uocn) call define_hist_field(n_vocn,"vocn","m/s",str2D_gov, cstr_gov, & - "ocean current (y)", & - "positive is y direction on U grid", c1, c0, & + "ocean current (y)", & + "positive is y direction on U grid", c1, c0, & ns1, f_vocn) call define_hist_field(n_ocnspd,"ocnspd","m/s",str2D_gou, cstr_gou, & - "ocean current speed", & - "vector magnitude", c1, c0, & + "ocean current speed", & + "vector magnitude", c1, c0, & ns1, f_ocnspd) call define_hist_field(n_ocndir,"ocndir","deg",str2D_gou, cstr_gou, & - "ocean current direction", & - "vector direction - going to", c1, c0, & + "ocean current direction", & + "vector direction - going to", c1, c0, & ns1, f_ocndir) 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, & + "freeze/melt potential", & + "if >0, new ice forms; if <0, ice melts", c1, c0, & ns1, f_frzmlt) call define_hist_field(n_fswfac,"scale_factor","1",tstr2D, tcstr, & - "shortwave scaling factor", & - "ratio of netsw new:old", c1, c0, & + "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, & @@ -930,194 +869,194 @@ subroutine init_hist (dt) 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, & + "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, & + "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", & - "averaged for coszen>0, weighted by aice", c100, c0, & + "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, & + "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, & + "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, & + "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, & + "near IR diffuse albedo", & + "scaled (divided) by aice", c100, c0, & ns1, f_alidf) call define_hist_field(n_alvdr_ai,"alvdr_ai","%",tstr2D, tcstr, & - "visible direct albedo", & - " ", c100, c0, & + "visible direct albedo", & + " ", c100, c0, & ns1, f_alvdr_ai) call define_hist_field(n_alidr_ai,"alidr_ai","%",tstr2D, tcstr, & - "near IR direct albedo", & - " ", c100, c0, & + "near IR direct albedo", & + " ", c100, c0, & ns1, f_alidr_ai) call define_hist_field(n_alvdf_ai,"alvdf_ai","%",tstr2D, tcstr, & - "visible diffuse albedo", & - " ", c100, c0, & + "visible diffuse albedo", & + " ", c100, c0, & ns1, f_alvdf_ai) call define_hist_field(n_alidf_ai,"alidf_ai","%",tstr2D, tcstr, & - "near IR diffuse albedo", & - " ", c100, c0, & + "near IR diffuse albedo", & + " ", c100, c0, & ns1, f_alidf_ai) call define_hist_field(n_albice,"albice","%",tstr2D, tcstr, & - "bare ice albedo", & - "averaged for coszen>0, weighted by aice", c100, c0, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "air temperature", & + "none", c1, -Tffresh, & ns1, f_Tair) call define_hist_field(n_Tref,"Tref","C",tstr2D, tcstr, & - "2m reference temperature", & - "none", c1, -Tffresh, & + "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, & + "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, & + "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, & + "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, & + "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", & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "salt flux ice to ocean", & + "weighted by ice area", c1, c0, & ns1, f_fsalt_ai) call define_hist_field(n_fbot,"fbot","W/m^2",tstr2D, tcstr, & @@ -1126,223 +1065,223 @@ subroutine init_hist (dt) ns1, f_fbot) 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, & + "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 (fhocn_ai)", & - "weighted by ice area", c1, c0, & + "heat flux ice to ocean (fhocn_ai)", & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "internal ice stress (y)", & + "positive is y direction on U grid", c1, c0, & ns1, f_strinty) - call define_hist_field(n_taubx,"taubx","N/m^2",ustr2D, ucstr, & - "seabed (basal) stress (x)", & - "positive is x direction on U grid", c1, c0, & + call define_hist_field(n_taubx,"taubx","N/m^2",ustr2D, ucstr, & + "seabed (basal) stress (x)", & + "positive is x direction on U grid", c1, c0, & ns1, f_taubx) - call define_hist_field(n_tauby,"tauby","N/m^2",ustr2D, ucstr, & - "seabed (basal) stress (y)", & - "positive is y direction on U grid", c1, c0, & + call define_hist_field(n_tauby,"tauby","N/m^2",ustr2D, ucstr, & + "seabed (basal) stress (y)", & + "positive is y direction on U grid", c1, c0, & ns1, f_tauby) call define_hist_field(n_strairxN,"strairxN","N/m^2",nstr2D, ncstr, & - "atm/ice stress (x)", & - "positive is x direction on N grid", c1, c0, & + "atm/ice stress (x)", & + "positive is x direction on N grid", c1, c0, & ns1, f_strairxN) call define_hist_field(n_strairyN,"strairyN","N/m^2",nstr2D, ncstr, & - "atm/ice stress (y)", & - "positive is y direction on N grid", c1, c0, & + "atm/ice stress (y)", & + "positive is y direction on N grid", c1, c0, & ns1, f_strairyN) call define_hist_field(n_strairxE,"strairxE","N/m^2",estr2D, ecstr, & - "atm/ice stress (x)", & - "positive is x direction on E grid", c1, c0, & + "atm/ice stress (x)", & + "positive is x direction on E grid", c1, c0, & ns1, f_strairxE) call define_hist_field(n_strairyE,"strairyE","N/m^2",estr2D, ecstr, & - "atm/ice stress (y)", & - "positive is y direction on E grid", c1, c0, & + "atm/ice stress (y)", & + "positive is y direction on E grid", c1, c0, & ns1, f_strairyE) call define_hist_field(n_strtltxN,"strtltxN","N/m^2",nstr2D, ncstr, & - "sea sfc tilt stress (x)", & - "positive is x direction on N grid", c1, c0, & + "sea sfc tilt stress (x)", & + "positive is x direction on N grid", c1, c0, & ns1, f_strtltxN) call define_hist_field(n_strtltyN,"strtltyN","N/m^2",nstr2D, ncstr, & - "sea sfc tilt stress (y)", & - "positive is y direction on N grid", c1, c0, & + "sea sfc tilt stress (y)", & + "positive is y direction on N grid", c1, c0, & ns1, f_strtltyN) call define_hist_field(n_strtltxE,"strtltxE","N/m^2",estr2D, ecstr, & - "sea sfc tilt stress (x)", & - "positive is x direction on E grid", c1, c0, & + "sea sfc tilt stress (x)", & + "positive is x direction on E grid", c1, c0, & ns1, f_strtltxE) call define_hist_field(n_strtltyE,"strtltyE","N/m^2",estr2D, ecstr, & - "sea sfc tilt stress (y)", & - "positive is y direction on E grid", c1, c0, & + "sea sfc tilt stress (y)", & + "positive is y direction on E grid", c1, c0, & ns1, f_strtltyE) call define_hist_field(n_strcorxN,"strcorxN","N/m^2",nstr2D, ncstr, & - "coriolis stress (x)", & - "positive is x direction on N grid", c1, c0, & + "coriolis stress (x)", & + "positive is x direction on N grid", c1, c0, & ns1, f_strcorxN) call define_hist_field(n_strcoryN,"strcoryN","N/m^2",nstr2D, ncstr, & - "coriolis stress (y)", & - "positive is y direction on N grid", c1, c0, & + "coriolis stress (y)", & + "positive is y direction on N grid", c1, c0, & ns1, f_strcoryN) call define_hist_field(n_strcorxE,"strcorxE","N/m^2",estr2D, ecstr, & - "coriolis stress (x)", & - "positive is x direction on E grid", c1, c0, & + "coriolis stress (x)", & + "positive is x direction on E grid", c1, c0, & ns1, f_strcorxE) call define_hist_field(n_strcoryE,"strcoryE","N/m^2",estr2D, ecstr, & - "coriolis stress (y)", & - "positive is y direction on E grid", c1, c0, & + "coriolis stress (y)", & + "positive is y direction on E grid", c1, c0, & ns1, f_strcoryE) call define_hist_field(n_strocnxN,"strocnxN","N/m^2",nstr2D, ncstr, & - "ocean/ice stress (x)", & - "positive is x direction on N grid", c1, c0, & + "ocean/ice stress (x)", & + "positive is x direction on N grid", c1, c0, & ns1, f_strocnxN) call define_hist_field(n_strocnyN,"strocnyN","N/m^2",nstr2D, ncstr, & - "ocean/ice stress (y)", & - "positive is y direction on N grid", c1, c0, & + "ocean/ice stress (y)", & + "positive is y direction on N grid", c1, c0, & ns1, f_strocnyN) call define_hist_field(n_strocnxE,"strocnxE","N/m^2",estr2D, ecstr, & - "ocean/ice stress (x)", & - "positive is x direction on E grid", c1, c0, & + "ocean/ice stress (x)", & + "positive is x direction on E grid", c1, c0, & ns1, f_strocnxE) call define_hist_field(n_strocnyE,"strocnyE","N/m^2",estr2D, ecstr, & - "ocean/ice stress (y)", & - "positive is y direction on E grid", c1, c0, & + "ocean/ice stress (y)", & + "positive is y direction on E grid", c1, c0, & ns1, f_strocnyE) call define_hist_field(n_strintxN,"strintxN","N/m^2",nstr2D, ncstr, & - "internal ice stress (x)", & - "positive is x direction on N grid", c1, c0, & + "internal ice stress (x)", & + "positive is x direction on N grid", c1, c0, & ns1, f_strintxN) call define_hist_field(n_strintyN,"strintyN","N/m^2",nstr2D, ncstr, & - "internal ice stress (y)", & - "positive is y direction on N grid", c1, c0, & + "internal ice stress (y)", & + "positive is y direction on N grid", c1, c0, & ns1, f_strintyN) call define_hist_field(n_strintxE,"strintxE","N/m^2",estr2D, ecstr, & - "internal ice stress (x)", & - "positive is x direction on E grid", c1, c0, & + "internal ice stress (x)", & + "positive is x direction on E grid", c1, c0, & ns1, f_strintxE) call define_hist_field(n_strintyE,"strintyE","N/m^2",estr2D, ecstr, & - "internal ice stress (y)", & - "positive is y direction on E grid", c1, c0, & + "internal ice stress (y)", & + "positive is y direction on E grid", c1, c0, & ns1, f_strintyE) - call define_hist_field(n_taubxN,"taubxN","N/m^2",nstr2D, ncstr, & + call define_hist_field(n_taubxN,"taubxN","N/m^2",nstr2D, ncstr, & "seabed (basal) stress (x)", & "positive is x direction on N grid", c1, c0, & ns1, f_taubxN) - call define_hist_field(n_taubyN,"taubyN","N/m^2",nstr2D, ncstr, & + call define_hist_field(n_taubyN,"taubyN","N/m^2",nstr2D, ncstr, & "seabed (basal) stress (y)", & "positive is y direction on N grid", c1, c0, & ns1, f_taubyN) - call define_hist_field(n_taubxE,"taubxE","N/m^2",estr2D, ecstr, & + call define_hist_field(n_taubxE,"taubxE","N/m^2",estr2D, ecstr, & "seabed (basal) stress (x)", & "positive is x direction on E grid", c1, c0, & ns1, f_taubxE) - call define_hist_field(n_taubyE,"taubyE","N/m^2",estr2D, ecstr, & + call define_hist_field(n_taubyE,"taubyE","N/m^2",estr2D, ecstr, & "seabed (basal) stress (y)", & "positive is y direction on E grid", c1, c0, & ns1, f_taubyE) call define_hist_field(n_strength,"strength","N/m",tstr2D, tcstr, & - "compressive ice strength", & - "none", c1, c0, & + "compressive ice strength", & + "none", c1, c0, & ns1, f_strength) call define_hist_field(n_divu,"divu","%/day",tstr2D, tcstr, & - "strain rate (divergence)", & - "divu is instantaneous, on T grid", secday*c100, c0, & + "strain rate (divergence)", & + "divu is instantaneous, on T grid", secday*c100, c0, & ns1, f_divu) call define_hist_field(n_shear,"shear","%/day",tstr2D, tcstr, & - "strain rate (shear)", & - "shear is instantaneous, on T grid", secday*c100, c0, & + "strain rate (shear)", & + "shear is instantaneous, on T grid", secday*c100, c0, & ns1, f_shear) call define_hist_field(n_vort,"vort","%/day",tstr2D, tcstr, & - "strain rate (vorticity)", & - "vort is instantaneous, on T grid", secday*c100, c0, & + "strain rate (vorticity)", & + "vort is instantaneous, on T grid", secday*c100, c0, & ns1, f_vort) select case (grid_ice) @@ -1357,68 +1296,68 @@ subroutine init_hist (dt) end select call define_hist_field(n_sig1,"sig1","1",gridstr2d, gridstr, & - "norm. principal stress 1", & - "sig1 is instantaneous" // trim(description), c1, c0, & + "norm. principal stress 1", & + "sig1 is instantaneous" // trim(description), c1, c0, & ns1, f_sig1) call define_hist_field(n_sig2,"sig2","1",gridstr2d, gridstr, & - "norm. principal stress 2", & - "sig2 is instantaneous" // trim(description), c1, c0, & + "norm. principal stress 2", & + "sig2 is instantaneous" // trim(description), c1, c0, & ns1, f_sig2) call define_hist_field(n_sigP,"sigP","N/m",gridstr2d, gridstr, & - "ice pressure", & - "sigP is instantaneous" // trim(description), c1, c0, & + "ice pressure", & + "sigP is instantaneous" // trim(description), c1, c0, & ns1, f_sigP) call define_hist_field(n_dvidtt,"dvidtt","cm/day",tstr2D, tcstr, & - "volume tendency thermo", & - "none", mps_to_cmpdy, c0, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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", & + 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", & + 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, & + "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, & + "ice area snapshot", & + "none", c1, c0, & ns1, f_aisnap) call define_hist_field(n_trsig,"trsig","N/m",tstr2D, tcstr, & @@ -1427,396 +1366,492 @@ subroutine init_hist (dt) 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, & + "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", & + 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", & + 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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "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, & + "first-year ice area", & + "weighted by ice area", c1, c0, & ns1, f_FY) - ! CMIP 2D variables + ! CMIP 2D variables (for "intensive" variables per Notz et al 2016 definition, + ! that is a weighted time average when ice is present) + ! Use avg_ice_present = 'init', 'final', 'pond', or 'ridge' to divide by + ! sum(aice), sum(apond), or sum(ardg) over time + ! aice is at the start of the timestep ('init') or the end of the timestep ('final') + ! avg_ice_present = 'none' produces a time average including zeroes when ice is not present + + call define_hist_field(n_siage,"siage","s",tstr2D, tcstr, & + "age of sea ice", & + "age of sea ice since its formation in open water", c1, c0, & + ns1, f_siage, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_sicompstren,"sicompstren","N m-1",tstr2D, tcstr, & + "compressive sea-ice strength", & + "computed strength of the ice pack", c1, c0, & + ns1, f_sicompstren, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_sidragbot,"sidragbot","1",tstr2D, tcstr, & + "ocean drag coefficient", & + "drag coefficient that is used to calculate the oceanic momentum drag on sea ice", & + c1, c0, & + ns1, f_sidragbot, avg_ice_present='init', mask_ice_free_points=.true.) + + call define_hist_field(n_sidragtop,"sidragtop","1",tstr2D, tcstr, & + "atmospheric drag coefficient", & + "drag coefficient that is used to calculate the atmospheric momentum drag on sea ice", & + c1, c0, & + ns1, f_sidragtop, avg_ice_present='init', mask_ice_free_points=.true.) - call define_hist_field(n_sithick,"sithick","m",tstr2D, tcstr, & - "sea ice thickness", & - "volume divided by area", c1, c0, & - ns1, f_sithick, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_siage,"siage","s",tstr2D, tcstr, & - "sea ice age", & - "none", c1, c0, & - ns1, f_siage, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sisnthick,"sisnthick","m",tstr2D, tcstr, & - "sea ice snow thickness", & - "snow volume divided by area", c1, c0, & - ns1, f_sisnthick, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sitemptop,"sitemptop","K",tstr2D, tcstr, & - "sea ice surface temperature", & - "none", c1, Tffresh, & - ns1, f_sitemptop, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sitempsnic,"sitempsnic","K",tstr2D, tcstr, & - "snow ice interface temperature", & - "surface temperature when no snow present", c1, Tffresh, & - ns1, f_sitempsnic, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sitempbot,"sitempbot","K",tstr2D, tcstr, & - "sea ice bottom temperature", & - "none", c1, Tffresh, & - ns1, f_sitempbot, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_siu,"siu","m/s",ustr2D, ucstr, & - "ice x velocity component", & - "none", c1, c0, & - ns1, f_siu, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_siv,"siv","m/s",ustr2D, ucstr, & - "ice y velocity component", & - "none", c1, c0, & - ns1, f_siv, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sidmasstranx,"sidmasstranx","kg/s",ustr2D, ucstr, & - "x component of snow and sea ice mass transport", & - "none", c1, c0, & - ns1, f_sidmasstranx) - - call define_hist_field(n_sidmasstrany,"sidmasstrany","kg/s",ustr2D, ucstr, & - "y component of snow and sea ice mass transport", & - "none", c1, c0, & - ns1, f_sidmasstrany) - - call define_hist_field(n_sistrxdtop,"sistrxdtop","N m-2",ustr2D, ucstr, & - "x component of atmospheric stress on sea ice", & - "none", c1, c0, & - ns1, f_sistrxdtop, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sistrydtop,"sistrydtop","N m-2",ustr2D, ucstr, & - "y component of atmospheric stress on sea ice", & - "none", c1, c0, & - ns1, f_sistrydtop, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sistrxubot,"sistrxubot","N m-2",ustr2D, ucstr, & - "x component of ocean stress on sea ice", & - "none", c1, c0, & - ns1, f_sistrxubot, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sistryubot,"sistryubot","N m-2",ustr2D, ucstr, & - "y component of ocean stress on sea ice", & - "none", c1, c0, & - ns1, f_sistryubot, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sicompstren,"sicompstren","N m-1",tstr2D, tcstr, & - "compressive sea ice strength", & - "none", c1, c0, & - ns1, f_sicompstren, avg_ice_present=.true., mask_ice_free_points=.true.) + call define_hist_field(n_sifb,"sifb","m",tstr2D, tcstr, & + "sea-ice freeboard", & + "mean height of sea-ice surface above sea level", & + c1, c0, & + ns1, f_sifb, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_siflcondbot,"siflcondbot","W/m2",tstr2D, tcstr, & + "net conductive heat flux in sea-ice at the base", & + "conductive heat flux from the centre of the lowermost vertical sea-ice grid box", & + c1, c0, & + ns1, f_siflcondbot, avg_ice_present='init', mask_ice_free_points=.true.) + + call define_hist_field(n_siflcondtop,"siflcondtop","W/m2",tstr2D, tcstr, & + "net conductive heat flux in sea-ice at the surface", & + "conductive heat flux from the centre of the uppermost vertical sea-ice grid box", & + c1, c0, & + ns1, f_siflcondtop, avg_ice_present='init', mask_ice_free_points=.true.) - call define_hist_field(n_sispeed,"sispeed","m/s",ustr2D, ucstr, & - "ice speed", & - "none", c1, c0, & - ns1, f_sispeed, avg_ice_present=.true., mask_ice_free_points=.true.) + call define_hist_field(n_siflfwbot,"siflfwbot","kg m-2 s-1",tstr2D, tcstr, & + "freshwater flux from sea-ice", & + "total flux of fresh water between ocean and sea ice", c1, c0, & + ns1, f_siflfwbot, avg_ice_present='final', mask_ice_free_points=.true.) - call define_hist_field(n_sidir,"sidir","deg",ustr2D, ucstr, & - "ice direction", & - "vector direction - going to", c1, c0, & - ns1, f_sidir) + call define_hist_field(n_siflfwdrain,"siflfwdrain","kg m-2 s-1",tstr2D, tcstr, & + "freshwater flux from sea-ice surface", & + "total flux of fresh water from sea-ice surface into underlying ocean", & + c1, c0, & + ns1, f_siflfwdrain, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_sifllattop,"sifllattop","W/m2",tstr2D, tcstr, & + "net latent heat flux over sea ice", & + "positive downward", c1, c0, & + ns1, f_sifllattop, avg_ice_present='init', mask_ice_free_points=.true.) + + call define_hist_field(n_sifllwdtop,"sifllwdtop","W/m2",tstr2D, tcstr, & + "downwelling longwave flux over sea ice", & + "downwelling longwave flux from the atmosphere to the sea-ice surface ", & + c1, c0, & + ns1, f_sifllwdtop, avg_ice_present='init', mask_ice_free_points=.true.) - call define_hist_field(n_sialb,"sialb","1",tstr2D, tcstr, & - "sea ice albedo", & - "none", c1, c0, & - ns1, f_sialb, avg_ice_present=.true., mask_ice_free_points=.true.) + call define_hist_field(n_sifllwutop,"sifllwutop","W/m2",tstr2D, tcstr, & + "upwelling longwave flux over sea ice", & + "upward longwave flux from the sea-ice surface to the atmosphere", & + c1, c0, & + ns1, f_sifllwutop, avg_ice_present='init', mask_ice_free_points=.true.) - call define_hist_field(n_sihc,"sihc","J m-2",tstr2D, tcstr, & - "sea ice heat content", & - "none", c1, c0, & - ns1, f_sihc) + call define_hist_field(n_siflsaltbot,"siflsaltbot","kg m-2 s-1",tstr2D, tcstr, & + "salt flux from sea ice", & + "total flux of salt between ocean and sea ice", c1, c0, & + ns1, f_siflsaltbot, avg_ice_present='final', mask_ice_free_points=.true.) - call define_hist_field(n_sisnhc,"sisnhc","J m-2",tstr2D, tcstr, & - "snow heat content", & - "none", c1, c0, & - ns1, f_sisnhc) - - call define_hist_field(n_sidconcth,"sidconcth","1/s",tstr2D, tcstr, & - "sea ice area change from thermodynamics", & - "none", c1, c0, & - ns1, f_sidconcth) - - call define_hist_field(n_sidconcdyn,"sidconcdyn","1/s",tstr2D, tcstr, & - "sea ice area change from dynamics", & - "none", c1, c0, & - ns1, f_sidconcdyn) - - call define_hist_field(n_sidmassth,"sidmassth","kg m-2 s-1",tstr2D, tcstr, & - "sea ice mass change from thermodynamics", & - "none", c1, c0, & - ns1, f_sidmassth) - - call define_hist_field(n_sidmassdyn,"sidmassdyn","kg m-2 s-1",tstr2D, tcstr, & - "sea ice mass change from dynamics", & - "none", c1, c0, & - ns1, f_sidmassdyn) - - call define_hist_field(n_sidmassgrowthwat,"sidmassgrowthwat","kg m-2 s-1",tstr2D, tcstr, & - "sea ice mass change from frazil", & - "none", c1, c0, & - ns1, f_sidmassgrowthwat) - - call define_hist_field(n_sidmassgrowthbot,"sidmassgrowthbot","kg m-2 s-1",tstr2D, tcstr, & - "sea ice mass change from basal growth", & - "none", c1, c0, & - ns1, f_sidmassgrowthbot) - - call define_hist_field(n_sidmasssi,"sidmasssi","kg m-2 s-1",tstr2D, tcstr, & - "sea ice mass change from snow-ice formation", & - "none", c1, c0, & - ns1, f_sidmasssi) - - call define_hist_field(n_sidmassevapsubl,"sidmassevapsubl","kg m-2 s-1",tstr2D, tcstr, & - "sea ice mass change from evaporation and sublimation", & - "none", c1, c0, & - ns1, f_sidmassevapsubl) - - call define_hist_field(n_sndmasssubl,"sndmasssubl","kg m-2 s-1",tstr2D, tcstr, & - "snow mass change from evaporation and sublimation", & - "none", c1, c0, & - ns1, f_sndmasssubl) - - call define_hist_field(n_sidmassmelttop,"sidmassmelttop","kg m-2 s-1",tstr2D, tcstr, & - "sea ice mass change top melt", & - "none", c1, c0, & - ns1, f_sidmassmelttop) - - call define_hist_field(n_sidmassmeltbot,"sidmassmeltbot","kg m-2 s-1",tstr2D, tcstr, & - "sea ice mass change bottom melt", & - "none", c1, c0, & - ns1, f_sidmassmeltbot) - - call define_hist_field(n_sidmasslat,"sidmasslat","kg m-2 s-1",tstr2D, tcstr, & - "sea ice mass change lateral melt", & - "none", c1, c0, & - ns1, f_sidmasslat) - - call define_hist_field(n_sndmasssnf,"sndmasssnf","kg m-2 s-1",tstr2D, tcstr, & - "snow mass change from snow fall", & - "none", c1, c0, & - ns1, f_sndmasssnf) - - call define_hist_field(n_sndmassmelt,"sndmassmelt","kg m-2 s-1",tstr2D, tcstr, & - "snow mass change from snow melt", & - "none", c1, c0, & - ns1, f_sndmassmelt) - - call define_hist_field(n_sndmassdyn,"sndmassdyn","kg m-2 s-1",tstr2D, tcstr, & - "snow mass change from dynamics ridging", & - "none", c1, c0, & - ns1, f_sndmassdyn) - - call define_hist_field(n_siflswdtop,"siflswdtop","W/m2",tstr2D, tcstr, & - "down shortwave flux over sea ice", & - "positive downward", c1, c0, & - ns1, f_siflswdtop, avg_ice_present=.true., mask_ice_free_points=.true.) + call define_hist_field(n_siflsensbot,"siflsensbot","W/m2",tstr2D, tcstr, & + "net upward sensible heat flux under sea ice", & + "net sensible heat flux under sea ice from or to the ocean", c1, c0, & + ns1, f_siflsensbot, avg_ice_present='init', mask_ice_free_points=.true.) - call define_hist_field(n_siflswutop,"siflswutop","W/m2",tstr2D, tcstr, & - "upward shortwave flux over sea ice", & - "positive downward", c1, c0, & - ns1, f_siflswutop, avg_ice_present=.true., mask_ice_free_points=.true.) + call define_hist_field(n_siflsenstop,"siflsenstop","W/m2",tstr2D, tcstr, & + "net downward sensible heat flux over sea ice", & + "positive downward", c1, c0, & + ns1, f_siflsenstop, avg_ice_present='init', mask_ice_free_points=.true.) call define_hist_field(n_siflswdbot,"siflswdbot","W/m2",tstr2D, tcstr, & - "down shortwave flux at bottom of ice", & - "positive downward", c1, c0, & - ns1, f_siflswdbot, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sifllwdtop,"sifllwdtop","W/m2",tstr2D, tcstr, & - "down longwave flux over sea ice", & - "positive downward", c1, c0, & - ns1, f_sifllwdtop, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sifllwutop,"sifllwutop","W/m2",tstr2D, tcstr, & - "upward longwave flux over sea ice", & - "positive downward", c1, c0, & - ns1, f_sifllwutop, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_siflsenstop,"siflsenstop","W/m2",tstr2D, tcstr, & - "sensible heat flux over sea ice", & - "positive downward", c1, c0, & - ns1, f_siflsenstop, avg_ice_present=.true., mask_ice_free_points=.true.) + "downwelling shortwave flux underneath sea ice", & + "amount of shortwave radiation that penetrates the sea ice", & + c1, c0, & + ns1, f_siflswdbot, avg_ice_present='init', mask_ice_free_points=.true.) + + call define_hist_field(n_siflswdtop,"siflswdtop","W/m2",tstr2D, tcstr, & + "downwelling shortwave flux over sea ice", & + "downwelling shortwave flux from the atmosphere to the sea-ice surface", & + c1, c0, & + ns1, f_siflswdtop, avg_ice_present='init', mask_ice_free_points=.true.) + + call define_hist_field(n_siflswutop,"siflswutop","W/m2",tstr2D, tcstr, & + "upwelling shortwave flux over sea ice", & + "upward shortwave flux from the sea-ice surface to the atmosphere", & + c1, c0, & + ns1, f_siflswutop, avg_ice_present='init', mask_ice_free_points=.true.) - call define_hist_field(n_siflsensupbot,"siflsensupbot","W/m2",tstr2D, tcstr, & - "sensible heat flux at bottom of sea ice", & - "positive downward", c1, c0, & - ns1, f_siflsensupbot, avg_ice_present=.true., mask_ice_free_points=.true.) + call define_hist_field(n_siforcecoriolx,"siforcecoriolx","N m-2",tstr2D, tcstr, & + "coriolis force term in force balance (x-component)", & + "x-component of the force on sea ice caused by the Coriolis force", & + c1, c0, & + ns1, f_siforcecoriolx, avg_ice_present='final', mask_ice_free_points=.true.) - call define_hist_field(n_sifllatstop,"sifllatstop","W/m2",tstr2D, tcstr, & - "latent heat flux over sea ice", & - "positive downward", c1, c0, & - ns1, f_sifllatstop, avg_ice_present=.true., mask_ice_free_points=.true.) + call define_hist_field(n_siforcecorioly,"siforcecorioly","N m-2",tstr2D, tcstr, & + "coriolis force term in force balance (y-component)", & + "y-component of the force on sea ice caused by the Coriolis force", & + c1, c0, & + ns1, f_siforcecorioly, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_siforceintstrx,"siforceintstrx","N m-2",tstr2D, tcstr, & + "internal stress force term in force balance (x-component)", & + "x-component of the force on sea ice caused by internal stress (divergence of sigma)", & + c1, c0, & + ns1, f_siforceintstrx, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_siforceintstry,"siforceintstry","N m-2",tstr2D, tcstr, & + "internal stress force term in force balance (y-component)", & + "y-component of the force on sea ice caused by internal stress (divergence of sigma)", & + c1, c0, & + ns1, f_siforceintstry, avg_ice_present='final', mask_ice_free_points=.true.) - call define_hist_field(n_siflcondtop,"siflcondtop","W/m2",tstr2D, tcstr, & - "conductive heat flux at top of sea ice", & - "positive downward", c1, c0, & - ns1, f_siflcondtop, avg_ice_present=.true., mask_ice_free_points=.true.) + call define_hist_field(n_siforcetiltx,"siforcetiltx","N m-2",tstr2D, tcstr, & + "sea-surface tilt term in force balance (x-component)", & + "x-component of the force on sea ice caused by sea-surface tilt", & + c1, c0, & + ns1, f_siforcetiltx, avg_ice_present='init', mask_ice_free_points=.true.) - call define_hist_field(n_siflcondbot,"siflcondbot","W/m2",tstr2D, tcstr, & - "conductive heat flux at bottom of sea ice", & - "positive downward", c1, c0, & - ns1, f_siflcondbot, avg_ice_present=.true., mask_ice_free_points=.true.) + call define_hist_field(n_siforcetilty,"siforcetilty","N m-2",tstr2D, tcstr, & + "sea surface tilt term in force balance (y-component)", & + "y-component of the force on sea ice caused by sea-surface tilt", & + c1, c0, & + ns1, f_siforcetilty, avg_ice_present='init', mask_ice_free_points=.true.) call define_hist_field(n_sipr,"sipr","kg m-2 s-1",tstr2D, tcstr, & - "rainfall over sea ice", & - "none", c1, c0, & - ns1, f_sipr, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sifb,"sifb","m",tstr2D, tcstr, & - "sea ice freeboard above sea level", & - "none", c1, c0, & - ns1, f_sifb, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_siflsaltbot,"siflsaltbot","kg m-2 s-1",tstr2D, tcstr, & - "salt flux from sea ice", & - "positive downward", c1, c0, & - ns1, f_siflsaltbot, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_siflfwbot,"siflfwbot","kg m-2 s-1",tstr2D, tcstr, & - "fresh water flux from sea ice", & - "positive downward", c1, c0, & - ns1, f_siflfwbot, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_siflfwdrain,"siflfwdrain","kg m-2 s-1",tstr2D, tcstr, & - "fresh water drainage through sea ice", & - "positive downward", c1, c0, & - ns1, f_siflfwdrain, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sidragtop,"sidragtop","1",tstr2D, tcstr, & - "atmospheric drag over sea ice", & - "none", c1, c0, & - ns1, f_sidragtop, avg_ice_present=.true., mask_ice_free_points=.true.) + "rainfall rate over sea ice", & + "mass of liquid precipitation falling onto sea ice", c1, c0, & + ns1, f_sipr, avg_ice_present='init', mask_ice_free_points=.true.) - call define_hist_field(n_sirdgthick,"sirdgthick","m",tstr2D, tcstr, & - "sea ice ridge thickness", & - "vrdg divided by ardg", c1, c0, & - ns1, f_sirdgthick, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_siforcetiltx,"siforcetiltx","N m-2",tstr2D, tcstr, & - "sea surface tilt term", & - "none", c1, c0, & - ns1, f_siforcetiltx, avg_ice_present=.true., mask_ice_free_points=.true.) + call define_hist_field(n_sisali,"sisali","ppt",tstr2D, tcstr, & + "sea-ice salinity", & + "mean sea-ice salinity of all sea ice in grid cell", & + c1, c0, & + ns1, f_sisali, avg_ice_present='final', mask_ice_free_points=.true.) - call define_hist_field(n_siforcetilty,"siforcetilty","N m-2",tstr2D, tcstr, & - "sea surface tile term", & - "none", c1, c0, & - ns1, f_siforcetilty, avg_ice_present=.true., mask_ice_free_points=.true.) + call define_hist_field(n_sispeed,"sispeed","m/s",ustr2D, ucstr, & + "sea-ice speed", & + "speed of ice (i.e. mean absolute velocity)", c1, c0, & + ns1, f_sispeed, avg_ice_present='final', mask_ice_free_points=.true.) + +! sidir is not actually in the CMIP7 table +! call define_hist_field(n_sidir,"sidir","deg",ustr2D, ucstr, & +! "ice direction", & +! "vector direction - going to", c1, c0, & +! ns1, f_sidir) + + call define_hist_field(n_sisnthick,"sisnthick","m",tstr2D, tcstr, & + "snow thickness", & + "snow volume divided by sea-ice area", c1, c0, & + ns1, f_sisnthick, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_sistrxdtop,"sistrxdtop","N m-2",ustr2D, ucstr, & + "x-component of atmospheric stress on sea ice", & + "x-component of the atmospheric stress on the surface of sea ice", & + c1, c0, & + ns1, f_sistrxdtop, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_sistrydtop,"sistrydtop","N m-2",ustr2D, ucstr, & + "y-component of atmospheric stress on sea ice", & + "y-component of the atmospheric stress on the surface of sea ice", & + c1, c0, & + ns1, f_sistrydtop, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_sistrxubot,"sistrxubot","N m-2",ustr2D, ucstr, & + "x-component of ocean stress on sea ice", & + "x-component of the ocean stress on the sea ice bottom ", c1, c0, & + ns1, f_sistrxubot, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_sistryubot,"sistryubot","N m-2",ustr2D, ucstr, & + "y-component of ocean stress on sea ice", & + "y-component of the ocean stress on the sea ice bottom ", c1, c0, & + ns1, f_sistryubot, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_sitempbot,"sitempbot","K",tstr2D, tcstr, & + "temperature at ice-ocean interface", & + "mean temperature at the base of the sea ice", c1, Tffresh, & + ns1, f_sitempbot, avg_ice_present='init', mask_ice_free_points=.true.) + + call define_hist_field(n_sitempsnic,"sitempsnic","K",tstr2D, tcstr, & + "temperature at snow-ice interface", & + "surface temperature when no snow present", c1, Tffresh, & + ns1, f_sitempsnic, avg_ice_present='init', mask_ice_free_points=.true.) + + call define_hist_field(n_sitemptop,"sitemptop","K",tstr2D, tcstr, & + "surface temperature of sea ice", & + "mean surface temperature of the sea-ice covered part of the grid cell", & + c1, Tffresh, & + ns1, f_sitemptop, avg_ice_present='final', mask_ice_free_points=.true.) - call define_hist_field(n_siforcecoriolx,"siforcecoriolx","N m-2",tstr2D, tcstr, & - "coriolis term", & - "none", c1, c0, & - ns1, f_siforcecoriolx, avg_ice_present=.true., mask_ice_free_points=.true.) + call define_hist_field(n_sithick,"sithick","m",tstr2D, tcstr, & + "sea-ice thickness", & + "volume divided by sea-ice area", c1, c0, & + ns1, f_sithick, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_siu,"siu","m/s",ustr2D, ucstr, & + "x-component of sea-ice velocity", & + "on native model grid", c1, c0, & + ns1, f_siu, avg_ice_present='final', mask_ice_free_points=.true.) + + call define_hist_field(n_siv,"siv","m/s",ustr2D, ucstr, & + "y-component of sea-ice velocity", & + "on native model grid", c1, c0, & + ns1, f_siv, avg_ice_present='final', mask_ice_free_points=.true.) + + ! CMIP 2D extensive variables + + call define_hist_field(n_siconc,"siconc","%",tstr2D, tcstr, & + "sea-ice area percentage (ocean grid)", & + "percentage of a given grid cell that is covered by sea ice on the ocean grid", & + c100, c0, & + ns1, f_siconc, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidconcdyn,"sidconcdyn","1/s",tstr2D, tcstr, & + "sea-ice area fraction tendency due to dynamics", & + "total rate of change in sea-ice area fraction through dynamics-related processes", & + c1, c0, & + ns1, f_sidconcdyn, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidconcth,"sidconcth","1/s",tstr2D, tcstr, & + "sea-ice area fraction tendency due to thermodynamics", & + "total rate of change in sea-ice area fraction through thermodynamic processes", & + c1, c0, & + ns1, f_sidconcth, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidmassdyn,"sidmassdyn","kg m-2 s-1",tstr2D, tcstr, & + "sea-ice mass change from dynamics", & + "total rate of change in sea-ice mass through dynamics-related processes", & + c1, c0, & + ns1, f_sidmassdyn, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidmassevapsubl,"sidmassevapsubl","kg m-2 s-1",tstr2D, tcstr, & + "sea-ice mass change through evaporation and sublimation", & + "rate of change of sea-ice mass through evaporation and sublimation", & + c1, c0, & + ns1, f_sidmassevapsubl, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidmassgrowthbot,"sidmassgrowthbot","kg m-2 s-1",tstr2D, tcstr, & + "sea-ice mass change through basal growth", & + "rate of change of sea-ice mass due to vertical growth of existing sea ice at its base", & + c1, c0, & + ns1, f_sidmassgrowthbot, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidmassgrowthsi,"sidmassgrowthsi","kg m-2 s-1",tstr2D, tcstr, & + "sea-ice mass change from snow-to-ice conversion", & + "rate of change of sea-ice mass due to transformation of snow to sea ice", & + c1, c0, & + ns1, f_sidmassgrowthsi, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidmassgrowthwat,"sidmassgrowthwat","kg m-2 s-1",tstr2D, tcstr, & + "sea-ice mass change through growth in supercooled open water (frazil)", & + "always positive or zero", c1, c0, & + ns1, f_sidmassgrowthwat, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidmassmeltbot,"sidmassmeltbot","kg m-2 s-1",tstr2D, tcstr, & + "sea-ice mass change through bottom melting", & + "rate of change of sea-ice mass through melting/dissolution at the ice bottom", & + c1, c0, & + ns1, f_sidmassmeltbot, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidmassmeltlat,"sidmassmeltlat","kg m-2 s-1",tstr2D, tcstr, & + "sea-ice mass change through lateral melting", & + "rate of change of sea-ice mass through lateral melting/dissolution", c1, c0, & + ns1, f_sidmassmeltlat, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidmassmelttop,"sidmassmelttop","kg m-2 s-1",tstr2D, tcstr, & + "sea-ice mass change through surface melting", & + "rate of change of sea-ice mass through melting at the ice surface", c1, c0, & + ns1, f_sidmassmelttop, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidmassth,"sidmassth","kg m-2 s-1",tstr2D, tcstr, & + "sea-ice mass change from thermodynamics", & + "total rate of change in sea-ice mass from thermodynamic processes", & + c1, c0, & + ns1, f_sidmassth, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidmasstranx,"sidmasstranx","kg/s",ustr2D, ucstr, & + "x-component of sea-ice mass transport", & + "x-component of the sea-ice drift-induced transport of snow and sea ice mass", & + c1, c0, & + ns1, f_sidmasstranx, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sidmasstrany,"sidmasstrany","kg/s",ustr2D, ucstr, & + "y-component of sea-ice mass transport", & + "y-component of the sea-ice drift-induced transport of snow and sea ice mass", & + c1, c0, & + ns1, f_sidmasstrany, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sihc,"sihc","J m-2",tstr2D, tcstr, & + "sea-ice heat content", & + "heat content of all ice in grid cell", c1, c0, & + ns1, f_sihc, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_simass,"simass","kg m-2",tstr2D, tcstr, & + "sea-ice mass per area", & + "total mass of sea ice", c1, c0, & + ns1, f_simass, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sisaltmass,"sisaltmass","kg m-2",tstr2D, tcstr, & + "mass of salt in sea-ice per area", & + "total mass of all salt in sea ice", c1, c0, & + ns1, f_sisaltmass, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sisnconc,"sisnconc","%",tstr2D, tcstr, & + "snow area percentage", & + "percentage of the sea-ice surface that is covered by snow", & + c100, c0, & + ns1, f_sisnconc, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sisndmassdyn,"sisndmassdyn","kg m-2 s-1",tstr2D, tcstr, & + "snow mass rate of change through advection by sea-ice dynamics", & + "rate of change of snow mass due to sea ice dynamics", c1, c0, & + ns1, f_sisndmassdyn, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sisndmassmelt,"sisndmassmelt","kg m-2 s-1",tstr2D, tcstr, & + "snow mass rate of change through melt", & + "always negative or zero", c1, c0, & + ns1, f_sisndmassmelt, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sisndmasssi,"sisndmasssi","kg m-2 s-1",tstr2D, tcstr, & + "snow mass rate of change through snow-to-ice conversion", & + "always negative or zero", c1, c0, & + ns1, f_sisndmasssi, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sisndmasssnf,"sisndmasssnf","kg m-2 s-1",tstr2D, tcstr, & + "snow mass change through snowfall", & + "rate of change of snow mass due to solid precipitation falling onto sea ice", & + c1, c0, & + ns1, f_sisndmasssnf, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sisndmasssubl,"sisndmasssubl","kg m-2 s-1",tstr2D, tcstr, & + "snow mass rate of change through evaporation or sublimation", & + "none", c1, c0, & + ns1, f_sisndmasssubl, avg_ice_present='none', mask_ice_free_points=.false.) - call define_hist_field(n_siforcecorioly,"siforcecorioly","N m-2",tstr2D, tcstr, & - "coriolis term", & - "none", c1, c0, & - ns1, f_siforcecorioly, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_siforceintstrx,"siforceintstrx","N m-2",tstr2D, tcstr, & - "internal stress term", & - "none", c1, c0, & - ns1, f_siforceintstrx, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_siforceintstry,"siforceintstry","N m-2",tstr2D, tcstr, & - "internal stress term", & - "none", c1, c0, & - ns1, f_siforceintstry, avg_ice_present=.true., mask_ice_free_points=.true.) - - call define_hist_field(n_sistreave,"sistreave","N m-1",ustr2D, ucstr, & - "average normal stress", & - "sistreave is instantaneous", c1, c0, & - ns1, f_sistreave) - - call define_hist_field(n_sistremax,"sistremax","N m-1",ustr2D, ucstr, & - "maximum shear stress", & - "sistremax is instantaneous", c1, c0, & - ns1, f_sistremax) + call define_hist_field(n_sisnhc,"sisnhc","J m-2",tstr2D, tcstr, & + "snow heat content", & + "heat content of all snow in grid cell", c1, c0, & + ns1, f_sisnhc, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sisnmass,"sisnmass","kg m-2",tstr2D, tcstr, & + "snow mass per area", & + "total mass of snow on sea ice", c1, c0, & + ns1, f_sisnmass, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sitimefrac,"sitimefrac","1",tstr2D, tcstr, & + "fraction of time steps with sea ice", & + "averaging period during which sea ice is present (siconc > 0) in a grid cell", & + c1, c0, & + ns1, f_icepresent, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sivol,"sivol","m",tstr2D, tcstr, & + "sea-ice volume per area", & + "total volume of sea ice divided by grid-cell area (equivalent thickness)", & + c1, c0, & + ns1, f_sivol, avg_ice_present='none', mask_ice_free_points=.false.) + + ! CMIP 2D instantaneous fields + + call define_hist_field(n_sidivvel,"sidivvel","s-1",ustr2D, ucstr, & + "divergence of the sea-ice velocity field", & + "sidivvel is instantaneous", c1, c0, & + ns1, f_sidivvel, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sishearvel,"sishearvel","s-1",ustr2D, ucstr, & + "maximum shear of the sea-ice velocity field", & + "sishearvel is instantaneous", c1, c0, & + ns1, f_sishearvel, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sistressave,"sistressave","N m-1",ustr2D, ucstr, & + "average normal stress in sea ice", & + "sistressave is instantaneous", c1, c0, & + ns1, f_sistressave, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_sistressmax,"sistressmax","N m-1",ustr2D, ucstr, & + "maximum shear stress in sea ice", & + "sistressmax is instantaneous", c1, c0, & + ns1, f_sistressmax, avg_ice_present='none', mask_ice_free_points=.false.) endif ! if (histfreq(ns1) /= 'x') then enddo ! ns1 @@ -1890,17 +1925,29 @@ subroutine init_hist (dt) ns1, f_keffn_top) ! CMIP 3D - call define_hist_field(n_siitdconc,"siitdconc","1",tstr3Dc, tcstr, & - "ice area, categories","none", c1, c0, & - ns1, f_siitdconc) - - call define_hist_field(n_siitdthick,"siitdthick","m",tstr3Dc, tcstr, & - "ice thickness, categories","none", c1, c0, & - ns1, f_siitdthick, avg_ice_present=.true.) + call define_hist_field(n_siitdconc,"siitdconc","%",tstr3Dc, tcstr, & + "sea-ice area percentages in ice thickness categories", & + "percentage of grid cell covered by each ice thickness category", & + c100, c0, & + ns1, f_siitdconc, avg_ice_present='none', mask_ice_free_points=.false.) + + call define_hist_field(n_siitdsnconc,"siitdsnconc","%",tstr3Dc, tcstr, & + "snow area percentages in ice thickness categories", & + "percentage of grid cell covered by snow in each ice thickness category", & + c100, c0, & + ns1, f_siitdsnconc, avg_ice_present='none', mask_ice_free_points=.false.) call define_hist_field(n_siitdsnthick,"siitdsnthick","m",tstr3Dc, tcstr, & - "snow thickness, categories","none", c1, c0, & - ns1, f_siitdsnthick, avg_ice_present=.true.) + "snow thickness in ice thickness categories", & + "actual thickness of snow in each ice thickness category", & + c1, c0, & + ns1, f_siitdsnthick, avg_ice_present='final', mask_ice_free_points=.false.) + + call define_hist_field(n_siitdthick,"siitdthick","m",tstr3Dc, tcstr, & + "sea-ice thickness in ice thickness categories", & + "actual (floe) thickness of sea ice in each thickness category", & + c1, c0, & + ns1, f_siitdthick, avg_ice_present='final', mask_ice_free_points=.false.) endif ! if (histfreq(ns1) /= 'x') then enddo ! ns1 @@ -2164,7 +2211,7 @@ subroutine accum_hist (dt) 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, scale_factor, fswabs, fswthru, alvdr, alvdf, alidr, alidf, & + frzmlt_init, scale_factor, fswabs, fswup, fswthru, alvdr, alvdf, alidr, alidf, & albice, albsno, albpnd, coszen, flat, fsens, flwout, evap, evaps, evapi, & Tair, Tref, Qref, congel, frazil, frazil_diag, snoice, dsnow, & melts, meltb, meltt, meltl, fresh, fsalt, fresh_ai, fsalt_ai, & @@ -2183,17 +2230,17 @@ subroutine accum_hist (dt) mlt_onset, frz_onset, dagedtt, dagedtd, fswint_ai, keffn_top, & snowfrac, alvdr_ai, alvdf_ai, alidr_ai, alidf_ai, update_ocn_f, & cpl_frazil - use ice_arrays_column, only: snowfracn, Cdn_atm + use ice_arrays_column, only: snowfracn, Cdn_atm, Cdn_ocn 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_mechred, only: n_alvl, n_ardg, f_ardg use ice_history_pond, only: accum_hist_pond + use ice_history_pond, only: n_apond_ai, f_apond_ai use ice_history_snow, only: accum_hist_snow, & f_rhos_cmp, f_rhos_cnt, n_rhos_cmp, n_rhos_cnt use ice_history_drag, only: accum_hist_drag - use icepack_intfc, only: icepack_mushy_density_brine, icepack_mushy_liquid_fraction - use icepack_intfc, only: icepack_mushy_temperature_mush use ice_history_fsd, only: accum_hist_fsd use ice_state ! almost everything use ice_timers, only: ice_timer_start, ice_timer_stop, timer_readwrite @@ -2219,19 +2266,19 @@ subroutine accum_hist (dt) sn ! temporary variable for salinity real (kind=dbl_kind), dimension (nx_block,ny_block) :: & - worka, workb, ravgip + worka, workb, ravgip, ravgip_init, ravgip_pond, ravgip_ridge, rho_ice, rho_ocn, sal_ice real (kind=dbl_kind), dimension (nx_block,ny_block,ncat_hist) :: & - ravgipn, worka3 + ravgipn real (kind=dbl_kind) :: awtvdr, awtidr, awtvdf, awtidf, puny, secday, rad_to_deg real (kind=dbl_kind) :: Tffresh, rhoi, rhos, rhow, ice_ref_salinity - real (kind=dbl_kind) :: rho_ice, rho_ocn, Tice, Sbr, phi, rhob, dfresh, dfsalt, sicen + real (kind=dbl_kind) :: dfresh, dfsalt, sicen logical (kind=log_kind) :: formdrag, skl_bgc - logical (kind=log_kind) :: tr_pond, tr_aero, tr_brine, tr_snow + logical (kind=log_kind) :: tr_pond, tr_aero, tr_brine, tr_snow, tr_pond_topo integer (kind=int_kind) :: ktherm integer (kind=int_kind) :: nt_sice, nt_qice, nt_qsno, nt_iage, nt_FY, nt_Tsfc, & - nt_alvl, nt_vlvl + nt_alvl, nt_apnd character (len=char_len) :: saltflux_option type (block) :: & @@ -2245,11 +2292,11 @@ subroutine accum_hist (dt) rhow_out=rhow, ice_ref_salinity_out=ice_ref_salinity) call icepack_query_parameters(formdrag_out=formdrag, skl_bgc_out=skl_bgc, ktherm_out=ktherm) call icepack_query_parameters(saltflux_option_out=saltflux_option) - call icepack_query_tracer_flags(tr_pond_out=tr_pond, tr_aero_out=tr_aero, & - tr_brine_out=tr_brine, tr_snow_out=tr_snow) + call icepack_query_tracer_flags(tr_aero_out=tr_aero, tr_brine_out=tr_brine, tr_snow_out=tr_snow, & + tr_pond_out=tr_pond, tr_pond_topo_out=tr_pond_topo) call icepack_query_tracer_indices(nt_sice_out=nt_sice, nt_qice_out=nt_qice, & nt_qsno_out=nt_qsno, nt_iage_out=nt_iage, nt_FY_out=nt_FY, nt_Tsfc_out=nt_Tsfc, & - nt_alvl_out=nt_alvl, nt_vlvl_out=nt_vlvl) + nt_alvl_out=nt_alvl, nt_apnd_out=nt_apnd) call icepack_warnings_flush(nu_diag) if (icepack_warnings_aborted()) call abort_ice(error_message=subname, & file=__FILE__, line=__LINE__) @@ -2330,8 +2377,8 @@ subroutine accum_hist (dt) #ifndef __INTEL_LLVM_COMPILER !$OMP PARALLEL DO PRIVATE(iblk,i,j,ilo,ihi,jlo,jhi,this_block, & - !$OMP k,n,qn,ns,sn,rho_ocn,rho_ice,Tice,Sbr,phi,rhob,dfresh,dfsalt,sicen, & - !$OMP worka,workb,worka3,Tinz4d,Sinz4d,Tsnz4d) + !$OMP k,n,qn,ns,sn,rho_ocn,rho_ice,sal_ice,dfresh,dfsalt,sicen, & + !$OMP worka,workb,Tinz4d,Sinz4d,Tsnz4d) #endif do iblk = 1, nblocks @@ -2357,6 +2404,8 @@ subroutine accum_hist (dt) 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_aice_init (1:1) /= 'x') & + call accum_hist_field(n_aice_init, iblk, aice_init(:,:,iblk), a2D) if (f_uvel (1:1) /= 'x') & call accum_hist_field(n_uvel, iblk, uvel(:,:,iblk), a2D) if (f_vvel (1:1) /= 'x') & @@ -2447,8 +2496,7 @@ subroutine accum_hist (dt) endif if (f_fswup(1:1) /= 'x') & - call accum_hist_field(n_fswup, iblk, & - (fsw(:,:,iblk)-fswabs(:,:,iblk)*workb(:,:)), a2D) + call accum_hist_field(n_fswup, iblk, fswup(:,:,iblk), a2D) if (f_fswdn (1:1) /= 'x') & call accum_hist_field(n_fswdn, iblk, fsw(:,:,iblk), a2D) if (f_flwdn (1:1) /= 'x') & @@ -2501,14 +2549,16 @@ subroutine accum_hist (dt) 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) + call accum_hist_field(n_fswabs_ai, iblk, fswabs(:,:,iblk)*aice(:,:,iblk), a2D) - if (f_albsni (1:1) /= 'x') & + if (f_albsni(1:1) /= 'x') then call accum_hist_field(n_albsni, iblk, & - (awtvdr*alvdr(:,:,iblk) & - + awtidr*alidr(:,:,iblk) & - + awtvdf*alvdf(:,:,iblk) & - + awtidf*alidf(:,:,iblk))*workb(:,:), a2D) + (awtvdr*alvdr_ai(:,:,iblk) & + + awtidr*alidr_ai(:,:,iblk) & + + awtvdf*alvdf_ai(:,:,iblk) & + + awtidf*alidf_ai(:,:,iblk)), a2D) + endif + if (f_alvdr (1:1) /= 'x') & call accum_hist_field(n_alvdr, iblk, alvdr(:,:,iblk), a2D) if (f_alidr (1:1) /= 'x') & @@ -2538,19 +2588,19 @@ subroutine accum_hist (dt) 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) + call accum_hist_field(n_flat_ai,iblk, flat(:,:,iblk)*aice(:,:,iblk), 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) + call accum_hist_field(n_fsens_ai,iblk, fsens(:,:,iblk)*aice(:,:,iblk), 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) + call accum_hist_field(n_flwup_ai,iblk, flwout(:,:,iblk)*aice(:,:,iblk), 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) + call accum_hist_field(n_evap_ai,iblk, evap(:,:,iblk)*aice(:,:,iblk), a2D) if (f_Tair (1:1) /= 'x') & call accum_hist_field(n_Tair, iblk, Tair(:,:,iblk), a2D) @@ -2701,10 +2751,10 @@ subroutine accum_hist (dt) 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) + call accum_hist_field(n_fsurf_ai,iblk, fsurf(:,:,iblk)*aice(:,:,iblk), a2D) if (f_fcondtop_ai(1:1)/= 'x') & call accum_hist_field(n_fcondtop_ai, iblk, & - fcondtop(:,:,iblk)*workb(:,:), a2D) + fcondtop(:,:,iblk)*aice(:,:,iblk), a2D) if (f_icepresent(1:1) /= 'x') then worka(:,:) = c0 @@ -2718,121 +2768,124 @@ subroutine accum_hist (dt) ! 2D CMIP fields - if (f_sithick(1:1) /= 'x') then + if (f_sitimefrac(1:1) /= 'x') then worka(:,:) = c0 do j = jlo, jhi do i = ilo, ihi - if (aice(i,j,iblk) > puny) worka(i,j) = vice(i,j,iblk) + if (aice(i,j,iblk) > puny) worka(i,j) = c1 enddo enddo - call accum_hist_field(n_sithick, iblk, worka(:,:), a2D) + call accum_hist_field(n_sitimefrac, iblk, worka(:,:), a2D) + endif + + if (f_sithick(1:1) /= 'x') then + call accum_hist_field(n_sithick, iblk, vice(:,:,iblk), a2D) + endif + + if (f_sivol(1:1) /= 'x') then + call accum_hist_field(n_sivol, iblk, vice(:,:,iblk), a2D) + endif + + rho_ice(:,:) = rhoi + rho_ocn(:,:) = rhow + sal_ice(:,:) = ice_ref_salinity + + if (ktherm == 2) then + do j = jlo, jhi + do i = ilo, ihi + call ice_brine_density(trcr(i,j,nt_qice:nt_qice+nzilyr-1,iblk),trcr(i,j,nt_sice:nt_sice+nzilyr-1,iblk), & + sss(i,j,iblk), rho_ice(i,j), rho_ocn(i,j), sal_ice(i,j)) + enddo + enddo + endif + + if (f_simass(1:1) /= 'x') then + call accum_hist_field(n_simass, iblk, rho_ice(:,:)*vice(:,:,iblk), a2D) + endif + + if (f_sisaltmass(1:1) /= 'x') then + call accum_hist_field(n_sisaltmass, iblk, rho_ice(:,:)*sal_ice(:,:)*vice(:,:,iblk), a2D) + endif + + if (f_sisali(1:1) /= 'x') then + call accum_hist_field(n_sisali, iblk, aice(:,:,iblk)*sal_ice(:,:), a2D) + endif + + if (f_siconc(1:1) /= 'x') then + call accum_hist_field(n_siconc, iblk, aice(:,:,iblk), a2D) + endif + + if (f_sisnconc(1:1) /= 'x') then + call accum_hist_field(n_sisnconc, iblk, snowfrac(:,:,iblk), a2D) + endif + + if (f_sisnmass(1:1) /= 'x') then + call accum_hist_field(n_sisnmass, iblk, rhos*vsno(:,:,iblk), a2D) endif if (f_siage(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) worka(i,j) = aice(i,j,iblk)*trcr(i,j,nt_iage,iblk) - enddo - enddo - call accum_hist_field(n_siage, iblk, worka(:,:), a2D) + call accum_hist_field(n_siage, iblk, aice(:,:,iblk)*trcr(:,:,nt_iage,iblk), a2D) endif if (f_sisnthick(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (vsno(i,j,iblk) > puny) & - worka(i,j) = vsno(i,j,iblk) - enddo - enddo - call accum_hist_field(n_sisnthick, iblk, worka(:,:), a2D) + call accum_hist_field(n_sisnthick, iblk, vsno(:,:,iblk), a2D) endif if (f_sitemptop(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - worka(i,j) = aice(i,j,iblk)*trcr(i,j,nt_Tsfc,iblk) - enddo - enddo - call accum_hist_field(n_sitemptop, iblk, worka(:,:), a2D) + call accum_hist_field(n_sitemptop, iblk, aice(:,:,iblk)*trcr(:,:,nt_Tsfc,iblk), a2D) endif ! Tsnice is already multiplied by aicen in icepack. if (f_sitempsnic(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - worka(i,j) = Tsnice(i,j,iblk) - enddo - enddo - call accum_hist_field(n_sitempsnic, iblk, worka(:,:), a2D) + call accum_hist_field(n_sitempsnic, iblk, Tsnice(:,:,iblk), a2D) endif if (f_sitempbot(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - worka(i,j) = aice(i,j,iblk)*Tbot(i,j,iblk) - enddo - enddo - call accum_hist_field(n_sitempbot, iblk, worka(:,:), a2D) + call accum_hist_field(n_sitempbot, iblk, aice_init(:,:,iblk)*Tbot(:,:,iblk), a2D) endif if (f_siu(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) worka(i,j) = aice(i,j,iblk)*uvel(i,j,iblk) - enddo - enddo - call accum_hist_field(n_siu, iblk, worka(:,:), a2D) + call accum_hist_field(n_siu, iblk, aice(:,:,iblk)*uvel(:,:,iblk), a2D) endif if (f_siv(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) worka(i,j) = aice(i,j,iblk)*vvel(i,j,iblk) - enddo - enddo - call accum_hist_field(n_siv, iblk, worka(:,:), a2D) + call accum_hist_field(n_siv, iblk, aice(:,:,iblk)*vvel(:,:,iblk), a2D) endif if (f_sispeed(1:1) /= 'x') then worka(:,:) = c0 do j = jlo, jhi do i = ilo, ihi - if (aice(i,j,iblk) > puny) worka(i,j) = aice(i,j,iblk) & + worka(i,j) = aice(i,j,iblk) & * sqrt(uvel(i,j,iblk)*uvel(i,j,iblk)+vvel(i,j,iblk)*vvel(i,j,iblk)) enddo enddo call accum_hist_field(n_sispeed, iblk, worka(:,:), a2D) endif - if (f_sidir(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (abs(uvel(i,j,iblk)) > puny .or. abs(vvel(i,j,iblk)) > puny) & - worka(i,j) = atan2(uvel(i,j,iblk),vvel(i,j,iblk))*rad_to_deg - if (worka(i,j) < 0.0 ) then - worka(i,j) = worka(i,j) + c360 - else - worka(i,j) = worka(i,j) * c1 - endif - enddo - enddo - call accum_hist_field(n_sidir, iblk, worka(:,:), a2D) - endif + +! if (f_sidir(1:1) /= 'x') then +! worka(:,:) = c0 +! do j = jlo, jhi +! do i = ilo, ihi +! if (abs(uvel(i,j,iblk)) > puny .or. abs(vvel(i,j,iblk)) > puny) & +! worka(i,j) = atan2(uvel(i,j,iblk),vvel(i,j,iblk))*rad_to_deg +! if (worka(i,j) < 0.0 ) then +! worka(i,j) = worka(i,j) + c360 +! else +! worka(i,j) = worka(i,j) * c1 +! endif +! enddo +! enddo +! call accum_hist_field(n_sidir, iblk, worka(:,:), a2D) +! endif + if (f_sidmasstranx(1:1) /= 'x') then worka(:,:) = c0 do j = jlo, jhi do i = ilo, ihi - if (aice(i,j,iblk) > puny) & - worka(i,j) = (rhoi*p5*(vice(i+1,j,iblk)+vice(i,j,iblk))*dyU(i,j,iblk) & - + rhos*p5*(vsno(i+1,j,iblk)+vsno(i,j,iblk))*dyU(i,j,iblk)) & - * p5*(uvel(i,j-1,iblk)+uvel(i,j,iblk)) + worka(i,j) = (rho_ice(i,j)*p5*(vice(i+1,j,iblk)+vice(i,j,iblk))*dyU(i,j,iblk) & + + rhos*p5*(vsno(i+1,j,iblk)+vsno(i,j,iblk))*dyU(i,j,iblk)) & + * p5*(uvel(i,j-1,iblk)+uvel(i,j,iblk)) enddo enddo call accum_hist_field(n_sidmasstranx, iblk, worka(:,:), a2D) @@ -2842,81 +2895,32 @@ subroutine accum_hist (dt) worka(:,:) = c0 do j = jlo, jhi do i = ilo, ihi - if (aice(i,j,iblk) > puny) & - worka(i,j) = (rhoi*p5*(vice(i,j+1,iblk)+vice(i,j,iblk))*dxU(i,j,iblk) & - + rhos*p5*(vsno(i,j+1,iblk)+vsno(i,j,iblk))*dxU(i,j,iblk)) & - * p5*(vvel(i-1,j,iblk)+vvel(i,j,iblk)) + worka(i,j) = (rho_ice(i,j)*p5*(vice(i,j+1,iblk)+vice(i,j,iblk))*dxU(i,j,iblk) & + + rhos*p5*(vsno(i,j+1,iblk)+vsno(i,j,iblk))*dxU(i,j,iblk)) & + * p5*(vvel(i-1,j,iblk)+vvel(i,j,iblk)) enddo enddo call accum_hist_field(n_sidmasstrany, iblk, worka(:,:), a2D) endif if (f_sistrxdtop(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice_init(i,j,iblk) > puny) & - worka(i,j) = aice(i,j,iblk)*(aice(i,j,iblk)*strairxU(i,j,iblk)/aice_init(i,j,iblk)) - enddo - enddo - call accum_hist_field(n_sistrxdtop, iblk, worka(:,:), a2D) + call accum_hist_field(n_sistrxdtop, iblk, aice(:,:,iblk)*strairxU(:,:,iblk), a2D) endif if (f_sistrydtop(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice_init(i,j,iblk) > puny) & - worka(i,j) = aice(i,j,iblk)*(aice(i,j,iblk)*strairyU(i,j,iblk)/aice_init(i,j,iblk)) - enddo - enddo - call accum_hist_field(n_sistrydtop, iblk, worka(:,:), a2D) + call accum_hist_field(n_sistrydtop, iblk, aice(:,:,iblk)*strairyU(:,:,iblk), a2D) endif if (f_sistrxubot(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) & - worka(i,j) = aice(i,j,iblk)*strocnxU(i,j,iblk) - enddo - enddo - call accum_hist_field(n_sistrxubot, iblk, worka(:,:), a2D) + call accum_hist_field(n_sistrxubot, iblk, aice(:,:,iblk)*strocnxU(:,:,iblk), a2D) endif if (f_sistryubot(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) & - worka(i,j) = aice(i,j,iblk)*strocnyU(i,j,iblk) - enddo - enddo - call accum_hist_field(n_sistryubot, iblk, worka(:,:), a2D) + call accum_hist_field(n_sistryubot, iblk, aice(:,:,iblk)*strocnyU(:,:,iblk), a2D) endif if (f_sicompstren(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) & - worka(i,j) = aice(i,j,iblk)*strength(i,j,iblk) - enddo - enddo - call accum_hist_field(n_sicompstren, iblk, worka(:,:), a2D) - endif - - if (f_sialb(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (fsw(i,j,iblk) > puny .and. aice_init(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*(fsw(i,j,iblk)-fswabs(i,j,iblk) & - * aice(i,j,iblk)/aice_init(i,j,iblk)) / fsw(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_sialb, iblk, worka(:,:), a2D) + call accum_hist_field(n_sicompstren, iblk, aice(:,:,iblk)*strength(:,:,iblk), a2D) endif if (f_sihc(1:1) /= 'x') then @@ -2944,342 +2948,127 @@ subroutine accum_hist (dt) endif if (f_sidconcth(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = daidtt(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_sidconcth, iblk, worka(:,:), a2D) + call accum_hist_field(n_sidconcth, iblk, daidtt(:,:,iblk), a2D) endif if (f_sidconcdyn(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = daidtd(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_sidconcdyn, iblk, worka(:,:), a2D) + call accum_hist_field(n_sidconcdyn, iblk, daidtd(:,:,iblk), a2D) endif if (f_sidmassth(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = dvidtt(i,j,iblk) * rhoi - endif - enddo - enddo - call accum_hist_field(n_sidmassth, iblk, worka(:,:), a2D) + call accum_hist_field(n_sidmassth, iblk, rho_ice(:,:)*dvidtt(:,:,iblk), a2D) endif if (f_sidmassdyn(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = dvidtd(i,j,iblk) * rhoi - endif - enddo - enddo - call accum_hist_field(n_sidmassdyn, iblk, worka(:,:), a2D) + call accum_hist_field(n_sidmassdyn, iblk, rho_ice(:,:)*dvidtd(:,:,iblk), a2D) endif if (f_sidmassgrowthwat(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice_init(i,j,iblk) > puny) then - worka(i,j) = frazil(i,j,iblk)*rhoi/dt - endif - enddo - enddo - call accum_hist_field(n_sidmassgrowthwat, iblk, worka(:,:), a2D) + call accum_hist_field(n_sidmassgrowthwat, iblk, rho_ice(:,:)*frazil(:,:,iblk)/dt, a2D) endif if (f_sidmassgrowthbot(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = congel(i,j,iblk)*rhoi/dt - endif - enddo - enddo - call accum_hist_field(n_sidmassgrowthbot, iblk, worka(:,:), a2D) + call accum_hist_field(n_sidmassgrowthbot, iblk, rho_ice(:,:)*congel(:,:,iblk)/dt, a2D) endif - if (f_sidmasssi(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = snoice(i,j,iblk)*rhoi/dt - endif - enddo - enddo - call accum_hist_field(n_sidmasssi, iblk, worka(:,:), a2D) + if (f_sidmassgrowthsi(1:1) /= 'x') then + call accum_hist_field(n_sidmassgrowthsi, iblk, rho_ice(:,:)*snoice(:,:,iblk)/dt, a2D) + endif + + if (f_sisndmasssi(1:1) /= 'x') then + call accum_hist_field(n_sisndmasssi, iblk, -snoice(:,:,iblk)*rhos/dt, a2D) endif if (f_sidmassevapsubl(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = evapi(i,j,iblk)*rhoi - endif - enddo - enddo - call accum_hist_field(n_sidmassevapsubl, iblk, worka(:,:), a2D) + call accum_hist_field(n_sidmassevapsubl, iblk, rho_ice(:,:)*evapi(:,:,iblk), a2D) endif if (f_sidmassmelttop(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = meltt(i,j,iblk)*rhoi/dt - endif - enddo - enddo - call accum_hist_field(n_sidmassmelttop, iblk, worka(:,:), a2D) + call accum_hist_field(n_sidmassmelttop, iblk, -rho_ice(:,:)*meltt(:,:,iblk)/dt, a2D) endif if (f_sidmassmeltbot(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = meltb(i,j,iblk)*rhoi/dt - endif - enddo - enddo - call accum_hist_field(n_sidmassmeltbot, iblk, worka(:,:), a2D) + call accum_hist_field(n_sidmassmeltbot, iblk, -rho_ice(:,:)*meltb(:,:,iblk)/dt, a2D) endif - if (f_sidmasslat(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = meltl(i,j,iblk)*rhoi/dt - endif - enddo - enddo - call accum_hist_field(n_sidmasslat, iblk, worka(:,:), a2D) + if (f_sidmassmeltlat(1:1) /= 'x') then + call accum_hist_field(n_sidmassmeltlat, iblk, -rho_ice(:,:)*meltl(:,:,iblk)/dt, a2D) endif - if (f_sndmasssubl(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = evaps(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_sndmasssubl, iblk, worka(:,:), a2D) + if (f_sisndmasssubl(1:1) /= 'x') then + call accum_hist_field(n_sisndmasssubl, iblk, evaps(:,:,iblk), a2D) endif - if (f_sndmasssnf(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*fsnow(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_sndmasssnf, iblk, worka(:,:), a2D) + if (f_sisndmasssnf(1:1) /= 'x') then + call accum_hist_field(n_sisndmasssnf, iblk, fsnow(:,:,iblk), a2D) endif - if (f_sndmassmelt(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = melts(i,j,iblk)*rhos/dt - endif - enddo - enddo - call accum_hist_field(n_sndmassmelt, iblk, worka(:,:), a2D) + if (f_sisndmassmelt(1:1) /= 'x') then + call accum_hist_field(n_sisndmassmelt, iblk, -melts(:,:,iblk)*rhos/dt, a2D) endif - if (f_sndmassdyn(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = dvsdtd(i,j,iblk)*rhos - endif - enddo - enddo - call accum_hist_field(n_sndmassdyn, iblk, worka(:,:), a2D) + if (f_sisndmassdyn(1:1) /= 'x') then + call accum_hist_field(n_sisndmassdyn, iblk, dvsdtd(:,:,iblk)*rhos, a2D) endif if (f_siflswdtop(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (fsw(i,j,iblk) > puny .and. aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*fsw(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_siflswdtop, iblk, worka(:,:), a2D) + call accum_hist_field(n_siflswdtop, iblk, aice_init(:,:,iblk)*fsw(:,:,iblk), a2D) endif if (f_siflswutop(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (fsw(i,j,iblk) > puny .and. aice_init(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*(fsw(i,j,iblk)-fswabs(i,j,iblk) & - * aice(i,j,iblk)/aice_init(i,j,iblk)) - endif - enddo - enddo - call accum_hist_field(n_siflswutop, iblk, worka(:,:), a2D) + call accum_hist_field(n_siflswutop, iblk, fswup(:,:,iblk), a2D) endif if (f_siflswdbot(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*fswthru(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_siflswdbot, iblk, worka(:,:), a2D) + call accum_hist_field(n_siflswdbot, iblk, aice_init(:,:,iblk)*fswthru(:,:,iblk), a2D) endif if (f_sifllwdtop(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*flw(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_sifllwdtop, iblk, worka(:,:), a2D) + call accum_hist_field(n_sifllwdtop, iblk, aice_init(:,:,iblk)*flw(:,:,iblk), a2D) endif if (f_sifllwutop(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*flwout(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_sifllwutop, iblk, worka(:,:), a2D) + call accum_hist_field(n_sifllwutop, iblk, aice(:,:,iblk)*flwout(:,:,iblk), a2D) endif if (f_siflsenstop(1:1) /= 'x') then worka(:,:) = c0 do j = jlo, jhi do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*fsens(i,j,iblk) - endif + worka(i,j) = aice(i,j,iblk)*fsens(i,j,iblk) enddo enddo - call accum_hist_field(n_siflsenstop, iblk, worka(:,:), a2D) + call accum_hist_field(n_siflsenstop, iblk, aice(:,:,iblk)*fsens(:,:,iblk), a2D) endif - if (f_siflsensupbot(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*fhocn(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_siflsensupbot, iblk, worka(:,:), a2D) + if (f_siflsensbot(1:1) /= 'x') then + call accum_hist_field(n_siflsensbot, iblk, aice(:,:,iblk)*fhocn(:,:,iblk), a2D) endif - if (f_sifllatstop(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*flat(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_sifllatstop, iblk, worka(:,:), a2D) + if (f_sifllattop(1:1) /= 'x') then + call accum_hist_field(n_sifllattop, iblk, aice(:,:,iblk)*flat(:,:,iblk), a2D) endif if (f_siflcondtop(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*fcondtop(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_siflcondtop, iblk, worka(:,:), a2D) + call accum_hist_field(n_siflcondtop, iblk, aice(:,:,iblk)*fcondtop(:,:,iblk), a2D) endif if (f_siflcondbot(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice_init(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*fcondbot(i,j,iblk)/aice_init(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_siflcondbot, iblk, worka(:,:), a2D) + call accum_hist_field(n_siflcondbot, iblk, aice(:,:,iblk)*fcondbot(:,:,iblk), a2D) endif if (f_sipr(1:1) /= 'x') then worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*frain(i,j,iblk) - endif - enddo - enddo + if (tr_pond_topo) worka(:,:) = aice_init(:,:,iblk)*frain(:,:,iblk) call accum_hist_field(n_sipr, iblk, worka(:,:), a2D) endif if (f_sifb(1:1) /= 'x') then - worka(:,:) = c0 - rho_ice = rhoi - rho_ocn = rhow do j = jlo, jhi do i = ilo, ihi if (aice(i,j,iblk) > puny) then - if (ktherm == 2) then - rho_ocn = icepack_mushy_density_brine(sss(i,j,iblk)) - rho_ice = c0 - do k = 1, nzilyr - Tice = icepack_mushy_temperature_mush(trcr(i,j,nt_qice+k-1,iblk),trcr(i,j,nt_sice+k-1,iblk)) - Sbr = trcr(i,j,nt_sice+k-1,iblk) - phi = icepack_mushy_liquid_fraction(Tice,Sbr) - rhob = icepack_mushy_density_brine(Sbr) - rho_ice = rho_ice + min(phi*rhob+(c1-phi)*rhoi,rho_ocn) - enddo - rho_ice = rho_ice / real(nzilyr,kind=dbl_kind) - endif - worka(i,j) = ((rho_ocn-rho_ice)*vice(i,j,iblk) - rhos*vsno(i,j,iblk))/rho_ocn -! if (worka(i,j) < c0) then -! write(nu_diag,*) 'negative fb',rho_ocn,rho_ice,rhos -! write(nu_diag,*) vice(i,j,iblk),vsno(i,j,iblk) -! endif + worka(i,j) = max(((rho_ocn(i,j)-rho_ice(i,j))*vice(i,j,iblk)-rhos*vsno(i,j,iblk)) & + / rho_ocn(i,j), c0) endif enddo enddo @@ -3340,70 +3129,33 @@ subroutine accum_hist (dt) worka(:,:) = c0 do j = jlo, jhi do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*(frain(i,j,iblk)+melts(i,j,iblk)+meltt(i,j,iblk)) - endif + worka(i,j) = aice(i,j,iblk)*(melts(i,j,iblk)+meltt(i,j,iblk)) enddo enddo call accum_hist_field(n_siflfwdrain, iblk, worka(:,:), a2D) endif if (f_sidragtop(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*Cdn_atm(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_sidragtop, iblk, worka(:,:), a2D) + call accum_hist_field(n_sidragtop, iblk, aice_init(:,:,iblk)*Cdn_atm(:,:,iblk), a2D) endif - if (f_sirdgthick(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk)*(c1 - trcr(i,j,nt_alvl,iblk)) > puny) then - worka(i,j) = vice(i,j,iblk) * (c1 - trcr(i,j,nt_vlvl,iblk)) & - / (aice(i,j,iblk) * (c1 - trcr(i,j,nt_alvl,iblk))) - endif - enddo - enddo - call accum_hist_field(n_sirdgthick, iblk, worka(:,:), a2D) + if (f_sidragbot(1:1) /= 'x') then + call accum_hist_field(n_sidragbot, iblk, aice_init(:,:,iblk)*Cdn_ocn(:,:,iblk), a2D) endif if (f_siforcetiltx(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*strtltxU(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_siforcetiltx, iblk, worka(:,:), a2D) + call accum_hist_field(n_siforcetiltx, iblk, aice_init(:,:,iblk)*strtltxU(:,:,iblk), a2D) endif if (f_siforcetilty(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*strtltyU(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_siforcetilty, iblk, worka(:,:), a2D) + call accum_hist_field(n_siforcetilty, iblk, aice_init(:,:,iblk)*strtltyU(:,:,iblk), a2D) endif if (f_siforcecoriolx(1:1) /= 'x') then worka(:,:) = c0 do j = jlo, jhi do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*fmU(i,j,iblk)*vvel(i,j,iblk) - endif + worka(i,j) = aice(i,j,iblk)*fmU(i,j,iblk)*vvel(i,j,iblk) enddo enddo call accum_hist_field(n_siforcecoriolx, iblk, worka(:,:), a2D) @@ -3413,36 +3165,18 @@ subroutine accum_hist (dt) worka(:,:) = c0 do j = jlo, jhi do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = -aice(i,j,iblk)*fmU(i,j,iblk)*uvel(i,j,iblk) - endif + worka(i,j) = -aice(i,j,iblk)*fmU(i,j,iblk)*uvel(i,j,iblk) enddo enddo call accum_hist_field(n_siforcecorioly, iblk, worka(:,:), a2D) endif if (f_siforceintstrx(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*strintxU(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_siforceintstrx, iblk, worka(:,:), a2D) + call accum_hist_field(n_siforceintstrx, iblk, aice(:,:,iblk)*strintxU(:,:,iblk), a2D) endif if (f_siforceintstry(1:1) /= 'x') then - worka(:,:) = c0 - do j = jlo, jhi - do i = ilo, ihi - if (aice(i,j,iblk) > puny) then - worka(i,j) = aice(i,j,iblk)*strintyU(i,j,iblk) - endif - enddo - enddo - call accum_hist_field(n_siforceintstry, iblk, worka(:,:), a2D) + call accum_hist_field(n_siforceintstry, iblk, aice(:,:,iblk)*strintyU(:,:,iblk), a2D) endif endif ! if (allocated(a2D)) @@ -3488,45 +3222,19 @@ subroutine accum_hist (dt) *aicen_init(:,:,1:ncat_hist,iblk), a3Dc) if (f_siitdconc (1:1) /= 'x') then - worka3(:,:,:) = c0 - do n = 1,ncat_hist - do j = jlo, jhi - do i = ilo, ihi - if (aicen(i,j,n,iblk) > puny) then - worka3(i,j,n) = aicen(i,j,n,iblk) - endif - enddo - enddo - enddo - call accum_hist_field(n_siitdconc-n2D, iblk, ncat_hist, worka3(:,:,:), a3Dc) + call accum_hist_field(n_siitdconc-n2D, iblk, ncat_hist, aicen(:,:,:,iblk), a3Dc) + endif + + if (f_siitdsnconc (1:1) /= 'x') then + call accum_hist_field(n_siitdsnconc-n2D, iblk, ncat_hist, snowfracn(:,:,:,iblk), a3Dc) endif if (f_siitdthick (1:1) /= 'x') then - worka3(:,:,:) = c0 - do n = 1,ncat_hist - do j = jlo, jhi - do i = ilo, ihi - if (aicen(i,j,n,iblk) > puny) then - worka3(i,j,n) = vicen(i,j,n,iblk) - endif - enddo - enddo - enddo - call accum_hist_field(n_siitdthick-n2D, iblk, ncat_hist, worka3(:,:,:), a3Dc) + call accum_hist_field(n_siitdthick-n2D, iblk, ncat_hist, vicen(:,:,:,iblk), a3Dc) endif if (f_siitdsnthick (1:1) /= 'x') then - worka3(:,:,:) = c0 - do n = 1,ncat_hist - do j = jlo, jhi - do i = ilo, ihi - if (aicen(i,j,n,iblk) > puny) then - worka3(i,j,n) = vsnon(i,j,n,iblk) - endif - enddo - enddo - enddo - call accum_hist_field(n_siitdsnthick-n2D, iblk, ncat_hist, worka3(:,:,:), a3Dc) + call accum_hist_field(n_siitdsnthick-n2D, iblk, ncat_hist, vsnon(:,:,:,iblk), a3Dc) endif endif ! if (allocated(a3Dc)) @@ -3661,7 +3369,8 @@ subroutine accum_hist (dt) ravgct = c1/avgct(ns) !$OMP PARALLEL DO PRIVATE(iblk,i,j,ilo,ihi,jlo,jhi,this_block, & - !$OMP n,nn,ravgctz,ravgip,ravgipn) + !$OMP n,nn,ravgctz,ravgip,ravgip_init,ravgip_pond, & + !$OMP ravgip_ridge,ravgipn) do iblk = 1, nblocks this_block = get_block(blocks_ice(iblk),iblk) ilo = this_block%ilo @@ -3672,36 +3381,79 @@ subroutine accum_hist (dt) ! Ice fraction really needs to be on one of the history ! streams, but in case it is not. - if (n_aice(ns) > 0) then - do j = jlo, jhi - do i = ilo, ihi - if (a2D(i,j,n_aice(ns),iblk) > puny) then - ravgip(i,j) = c1/(a2D(i,j,n_aice(ns),iblk)) - else - ravgip(i,j) = c0 - endif - enddo ! i - enddo ! j + if (n_aice(ns) > 0 .and. any(avail_hist_fields(:)%avg_ice_present == 'final')) then + do j = jlo, jhi + do i = ilo, ihi + if (a2D(i,j,n_aice(ns),iblk) > puny) then + ravgip(i,j) = c1/(a2D(i,j,n_aice(ns),iblk)) + else + ravgip(i,j) = c0 + endif + enddo ! i + enddo ! j + elseif (f_aice(1:1) == 'x' .and. any(avail_hist_fields(:)%avg_ice_present == 'final')) then + call abort_ice(subname//' ERROR: f_aice must be defined', file=__FILE__, line=__LINE__) endif - if (n_aicen(ns) > n2D) then - do k=1,ncat_hist - do j = jlo, jhi - do i = ilo, ihi - if (a3Dc(i,j,k,n_aicen(ns)-n2D,iblk) > puny) then - ravgipn(i,j,k) = c1/(a3Dc(i,j,k,n_aicen(ns)-n2D,iblk)) - else - ravgipn(i,j,k) = c0 - endif - enddo ! i - enddo ! j - enddo ! k + if (n_aice_init(ns) > 0 .and. any(avail_hist_fields(:)%avg_ice_present == 'init')) then + do j = jlo, jhi + do i = ilo, ihi + if (a2D(i,j,n_aice_init(ns),iblk) > puny) then + ravgip_init(i,j) = c1/(a2D(i,j,n_aice_init(ns),iblk)) + else + ravgip_init(i,j) = c0 + endif + enddo ! i + enddo ! j + elseif (f_aice_init(1:1) == 'x' .and. any(avail_hist_fields(:)%avg_ice_present == 'init')) then + call abort_ice(subname//' ERROR: f_aice_init must be defined', file=__FILE__, line=__LINE__) + endif + if (n_apond_ai(ns) > 0 .and. any(avail_hist_fields(:)%avg_ice_present == 'pond')) then + do j = jlo, jhi + do i = ilo, ihi + if (a2D(i,j,n_apond_ai(ns),iblk) > puny) then + ravgip_pond(i,j) = c1/a2D(i,j,n_apond_ai(ns),iblk) + else + ravgip_pond(i,j) = c0 + endif + enddo ! i + enddo ! j + elseif (f_apond_ai(1:1) == 'x' .and. any(avail_hist_fields(:)%avg_ice_present == 'pond')) then + call abort_ice(subname//' ERROR: f_apond_ai must be defined', file=__FILE__, line=__LINE__) + endif + if (n_ardg(ns) > 0 .and. any(avail_hist_fields(:)%avg_ice_present == 'ridge')) then + do j = jlo, jhi + do i = ilo, ihi + if (a2D(i,j,n_ardg(ns),iblk) > puny) then + ravgip_ridge(i,j) = c1/a2D(i,j,n_ardg(ns),iblk) + else + ravgip_ridge(i,j) = c0 + endif + enddo ! i + enddo ! j + elseif (f_ardg(1:1) == 'x' .and. any(avail_hist_fields(:)%avg_ice_present == 'ridge')) then + call abort_ice(subname//' ERROR: f_ardg must be defined', file=__FILE__, line=__LINE__) + endif + if (n_aicen(ns) > n2D .and. any(avail_hist_fields(:)%avg_ice_present == 'final')) then + do k=1,ncat_hist + do j = jlo, jhi + do i = ilo, ihi + if (a3Dc(i,j,k,n_aicen(ns)-n2D,iblk) > puny) then + ravgipn(i,j,k) = c1/(a3Dc(i,j,k,n_aicen(ns)-n2D,iblk)) + else + ravgipn(i,j,k) = c0 + endif + enddo ! i + enddo ! j + enddo ! k + elseif (f_aicen(1:1) == 'x' .and. any(avail_hist_fields(:)%avg_ice_present == 'final')) then + call abort_ice(subname//' ERROR: f_aicen must be defined', file=__FILE__, line=__LINE__) endif do n = 1, num_avail_hist_fields_2D if (avail_hist_fields(n)%vhistfreq == histfreq(ns)) then ! Only average when/where ice present - if (avail_hist_fields(n)%avg_ice_present) then + if (trim(avail_hist_fields(n)%avg_ice_present) == 'final') then do j = jlo, jhi do i = ilo, ihi if (.not. tmask(i,j,iblk)) then @@ -3712,6 +3464,39 @@ subroutine accum_hist (dt) endif enddo ! i enddo ! j + elseif (trim(avail_hist_fields(n)%avg_ice_present) == 'init') then + do j = jlo, jhi + do i = ilo, ihi + if (.not. tmask(i,j,iblk)) then + a2D(i,j,n,iblk) = spval_dbl + else ! convert units + a2D(i,j,n,iblk) = avail_hist_fields(n)%cona*a2D(i,j,n,iblk) & + * ravgip_init(i,j) + avail_hist_fields(n)%conb + endif + enddo ! i + enddo ! j + elseif (trim(avail_hist_fields(n)%avg_ice_present) == 'pond') then + do j = jlo, jhi + do i = ilo, ihi + if (.not. tmask(i,j,iblk)) then + a2D(i,j,n,iblk) = spval_dbl + else ! convert units + a2D(i,j,n,iblk) = avail_hist_fields(n)%cona*a2D(i,j,n,iblk) & + * ravgip_pond(i,j) + avail_hist_fields(n)%conb + endif + enddo ! i + enddo ! j + elseif (trim(avail_hist_fields(n)%avg_ice_present) == 'ridge') then + do j = jlo, jhi + do i = ilo, ihi + if (.not. tmask(i,j,iblk)) then + a2D(i,j,n,iblk) = spval_dbl + else ! convert units + a2D(i,j,n,iblk) = avail_hist_fields(n)%cona*a2D(i,j,n,iblk) & + * ravgip_ridge(i,j) + avail_hist_fields(n)%conb + endif + enddo ! i + enddo ! j else do j = jlo, jhi do i = ilo, ihi @@ -3726,7 +3511,8 @@ subroutine accum_hist (dt) endif ! Mask ice-free points - if (avail_hist_fields(n)%mask_ice_free_points) then + if (avail_hist_fields(n)%mask_ice_free_points .and. & + trim(avail_hist_fields(n)%avg_ice_present) == 'final') then do j = jlo, jhi do i = ilo, ihi if (ravgip(i,j) == c0) a2D(i,j,n,iblk) = spval_dbl @@ -3734,73 +3520,74 @@ subroutine accum_hist (dt) enddo ! j endif - ! CMIP albedo: also mask points below horizon - if (index(avail_hist_fields(n)%vname,'sialb') /= 0) then - do j = jlo, jhi - do i = ilo, ihi - if (albcnt(i,j,iblk,ns) <= puny) a2D(i,j,n,iblk) = spval_dbl - enddo ! i - enddo ! j + if (avail_hist_fields(n)%mask_ice_free_points .and. & + trim(avail_hist_fields(n)%avg_ice_present) == 'init') then + do j = jlo, jhi + do i = ilo, ihi + if (ravgip_init(i,j) == c0) a2D(i,j,n,iblk) = spval_dbl + enddo ! i + enddo ! j + endif endif ! 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 - if (avail_hist_fields(n)%vname(1:8) == 'alvdr_ai') 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_alvdr_ai (1:1) /= 'x' .and. n_alvdr_ai(ns) /= 0) & - a2D(i,j,n_alvdr_ai(ns),iblk) = & - a2D(i,j,n_alvdr_ai(ns),iblk)*avgct(ns)*ravgctz - if (f_alvdf_ai (1:1) /= 'x' .and. n_alvdf_ai(ns) /= 0) & - a2D(i,j,n_alvdf_ai(ns),iblk) = & - a2D(i,j,n_alvdf_ai(ns),iblk)*avgct(ns)*ravgctz - if (f_alidr_ai (1:1) /= 'x' .and. n_alidr_ai(ns) /= 0) & - a2D(i,j,n_alidr_ai(ns),iblk) = & - a2D(i,j,n_alidr_ai(ns),iblk)*avgct(ns)*ravgctz - if (f_alidf_ai (1:1) /= 'x' .and. n_alidf_ai(ns) /= 0) & - a2D(i,j,n_alidf_ai(ns),iblk) = & - a2D(i,j,n_alidf_ai(ns),iblk)*avgct(ns)*ravgctz - endif - enddo ! i - enddo ! j - endif +! 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 +! if (avail_hist_fields(n)%vname(1:8) == 'alvdr_ai') 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_alvdr_ai (1:1) /= 'x' .and. n_alvdr_ai(ns) /= 0) & +! a2D(i,j,n_alvdr_ai(ns),iblk) = & +! a2D(i,j,n_alvdr_ai(ns),iblk)*avgct(ns)*ravgctz +! if (f_alvdf_ai (1:1) /= 'x' .and. n_alvdf_ai(ns) /= 0) & +! a2D(i,j,n_alvdf_ai(ns),iblk) = & +! a2D(i,j,n_alvdf_ai(ns),iblk)*avgct(ns)*ravgctz +! if (f_alidr_ai (1:1) /= 'x' .and. n_alidr_ai(ns) /= 0) & +! a2D(i,j,n_alidr_ai(ns),iblk) = & +! a2D(i,j,n_alidr_ai(ns),iblk)*avgct(ns)*ravgctz +! if (f_alidf_ai (1:1) /= 'x' .and. n_alidf_ai(ns) /= 0) & +! a2D(i,j,n_alidf_ai(ns),iblk) = & +! a2D(i,j,n_alidf_ai(ns),iblk)*avgct(ns)*ravgctz +! endif +! enddo ! i +! enddo ! j +! endif ! snwcnt averaging is not working correctly ! for now, these history fields will have zeroes includes in the averages @@ -3833,14 +3620,13 @@ subroutine accum_hist (dt) ! enddo ! j ! endif - endif ! avail_hist_fields(n)%vhistfreq == histfreq(ns) enddo ! n do n = 1, num_avail_hist_fields_3Dc nn = n2D + n if (avail_hist_fields(nn)%vhistfreq == histfreq(ns)) then - if (avail_hist_fields(nn)%avg_ice_present) then + if (trim(avail_hist_fields(nn)%avg_ice_present) /= 'none') then do k = 1, ncat_hist do j = jlo, jhi do i = ilo, ihi @@ -4040,8 +3826,10 @@ subroutine accum_hist (dt) if (n_sig1 (ns) /= 0) a2D(i,j,n_sig1(ns), iblk) = spval_dbl if (n_sig2 (ns) /= 0) a2D(i,j,n_sig2(ns), iblk) = spval_dbl if (n_sigP (ns) /= 0) a2D(i,j,n_sigP(ns), iblk) = spval_dbl - if (n_sistreave(ns) /= 0) a2D(i,j,n_sistreave(ns),iblk) = spval_dbl - if (n_sistremax(ns) /= 0) a2D(i,j,n_sistremax(ns),iblk) = spval_dbl + if (n_sistressave(ns) /= 0) a2D(i,j,n_sistressave(ns),iblk) = spval_dbl + if (n_sistressmax(ns) /= 0) a2D(i,j,n_sistressmax(ns),iblk) = spval_dbl + if (n_sidivvel(ns) /= 0) a2D(i,j,n_sidivvel(ns),iblk) = spval_dbl + if (n_sishearvel(ns) /= 0) a2D(i,j,n_sishearvel(ns),iblk) = spval_dbl if (n_mlt_onset(ns) /= 0) a2D(i,j,n_mlt_onset(ns),iblk) = spval_dbl if (n_frz_onset(ns) /= 0) a2D(i,j,n_frz_onset(ns),iblk) = spval_dbl if (n_hisnap (ns) /= 0) a2D(i,j,n_hisnap(ns), iblk) = spval_dbl @@ -4074,10 +3862,14 @@ subroutine accum_hist (dt) sig2 (i,j,iblk)*avail_hist_fields(n_sig2(ns))%cona if (n_sigP (ns) /= 0) a2D(i,j,n_sigP(ns),iblk) = & sigP (i,j,iblk)*avail_hist_fields(n_sigP(ns))%cona - if (n_sistreave(ns) /= 0) a2D(i,j,n_sistreave(ns),iblk) = & - p5*(sig1(i,j,iblk)+sig2(i,j,iblk))*avail_hist_fields(n_sistreave(ns))%cona - if (n_sistremax(ns) /= 0) a2D(i,j,n_sistremax(ns),iblk) = & - p5*(sig1(i,j,iblk)-sig2(i,j,iblk))*avail_hist_fields(n_sistremax(ns))%cona + if (n_sistressave(ns) /= 0) a2D(i,j,n_sistressave(ns),iblk) = & + p5*(sig1(i,j,iblk)+sig2(i,j,iblk))*avail_hist_fields(n_sistressave(ns))%cona + if (n_sistressmax(ns) /= 0) a2D(i,j,n_sistressmax(ns),iblk) = & + p5*(sig1(i,j,iblk)-sig2(i,j,iblk))*avail_hist_fields(n_sistressmax(ns))%cona + if (n_sidivvel(ns) /= 0) a2D(i,j,n_sidivvel(ns),iblk) = & + divu (i,j,iblk)*avail_hist_fields(n_sidivvel(ns))%cona + if (n_sishearvel(ns) /= 0) a2D(i,j,n_sishearvel(ns),iblk) = & + shear (i,j,iblk)*avail_hist_fields(n_sishearvel(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) = & diff --git a/cicecore/cicedyn/analysis/ice_history_mechred.F90 b/cicecore/cicedyn/analysis/ice_history_mechred.F90 index e0d15fcf2..f3707f6f6 100644 --- a/cicecore/cicedyn/analysis/ice_history_mechred.F90 +++ b/cicecore/cicedyn/analysis/ice_history_mechred.F90 @@ -37,6 +37,11 @@ module ice_history_mechred f_aredistn = 'x', f_vredistn = 'x', & f_araftn = 'x', f_vraftn = 'x' + ! CMIP ridging variables. + + character (len=max_nstrm), public :: & + f_sirdgconc = 'm', f_sirdgthick = 'm' + !--------------------------------------------------------------- ! namelist variables !--------------------------------------------------------------- @@ -51,13 +56,14 @@ module ice_history_mechred f_dvirdgndt, & f_aparticn, f_krdgn , & f_aredistn, f_vredistn , & - f_araftn, f_vraftn + f_araftn, f_vraftn , & + f_sirdgconc, f_sirdgthick !--------------------------------------------------------------- ! field indices !--------------------------------------------------------------- - integer (kind=int_kind), dimension(max_nstrm) :: & + integer (kind=int_kind), dimension(max_nstrm), public :: & n_ardg , n_vrdg , & n_alvl , n_vlvl , & n_dardg1dt , n_dardg2dt , & @@ -69,6 +75,9 @@ module ice_history_mechred n_aredistn , n_vredistn , & n_araftn , n_vraftn + integer (kind=int_kind), dimension(max_nstrm) :: & + n_sirdgconc, n_sirdgthick + !======================================================================= contains @@ -152,6 +161,8 @@ subroutine init_hist_mechred_2D f_vrdgn = 'x' f_araftn = 'x' f_vraftn = 'x' + f_sirdgconc = 'x' + f_sirdgthick = 'x' endif if (f_araftn /= 'x' .or. f_vraftn /= 'x') f_ardgn = f_araftn @@ -174,6 +185,8 @@ subroutine init_hist_mechred_2D call broadcast_scalar (f_vredistn, master_task) call broadcast_scalar (f_araftn, master_task) call broadcast_scalar (f_vraftn, master_task) + call broadcast_scalar (f_sirdgconc, master_task) + call broadcast_scalar (f_sirdgthick, master_task) ! 2D variables @@ -182,49 +195,63 @@ subroutine init_hist_mechred_2D if (f_alvl(1:1) /= 'x') & call define_hist_field(n_alvl,"alvl","1",tstr2D, tcstr, & - "level ice area fraction", & - "none", c1, c0, & + "level ice area fraction", & + "none", c1, c0, & ns, f_alvl) if (f_vlvl(1:1) /= 'x') & call define_hist_field(n_vlvl,"vlvl","m",tstr2D, tcstr, & - "level ice volume", & - "grid cell mean level ice thickness", c1, c0, & + "level ice volume", & + "grid cell mean level ice thickness", c1, c0, & ns, f_vlvl) if (f_ardg(1:1) /= 'x') & call define_hist_field(n_ardg,"ardg","1",tstr2D, tcstr, & - "ridged ice area fraction", & - "none", c1, c0, & + "ridged ice area fraction", & + "none", c1, c0, & ns, f_ardg) if (f_vrdg(1:1) /= 'x') & call define_hist_field(n_vrdg,"vrdg","m",tstr2D, tcstr, & - "ridged ice volume", & - "grid cell mean level ridged thickness", c1, c0, & + "ridged ice volume", & + "grid cell mean level ridged thickness", c1, c0, & ns, f_vrdg) if (f_dardg1dt(1:1) /= 'x') & call define_hist_field(n_dardg1dt,"dardg1dt","%/day",tstr2D, tcstr, & - "ice area ridging rate", & - "none", secday*c100, c0, & + "ice area ridging rate", & + "none", secday*c100, c0, & ns, f_dardg1dt) if (f_dardg2dt(1:1) /= 'x') & call define_hist_field(n_dardg2dt,"dardg2dt","%/day",tstr2D, tcstr, & - "ridge area formation rate", & - "none", secday*c100, c0, & + "ridge area formation rate", & + "none", secday*c100, c0, & ns, f_dardg2dt) if (f_dvirdgdt(1:1) /= 'x') & call define_hist_field(n_dvirdgdt,"dvirdgdt","cm/day",tstr2D, tcstr, & - "ice volume ridging rate", & - "none", mps_to_cmpdy, c0, & + "ice volume ridging rate", & + "none", mps_to_cmpdy, c0, & ns, f_dvirdgdt) if (f_opening(1:1) /= 'x') & call define_hist_field(n_opening,"opening","%/day",tstr2D, tcstr, & - "lead area opening rate", & - "none", secday*c100, c0, & + "lead area opening rate", & + "none", secday*c100, c0, & ns, f_opening) + if (f_sirdgconc(1:1) /= 'x') & + call define_hist_field(n_sirdgconc,"sirdgconc","%",tstr2D, tcstr, & + "percentage of ridged sea ice", & + "area percentage of sea ice surface that is ridged sea ice", & + c100, c0, & + ns, f_sirdgconc, avg_ice_present='none', mask_ice_free_points=.false.) + + if (f_sirdgthick(1:1) /= 'x') & + call define_hist_field(n_sirdgthick,"sirdgthick","m",tstr2D, tcstr, & + "ridged ice thickness", & + "total volume of ridged sea ice divided by area of ridges", & + c1, c0, & + ns, f_sirdgthick, avg_ice_present='ridge', mask_ice_free_points=.true.) + endif ! histfreq(ns) /= 'x' enddo ! nstreams @@ -331,6 +358,7 @@ end subroutine init_hist_mechred_3Dc subroutine accum_hist_mechred (iblk) + use ice_blocks, only: nx_block, ny_block use ice_history_shared, only: n2D, a2D, a3Dc, ncat_hist, & accum_hist_field use ice_state, only: aice, vice, trcr, aicen, vicen, trcrn @@ -380,6 +408,15 @@ subroutine accum_hist_mechred (iblk) if (f_opening(1:1) /= 'x') & call accum_hist_field(n_opening, iblk, opening(:,:,iblk), a2D) + if (f_sirdgconc(1:1)/= 'x') & + call accum_hist_field(n_sirdgconc, iblk, & + aice(:,:,iblk) * (c1 - trcr(:,:,nt_alvl,iblk)), a2D) + + if (f_sirdgthick(1:1)/= 'x') then + call accum_hist_field(n_sirdgthick, iblk, & + vice(:,:,iblk) * (c1 - trcr(:,:,nt_vlvl,iblk)), a2D) + endif + endif ! allocated(a2D) ! 3D category fields diff --git a/cicecore/cicedyn/analysis/ice_history_pond.F90 b/cicecore/cicedyn/analysis/ice_history_pond.F90 index 88b5fa899..17d1c50b3 100644 --- a/cicecore/cicedyn/analysis/ice_history_pond.F90 +++ b/cicecore/cicedyn/analysis/ice_history_pond.F90 @@ -8,7 +8,7 @@ module ice_history_pond use ice_kinds_mod use ice_domain_size, only: max_nstrm - use ice_constants, only: c0, c1 + use ice_constants, only: c0, c1, c100 use ice_fileunits, only: nu_nml, nml_filename, & get_fileunit, release_fileunit use ice_fileunits, only: nu_diag @@ -40,6 +40,11 @@ module ice_history_pond f_dpnd_freebdn= 'x', f_dpnd_initialn= 'x', & f_dpnd_dlidn = 'x' + ! CMIP related pond variables + character (len=max_nstrm), public :: & + f_simpconc = 'm', f_simpeffconc = 'm', & + f_simpthick = 'm', f_simprefrozen = 'm' + !--------------------------------------------------------------- ! namelist variables !--------------------------------------------------------------- @@ -62,13 +67,17 @@ module ice_history_pond f_dpnd_exponn , & f_dpnd_freebdn , & f_dpnd_initialn , & - f_dpnd_dlidn + f_dpnd_dlidn , & + f_simpconc , & + f_simpeffconc , & + f_simpthick , & + f_simprefrozen !--------------------------------------------------------------- ! field indices !--------------------------------------------------------------- - integer (kind=int_kind), dimension(max_nstrm) :: & + integer (kind=int_kind), dimension(max_nstrm), public :: & n_apondn , n_apeffn , & n_hpondn , & n_apond , n_apond_ai, & @@ -83,6 +92,11 @@ module ice_history_pond n_dpnd_freebdn, n_dpnd_initialn, & n_dpnd_dlidn + ! CMIP related melt pond variables + integer (kind=int_kind), dimension(max_nstrm) :: & + n_simpconc , n_simpeffconc, & + n_simpthick , n_simprefrozen + !======================================================================= contains @@ -177,6 +191,10 @@ subroutine init_hist_pond_2D f_dpnd_freebdn = 'x' f_dpnd_initialn= 'x' f_dpnd_dlidn = 'x' + f_simpconc = 'x' + f_simpeffconc = 'x' + f_simpthick = 'x' + f_simprefrozen = 'x' endif if (tr_pond_topo) then @@ -203,6 +221,10 @@ subroutine init_hist_pond_2D call broadcast_scalar (f_hpond, master_task) call broadcast_scalar (f_ipond, master_task) call broadcast_scalar (f_apeff, master_task) + call broadcast_scalar (f_simpconc, master_task) + call broadcast_scalar (f_simpeffconc, master_task) + call broadcast_scalar (f_simpthick, master_task) + call broadcast_scalar (f_simprefrozen, master_task) call broadcast_scalar (f_apond_ai, master_task) call broadcast_scalar (f_hpond_ai, master_task) call broadcast_scalar (f_ipond_ai, master_task) @@ -234,8 +256,8 @@ subroutine init_hist_pond_2D if (f_apond_ai(1:1) /= 'x') & call define_hist_field(n_apond_ai,"apond_ai","1",tstr2D, tcstr, & - "melt pond fraction of grid cell", & - "weighted by ice area", c1, c0, & + "melt pond fraction of grid cell", & + "weighted by ice area", c1, c0, & ns, f_apond_ai) if (f_hpond(1:1) /= 'x') & @@ -246,8 +268,8 @@ subroutine init_hist_pond_2D if (f_hpond_ai(1:1) /= 'x') & call define_hist_field(n_hpond_ai,"hpond_ai","m",tstr2D, tcstr, & - "mean melt pond depth over grid cell", & - "weighted by ice area", c1, c0, & + "mean melt pond depth over grid cell", & + "weighted by ice area", c1, c0, & ns, f_hpond) if (f_ipond(1:1) /= 'x') & @@ -258,64 +280,93 @@ subroutine init_hist_pond_2D if (f_ipond_ai(1:1) /= 'x') & call define_hist_field(n_ipond_ai,"ipond_ai","m",tstr2D, tcstr, & - "mean pond ice thickness over grid cell", & - "weighted by ice area", c1, c0, & + "mean pond ice thickness over grid cell", & + "weighted by ice area", c1, c0, & ns, f_ipond_ai) if (f_apeff(1:1) /= 'x') & call define_hist_field(n_apeff,"apeff","1",tstr2D, tcstr, & "radiation-effective pond area fraction of sea ice", & - "none", c1, c0, & + "none", c1, c0, & ns, f_apeff) if (f_apeff_ai(1:1) /= 'x') & call define_hist_field(n_apeff_ai,"apeff_ai","1",tstr2D, tcstr, & - "radiation-effective pond area fraction over grid cell", & - "weighted by ice area", c1, c0, & + "radiation-effective pond area fraction over grid cell", & + "weighted by ice area", c1, c0, & ns, f_apeff_ai) if (f_dpnd_flush(1:1) /= 'x') & call define_hist_field(n_dpnd_flush,"dpnd_flush","m/s",tstr2D, tcstr, & - "pond flushing rate due to ice permeability", & - "none", c1, c0, & + "pond flushing rate due to ice permeability", & + "none", c1, c0, & ns, f_dpnd_flush) if (f_dpnd_expon(1:1) /= 'x') & call define_hist_field(n_dpnd_expon,"dpnd_expon","m/s",tstr2D, tcstr, & - "exponential pond drainage rate", & - "none", c1, c0, & + "exponential pond drainage rate", & + "none", c1, c0, & ns, f_dpnd_expon) if (f_dpnd_freebd(1:1) /= 'x') & call define_hist_field(n_dpnd_freebd,"dpnd_freebd","m/s",tstr2D, tcstr, & - "pond drainage rate due to freeboard constraint", & - "none", c1, c0, & + "pond drainage rate due to freeboard constraint", & + "none", c1, c0, & ns, f_dpnd_freebd) if (f_dpnd_initial(1:1) /= 'x') & call define_hist_field(n_dpnd_initial,"dpnd_initial","m/s",tstr2D, tcstr, & - "runoff rate due to rfrac", & - "none", c1, c0, & + "runoff rate due to rfrac", & + "none", c1, c0, & ns, f_dpnd_initial) if (f_dpnd_dlid(1:1) /= 'x') & call define_hist_field(n_dpnd_dlid,"dpnd_dlid","m/s",tstr2D, tcstr, & - "pond loss / gain to ice lid freezing / melting", & - "none", c1, c0, & + "pond loss / gain to ice lid freezing / melting", & + "none", c1, c0, & ns, f_dpnd_dlid) if (f_dpnd_melt(1:1) /= 'x') & call define_hist_field(n_dpnd_melt,"dpnd_melt","m/s",tstr2D, tcstr, & - "pond drainage due to ice melting", & - "none", c1, c0, & + "pond drainage due to ice melting", & + "none", c1, c0, & ns, f_dpnd_melt) if (f_dpnd_ridge(1:1) /= 'x') & call define_hist_field(n_dpnd_ridge,"dpnd_ridge","m",tstr2D, tcstr, & - "pond drainage due to ridging", & - "none", c1, c0, & + "pond drainage due to ridging", & + "none", c1, c0, & ns, f_dpnd_ridge) + ! CMIP melt pond variables + if (f_simpconc(1:1) /= 'x') & + call define_hist_field(n_simpconc,"simpconc","%",tstr2D, tcstr, & + "percentage of sea ice covered by melt ponds", & + "area percentage of sea-ice surface that is covered by melt ponds", & + c100, c0, & + ns, f_simpconc, avg_ice_present='none', mask_ice_free_points=.false.) + + if (f_simpeffconc(1:1) /= 'x') & + call define_hist_field(n_simpeffconc,"simpeffconc","%",tstr2D, tcstr, & + "percentage of sea ice covered by effective melt ponds", & + "area percentage of sea-ice surface that is covered by open melt ponds", & + c100, c0, & + ns, f_simpeffconc, avg_ice_present='none', mask_ice_free_points=.false.) + + if (f_simprefrozen(1:1) /= 'x') & + call define_hist_field(n_simprefrozen,"simprefrozen","m",tstr2D, tcstr, & + "thickness of refrozen ice on melt ponds", & + "volume of refrozen ice on melt ponds divided by melt pond covered area", & + c1, c0, & + ns, f_simprefrozen, avg_ice_present='pond', mask_ice_free_points=.true.) + + if (f_simpthick(1:1) /= 'x') & + call define_hist_field(n_simpthick,"simpthick","m",tstr2D, tcstr, & + "melt pond depth", & + "average depth of melt ponds on sea ice, that is melt pond volume divided by melt pond area", & + c1, c0, & + ns, f_simpthick, avg_ice_present='pond', mask_ice_free_points=.true.) + endif ! histfreq(ns) /= 'x' enddo ! nstreams @@ -413,7 +464,7 @@ subroutine accum_hist_pond (iblk) use ice_flux, only: dpnd_flushn, dpnd_exponn, dpnd_freebdn, dpnd_initialn use ice_history_shared, only: n2D, a2D, a3Dc, ncat_hist, & accum_hist_field - use ice_state, only: aice, trcr, trcrn + use ice_state, only: aice, aice_init, trcr, trcrn integer (kind=int_kind), intent(in) :: & iblk ! block index @@ -455,59 +506,38 @@ subroutine accum_hist_pond (iblk) if (allocated(a2D)) then + worka(:,:) = c0 if (tr_pond_lvl) then + worka(:,:) = trcr(:,:,nt_alvl,iblk) * trcr(:,:,nt_apnd,iblk) + elseif (tr_pond_sealvl .or. tr_pond_topo) then + worka(:,:) = trcr(:,:,nt_apnd,iblk) + endif if (f_apond(1:1)/= 'x') & - call accum_hist_field(n_apond, iblk, & - trcr(:,:,nt_alvl,iblk) * trcr(:,:,nt_apnd,iblk), a2D) + call accum_hist_field(n_apond, iblk, worka(:,:), a2D) if (f_apond_ai(1:1)/= 'x') & - call accum_hist_field(n_apond_ai, iblk, & - aice(:,:,iblk) & - * trcr(:,:,nt_alvl,iblk) * trcr(:,:,nt_apnd,iblk), a2D) + call accum_hist_field(n_apond_ai, iblk, aice_init(:,:,iblk)*worka(:,:), a2D) if (f_hpond(1:1)/= 'x') & - call accum_hist_field(n_hpond, iblk, & - trcr(:,:,nt_alvl,iblk) * trcr(:,:,nt_apnd,iblk) & - * trcr(:,:,nt_hpnd,iblk), a2D) + call accum_hist_field(n_hpond, iblk, worka(:,:)*trcr(:,:,nt_hpnd,iblk), a2D) if (f_hpond_ai(1:1)/= 'x') & - call accum_hist_field(n_hpond_ai, iblk, & - aice(:,:,iblk) & - * trcr(:,:,nt_alvl,iblk) * trcr(:,:,nt_apnd,iblk) & - * trcr(:,:,nt_hpnd,iblk), a2D) + call accum_hist_field(n_hpond_ai, iblk, aice_init(:,:,iblk)*worka(:,:)*trcr(:,:,nt_hpnd,iblk), a2D) if (f_ipond(1:1)/= 'x') & - call accum_hist_field(n_ipond, iblk, & - trcr(:,:,nt_alvl,iblk) * trcr(:,:,nt_apnd,iblk) & - * trcr(:,:,nt_ipnd,iblk), a2D) + call accum_hist_field(n_ipond, iblk, worka(:,:)*trcr(:,:,nt_ipnd,iblk), a2D) if (f_ipond_ai(1:1)/= 'x') & - call accum_hist_field(n_ipond_ai, iblk, & - aice(:,:,iblk) & - * trcr(:,:,nt_alvl,iblk) * trcr(:,:,nt_apnd,iblk) & - * trcr(:,:,nt_ipnd,iblk), a2D) + call accum_hist_field(n_ipond_ai, iblk, aice_init(:,:,iblk)*worka(:,:)*trcr(:,:,nt_ipnd,iblk), a2D) - elseif (tr_pond_topo .or. tr_pond_sealvl) then + ! CMIP pond related variables + if (f_simpeffconc (1:1) /= 'x') & + call accum_hist_field(n_simpeffconc, iblk, apeff_ai(:,:,iblk), a2D) - if (f_apond(1:1)/= 'x') & - call accum_hist_field(n_apond, iblk, & - trcr(:,:,nt_apnd,iblk), a2D) - if (f_apond_ai(1:1)/= 'x') & - call accum_hist_field(n_apond_ai, iblk, & - aice(:,:,iblk) * trcr(:,:,nt_apnd,iblk), a2D) - if (f_hpond(1:1)/= 'x') & - call accum_hist_field(n_hpond, iblk, & - trcr(:,:,nt_apnd,iblk) & - * trcr(:,:,nt_hpnd,iblk), a2D) - if (f_hpond_ai(1:1)/= 'x') & - call accum_hist_field(n_hpond_ai, iblk, & - aice(:,:,iblk) * trcr(:,:,nt_apnd,iblk) & - * trcr(:,:,nt_hpnd,iblk), a2D) - if (f_ipond(1:1)/= 'x') & - call accum_hist_field(n_ipond, iblk, & - trcr(:,:,nt_apnd,iblk) & - * trcr(:,:,nt_ipnd,iblk), a2D) - if (f_ipond_ai(1:1)/= 'x') & - call accum_hist_field(n_ipond_ai, iblk, & - aice(:,:,iblk) * trcr(:,:,nt_apnd,iblk) & - * trcr(:,:,nt_ipnd,iblk), a2D) - endif ! ponds + if (f_simpconc(1:1)/= 'x') & + call accum_hist_field(n_simpconc, iblk, worka(:,:), a2D) + + if (f_simpthick(1:1)/= 'x') & + call accum_hist_field(n_simpthick, iblk, aice_init(:,:,iblk)*worka(:,:)*trcr(:,:,nt_hpnd,iblk), a2D) + + if (f_simprefrozen(1:1)/= 'x') & + call accum_hist_field(n_simprefrozen, iblk, aice_init(:,:,iblk)*worka(:,:)*trcr(:,:,nt_ipnd,iblk), a2D) this_block = get_block(blocks_ice(iblk),iblk) ilo = this_block%ilo diff --git a/cicecore/cicedyn/analysis/ice_history_shared.F90 b/cicecore/cicedyn/analysis/ice_history_shared.F90 index 83b8d8eca..2d2c500e8 100644 --- a/cicecore/cicedyn/analysis/ice_history_shared.F90 +++ b/cicecore/cicedyn/analysis/ice_history_shared.F90 @@ -1,4 +1,4 @@ -!======================================================================= +! ======================================================================= ! ! Output files: netCDF or binary data, Fortran unformatted dumps ! @@ -34,7 +34,7 @@ module ice_history_shared implicit none private - public :: define_hist_field, accum_hist_field, icefields_nml, construct_filename + public :: define_hist_field, accum_hist_field, icefields_nml, construct_filename, ice_brine_density integer (kind=int_kind), public :: history_precision @@ -80,22 +80,22 @@ module ice_history_shared !--------------------------------------------------------------- type, public :: ice_hist_field - character (len=16) :: vname ! variable name + character (len=20) :: vname ! variable name character (len=16) :: vunit ! variable units character (len=25) :: vcoord ! variable coordinates character (len=16) :: vcellmeas ! variable cell measures - character (len=55) :: vdesc ! variable description - character (len=55) :: vcomment ! variable description + character (len=80) :: vdesc ! variable description + character (len=80) :: vcomment ! variable description real (kind=dbl_kind) :: cona ! multiplicative conversion factor real (kind=dbl_kind) :: conb ! additive conversion factor character (len=1) :: vhistfreq ! frequency of history output integer (kind=int_kind) :: vhistfreq_n ! number of vhistfreq intervals - logical (kind=log_kind) :: avg_ice_present ! only average where ice is present + character (len=16) :: avg_ice_present ! only average where ice is present 'init', 'final', 'none' logical (kind=log_kind) :: mask_ice_free_points ! mask ice-free points end type integer (kind=int_kind), parameter, public :: & - max_avail_hist_fields = 800 ! Max number of history fields + max_avail_hist_fields = 1200 ! Max number of history fields integer (kind=int_kind), public :: & num_avail_hist_fields_tot = 0, & ! Current, total number of defined fields @@ -218,9 +218,11 @@ module ice_history_shared character (len=max_nstrm), public :: & ! f_example = 'md', & + f_CMIP = 'x', & f_hi = 'm', f_hs = 'm', & f_snowfrac = 'x', f_snowfracn = 'x', & f_Tsfc = 'm', f_aice = 'm', & + f_aice_init = 'x', & f_uvel = 'm', f_vvel = 'm', & f_icespd = 'm', f_icedir = 'm', & f_uvelE = 'x', f_vvelE = 'x', & @@ -289,40 +291,43 @@ module ice_history_shared f_mlt_onset = 'm', f_frz_onset = 'm', & f_iage = 'm', f_FY = 'm', & f_hisnap = 'm', f_aisnap = 'm', & - f_CMIP = 'x' , & f_sithick = 'x', f_sisnthick = 'x', & - f_siage = 'x', & + f_siage = 'x', f_siconc = 'x', & + f_sisnconc = 'x', f_sisnmass = 'x', & f_sitemptop = 'x', f_sitempsnic = 'x', & - f_sitempbot = 'x', & + f_sitempbot = 'x', f_sivol = 'x', & f_sispeed = 'x', f_sidir = 'x', & f_siu = 'x', f_siv = 'x', & f_sidmasstranx = 'x', f_sidmasstrany = 'x', & f_sistrxdtop = 'x', f_sistrydtop = 'x', & f_sistrxubot = 'x', f_sistryubot = 'x', & f_sicompstren = 'x', & - f_sialb = 'x', & + f_sisali = 'x', & f_sihc = 'x', f_sisnhc = 'x', & f_sidconcth = 'x', f_sidconcdyn = 'x', & f_sidmassth = 'x', f_sidmassdyn = 'x', & f_sidmassgrowthwat = 'x', & f_sidmassgrowthbot = 'x', & - f_sidmasssi = 'x', & + f_simass = 'x', & + f_sisaltmass = 'x', & + f_sidmassgrowthsi = 'x', & f_sidmassevapsubl = 'x', & - f_sndmasssubl = 'x', & + f_sisndmasssubl = 'x', & f_sidmassmelttop = 'x', & f_sidmassmeltbot = 'x', & - f_sidmasslat = 'x', & - f_sndmasssnf = 'x', & - f_sndmassmelt = 'x', & - f_sndmassdyn = 'x', & + f_sidmassmeltlat = 'x', & + f_sisndmasssnf = 'x', & + f_sisndmassmelt = 'x', & + f_sisndmassdyn = 'x', & + f_sisndmasssi = 'x', & f_siflswdtop = 'x', & f_siflswutop = 'x', & f_siflswdbot = 'x', & f_sifllwdtop = 'x', & f_sifllwutop = 'x', & f_siflsenstop = 'x', & - f_siflsensupbot = 'x', & - f_sifllatstop = 'x', & + f_siflsensbot = 'x', & + f_sifllattop = 'x', & f_siflcondtop = 'x', & f_siflcondbot = 'x', & f_sipr = 'x', & @@ -337,12 +342,16 @@ module ice_history_shared f_siforceintstrx = 'x', & f_siforceintstry = 'x', & f_siitdconc = 'x', & + f_siitdsnconc = 'x', & f_siitdthick = 'x', & f_siitdsnthick = 'x', & + f_sidragbot = 'x', & f_sidragtop = 'x', & - f_sirdgthick = 'x', & - f_sistreave = 'x', & - f_sistremax = 'x', & + f_sistressave = 'x', & + f_sistressmax = 'x', & + f_sidivvel = 'x', & + f_sishearvel = 'x', & + f_sitimefrac = 'x', & f_aicen = 'x', f_vicen = 'x', & f_vsnon = 'x', & f_trsig = 'm', f_icepresent = 'm', & @@ -389,9 +398,11 @@ module ice_history_shared f_VGRDb , f_VGRDa , & f_NFSD , & ! f_example , & + f_CMIP , & f_hi, f_hs , & f_snowfrac, f_snowfracn, & f_Tsfc, f_aice , & + f_aice_init, & f_uvel, f_vvel , & f_icespd, f_icedir , & ! For now, C and CD grid quantities are controlled by the generic (originally B-grid) namelist flag @@ -461,40 +472,43 @@ module ice_history_shared f_mlt_onset, f_frz_onset, & f_iage, f_FY , & f_hisnap, f_aisnap , & - f_CMIP, & f_sithick, f_sisnthick, & - f_siage, & + f_siage, f_siconc , & + f_sisnconc, f_sisnmass , & f_sitemptop, f_sitempsnic,& - f_sitempbot, & + f_sitempbot, f_sivol, & f_sispeed, f_sidir, & f_siu, f_siv, & f_sidmasstranx, f_sidmasstrany, & f_sistrxdtop, f_sistrydtop, & f_sistrxubot, f_sistryubot, & f_sicompstren, & - f_sialb, & + f_sisali, & f_sihc, f_sisnhc, & f_sidconcth, f_sidconcdyn,& f_sidmassth, f_sidmassdyn,& f_sidmassgrowthwat, & f_sidmassgrowthbot, & - f_sidmasssi, & + f_simass, & + f_sisaltmass, & + f_sidmassgrowthsi, & f_sidmassevapsubl, & - f_sndmasssubl, & + f_sisndmasssubl, & f_sidmassmelttop, & f_sidmassmeltbot, & - f_sidmasslat, & - f_sndmasssnf, & - f_sndmassmelt, & - f_sndmassdyn, & + f_sidmassmeltlat, & + f_sisndmasssnf, & + f_sisndmassmelt, & + f_sisndmassdyn, & + f_sisndmasssi, & f_siflswdtop, & f_siflswutop, & f_siflswdbot, & f_sifllwdtop, & f_sifllwutop, & f_siflsenstop, & - f_siflsensupbot, & - f_sifllatstop, & + f_siflsensbot, & + f_sifllattop, & f_siflcondtop, & f_siflcondbot, & f_sipr, & @@ -509,12 +523,16 @@ module ice_history_shared f_siforceintstrx, & f_siforceintstry, & f_siitdconc, & + f_siitdsnconc, & f_siitdthick, & f_siitdsnthick, & + f_sidragbot, & f_sidragtop, & - f_sirdgthick, & - f_sistreave, & - f_sistremax, & + f_sistressave, & + f_sistressmax, & + f_sidivvel, & + f_sishearvel, & + f_sitimefrac, & f_aicen, f_vicen , & f_vsnon, & f_trsig, f_icepresent,& @@ -593,6 +611,7 @@ module ice_history_shared n_hi , n_hs , & n_snowfrac , n_snowfracn , & n_Tsfc , n_aice , & + n_aice_init , & n_uvel , n_vvel , & n_icespd , n_icedir , & n_uvelE , n_vvelE , & @@ -663,38 +682,42 @@ module ice_history_shared n_mlt_onset , n_frz_onset , & n_hisnap , n_aisnap , & n_sithick , n_sisnthick , & - n_siage, & + n_siage , n_siconc , & + n_sisnconc , n_sisnmass , & n_sitemptop , n_sitempsnic , & - n_sitempbot , & + n_sitempbot , n_sivol , & n_sispeed , n_sidir , & - n_siu, n_siv, & + n_siu , n_siv , & n_sidmasstranx, n_sidmasstrany, & n_sistrxdtop, n_sistrydtop, & n_sistrxubot, n_sistryubot, & n_sicompstren, & - n_sialb, & + n_sisali, & n_sihc , n_sisnhc, & n_sidconcth , n_sidconcdyn, & n_sidmassth , n_sidmassdyn, & n_sidmassgrowthwat, & n_sidmassgrowthbot, & - n_sidmasssi, & + n_simass, & + n_sisaltmass, & + n_sidmassgrowthsi, & n_sidmassevapsubl, & - n_sndmasssubl, & + n_sisndmasssubl, & n_sidmassmelttop, & n_sidmassmeltbot, & - n_sidmasslat, & - n_sndmasssnf, & - n_sndmassmelt, & - n_sndmassdyn, & + n_sidmassmeltlat, & + n_sisndmasssnf, & + n_sisndmassmelt, & + n_sisndmassdyn, & + n_sisndmasssi, & n_siflswdtop, & n_siflswutop, & n_siflswdbot, & n_sifllwdtop, & n_sifllwutop, & n_siflsenstop, & - n_siflsensupbot, & - n_sifllatstop, & + n_siflsensbot, & + n_sifllattop, & n_siflcondtop, & n_siflcondbot, & n_sipr, & @@ -709,12 +732,16 @@ module ice_history_shared n_siforceintstrx, & n_siforceintstry, & n_siitdconc, & + n_siitdsnconc, & n_siitdthick, & n_siitdsnthick, & + n_sidragbot, & n_sidragtop, & - n_sirdgthick, & - n_sistreave, & - n_sistremax, & + n_sistressave, & + n_sistressmax, & + n_sidivvel, & + n_sishearvel, & + n_sitimefrac, & n_trsig , n_icepresent , & n_iage , n_FY , & n_fsurf_ai , & @@ -888,6 +915,9 @@ subroutine define_hist_field(id, vname, vunit, vcoord, vcellmeas, & vdesc , & ! variable descriptions vcomment ! variable comments + character (len=*), optional, intent(in) :: & + avg_ice_present ! compute average only when ice is present + real (kind=dbl_kind), intent(in) :: & cona , & ! multiplicative conversion factor conb ! additive conversion factor @@ -899,7 +929,6 @@ subroutine define_hist_field(id, vname, vunit, vcoord, vcellmeas, & ns ! history file stream index logical (kind=log_kind), optional, intent(in) :: & - avg_ice_present , & ! compute average only when ice is present mask_ice_free_points ! mask ice-free points integer (kind=int_kind) :: & @@ -908,13 +937,14 @@ subroutine define_hist_field(id, vname, vunit, vcoord, vcellmeas, & character (len=40) :: stmp + character (len=16) :: l_avg_ice_present ! compute average only when ice is present + logical (kind=log_kind) :: & - l_avg_ice_present , & ! compute average only when ice is present l_mask_ice_free_points ! mask ice-free points character(len=*), parameter :: subname = '(define_hist_field)' - l_avg_ice_present = .false. + l_avg_ice_present = 'none' l_mask_ice_free_points = .false. if(present(avg_ice_present)) l_avg_ice_present = avg_ice_present @@ -1199,6 +1229,55 @@ end subroutine accum_hist_field_4D !======================================================================= +! Computes total density of brine plus fresh ice. Used for mass +! related CMIP variables. +! +! 2025 Created by D. Bailey + + subroutine ice_brine_density (qice,sice,sss,rho_ice,rho_ocn,salt_ice) + + use ice_constants, only: c0, c1 + use icepack_intfc, only: icepack_mushy_density_brine, icepack_mushy_liquid_fraction + use icepack_intfc, only: icepack_mushy_temperature_mush, icepack_query_parameters + + real (kind=dbl_kind), intent(in), dimension(:) :: & + qice, & ! sea ice enthalpy of each layer (J m-3) + sice ! sea ice salinity in each layer (psu) + + real (kind=dbl_kind), intent(in) :: & + sss ! sea surface (ocean) salinity (psu) + + real (kind=dbl_kind), intent(out) :: & + rho_ice, & ! combined brine + ice density (kg m-3) + rho_ocn, & ! ocean density from sss (kg m-3) + salt_ice ! bulk salinity of brine + ice (psu) + + integer (kind=int_kind) :: k + real (kind=dbl_kind) :: rhoi ! constant fresh ice density (kg m-3) + real (kind=dbl_kind) :: & + Tice, & ! sea ice temperature in each layer (C) + Sbr, & ! salinity of brine in each layer (psu) + phi, & ! brine fraction in each layer + rhob ! density of brine (kg m-3) + + call icepack_query_parameters(rhoi_out=rhoi) + + rho_ocn = icepack_mushy_density_brine(sss) + rho_ice = c0 + salt_ice = c0 + do k = 1, nzilyr + Sbr = sice(k) + Tice = icepack_mushy_temperature_mush(qice(k),Sbr) + salt_ice = salt_ice + Sbr / real(nzilyr,kind=dbl_kind) + phi = icepack_mushy_liquid_fraction(Tice,Sbr) + rhob = icepack_mushy_density_brine(Sbr) + rho_ice = rho_ice + min(phi*rhob+(c1-phi)*rhoi,rho_ocn) + enddo + rho_ice = rho_ice / real(nzilyr,kind=dbl_kind) + + end subroutine ice_brine_density +!======================================================================= + end module ice_history_shared !======================================================================= diff --git a/cicecore/cicedyn/analysis/ice_history_snow.F90 b/cicecore/cicedyn/analysis/ice_history_snow.F90 index 19722b014..bb702defa 100644 --- a/cicecore/cicedyn/analysis/ice_history_snow.F90 +++ b/cicecore/cicedyn/analysis/ice_history_snow.F90 @@ -31,6 +31,9 @@ module ice_history_snow f_rsnw = 'm', f_rsnwn = 'x', & f_meltsliq = 'm', f_fsloss = 'x' + character (len=max_nstrm), public :: & + f_sisndmasswind = 'm' + !--------------------------------------------------------------- ! namelist variables !--------------------------------------------------------------- @@ -41,7 +44,8 @@ module ice_history_snow f_rhos_cmp, f_rhos_cmpn, & f_rhos_cnt, f_rhos_cntn, & f_rsnw, f_rsnwn, & - f_meltsliq, f_fsloss + f_meltsliq, f_fsloss, & + f_sisndmasswind !--------------------------------------------------------------- ! field indices @@ -55,6 +59,9 @@ module ice_history_snow n_rsnw, n_rsnwn, & n_meltsliq, n_fsloss + integer (kind=int_kind), dimension(max_nstrm), public :: & + n_sisndmasswind + !======================================================================= contains @@ -147,6 +154,7 @@ subroutine init_hist_snow_2D (dt) f_rsnwn = 'x' f_meltsliq = 'x' f_fsloss = 'x' + f_sisndmasswind = 'x' endif call broadcast_scalar (f_smassice, master_task) @@ -161,6 +169,7 @@ subroutine init_hist_snow_2D (dt) call broadcast_scalar (f_rsnwn, master_task) call broadcast_scalar (f_meltsliq, master_task) call broadcast_scalar (f_fsloss, master_task) + call broadcast_scalar (f_sisndmasswind, master_task) if (tr_snow) then @@ -210,6 +219,13 @@ subroutine init_hist_snow_2D (dt) "none", c1, c0, & ns, f_fsloss) + if (f_sisndmasswind(1:1) /= 'x') & + call define_hist_field(n_sisndmasswind,"sisndmasswind","kg/m^2/s",tstr2D, tcstr, & + "snow mass rate of change through wind drift of snow", & + "rate of change of snow mass due to wind-driven transport into the ocean", & + c1, c0, & + ns, f_sisndmasswind, avg_ice_present='none', mask_ice_free_points=.false.) + endif ! histfreq(ns) /= 'x' enddo ! nstreams endif ! tr_snow @@ -376,6 +392,10 @@ subroutine accum_hist_snow (iblk) call accum_hist_field(n_fsloss, iblk, & fsloss(:,:,iblk), a2D) + if (f_sisndmasswind(1:1)/= 'x') & + call accum_hist_field(n_sisndmasswind, iblk, & + fsloss(:,:,iblk), a2D) + endif ! allocated(a2D) ! 3D category fields diff --git a/cicecore/cicedyn/general/ice_flux.F90 b/cicecore/cicedyn/general/ice_flux.F90 index 16d3a6edc..e6f38acda 100644 --- a/cicecore/cicedyn/general/ice_flux.F90 +++ b/cicecore/cicedyn/general/ice_flux.F90 @@ -400,6 +400,7 @@ module ice_flux vatmT , & ! vatm on T grid (m/s) wlat , & ! lateral heat rate (m/s) fsw , & ! incoming shortwave radiation (W/m^2) + fswup , & ! outgoing shortwave radiation (W/m^2) coszen , & ! cosine solar zenith angle, < 0 for sun below horizon rdg_conv, & ! convergence term for ridging (1/s) rdg_shear ! shear term for ridging (1/s) @@ -585,6 +586,7 @@ subroutine alloc_flux vatmT (nx_block,ny_block,max_blocks), & ! vatm on T grid wlat (nx_block,ny_block,max_blocks), & ! lateral melt rate (m/s) fsw (nx_block,ny_block,max_blocks), & ! incoming shortwave radiation (W/m^2) + fswup (nx_block,ny_block,max_blocks), & ! outgoing shortwave radiation (W/m^2) coszen (nx_block,ny_block,max_blocks), & ! cosine solar zenith angle, < 0 for sun below horizon rdg_conv (nx_block,ny_block,max_blocks), & ! convergence term for ridging (1/s) rdg_shear (nx_block,ny_block,max_blocks), & ! shear term for ridging (1/s) @@ -785,6 +787,7 @@ subroutine alloc_flux vatmT = c0 wlat = c0 fsw = c0 + fswup = c0 coszen = c0 rdg_conv = c0 rdg_shear = c0 @@ -1070,8 +1073,6 @@ subroutine init_coupler_flux flat (:,:,:) = c0 fswabs (:,:,:) = c0 fswint_ai(:,:,:) = c0 - flwout (:,:,:) = -stefan_boltzmann*Tffresh**4 - ! in case atm model diagnoses Tsfc from flwout evap (:,:,:) = c0 evaps (:,:,:) = c0 evapi (:,:,:) = c0 @@ -1128,6 +1129,7 @@ subroutine init_coupler_flux coszen (:,:,:) = c0 ! Cosine of the zenith angle fsw (:,:,:) = c0 ! shortwave radiation (W/m^2) + fswup (:,:,:) = c0 ! shortwave radiation (W/m^2) scale_factor(:,:,:) = c1 ! shortwave scaling factor wind (:,:,:) = sqrt(uatm(:,:,:)**2 & + vatm(:,:,:)**2) ! wind speed, (m/s) @@ -1170,6 +1172,8 @@ subroutine init_flux_atm fswabs (:,:,:) = c0 flwout (:,:,:) = c0 evap (:,:,:) = c0 + evaps (:,:,:) = c0 + evapi (:,:,:) = c0 Tref (:,:,:) = c0 Qref (:,:,:) = c0 Uref (:,:,:) = c0 @@ -1563,9 +1567,6 @@ subroutine scale_fluxes (nx_block, ny_block, & fsens (i,j) = fsens (i,j) * ar flat (i,j) = flat (i,j) * ar fswabs (i,j) = fswabs (i,j) * ar - ! Special case where aice_init was zero and aice > 0. - if (flwout(i,j) > -puny) & - flwout (i,j) = -stefan_boltzmann *(Tf(i,j) + Tffresh)**4 flwout (i,j) = flwout (i,j) * ar evap (i,j) = evap (i,j) * ar Tref (i,j) = Tref (i,j) * ar diff --git a/cicecore/cicedyn/infrastructure/io/io_binary/ice_history_write.F90 b/cicecore/cicedyn/infrastructure/io/io_binary/ice_history_write.F90 index af4a88007..030977ccc 100644 --- a/cicecore/cicedyn/infrastructure/io/io_binary/ice_history_write.F90 +++ b/cicecore/cicedyn/infrastructure/io/io_binary/ice_history_write.F90 @@ -163,7 +163,8 @@ subroutine ice_write_hist(ns) .or. n==n_vort(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_sidivvel(ns) .or. n==n_sishearvel(ns) & + .or. n==n_sistressave(ns) .or. n==n_sistressmax(ns) & .or. n==n_mlt_onset(ns) .or. n==n_frz_onset(ns) & .or. n==n_hisnap(ns) .or. n==n_aisnap(ns)) then write (nu_hdr, 996) nrec,trim(avail_hist_fields(n)%vname), & diff --git a/cicecore/cicedyn/infrastructure/io/io_netcdf/ice_history_write.F90 b/cicecore/cicedyn/infrastructure/io/io_netcdf/ice_history_write.F90 index 2d6c5915a..668badb14 100644 --- a/cicecore/cicedyn/infrastructure/io/io_netcdf/ice_history_write.F90 +++ b/cicecore/cicedyn/infrastructure/io/io_netcdf/ice_history_write.F90 @@ -270,14 +270,15 @@ subroutine ice_write_hist (ns) endif ! Define coord time_bounds if hist_avg is true - ! bounds inherit attributes if (hist_avg(ns) .and. .not. write_ic) then - time_coord = coord_attributes('time_bounds', 'undefined', 'undefined', 'undefined') + time_coord = coord_attributes('time_bounds', 'time interval endpoints', trim(cal_units), 'undefined') dimid(1) = boundid dimid(2) = timid call ice_hist_coord_def(ncid, time_coord, nf90_double, dimid(1:2), varid) + status = nf90_put_att(ncid,varid,'calendar',cal_att) !extra attribute + call ice_check_nc(status, subname//' ERROR: defining att calendar: '//cal_att,file=__FILE__,line=__LINE__) endif endif ! histfreq(ns)/='g' @@ -406,15 +407,22 @@ subroutine ice_write_hist (ns) ! bounds fields are required for CF compliance ! dimensions (nx,ny,nverts) - ! bounds inherit attributes - var_nverts(n_lont_bnds) = coord_attributes('lont_bounds','und','und','und') - var_nverts(n_latt_bnds) = coord_attributes('latt_bounds','und','und','und') - var_nverts(n_lonu_bnds) = coord_attributes('lonu_bounds','und','und','und') - var_nverts(n_latu_bnds) = coord_attributes('latu_bounds','und','und','und') - var_nverts(n_lonn_bnds) = coord_attributes('lonn_bounds','und','und','und') - var_nverts(n_latn_bnds) = coord_attributes('latn_bounds','und','und','und') - var_nverts(n_lone_bnds) = coord_attributes('lone_bounds','und','und','und') - var_nverts(n_late_bnds) = coord_attributes('late_bounds','und','und','und') + var_nverts(n_lont_bnds) = coord_attributes('lont_bounds', & + 'longitude of gridbox corners for T points','degrees_east','und') + var_nverts(n_latt_bnds) = coord_attributes('latt_bounds', & + 'latitude of gridbox corners for T points','degrees_north','und') + var_nverts(n_lonu_bnds) = coord_attributes('lonu_bounds', & + 'longitude of gridbox corners for U points','degrees_east','und') + var_nverts(n_latu_bnds) = coord_attributes('latu_bounds', & + 'latitude of gridbox corners for U points','degrees_north','und') + var_nverts(n_lonn_bnds) = coord_attributes('lonn_bounds', & + 'longitude of gridbox corners for N points','degrees_east','und') + var_nverts(n_latn_bnds) = coord_attributes('latn_bounds', & + 'latitude of gridbox corners for N points','degrees_north','und') + var_nverts(n_lone_bnds) = coord_attributes('lone_bounds', & + 'longitude of gridbox corners for E points','degrees_east','und') + var_nverts(n_late_bnds) = coord_attributes('late_bounds', & + 'latitude of gridbox corners for E points','degrees_north','und') !----------------------------------------------------------------- ! define attributes for time-invariant variables @@ -473,7 +481,6 @@ subroutine ice_write_hist (ns) enddo ! bounds fields with dimensions (nverts,nx,ny) - ! bounds inherits attributes dimid_nverts(1) = nvertexid dimid_nverts(2) = imtid dimid_nverts(3) = jmtid @@ -1210,7 +1217,7 @@ end subroutine ice_write_hist !======================================================================= ! Defines a (time-dependent) history var in the history file -! variables have short_name, long_name and units, coordiantes and cell_measures attributes, +! variables have short_name, long_name and units, coordinates and cell_measures attributes, ! and are compressed and chunked for 'hdf5' subroutine ice_hist_field_def(ncid, hfield, lprecision, dimids, ns) @@ -1284,10 +1291,20 @@ subroutine ice_hist_field_def(ncid, hfield, lprecision, dimids, ns) if (hist_avg(ns) .and. .not. write_ic) then if (TRIM(hfield%vname(1:4))/='sig1' & .and.TRIM(hfield%vname(1:4))/='sig2' & - .and.TRIM(hfield%vname(1:9))/='sistreave' & - .and.TRIM(hfield%vname(1:9))/='sistremax' & + .and.TRIM(hfield%vname(1:5))/='trsig' & + .and.TRIM(hfield%vname(1:4))/='divu' & + .and.TRIM(hfield%vname(1:5))/='shear' & + .and.TRIM(hfield%vname(1:4))/='vort' & + .and.TRIM(hfield%vname(1:9))/='frz_onset' & + .and.TRIM(hfield%vname(1:9))/='mlt_onset' & + .and.TRIM(hfield%vname(1:6))/='aisnap' & + .and.TRIM(hfield%vname(1:6))/='hisnap' & + .and.TRIM(hfield%vname(1:8))/='sidivvel' & + .and.TRIM(hfield%vname(1:10))/='sishearvel' & + .and.TRIM(hfield%vname(1:11))/='sistressave' & + .and.TRIM(hfield%vname(1:11))/='sistressmax' & .and.TRIM(hfield%vname(1:4))/='sigP') then - if (hfield%avg_ice_present) then + if (trim(hfield%avg_ice_present) /= 'none') then status = nf90_put_att(ncid,varid,'cell_methods','area: time: mean where sea ice (mask=siconc)') call ice_check_nc(status, subname// ' ERROR: defining cell methods for '//hfield%vname, & file=__FILE__, line=__LINE__) @@ -1309,8 +1326,10 @@ subroutine ice_hist_field_def(ncid, hfield, lprecision, dimids, ns) .or.TRIM(hfield%vname(1:4))=='sig2' & .or.TRIM(hfield%vname(1:4))=='sigP' & .or.TRIM(hfield%vname(1:5))=='trsig' & - .or.TRIM(hfield%vname(1:9))=='sistreave' & - .or.TRIM(hfield%vname(1:9))=='sistremax' & + .or.TRIM(hfield%vname(1:8))=='sidivvel' & + .or.TRIM(hfield%vname(1:10))=='sishearvel' & + .or.TRIM(hfield%vname(1:11))=='sistressave' & + .or.TRIM(hfield%vname(1:11))=='sistressmax' & .or.TRIM(hfield%vname(1:9))=='mlt_onset' & .or.TRIM(hfield%vname(1:9))=='frz_onset' & .or.TRIM(hfield%vname(1:6))=='hisnap' & diff --git a/cicecore/cicedyn/infrastructure/io/io_pio2/ice_history_write.F90 b/cicecore/cicedyn/infrastructure/io/io_pio2/ice_history_write.F90 index f4c6e51db..5c41f0eda 100644 --- a/cicecore/cicedyn/infrastructure/io/io_pio2/ice_history_write.F90 +++ b/cicecore/cicedyn/infrastructure/io/io_pio2/ice_history_write.F90 @@ -268,14 +268,15 @@ subroutine ice_write_hist (ns) endif ! Define coord time_bounds if hist_avg is true - ! bounds inherit attributes if (hist_avg(ns) .and. .not. write_ic) then - time_coord = coord_attributes('time_bounds', 'undefined', 'undefined', 'undefined') + time_coord = coord_attributes('time_bounds', 'time interval bounds', trim(cal_units), 'undefined') dimid2(1) = boundid dimid2(2) = timid call ice_hist_coord_def(File, time_coord, pio_double, dimid2, varid) + call ice_pio_check(pio_put_att(File,varid,'calendar',cal_att), & + subname//' ERROR: defining att calendar: '//cal_att,file=__FILE__,line=__LINE__) endif endif ! histfreq(ns)/='g' @@ -404,15 +405,14 @@ subroutine ice_write_hist (ns) ! bounds fields are required for CF compliance ! dimensions (nx,ny,nverts) - ! bounds inherit attributes - var_nverts(n_lont_bnds) = coord_attributes('lont_bounds','und','und','und') - var_nverts(n_latt_bnds) = coord_attributes('latt_bounds','und','und','und') - var_nverts(n_lonu_bnds) = coord_attributes('lonu_bounds','und','und','und') - var_nverts(n_latu_bnds) = coord_attributes('latu_bounds','und','und','und') - var_nverts(n_lonn_bnds) = coord_attributes('lonn_bounds','und','und','und') - var_nverts(n_latn_bnds) = coord_attributes('latn_bounds','und','und','und') - var_nverts(n_lone_bnds) = coord_attributes('lone_bounds','und','und','und') - var_nverts(n_late_bnds) = coord_attributes('late_bounds','und','und','und') + var_nverts(n_lont_bnds) = coord_attributes('lont_bounds','longitude bounds (T-cell)','degrees_east','und') + var_nverts(n_latt_bnds) = coord_attributes('latt_bounds','latitude bounds (T-cell)','degrees_north','und') + var_nverts(n_lonu_bnds) = coord_attributes('lonu_bounds','longitude bounds (U-cell)','degrees_east','und') + var_nverts(n_latu_bnds) = coord_attributes('latu_bounds','latitude bounds (U-cell)','degrees_north','und') + var_nverts(n_lonn_bnds) = coord_attributes('lonn_bounds','longitude bounds (N-cell)','degrees_east','und') + var_nverts(n_latn_bnds) = coord_attributes('latn_bounds','latitude bounds (N-cell)','degrees_north','und') + var_nverts(n_lone_bnds) = coord_attributes('lone_bounds','longitude bounds (E-cell)','degrees_east','und') + var_nverts(n_late_bnds) = coord_attributes('late_bounds','latitude bounds (E-cell)','degrees_north','und') !----------------------------------------------------------------- ! define attributes for time-invariant variables @@ -467,7 +467,6 @@ subroutine ice_write_hist (ns) enddo ! bounds fields with dimensions (nverts,nx,ny) - ! bounds inherit attributes dimid_nverts(1) = nvertexid dimid_nverts(2) = imtid dimid_nverts(3) = jmtid @@ -1351,7 +1350,7 @@ end subroutine ice_hist_coord_def !======================================================================= ! Defines a (time-dependent) history var in the history file -! variables have short_name, long_name and units, coordiantes and cell_measures attributes, +! variables have short_name, long_name and units, coordinates and cell_measures attributes, ! and are compressed and chunked for 'hdf5' subroutine ice_hist_field_def(File, hfield,lprecision, dimids, ns) @@ -1425,10 +1424,20 @@ subroutine ice_hist_field_def(File, hfield,lprecision, dimids, ns) if (hist_avg(ns) .and. .not. write_ic) then if (TRIM(hfield%vname(1:4))/='sig1' & .and.TRIM(hfield%vname(1:4))/='sig2' & - .and.TRIM(hfield%vname(1:9))/='sistreave' & - .and.TRIM(hfield%vname(1:9))/='sistremax' & + .and.TRIM(hfield%vname(1:5))/='trsig' & + .and.TRIM(hfield%vname(1:4))/='divu' & + .and.TRIM(hfield%vname(1:5))/='shear' & + .and.TRIM(hfield%vname(1:4))/='vort' & + .and.TRIM(hfield%vname(1:9))/='frz_onset' & + .and.TRIM(hfield%vname(1:9))/='mlt_onset' & + .and.TRIM(hfield%vname(1:6))/='aisnap' & + .and.TRIM(hfield%vname(1:6))/='hisnap' & + .and.TRIM(hfield%vname(1:8))/='sidivvel' & + .and.TRIM(hfield%vname(1:10))/='sishearvel' & + .and.TRIM(hfield%vname(1:11))/='sistressave' & + .and.TRIM(hfield%vname(1:11))/='sistressmax' & .and.TRIM(hfield%vname(1:4))/='sigP') then - if (hfield%avg_ice_present) then + if (trim(hfield%avg_ice_present) /= 'none') then call ice_pio_check(pio_put_att(File,varid,'cell_methods', & 'area: time: mean where sea ice (mask=siconc)'), & subname//' ERROR: defining att cell_methods',file=__FILE__,line=__LINE__) @@ -1449,8 +1458,10 @@ subroutine ice_hist_field_def(File, hfield,lprecision, dimids, ns) .or.TRIM(hfield%vname(1:4))=='sig2' & .or.TRIM(hfield%vname(1:4))=='sigP' & .or.TRIM(hfield%vname(1:5))=='trsig' & - .or.TRIM(hfield%vname(1:9))=='sistreave' & - .or.TRIM(hfield%vname(1:9))=='sistremax' & + .or.TRIM(hfield%vname(1:8))=='sidivvel' & + .or.TRIM(hfield%vname(1:10))=='sishearvel' & + .or.TRIM(hfield%vname(1:11))=='sistressave' & + .or.TRIM(hfield%vname(1:11))=='sistressmax' & .or.TRIM(hfield%vname(1:9))=='mlt_onset' & .or.TRIM(hfield%vname(1:9))=='frz_onset' & .or.TRIM(hfield%vname(1:6))=='hisnap' & diff --git a/cicecore/drivers/direct/hadgem3/CICE_RunMod.F90 b/cicecore/drivers/direct/hadgem3/CICE_RunMod.F90 index ca0099680..874b382ea 100644 --- a/cicecore/drivers/direct/hadgem3/CICE_RunMod.F90 +++ b/cicecore/drivers/direct/hadgem3/CICE_RunMod.F90 @@ -358,7 +358,7 @@ subroutine coupling_prep (iblk) fswthru_ai, fhocn, fswthru, scale_factor, snowfrac, & fswthru_vdr, fswthru_vdf, fswthru_idr, fswthru_idf, & swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyt, & - fsens, flat, fswabs, flwout, evap, Tref, Qref, & + fsens, flat, fswabs, fsw, fswup, flwout, evap, Tref, Qref, & fsurfn_f, flatn_f, scale_fluxes, frzmlt_init, frzmlt use ice_flux_bgc, only: faero_ocn, flux_bio, flux_bio_ai use ice_grid, only: tmask @@ -507,6 +507,8 @@ subroutine coupling_prep (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) + fswup (i,j,iblk) = aice_init(i,j,iblk) & + * fsw (i,j,iblk) - fswabs(i,j,iblk) if (nbtrcr > 0) then do k = 1, nbtrcr diff --git a/cicecore/drivers/direct/nemo_concepts/CICE_RunMod.F90 b/cicecore/drivers/direct/nemo_concepts/CICE_RunMod.F90 index 23eb990a0..0fa89f984 100644 --- a/cicecore/drivers/direct/nemo_concepts/CICE_RunMod.F90 +++ b/cicecore/drivers/direct/nemo_concepts/CICE_RunMod.F90 @@ -358,7 +358,7 @@ subroutine coupling_prep (iblk) fswthru_ai, fhocn, fswthru, scale_factor, snowfrac, & fswthru_vdr, fswthru_vdf, fswthru_idr, fswthru_idf, & swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyt, & - fsens, flat, fswabs, flwout, evap, Tref, Qref, & + fsens, flat, fswabs, fsw, fswup, flwout, evap, Tref, Qref, & fsurfn_f, flatn_f, scale_fluxes, frzmlt_init, frzmlt use ice_flux_bgc, only: faero_ocn, flux_bio, flux_bio_ai use ice_grid, only: tmask @@ -507,6 +507,8 @@ subroutine coupling_prep (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) + fswup (i,j,iblk) = aice_init(i,j,iblk) & + * fsw (i,j,iblk) - fswabs(i,j,iblk) if (nbtrcr > 0) then do k = 1, nbtrcr diff --git a/cicecore/drivers/mapl/geos/CICE_RunMod.F90 b/cicecore/drivers/mapl/geos/CICE_RunMod.F90 index ca5cf8739..04c7fcc97 100644 --- a/cicecore/drivers/mapl/geos/CICE_RunMod.F90 +++ b/cicecore/drivers/mapl/geos/CICE_RunMod.F90 @@ -370,7 +370,7 @@ subroutine coupling_prep (iblk) fswthru_vdr, fswthru_vdf, fswthru_idr, fswthru_idf, & fswthru_uvrdr, fswthru_uvrdf, fswthru_pardr, fswthru_pardf, & swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyT, & - fsens, flat, fswabs, flwout, evap, Tref, Qref, & + fsens, flat, fswabs, fsw, fswup, flwout, evap, Tref, Qref, & scale_fluxes, frzmlt_init, frzmlt, Uref, wind use ice_flux_bgc, only: faero_ocn, fiso_ocn, Qref_iso, fiso_evap, & flux_bio, flux_bio_ai, & @@ -556,6 +556,8 @@ subroutine coupling_prep (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) + fswup (i,j,iblk) = aice_init(i,j,iblk) & + * fsw (i,j,iblk) - fswabs(i,j,iblk) if (nbtrcr > 0) then do k = 1, nbtrcr diff --git a/cicecore/drivers/mct/cesm1/CICE_RunMod.F90 b/cicecore/drivers/mct/cesm1/CICE_RunMod.F90 index d22570ae1..f951e762d 100644 --- a/cicecore/drivers/mct/cesm1/CICE_RunMod.F90 +++ b/cicecore/drivers/mct/cesm1/CICE_RunMod.F90 @@ -429,7 +429,7 @@ subroutine coupling_prep (iblk) fswthru_ai, fhocn, scale_factor, snowfrac, & fswthru, fswthru_vdr, fswthru_vdf, fswthru_idr, fswthru_idf, & swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyT, & - fsens, flat, fswabs, flwout, evap, Tref, Qref, & + fsens, flat, fswabs, fsw, fswup, flwout, evap, Tref, Qref, & scale_fluxes, frzmlt_init, frzmlt, Uref, wind use ice_flux_bgc, only: faero_ocn, fiso_ocn, Qref_iso, fiso_evap, & flux_bio, flux_bio_ai, & @@ -585,6 +585,8 @@ subroutine coupling_prep (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) + fswup (i,j,iblk) = aice_init(i,j,iblk) & + * fsw (i,j,iblk) - fswabs(i,j,iblk) if (nbtrcr > 0) then do k = 1, nbtrcr diff --git a/cicecore/drivers/nuopc/cmeps/CICE_RunMod.F90 b/cicecore/drivers/nuopc/cmeps/CICE_RunMod.F90 index 876f10512..9d487ea60 100644 --- a/cicecore/drivers/nuopc/cmeps/CICE_RunMod.F90 +++ b/cicecore/drivers/nuopc/cmeps/CICE_RunMod.F90 @@ -420,7 +420,7 @@ subroutine coupling_prep (iblk) fswthru_ai, fhocn, fswthru, scale_factor, snowfrac, & fswthru_vdr, fswthru_vdf, fswthru_idr, fswthru_idf, & swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyT, & - fsens, flat, fswabs, flwout, evap, Tref, Qref, & + fsens, flat, fswabs, fsw, fswup, flwout, evap, Tref, Qref, & scale_fluxes, frzmlt_init, frzmlt, Uref, wind use ice_flux_bgc, only: faero_ocn, fiso_ocn, Qref_iso, fiso_evap, & flux_bio, flux_bio_ai, fnit, fsil, famm, fdmsp, fdms, fhum, & @@ -577,6 +577,8 @@ subroutine coupling_prep (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) + fswup (i,j,iblk) = aice_init(i,j,iblk) & + * fsw (i,j,iblk) - fswabs(i,j,iblk) if (nbtrcr > 0) then do k = 1, nbtrcr diff --git a/cicecore/drivers/nuopc/cmeps/ice_import_export.F90 b/cicecore/drivers/nuopc/cmeps/ice_import_export.F90 index 412b2cce5..91236c11e 100644 --- a/cicecore/drivers/nuopc/cmeps/ice_import_export.F90 +++ b/cicecore/drivers/nuopc/cmeps/ice_import_export.F90 @@ -1073,6 +1073,7 @@ subroutine ice_export( exportState, rc ) call state_setexport(exportState, 'Si_imask', input=ocn_gridcell_frac, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return else + tempfld(:,:,:) = c0 do iblk = 1, nblocks this_block = get_block(blocks_ice(iblk),iblk) ilo = this_block%ilo @@ -1134,6 +1135,7 @@ subroutine ice_export( exportState, rc ) if (ChkErr(rc,__LINE__,u_FILE_u)) return ! Snow height + tempfld(:,:,:) = c0 do iblk = 1, nblocks this_block = get_block(blocks_ice(iblk),iblk) ilo = this_block%ilo @@ -1191,11 +1193,29 @@ subroutine ice_export( exportState, rc ) areacor=mod2med_areacor, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return + ! Fix outgoing longwave if aice_init = 0, but aice > 0. + tempfld(:,:,:) = flwout(:,:,:) + do iblk = 1, nblocks + this_block = get_block(blocks_ice(iblk),iblk) + ilo = this_block%ilo + ihi = this_block%ihi + jlo = this_block%jlo + jhi = this_block%jhi + do j = jlo, jhi + do i = ilo, ihi + if ( tmask(i,j,iblk) .and. ailohi(i,j,iblk) > c0 .and. flwout(i,j,iblk) > -puny) then + tempfld(i,j,iblk) = (-stefan_boltzmann *(Tf(i,j) + Tffresh)**4) / ailohi(i,j,iblk) + end if + end do + end do + end do ! longwave outgoing (upward), average over ice fraction only - call state_setexport(exportState, 'Faii_lwup' , input=flwout, lmask=tmask, ifrac=ailohi, & + call state_setexport(exportState, 'Faii_lwup' , input=tempfld, lmask=tmask, ifrac=ailohi, & areacor=mod2med_areacor, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return + deallocate(tempfld) + ! Evaporative water flux (kg/m^2/s) call state_setexport(exportState, 'Faii_evap' , input=evap, lmask=tmask, ifrac=ailohi, & areacor=mod2med_areacor, rc=rc) diff --git a/cicecore/drivers/nuopc/dmi/CICE_RunMod.F90 b/cicecore/drivers/nuopc/dmi/CICE_RunMod.F90 index 5f8fb52bc..2de15f3a4 100644 --- a/cicecore/drivers/nuopc/dmi/CICE_RunMod.F90 +++ b/cicecore/drivers/nuopc/dmi/CICE_RunMod.F90 @@ -442,14 +442,13 @@ subroutine coupling_prep (iblk) fswthru_ai, fhocn, scale_factor, snowfrac, & fswthru, fswthru_vdr, fswthru_vdf, fswthru_idr, fswthru_idf, & swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyT, & - fsens, flat, fswabs, flwout, evap, Tref, Qref, & + fsens, flat, fswabs, fsw, fswup, flwout, evap, Tref, Qref, & scale_fluxes, frzmlt_init, frzmlt use ice_flux_bgc, only: faero_ocn, fiso_ocn, Qref_iso, fiso_evap, & flux_bio, flux_bio_ai use ice_grid, only: tmask - use ice_state, only: aicen, aice + use ice_state, only: aicen, aice, aice_init #ifdef CICE_IN_NEMO - use ice_state, only: aice_init use ice_flux, only: flatn_f, fsurfn_f #endif use ice_step_mod, only: ocean_mixed_layer @@ -596,6 +595,8 @@ subroutine coupling_prep (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) + fswup (i,j,iblk) = aice_init(i,j,iblk) & + * fsw (i,j,iblk) - fswabs(i,j,iblk) if (nbtrcr > 0) then do k = 1, nbtrcr diff --git a/cicecore/drivers/standalone/cice/CICE_RunMod.F90 b/cicecore/drivers/standalone/cice/CICE_RunMod.F90 index 53476776b..335fc1e50 100644 --- a/cicecore/drivers/standalone/cice/CICE_RunMod.F90 +++ b/cicecore/drivers/standalone/cice/CICE_RunMod.F90 @@ -435,14 +435,13 @@ subroutine coupling_prep (iblk) fswthru_ai, fhocn, scale_factor, snowfrac, & fswthru, fswthru_vdr, fswthru_vdf, fswthru_idr, fswthru_idf, & swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyT, & - fsens, flat, fswabs, flwout, evap, Tref, Qref, & + fsens, flat, fswabs, fsw, fswup, flwout, evap, Tref, Qref, & scale_fluxes, frzmlt_init, frzmlt use ice_flux_bgc, only: faero_ocn, fiso_ocn, Qref_iso, fiso_evap, & flux_bio, flux_bio_ai use ice_grid, only: tmask - use ice_state, only: aicen, aice + use ice_state, only: aicen, aice, aice_init #ifdef CICE_IN_NEMO - use ice_state, only: aice_init use ice_flux, only: flatn_f, fsurfn_f #endif use ice_step_mod, only: ocean_mixed_layer @@ -585,10 +584,18 @@ subroutine coupling_prep (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) + + + !---------------------------------------------------------------- + ! Store fluxes before scaling by aice + !---------------------------------------------------------------- + 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) + fswup (i,j,iblk) = aice_init(i,j,iblk) & + * fsw (i,j,iblk) - fswabs(i,j,iblk) if (nbtrcr > 0) then do k = 1, nbtrcr diff --git a/cicecore/drivers/unittest/opticep/CICE_RunMod.F90 b/cicecore/drivers/unittest/opticep/CICE_RunMod.F90 index 53476776b..4ee0e2bc7 100644 --- a/cicecore/drivers/unittest/opticep/CICE_RunMod.F90 +++ b/cicecore/drivers/unittest/opticep/CICE_RunMod.F90 @@ -435,14 +435,13 @@ subroutine coupling_prep (iblk) fswthru_ai, fhocn, scale_factor, snowfrac, & fswthru, fswthru_vdr, fswthru_vdf, fswthru_idr, fswthru_idf, & swvdr, swidr, swvdf, swidf, Tf, Tair, Qa, strairxT, strairyT, & - fsens, flat, fswabs, flwout, evap, Tref, Qref, & + fsens, flat, fswabs, fsw, fswup, flwout, evap, Tref, Qref, & scale_fluxes, frzmlt_init, frzmlt use ice_flux_bgc, only: faero_ocn, fiso_ocn, Qref_iso, fiso_evap, & flux_bio, flux_bio_ai use ice_grid, only: tmask - use ice_state, only: aicen, aice + use ice_state, only: aicen, aice, aice_init #ifdef CICE_IN_NEMO - use ice_state, only: aice_init use ice_flux, only: flatn_f, fsurfn_f #endif use ice_step_mod, only: ocean_mixed_layer @@ -589,6 +588,8 @@ subroutine coupling_prep (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) + fswup (i,j,iblk) = aice_init(i,j,iblk) & + * fsw (i,j,iblk) - fswabs(i,j,iblk) if (nbtrcr > 0) then do k = 1, nbtrcr diff --git a/cicecore/shared/ice_calendar.F90 b/cicecore/shared/ice_calendar.F90 index 829e78218..8f0e793aa 100644 --- a/cicecore/shared/ice_calendar.F90 +++ b/cicecore/shared/ice_calendar.F90 @@ -443,6 +443,11 @@ subroutine calendar() if (mod(elapsed_days,histfreq_n(ns))==0) & write_history(ns) = .true. endif + case ("n", "N") + if (new_day .and. histfreq_n(ns)/=0) then + if (mday == histfreq_n(ns)) & + write_history(ns) = .true. + endif case ("h", "H") if (new_hour .and. histfreq_n(ns)/=0) then if (mod(elapsed_hours,histfreq_n(ns))==0) & diff --git a/configuration/scripts/ice_in b/configuration/scripts/ice_in index 04f7c0d41..5c2bb87ee 100644 --- a/configuration/scripts/ice_in +++ b/configuration/scripts/ice_in @@ -513,12 +513,15 @@ f_VGRDa = .true. f_bounds = .false. f_aice = 'm' + f_aice_init = 'x' f_hi = 'm' f_hs = 'm' f_Tsfc = 'm' f_sice = 'm' f_uvel = 'm' f_vvel = 'm' + f_icedir = 'x' + f_icespd = 'x' f_uatm = 'm' f_vatm = 'm' f_fswdn = 'm' @@ -628,7 +631,78 @@ f_fmelttn_ai = 'x' f_flatn_ai = 'x' f_fsensn_ai = 'x' - f_CMIP = 'x' + f_siage = 'x' + f_siconc = 'x' + f_sifb = 'x' + f_siflcondbot = 'x' + f_siflcondtop = 'x' + f_sifllattop = 'x' + f_sifllwdtop = 'x' + f_sifllwutop = 'x' + f_siflsensbot = 'x' + f_siflsenstop = 'x' + f_siflswdtop = 'x' + f_siflswutop = 'x' + f_sihc = 'x' + f_siitdconc = 'x' + f_siitdsnconc = 'x' + f_siitdsnthick = 'x' + f_siitdthick = 'x' + f_sisali = 'x' + f_sisaltmass = 'x' + f_sisnhc = 'x' + f_sisnthick = 'x' + f_sitempbot = 'x' + f_sitempsnic = 'x' + f_sitemptop = 'x' + f_sithick = 'x' + f_sitimefrac = 'x' + f_siu = 'x' + f_siv = 'x' + f_sicompstren = 'x' + f_sidconcdyn = 'x' + f_sidconcth = 'x' + f_sidivvel = 'x' + f_sidmassdyn = 'x' + f_sidmassevapsubl = 'x' + f_sidmassgrowthbot = 'x' + f_sidmassgrowthsi = 'x' + f_sidmassgrowthwat = 'x' + f_sidmassmeltbot = 'x' + f_sidmassmeltlat = 'x' + f_sidmassmelttop = 'x' + f_sidmassth = 'x' + f_sidmasstranx = 'x' + f_sidmasstrany = 'x' + f_sidragbot = 'x' + f_sidragtop = 'x' + f_siflfwbot = 'x' + f_siflfwdrain = 'x' + f_siflsaltbot = 'x' + f_siforcecoriolx = 'x' + f_siforcecorioly = 'x' + f_siforceintstrx = 'x' + f_siforceintstry = 'x' + f_siforcetiltx = 'x' + f_siforcetilty = 'x' + f_simass = 'x' + f_sipr = 'x' + f_sishearvel = 'x' + f_sisnconc = 'x' + f_sisndmassdyn = 'x' + f_sisndmassmelt = 'x' + f_sisndmasssi = 'x' + f_sisndmasssnf = 'x' + f_sisndmasssubl = 'x' + f_sisnmass = 'x' + f_sispeed = 'x' + f_sistressave = 'x' + f_sistressmax = 'x' + f_sistrxdtop = 'x' + f_sistrydtop = 'x' + f_sistrxubot = 'x' + f_sistryubot = 'x' + f_sivol = 'x' / &icefields_mechred_nml @@ -651,6 +725,8 @@ f_vredistn = 'x' f_araftn = 'x' f_vraftn = 'x' + f_sirdgconc = 'x' + f_sirdgthick = 'x' / &icefields_pond_nml @@ -677,6 +753,10 @@ f_hpond_ai = 'm' f_ipond_ai = 'm' f_apeff_ai = 'm' + f_simpconc = 'x' + f_simpeffconc = 'x' + f_simprefrozen = 'x' + f_simpthick = 'x' / &icefields_snow_nml @@ -692,6 +772,7 @@ f_rsnw = 'm' f_meltsliq = 'm' f_fsloss = 'm' + f_sisndmasswind = 'x' / &icefields_bgc_nml diff --git a/configuration/scripts/options/set_nml.cmip b/configuration/scripts/options/set_nml.cmip new file mode 100644 index 000000000..80d30a755 --- /dev/null +++ b/configuration/scripts/options/set_nml.cmip @@ -0,0 +1,143 @@ + histfreq = 'm','d','n','x','x' + histfreq_n = 1,1,15,1,1 + hist_avg = .true.,.true.,.false.,.false,.false. + f_aice = 'md' + f_aice_init = 'md' + f_aicen = 'md' + f_siage = 'md' + f_siconc = 'md' + f_sifb = 'md' + f_siflcondbot = 'md' + f_siflcondtop = 'md' + f_sifllattop = 'md' + f_sifllwdtop = 'md' + f_sifllwutop = 'md' + f_siflsensbot = 'md' + f_siflsenstop = 'md' + f_siflswdtop = 'md' + f_siflswutop = 'md' + f_sihc = 'md' + f_siitdconc = 'md' + f_siitdsnconc = 'md' + f_siitdsnthick = 'md' + f_siitdthick = 'md' + f_simpconc = 'md' + f_simpeffconc = 'md' + f_simprefrozen = 'md' + f_simpthick = 'md' + f_sirdgconc = 'md' + f_sirdgthick = 'md' + f_sisali = 'md' + f_sisaltmass = 'md' + f_sisnhc = 'md' + f_sisnthick = 'md' + f_sitempbot = 'md' + f_sitempsnic = 'md' + f_sitemptop = 'md' + f_sithick = 'md' + f_sitimefrac = 'md' + f_siu = 'md' + f_siv = 'md' + f_sicompstren = 'm' + f_sidconcdyn = 'm' + f_sidconcth = 'm' + f_sidivvel = 'n' + f_sidmassdyn = 'm' + f_sidmassevapsubl = 'm' + f_sidmassgrowthbot = 'm' + f_sidmassgrowthsi = 'm' + f_sidmassgrowthwat = 'm' + f_sidmassmeltbot = 'm' + f_sidmassmeltlat = 'm' + f_sidmassmelttop = 'm' + f_sidmassth = 'm' + f_sidmasstranx = 'm' + f_sidmasstrany = 'm' + f_sidragbot = 'm' + f_sidragtop = 'm' + f_siflfwbot = 'm' + f_siflfwdrain = 'm' + f_siflsaltbot = 'm' + f_siforcecoriolx = 'm' + f_siforcecorioly = 'm' + f_siforceintstrx = 'm' + f_siforceintstry = 'm' + f_siforcetiltx = 'm' + f_siforcetilty = 'm' + f_simass = 'm' + f_sipr = 'm' + f_sishearvel = 'n' + f_sisnconc = 'm' + f_sisndmassdyn = 'm' + f_sisndmassmelt = 'm' + f_sisndmasssi = 'm' + f_sisndmasssnf = 'm' + f_sisndmasssubl = 'm' + f_sisndmasswind = 'm' + f_sisnmass = 'm' + f_sispeed = 'm' + f_sistressave = 'n' + f_sistressmax = 'n' + f_sistrxdtop = 'm' + f_sistrydtop = 'm' + f_sistrxubot = 'm' + f_sistryubot = 'm' + f_sivol = 'm' + f_icepresent = 'x' + f_hi = 'x' + f_hs = 'x' + f_Tsfc = 'x' + f_iage = 'x' + f_fswdn = 'x' + f_fswup = 'x' + f_fswthru = 'x' + f_flwdn = 'x' + f_flwup = 'x' + f_fsens = 'x' + f_flat = 'x' + f_fhocn = 'x' + f_rain = 'x' + f_snow = 'x' + f_evap = 'x' + f_fsalt = 'x' + f_fresh = 'x' + f_meltl = 'x' + f_melts = 'x' + f_meltt = 'x' + f_meltb = 'x' + f_strairx = 'x' + f_strairy = 'x' + f_strocnx = 'x' + f_strocny = 'x' + f_strtltx = 'x' + f_strtlty = 'x' + f_strintx = 'x' + f_strinty = 'x' + f_strcorx = 'x' + f_strcory = 'x' + f_strength = 'x' + f_daidtt = 'x' + f_daidtd = 'x' + f_dvidtt = 'x' + f_dvidtd = 'x' + f_congel = 'x' + f_frazil = 'x' + f_snoice = 'x' + f_vicen = 'x' + f_vsnon = 'x' + f_divu = 'x' + f_shear = 'x' + f_snowfrac = 'x' + f_snowfracn = 'x' + f_uvel = 'x' + f_vvel = 'x' + f_icespd = 'x' + f_apeff = 'x' + f_apond = 'md' + f_hpond = 'x' + f_ipond = 'x' + f_fsloss = 'x' + f_alvl = 'md' + f_vlvl = 'md' + f_ardg = 'md' + f_vrdg = 'md' diff --git a/configuration/scripts/options/set_nml.histall b/configuration/scripts/options/set_nml.histall index 2b9810bdd..5423e0216 100644 --- a/configuration/scripts/options/set_nml.histall +++ b/configuration/scripts/options/set_nml.histall @@ -21,7 +21,6 @@ f_VGRDb = .true. f_VGRDa = .true. f_bounds = .true. - f_CMIP = 'm' f_aice = 'md1h' f_hi = 'h1dm' f_hs = 'd1m' @@ -29,6 +28,8 @@ f_sice = 'md' f_uvel = 'md' f_vvel = 'dm' + f_icespd = 'md' + f_icedir = 'md' f_uatm = 'dm' f_vatm = 'dm' f_fswdn = 'dm' @@ -237,3 +238,83 @@ f_aice_ww = 'md' f_diam_ww = 'md' f_hice_ww = 'md' + f_aice_init = 'md' + f_siage = 'md' + f_siconc = 'md' + f_sifb = 'md' + f_siflcondbot = 'md' + f_siflcondtop = 'md' + f_sifllattop = 'md' + f_sifllwdtop = 'md' + f_sifllwutop = 'md' + f_siflsensbot = 'md' + f_siflsenstop = 'md' + f_siflswdtop = 'md' + f_siflswutop = 'md' + f_sihc = 'md' + f_siitdconc = 'md' + f_siitdsnconc = 'md' + f_siitdsnthick = 'md' + f_siitdthick = 'md' + f_simpconc = 'md' + f_simpeffconc = 'md' + f_simprefrozen = 'md' + f_simpthick = 'md' + f_sirdgconc = 'md' + f_sirdgthick = 'md' + f_sisali = 'md' + f_sisaltmass = 'md' + f_sisnhc = 'md' + f_sisnthick = 'md' + f_sitempbot = 'md' + f_sitempsnic = 'md' + f_sitemptop = 'md' + f_sithick = 'md' + f_sitimefrac = 'md' + f_siu = 'md' + f_siv = 'md' + f_sicompstren = 'md' + f_sidconcdyn = 'md' + f_sidconcth = 'md' + f_sidivvel = 'md' + f_sidmassdyn = 'md' + f_sidmassevapsubl = 'md' + f_sidmassgrowthbot = 'md' + f_sidmassgrowthsi = 'md' + f_sidmassgrowthwat = 'md' + f_sidmassmeltbot = 'md' + f_sidmassmeltlat = 'md' + f_sidmassmelttop = 'md' + f_sidmassth = 'md' + f_sidmasstranx = 'md' + f_sidmasstrany = 'md' + f_sidragbot = 'md' + f_sidragtop = 'md' + f_siflfwbot = 'md' + f_siflfwdrain = 'md' + f_siflsaltbot = 'md' + f_siforcecoriolx = 'md' + f_siforcecorioly = 'md' + f_siforceintstrx = 'md' + f_siforceintstry = 'md' + f_siforcetiltx = 'md' + f_siforcetilty = 'md' + f_simass = 'md' + f_sipr = 'md' + f_sishearvel = 'md' + f_sisnconc = 'md' + f_sisndmassdyn = 'md' + f_sisndmassmelt = 'md' + f_sisndmasssi = 'md' + f_sisndmasssnf = 'md' + f_sisndmasssubl = 'md' + f_sisndmasswind = 'md' + f_sisnmass = 'md' + f_sispeed = 'md' + f_sistressave = 'md' + f_sistressmax = 'md' + f_sistrxdtop = 'md' + f_sistrydtop = 'md' + f_sistrxubot = 'md' + f_sistryubot = 'md' + f_sivol = 'md' diff --git a/configuration/scripts/options/set_nml.histdbg b/configuration/scripts/options/set_nml.histdbg index a70e734e5..23a7685c7 100644 --- a/configuration/scripts/options/set_nml.histdbg +++ b/configuration/scripts/options/set_nml.histdbg @@ -20,7 +20,6 @@ f_VGRDb = .true. f_VGRDa = .true. f_bounds = .true. - f_CMIP = 'm' f_aice = 'md1h' f_hi = 'h1dm' f_hs = 'd1m' diff --git a/configuration/scripts/tests/io_suite.ts b/configuration/scripts/tests/io_suite.ts index 71bcd00a5..02e9e015e 100644 --- a/configuration/scripts/tests/io_suite.ts +++ b/configuration/scripts/tests/io_suite.ts @@ -68,3 +68,4 @@ restart gx3 8x4 isotope,histall,iopio2,iohdf5,precision8 restart gx3 12x2 fsd12,histall,iopio2,iocdf1,precision8 restart gx3 16x2 debug,histall,iopio2,iocdf2,histinst,precision8 +restart gx3 12x2 cmip,ionetcdf,iocdf2 diff --git a/doc/source/user_guide/ug_case_settings.rst b/doc/source/user_guide/ug_case_settings.rst index 486dfdfc1..2f00fed85 100644 --- a/doc/source/user_guide/ug_case_settings.rst +++ b/doc/source/user_guide/ug_case_settings.rst @@ -189,12 +189,13 @@ setup_nml "", "``zero``", "restart output frequency relative to year-month-day of 0000-01-01", "" "``dumpfreq_n``", "integer array", "write restart frequency with ``dumpfreq``", "1,1,1,1,1" "``dump_last``", "logical", "write restart on last time step of simulation", "``.false.``" - "``histfreq``", "``d``", "write history every ``histfreq_n`` days", "'1','h','d','m','y'" + "``histfreq``", "``d``", "write history every ``histfreq_n`` days", "'1','h','d','m','y','n'" "", "``h``", "write history every ``histfreq_n`` hours", "" "", "``m``", "write history every ``histfreq_n`` months", "" "", "``x``", "unused frequency stream (not written)", "" "", "``y``", "write history every ``histfreq_n`` years", "" "", "``1``", "write history every ``histfreq_n`` time step", "" + "", "``n``", "write history on day ``histfreq_n``", "" "``histfreq_base``", "``init``", "history output frequency relative to year_init, month_init, day_init", "'zero','zero','zero','zero','zero'" "", "``zero``", "history output frequency relative to year-month-day of 0000-01-01", "" "``histfreq_n``", "integer array", "frequency history output is written with ``histfreq``", "1,1,1,1,1" @@ -948,6 +949,7 @@ source code for a full list of supported output fields. "", "``x``", "do not write var to history", "" "", "``y``", "write field var every ``histfreq_n`` years", "" "", "``1``", "write field var every time step", "" + "", "``n``", "write field var on day ``histfreq_n`` instantaneously", "" "", "``md``", "*e.g.,* write both monthly and daily files", "" "``f__ai``", "``d``", "write field cell average var every ``histfreq_n`` days", "" "", "``h``", "write field cell average var every ``histfreq_n`` hours", "" diff --git a/doc/source/user_guide/ug_implementation.rst b/doc/source/user_guide/ug_implementation.rst index 17b1f618e..a79e3bc1e 100644 --- a/doc/source/user_guide/ug_implementation.rst +++ b/doc/source/user_guide/ug_implementation.rst @@ -1296,8 +1296,8 @@ collected in their own history modules (**ice_history_bgc.F90**, **ice_history_drag.F90**, **ice_history_mechred.F90**, **ice_history_pond.F90**). -The history modules allow output at different frequencies. Five output -options (``1``, ``h``, ``d``, ``m``, ``y``) are available simultaneously for ``histfreq`` +The history modules allow output at different frequencies. Six output +options (``1``, ``h``, ``d``, ``m``, ``y``, ``n``) are available simultaneously for ``histfreq`` during a run, and each stream must have a unique value for ``histfreq``. In other words, ``d`` cannot be used by two different streams. Each stream has an associated frequency set by ``histfreq_n``. The frequency is @@ -1396,17 +1396,26 @@ subroutine **define_hist_field**. ``cona`` and ``conb`` are multiplicative and terms respectively that are hardwired into the source code to convert model units to history units. -Beginning with CICE v6, history variables requested by the Sea Ice Model Intercomparison -Project (SIMIP) :cite:`Notz16` have been added as possible history output variables (e.g. -``f_sithick``, ``f_sidmassgrowthbottom``, etc.). The lists of -`monthly `_ and -`daily `_ -requested SIMIP variables provide the names of possible history fields in CICE. -However, each of the additional variables can be output at any temporal frequency -specified in the **icefields_nml** section of **ice_in** as detailed above. -Additionally, a new history output variable, ``f_CMIP``, has been added. When ``f_CMIP`` -is added to the **icefields_nml** section of **ice_in** then all SIMIP variables -will be turned on for output at the frequency specified by ``f_CMIP``. +Beginning with CICE v6, history variables requested by the Sea Ice Model Intercomparison +Project (SIMIP) :cite:`Notz16` are available as history output variables +(e.g. ``f_sithick``, ``f_sidmassgrowthbottom``, etc.). The lists of +`monthly `_ and +`daily `_ +requested SIMIP variables provide their history field names in CICE. +These variables have been updated for the +`CMIP7 data request `_. + +The ``f_CMIP`` flag has been removed. This is now a ``set_nml.cmip`` namelist option +which can be invoked with the ``-s cmip`` option during cice.setup. This optional +namelist setting will turn on the CMIP data request and turn off CICE duplicates of +SIMIP variables. However, these can be changed by the user in their case ``ice_in`` file. +Note that all SIMIP variables have been updated to correspond to the new +`CMIP7 data request `_. + +Note that some SIMIP variables require division by ice or sub-ice areas, which can be extremely +small and cause the output variables to appear unphysically large. Please interpret these +quantities (such as ``sithick``) very carefully. A future release will have an option to mask +these regions. It may also be helpful for debugging to increase the precision of the history file output from 4 bytes to 8 bytes. This is changed through the ``history_precision``