subroutine globalsum (ind) #if defined O_global_sums || defined O_co2emit_diag !======================================================================= ! calculate the global sums of heat and fresh water ! input: ! ind = index for calculation (1=initial, 2=current, 3=final) !======================================================================= implicit none character(120) :: fname character(32) :: nstamp integer i, iou, id_time, ind, j, jrow, k, l, n, nc, ntrec integer nyear, nmonth, nday, nhour, nmin, nsec, it(10) integer layer logical defined real area, cstdyt, dtoih1, dtoih2, dtoic1, dtoic2 real dtah1, dtah2, dtaf1, dtaf2, dtac1, dtac2 real dtsh1, dtsh2, dtsf1, dtsf2, dtsc1, dtsc2 real dtih1, dtih2, dtif1, dtif2, dtic1, dtic2 real dtlh1, dtlh2, dtlf1, dtlf2, dtlc1, dtlc2 real dtoh1, dtof1, dtoh2, dtof2, dtoc1, dtoc2 real dth1, dth2, dtf1, dtf2, dtc1, dtc2 real tah(3), taf(3), tac(3), tsh(3), tsf(3), tsc(3), tih(3) real tif(3), tic(3), tlh(3), tlf(3), tlc(3), toh(3), tof(3) real toc(3), th, tf, tc, RFTIME, time, tmp, c1e20, vol, avgper save tah, taf, tac, tsh, tsf, tsc, tih, tif, tic, tlh, tlf, tlc save toh, tof, toc, dtoih1, dtoic1 include "size.h" include "param.h" include "pconst.h" include "stdunits.h" include "cembm.h" # if defined O_embm include "atm.h" # endif # if defined O_ice && defined O_embm # if defined O_ice_cpts include "cpts.h" # endif include "ice.h" # endif include "calendar.h" include "coord.h" include "grdvar.h" include "levind.h" include "csbc.h" include "iounit.h" include "switch.h" # if defined O_mom include "mw.h" # if defined O_npzd include "npzd.h" # endif # endif # if defined O_mtlm include "mtlm.h" # endif include "tmngr.h" !----------------------------------------------------------------------- ! sum heat, water, carbon in atmosphere, snow, ice, land, ocean !----------------------------------------------------------------------- c1e20 = 1.e20 if (ind .eq. 1) then dtoih = 0.0 dtoih1 = 0.0 dtoic = 0. dtoic1 = 0.0 endif ! zero accumulators do n=min(3,ind*2-1),3 tah(n) = 0. taf(n) = 0. tac(n) = 0. tsh(n) = 0. tsf(n) = 0. tsc(n) = 0. tih(n) = 0. tif(n) = 0. tic(n) = 0. tlh(n) = 0. tlf(n) = 0. tlc(n) = 0. toh(n) = 0. tof(n) = 0. toc(n) = 0. enddo # if defined O_mtlm # if defined O_mtlm_segday RFTIME = DAY_TRIF*segtim/DAY_YEAR # else RFTIME = DAY_TRIF/DAY_YEAR # endif # endif ! calculate totals L = 0 do jrow=2,jmtm1 # if defined O_mom if (wide_open_mw) then j = jrow else j = jmw call getrow (latdisk(taup1disk), nslab, jrow &, u(1,1,j,1,taup1), t(1,1,j,1,taup1)) endif # endif cstdyt = cst(jrow)*dyt(jrow) do i=2,imtm1 area = cstdyt*dxt(i) if (kmt(i,jrow) .eq. 0) L = L + 1 ! atmosphere # if defined O_embm tah(3) = tah(3) + at(i,jrow,2,isat)*area taf(3) = taf(3) + at(i,jrow,2,ishum)*area # if defined O_carbon_co2_2d tac(3) = tac(3) + at(i,jrow,2,ico2)*area # else tac(3) = tac(3) + co2ccn*area # endif # endif ! snow and ice # if defined O_ice && defined O_embm tsf(3) = tsf(3) + hsno(i,jrow,2)*area tif(3) = tif(3) + hice(i,jrow,2)*area # if defined O_ice_cpts do nc=1,ncat tsf(3) = tsf(3) + hseff(i,jrow,2,nc)*area tif(3) = tif(3) + heff(i,jrow,2,nc)*area enddo do layer=1,ntilay tih(3) = tih(3) + E(i,jrow,2,layer)*area & + hice(i,jrow,2)*area enddo # endif # endif # if defined O_mtlm && defined O_embm if (kmt(i,jrow) .eq. 0) & tsf(3) = tsf(3) + LYING_SNOW(L)*area*0.1/rhosno # endif # if defined O_embm ! land tlf(3) = tlf(3) + soilm(i,jrow,2)*area*1.0e-3 # if defined O_mtlm if (kmt(i,jrow) .eq. 0) then tlh(3) = tlh(3) + (TS1(L)*HCAP_SOIL*ROOTDEP ! add in most recent fluxes not felt by atmosphere yet & + (sbc(i,jrow,ievap)*vlocn + sbc(i,jrow,isens) & + sbc(i,jrow,ilwr))*segtim*SEC_DAY*1.e-3)*area*1.e-4 tlf(3) = tlf(3) + (M(L) + MNEG(L) ! add back in most recent fluxes not felt by atmosphere yet & + sbc(i,jrow,ievap)*segtim*SEC_DAY*10.)*area*1.e-4 tlc(3) = tlc(3) + (CV(L) + CS(L) # if defined O_carbon ! add in most recent fluxes not felt by atmosphere yet # if defined O_mtlm_segday & + (sbc(i,jrow,isr) - sbc(i,jrow,inpp) & + sbc(i,jrow,iburn)*(DAY_TRIF-ntlbc))*segtim*SEC_DAY # else & + ((sbc(i,jrow,isr) - sbc(i,jrow,inpp))*segtim & + sbc(i,jrow,iburn)*(DAY_TRIF-ntlbc*segtim))*SEC_DAY # endif # endif ! add in driving fluxes not felt by triffid yet & - RESP_S_DR(L)*RFTIME & )*area*1.e-4 do N=1,NPFT tlc(3) = tlc(3) + NPP_DR(L,N)*FRAC(L,N)*RFTIME*area*1.0e-4 enddo endif # endif # endif ! ocean # if defined O_mom do k=1,km vol = dzt(k)*area toh(3) = toh(3) + t(i,k,j,itemp,taup1)*vol tof(3) = tof(3) + t(i,k,j,isalt,taup1)*vol # if defined O_carbon toc(3) = toc(3) + t(i,k,j,idic,taup1)*vol # if defined O_npzd toc(3) = toc(3) + t(i,k,j,iphyt,taup1)*vol*redctn toc(3) = toc(3) + t(i,k,j,izoop,taup1)*vol*redctn toc(3) = toc(3) + t(i,k,j,idetr,taup1)*vol*redctn # if defined O_npzd_nitrogen toc(3) = toc(3) + t(i,k,j,idiaz,taup1)*vol*redctn # endif # endif # endif enddo # if defined O_sealev_data && defined O_sealev_salinity if (kmt(i,jrow) .ne. 0) & tof(3) = tof(3) + sealev*area*socn # endif # endif enddo enddo ! convert units to Joules and kilograms # if defined O_embm ! atmosphere tah(3) = (taf(3)*rhoatm*shq*vlocn+tah(3)*cpatm*rhoatm*sht)*1.0e-7 taf(3) = taf(3)*rhoatm*shq*1.0e-3 ! 4.138e-7 => 12e-6 g/umol carbon / 29 g/mol air tac(3) = tac(3)*4.138e-7*rhoatm*shc*1.0e-3 # endif # if defined O_ice && defined O_embm ! snow and ice tsh(3) = -tsf(3)*rhosno*flice*1.0e-7 tsf(3) = tsf(3)*rhosno*1.0e-3 # if defined O_ice_cpts tih(3) = tih(3)*1.0e-7 ! cpts ice energy of melting is negative # else tih(3) = -tif(3)*rhoice*flice*1.0e-7 # endif tif(3) = tif(3)*rhoice*1.0e-3 # endif ! land tlh(3) = tlh(3) tlf(3) = tlf(3) tlc(3) = tlc(3) # if defined O_mom ! ocean toh(3) = toh(3)/0.2389 tof(3) = -tof(3)/(socn*1.e3) ! 12e-9 kg/umol toc(3) = toc(3)*12.e-9 # endif ! dtoic includes any sediment burial, weathering or emissions ! 12e-9 kg/umol dtoic = dtoic*12e-9 ! calculate differences from start of run dtah1 = tah(3) - tah(1) dtaf1 = taf(3) - taf(1) dtac1 = tac(3) - tac(1) dtsh1 = tsh(3) - tsh(1) dtsf1 = tsf(3) - tsf(1) dtsc1 = tsc(3) - tsc(1) dtih1 = tih(3) - tih(1) dtif1 = tif(3) - tif(1) dtic1 = tic(3) - tic(1) dtlh1 = tlh(3) - tlh(1) dtlf1 = tlf(3) - tlf(1) dtlc1 = tlc(3) - tlc(1) dtoh1 = toh(3) - toh(1) dtof1 = tof(3) - tof(1) dtoc1 = toc(3) - toc(1) dtoih1 = dtoih1 + dtoih dtoic1 = dtoic1 + dtoic ! calculate differences from last time step dtah2 = tah(3) - tah(2) dtaf2 = taf(3) - taf(2) dtac2 = tac(3) - tac(2) dtsh2 = tsh(3) - tsh(2) dtsf2 = tsf(3) - tsf(2) dtsc2 = tsc(3) - tsc(2) dtih2 = tih(3) - tih(2) dtif2 = tif(3) - tif(2) dtic2 = tic(3) - tic(2) dtlh2 = tlh(3) - tlh(2) dtlf2 = tlf(3) - tlf(2) dtlc2 = tlc(3) - tlc(2) dtoh2 = toh(3) - toh(2) dtof2 = tof(3) - tof(2) dtoc2 = toc(3) - toc(2) dtoih2 = dtoih dtoic2 = dtoic dtoih = 0.0 dtoic = 0.0 do n=ind,2 tah(n) = tah(3) taf(n) = taf(3) tac(n) = tac(3) tsh(n) = tsh(3) tsf(n) = tsf(3) tsc(n) = tsc(3) tih(n) = tih(3) tif(n) = tif(3) tic(n) = tic(3) tlh(n) = tlh(3) tlf(n) = tlf(3) tlc(n) = tlc(3) toh(n) = toh(3) tof(n) = tof(3) toc(n) = toc(3) enddo ! write differences or totals if (ind .eq. 2) then dth1 = dtah1 + dtsh1 + dtih1 + dtlh1 + dtoh1 + dtoih1 dtf1 = dtaf1 + dtsf1 + dtif1 + dtlf1 + dtof1 dtc1 = dtac1 + dtsc1 + dtic1 + dtlc1 + dtoc1 + dtoic1 dth2 = dtah2 + dtsh2 + dtih2 + dtlh2 + dtoh2 + dtoih2 dtf2 = dtaf2 + dtsf2 + dtif2 + dtlf2 + dtof2 dtc2 = dtac2 + dtsc2 + dtic2 + dtlc2 + dtoc2 + dtoic2 if (tsits) then time = year0 + accel_yr0 + (relyr - accel_yr0)*accel call rdstmp (stamp, nyear, nmonth, nday, nhour, nmin, nsec) nyear = time call mkstmp (nstamp, nyear, nmonth, nday, nhour, nmin, nsec) call def_tsi call def_tsi_gsums (fname) # if !defined O_save_time_relyear0 ! make output time relative to year 1 time = time - 1. # endif avgper = tsiper*accel if (avgper .le. 1e-6) avgper = 0. # if defined O_save_time_endper tmp = 0. # elif defined O_save_time_startper tmp = 1. # else tmp = 0.5 # endif # if defined O_units_time_years # if defined O_calendar_360_day time = time - tmp*avgper/360. # elif defined O_calendar_gregorian time = time - tmp*avgper/365.25 # else time = time - tmp*avgper/365. # endif # else # if defined O_calendar_360_day time = time*360. - tmp*avgper # elif defined O_calendar_gregorian time = time*365.25 - tmp*avgper # else time = time*365. - tmp*avgper # endif # endif call opennext (fname, time, ntrec, iou) if (ntrec .le. 0) ntrec = 1 call putvars ('time', iou, ntrec, time, c1, c0) call rdstmp (stamp, nyear, nmonth, nday, nhour, nmin, nsec) call putvars ('T_avgper', iou, ntrec, avgper, c1, c0) # if defined O_co2emit_diag call putvars ('F_co2diag', iou, ntrec # if defined O_save_carbon_totals &, dtc2, tsiint*1.e12/yrlen, c0) # else &, dtc2, tsiint*daylen, c0) # endif carbemit = carbemit + dtc2*1e-12 dtoic1 = dtoic1 - dtc1 dtoic2 = dtoic2 - dtc2 dtc1 = 0. dtc2 = 0. # endif # if defined O_global_sums call putvars ('G_dtheat', iou, ntrec, dth1, c1, c0) call putvars ('G_dtwater', iou, ntrec, dtf1, c1, c0) call putvars ('G_dtcarb', iou, ntrec, dtc1, c1, c0) # endif # if defined O_global_sums if (iotsi .eq. stdout .or. iotsi .lt. 0) then write (*,'(a,a)') 'Changes in heat, fresh water and carbon ' &, 'from start of run and last calculation' write (*,'(a,e22.14,a,e22.14,a,e22.14,a &, e22.14,a,e22.14,a,e22.14,a)') & ' d atm ',dtah1, ' J ',dtaf1,' kg ',dtac1, ' kg' &, dtah2, ' J ',dtaf2,' kg ',dtac2, ' kg' &, ' d snow ',dtsh1, ' J ',dtsf1,' kg ',dtsc1, ' kg' &, dtsh2, ' J ',dtsf2,' kg ',dtsc2, ' kg' &, ' d ice ',dtih1, ' J ',dtif1,' kg ',dtic1, ' kg' &, dtih2, ' J ',dtif2,' kg ',dtic2, ' kg' &, ' d lnd ',dtlh1, ' J ',dtlf1,' kg ',dtlc1, ' kg' &, dtlh2, ' J ',dtlf2,' kg ',dtlc2, ' kg' &, ' d ocn ',dtoh1, ' J ',dtof1,' kg ',dtoc1, ' kg' &, dtoh2, ' J ',dtof2,' kg ',dtoc2, ' kg' &, ' d out-in ',dtoih1,' J ',0.0, ' kg ',dtoic1,' kg' &, dtoih2,' J ',0.0, ' kg ',dtoic2,' kg' &, ' d total ',dth1, ' J ',dtf1, ' kg ',dtc1, ' kg' &, dth2, ' J ',dtf2, ' kg ',dtc2, ' kg' endif endif else th = tah(ind) + tsh(ind) + tih(ind) + tlh(ind) + toh(ind) tf = taf(ind) + tsf(ind) + tif(ind) + tlf(ind) + tof(ind) tc = tac(ind) + tsc(ind) + tic(ind) + tlc(ind) + toc(ind) write (*,'(/,a,a)') 'Total heat (in Joules referenced to 0 C' &, ' and no ice or snow) and fresh water (in kg)' # if defined O_embm write (*,'(a,a)') 'Total ocean fresh water is the equivalent' &, ' difference from the ocean volume referenced to socn' # else write (*,'(a,a)') 'Total ocean fresh water is the equivalent' &, ' difference from the ocean volume referenced to 34.9 ppt' # endif write (*,'(a,e22.14,a,e22.14,a,e22.14,a)') & ' t atm ', tah(ind),' J ',taf(ind),' kg ',tac(ind),' kg' &, ' t snow ', tsh(ind),' J ',tsf(ind),' kg ',tsc(ind),' kg' &, ' t ice ', tih(ind),' J ',tif(ind),' kg ',tic(ind),' kg' &, ' t lnd ', tlh(ind),' J ',tlf(ind),' kg ',tlc(ind),' kg' &, ' t ocn ', toh(ind),' J ',tof(ind),' kg ',toc(ind),' kg' &, ' t total ', th, ' J ',tf, ' kg ',tc, ' kg' # else endif # endif endif return end subroutine gsums_tsi_def (fname, calendar, expnam, runstamp) !======================================================================= ! output routine for atmospheric global sums step integrals ! inputs: ! fname = file name ! calendar = calendar ! expnam = experiment name ! runstamp = run stamp !======================================================================= implicit none character(*) :: fname, calendar, expnam, runstamp integer id(1), id_time, iou real c0, c1, c1e3, c1e20 c0 = 0. c1 = 1. c1e3 = 1.e3 c1e20 = 1.e20 !----------------------------------------------------------------------- ! open file !----------------------------------------------------------------------- call openfile (fname, iou) !----------------------------------------------------------------------- ! start definitions !----------------------------------------------------------------------- call redef (iou) !----------------------------------------------------------------------- ! write global attributes !----------------------------------------------------------------------- call putatttext (iou, 'global', 'Conventions', 'CF-1.0') call putatttext (iou, 'global', 'experiment_name', expnam) call putatttext (iou, 'global', 'run_stamp', runstamp) !----------------------------------------------------------------------- ! define dimensions !----------------------------------------------------------------------- call defdim ('time', iou, 0, id_time) id(1) = id_time !----------------------------------------------------------------------- ! define data !----------------------------------------------------------------------- call defvar ('time', iou, 1, id, c0, c0, 'T', 'D' # if defined O_units_time_years # if !defined O_save_time_relyear0 &, 'time', 'time', 'years since 1-1-1') # else &, 'time', 'time', 'years since 0-1-1') # endif # else # if !defined O_save_time_relyear0 &, 'time', 'time', 'days since 1-1-1') # else &, 'time', 'time', 'days since 0-1-1') # endif # endif call putatttext (iou, 'time', 'calendar', calendar) call defvar ('T_avgper', iou, 1, id, c0, c0, ' ', 'F' &, 'averaging period', ' ','day') # if defined O_co2emit_diag call defvar ('F_co2diag', iou, 1, id, -c1e20, c1e20, ' ' # if defined O_save_carbon_totals &, 'F', 'diagnosed carbon emissions', ' ', 'Pg yr-1') # else &, 'F', 'diagnosed carbon emissions', ' ', 'kg s-1') # endif # endif # if defined O_global_sums call defvar ('G_dtheat', iou, 1, id, -c1e20 &, c1e20, ' ', 'F', 'heat conservation error', ' ', 'J') call defvar ('G_dtwater', iou, 1, id, -c1e20 &, c1e20, ' ', 'F', 'water conservation error', ' ', 'kg') call defvar ('G_dtcarb', iou, 1, id, -c1e20 &, c1e20, ' ', 'F', 'carbon conservation error', ' ', 'kg') # endif !----------------------------------------------------------------------- ! end definitions !----------------------------------------------------------------------- call enddef (iou) #endif return end