diff --git a/cpl/gcam_cpl_indices.F90 b/cpl/gcam_cpl_indices.F90 index c4f8f9d..8aa2ddd 100644 --- a/cpl/gcam_cpl_indices.F90 +++ b/cpl/gcam_cpl_indices.F90 @@ -43,6 +43,9 @@ module gcam_cpl_indices integer, pointer, public ::index_x2z_Sl_hr(:) ! total heterotrophic respiration integer, pointer, public ::index_x2z_Sl_npp(:) ! net primary production integer, pointer, public ::index_x2z_Sl_pftwgt(:) ! pft weights for each cell + integer, public ::index_x2z_Sl_forc_hdm = 0 ! human population density + integer, public ::index_x2z_Sl_hdd = 0 ! heating degree days + integer, public ::index_x2z_Sl_cdd = 0 ! cooling degree days integer, public ::nflds_x2z = 0 !----------------------------------------------------------------------- @@ -78,8 +81,10 @@ subroutine gcam_cpl_indices_init( ) gcam_config, base_gcam_co2_file, base_gcam_lu_wh_file, & base_co2_surface_file, base_co2_shipment_file, base_co2_aircraft_file, & base_npp_file, base_hr_file, base_pft_file, & + base_hdd_file, base_cdd_file, & gcam2elm_co2_mapping_file, gcam2elm_luc_mapping_file, & gcam2elm_woodharvest_mapping_file, gcam2elm_cdensity_mapping_file, & + gcam2elm_degdays_mapping_file, & gcam_gridfile, elm2gcam_mapping_file, & gcam2glm_glumap, gcam2glm_baselu, gcam2glm_basebiomass, & country2grid_map, country2region_map, pop_iiasa_file, gdp_iiasa_file, & @@ -89,8 +94,8 @@ subroutine gcam_cpl_indices_init( ) pasture_addtreeonly, pasture_setherbfracrem, pasture_setavailtreefracrem, & fdyndat_ehc, & read_scalars, scalar_source_dir, & - write_scalars, write_co2, & - elm_ehc_agyield_scaling, elm_ehc_carbon_scaling, ehc_eam_co2_emissions,& + write_scalars, read_hdd_cdd, write_hdd_cdd, write_co2, & + elm_ehc_agyield_scaling, elm_ehc_carbon_scaling, elm_ehc_hdd_cdd, ehc_eam_co2_emissions,& gcam_spinup, run_gcam nlfilename_iac = "gcam_in" @@ -213,6 +218,11 @@ subroutine gcam_cpl_indices_set( ) end do + ! Scalar per-gridcell lnd->iac field (no PFT loop needed) + index_x2z_Sl_forc_hdm = mct_avect_indexra(x2z, 'Sl_forc_hdm') + index_x2z_Sl_hdd = mct_avect_indexra(x2z, 'Sl_hdd') + index_x2z_Sl_cdd = mct_avect_indexra(x2z, 'Sl_cdd') + ! iac -> atm ! Monthly sfc, low alt air, high alt air do m=1,12 diff --git a/cpl/iac_import_export.F90 b/cpl/iac_import_export.F90 index 1a946c3..7c99d2c 100644 --- a/cpl/iac_import_export.F90 +++ b/cpl/iac_import_export.F90 @@ -33,6 +33,8 @@ subroutine iac_import(x2z, lnd2iac_vars) integer :: n,n1,p,g,i,j integer :: begg, endg character(len=32), parameter :: sub = '(iac_import)' + real(r8), parameter :: dayspy = 365._r8 ! days per year + real(r8), parameter :: k2f = 9.0_r8/5.0_r8 ! Kelvin-day to Fahrenheit-day conversion ! Gcam expects things in npp_m[lon][lat][pft] format, so we need ! to extract from the flattened column representation. @@ -71,6 +73,15 @@ subroutine iac_import(x2z, lnd2iac_vars) end do ! global index g end do ! pft index p + ! Scalar per-gridcell fields + do g=iac_ctl%begg,iac_ctl%endg + i=iac_ctl%ilon(g) + j=iac_ctl%jlat(g) + lnd2iac_vars%forc_hdm(i,j) = 0._r8; if (index_x2z_Sl_forc_hdm > 0) lnd2iac_vars%forc_hdm(i,j) = x2z(index_x2z_Sl_forc_hdm, g) + lnd2iac_vars%hdd(i,j) = 0._r8; if (index_x2z_Sl_hdd > 0) lnd2iac_vars%hdd(i,j) = x2z(index_x2z_Sl_hdd,g) * dayspy * k2f + lnd2iac_vars%cdd(i,j) = 0._r8; if (index_x2z_Sl_cdd > 0) lnd2iac_vars%cdd(i,j) = x2z(index_x2z_Sl_cdd,g) * dayspy * k2f + end do ! global index g + end subroutine iac_import !=============================================================================== diff --git a/iac/coupling/gcam_comp_mod.F90 b/iac/coupling/gcam_comp_mod.F90 index 1a54fff..237540e 100644 --- a/iac/coupling/gcam_comp_mod.F90 +++ b/iac/coupling/gcam_comp_mod.F90 @@ -236,6 +236,8 @@ subroutine gcam_init_mod(gcamo, gcamoemis,gcamoco2sfcjan, gcamoco2sfcfeb, & gcam2elm_woodharvest_mapping_file(char_len+1:char_len+1) = c_null_char char_len = len_trim(gcam2elm_cdensity_mapping_file) gcam2elm_cdensity_mapping_file(char_len+1:char_len+1) = c_null_char + char_len = len_trim(gcam2elm_degdays_mapping_file) + gcam2elm_degdays_mapping_file(char_len+1:char_len+1) = c_null_char char_len = len_trim(base_co2_surface_file) base_co2_surface_file(char_len+1:char_len+1) = c_null_char char_len = len_trim(base_co2_shipment_file) @@ -281,6 +283,7 @@ subroutine gcam_init_mod(gcamo, gcamoemis,gcamoco2sfcjan, gcamoco2sfcfeb, & trim(gcam2elm_luc_mapping_file),& trim(gcam2elm_woodharvest_mapping_file),& trim(gcam2elm_cdensity_mapping_file),& + trim(gcam2elm_degdays_mapping_file),& trim(base_gcam_co2_file), trim(base_co2_surface_file),& trim(base_co2_shipment_file), trim(base_co2_aircraft_file),& iac_ctl%area, iac_ctl%nlon, iac_ctl%nlat, num_gcam_energy_regions, num_emiss_sectors, rr) @@ -360,12 +363,14 @@ subroutine gcam_run_mod(gcamo, gcamoemis, gcamoco2sfcjan, gcamoco2sfcfeb, & ! !LOCAL VARIABLES: integer :: ymd, tod, dt - integer :: i,j,wc,gs,cs,rs,ws,rr,ays + integer :: i,j,wc,gs,cs,rs,ws,rdd,wdd,rr,ays,dd character(len=256) :: scalar_source_dir_loc character(len=256) :: elm2gcam_mapping_file_loc character(len=256) :: base_npp_file_loc character(len=256) :: base_hr_file_loc character(len=256) :: base_pft_file_loc + character(len=256) :: base_hdd_file_loc + character(len=256) :: base_cdd_file_loc character(len=*),parameter :: subname='(gcam_run_mod)' @@ -412,6 +417,18 @@ subroutine gcam_run_mod(gcamo, gcamoemis, gcamoco2sfcjan, gcamoco2sfcfeb, & ws = 0 end if + ! for reading/writing hdd and cdd from/to diagnostic files + if ( read_hdd_cdd ) then + rdd = 1 + else + rdd = 0 + end if + if ( write_hdd_cdd ) then + wdd = 1 + else + wdd = 0 + end if + ! for ag yield scaling if ( elm_ehc_agyield_scaling ) then ays = 1 @@ -426,6 +443,13 @@ subroutine gcam_run_mod(gcamo, gcamoemis, gcamoco2sfcjan, gcamoco2sfcfeb, & cs = 0 end if + ! for passing hdd and cdd from elm to ehc + if ( elm_ehc_hdd_cdd ) then + dd = 1 + else + dd = 0 + end if + ! get some file names for scalars ! use local variables to avoid adding multiple null characters to the orig scalar_source_dir_loc=trim(scalar_source_dir)//c_null_char @@ -433,14 +457,16 @@ subroutine gcam_run_mod(gcamo, gcamoemis, gcamoco2sfcjan, gcamoco2sfcfeb, & base_npp_file_loc=trim(base_npp_file)//c_null_char base_hr_file_loc=trim(base_hr_file)//c_null_char base_pft_file_loc=trim(base_pft_file)//c_null_char + base_hdd_file_loc=trim(base_hdd_file)//c_null_char + base_cdd_file_loc=trim(base_cdd_file)//c_null_char ! Call runcGCAM method of E3SM Interface ! The yields and carbon density scalars are set within this function also call runcGCAM(ymd, gcamo, gcamoemis, trim(base_gcam_lu_wh_file), trim(base_gcam_co2_file), gs, & - iac_ctl%area, lnd2iac_vars%pftwgt, lnd2iac_vars%npp, lnd2iac_vars%hr, & - iac_ctl%nlon, iac_ctl%nlat, iac_ctl%npft, num_gcam_energy_regions, num_emiss_ctys, num_emiss_sectors, num_periods,& - elm2gcam_mapping_file_loc, iac_first_coupled_year, rs, scalar_source_dir_loc, ws, ays, cs,& - base_npp_file_loc, base_hr_file_loc, base_pft_file_loc, rr) + iac_ctl%area, lnd2iac_vars%pftwgt, lnd2iac_vars%npp, lnd2iac_vars%hr, lnd2iac_vars%hdd, lnd2iac_vars%cdd, lnd2iac_vars%forc_hdm, & + iac_ctl%landfrac, iac_ctl%nlon, iac_ctl%nlat, iac_ctl%npft, num_gcam_energy_regions, num_emiss_ctys, num_emiss_sectors, num_periods,& + elm2gcam_mapping_file_loc, iac_first_coupled_year, rs, scalar_source_dir_loc, ws, rdd, wdd, ays, cs, dd, & + base_npp_file_loc, base_hr_file_loc, base_pft_file_loc, base_hdd_file_loc, base_cdd_file_loc, rr) ! If co2 emissions need to be passed from GCAM to EAM, then call downscale CO2 if ( ehc_eam_co2_emissions ) then diff --git a/iac/coupling/gcam_var_mod.F90 b/iac/coupling/gcam_var_mod.F90 index b226033..7e7e1a2 100644 --- a/iac/coupling/gcam_var_mod.F90 +++ b/iac/coupling/gcam_var_mod.F90 @@ -127,10 +127,13 @@ module gcam_var_mod character(len=256), public :: base_npp_file ! by pft, lat, lon, for land scalars character(len=256), public :: base_hr_file ! by pft, lat, lon, for land scalars character(len=256), public :: base_pft_file ! by pft, lat, lon, for land scalars + character(len=256), public :: base_hdd_file ! by lat, lon, for baseline hdd + character(len=256), public :: base_cdd_file ! by lat, lon, for baseline cdd character(len=256), public :: gcam2elm_co2_mapping_file ! define gcam co2 out character(len=256), public :: gcam2elm_luc_mapping_file ! def gcam lu out character(len=256), public :: gcam2elm_woodharvest_mapping_file ! def gcam wh out character(len=256), public :: gcam2elm_cdensity_mapping_file ! def gcam cdensity in + character(len=256), public :: gcam2elm_degdays_mapping_file ! def gcam degday in ! Grid and mapping and initialization files character(len=256), public :: gcam_gridfile ! definition of iac grid @@ -165,6 +168,8 @@ module gcam_var_mod character(len=256), public :: scalar_source_dir logical, public :: write_scalars ! scalars will be written to a file. ! hr, area, pft weight) are passed from e3sm. + logical, public :: read_hdd_cdd ! if .FALSE., HDD/CDD are calculated; if .TRUE. HDD/CDD are read from file + logical, public :: write_hdd_cdd ! hdd cdd will be written to a file. logical, public :: write_co2 ! gridded co2 emissions will be ! written to a file (in addition to passed in code). @@ -174,6 +179,8 @@ module gcam_var_mod ! land productivity from elm scale ag yield in gcam. logical, public :: elm_ehc_carbon_scaling ! if .TRUE., changes in ! land productivity from elm scale carbon density in gcam. + logical, public :: elm_ehc_hdd_cdd ! if .TRUE., HDD and CDD + ! passed from elm to gcam. logical, public :: ehc_eam_co2_emissions ! if .TRUE., energy system ! co2 is passed from gcam to eam. diff --git a/iac/coupling/iac_data_mod.F90 b/iac/coupling/iac_data_mod.F90 index 2ab5619..4f99aeb 100644 --- a/iac/coupling/iac_data_mod.F90 +++ b/iac/coupling/iac_data_mod.F90 @@ -61,6 +61,9 @@ module iac_data_mod real(r8), allocatable :: npp(:,:,:) real(r8), allocatable :: hr(:,:,:) real(r8), allocatable :: pftwgt(:,:,:) + real(r8), allocatable :: forc_hdm(:,:) ! population density (lon,lat) + real(r8), allocatable :: hdd(:,:) + real(r8), allocatable :: cdd(:,:) end type lnd2iac_type type, public :: iac2lnd_type diff --git a/iac/coupling/iac_init_mod.F90 b/iac/coupling/iac_init_mod.F90 index 7cebe43..d52ef64 100644 --- a/iac/coupling/iac_init_mod.F90 +++ b/iac/coupling/iac_init_mod.F90 @@ -59,8 +59,10 @@ subroutine iac_init(EClock) gcam_config, base_gcam_co2_file, base_gcam_lu_wh_file, & base_co2_surface_file, base_co2_shipment_file, base_co2_aircraft_file, & base_npp_file, base_hr_file, base_pft_file, & + base_hdd_file, base_cdd_file, & gcam2elm_co2_mapping_file, gcam2elm_luc_mapping_file, & gcam2elm_woodharvest_mapping_file, gcam2elm_cdensity_mapping_file, & + gcam2elm_degdays_mapping_file, & gcam_gridfile, elm2gcam_mapping_file, & gcam2glm_glumap, gcam2glm_baselu, gcam2glm_basebiomass, & country2grid_map, country2region_map, pop_iiasa_file, gdp_iiasa_file, & @@ -70,8 +72,8 @@ subroutine iac_init(EClock) pasture_addtreeonly, pasture_setherbfracrem, pasture_setavailtreefracrem, & fdyndat_ehc, & read_scalars, scalar_source_dir, & - write_scalars, write_co2, & - elm_ehc_agyield_scaling, elm_ehc_carbon_scaling, ehc_eam_co2_emissions, & + write_scalars, read_hdd_cdd, write_hdd_cdd, write_co2, & + elm_ehc_agyield_scaling, elm_ehc_carbon_scaling, elm_ehc_hdd_cdd, ehc_eam_co2_emissions, & gcam_spinup, run_gcam nlfilename_iac = "gcam_in" @@ -135,10 +137,13 @@ subroutine iac_init(EClock) write(iulog, '(A,A)') "base_npp_file = ", trim(base_npp_file ) write(iulog, '(A,A)') "base_hr_file = ", trim(base_hr_file) write(iulog, '(A,A)') "base_pft_file = ", trim(base_pft_file ) + write(iulog, '(A,A)') "base_hdd_file = ", trim(base_hdd_file ) + write(iulog, '(A,A)') "base_cdd_file = ", trim(base_cdd_file ) write(iulog, '(A,A)') "gcam2elm_co2_mapping_file = ", trim(gcam2elm_co2_mapping_file ) write(iulog, '(A,A)') "gcam2elm_luc_mapping_file = ", trim(gcam2elm_luc_mapping_file) write(iulog, '(A,A)') "gcam2elm_woodharvest_mapping_file = ", trim(gcam2elm_woodharvest_mapping_file) write(iulog, '(A,A)') "gcam2elm_cdensity_mapping_file = ", trim(gcam2elm_cdensity_mapping_file) + write(iulog, '(A,A)') "gcam2elm_degdays_mapping_file = ", trim(gcam2elm_degdays_mapping_file) write(iulog,*) 'grid mapping and initialization files:' write(iulog, '(A,A)') "gcam_gridfile = ", trim(gcam_gridfile) @@ -173,12 +178,27 @@ subroutine iac_init(EClock) write(iulog, '(A,L)') "read_scalars = ",read_scalars write(iulog, '(A,A)') "scalar_source_dir = ", trim(scalar_source_dir) write(iulog, '(A,L10)') "write_scalars = ",write_scalars + write(iulog, '(A,L10)') "read_hdd_cdd = ",read_hdd_cdd + write(iulog, '(A,L10)') "write_hdd_cdd = ",write_hdd_cdd write(iulog, '(A,L10)') "write_co2 = ",write_co2 write(iulog, '(A,L10)') "elm_ehc_agyield_scaling = ", elm_ehc_agyield_scaling write(iulog, '(A,L10)') "elm_ehc_carbon_scaling = ", elm_ehc_carbon_scaling + write(iulog, '(A,L10)') "elm_ehc_hdd_cdd = ", elm_ehc_hdd_cdd write(iulog, '(A,L10)') "ehc_eam_co2_emissions = ", ehc_eam_co2_emissions write(iulog, '(A,L10)') "gcam_spinup = ",gcam_spinup write(iulog, '(A,L10)') "run_gcam = ",run_gcam + if (read_scalars .and. write_scalars) then + write(iulog,*) '('//trim(subname)//') WARNING: read_scalars and write_scalars are both true.' + write(iulog,*) '('//trim(subname)//') read_scalars takes precedence, so write_scalars is effectively false:' + write(iulog,*) '('//trim(subname)//') no scalar calculations or writes will be performed.' + write_scalars = .false. + end if + if (read_hdd_cdd .and. write_hdd_cdd) then + write(iulog,*) '('//trim(subname)//') WARNING: read_hdd_cdd and write_hdd_cdd are both true.' + write(iulog,*) '('//trim(subname)//') read_hdd_cdd takes precedence, so write_hdd_cdd is effectively false:' + write(iulog,*) '('//trim(subname)//') no HDD/CDD calculations or writes will be performed.' + write_hdd_cdd = .false. + end if !if (nsrest == nsrStartup .and. finidat_rtm /= ' ') then ! write(iulog,*) ' MOSART initial data = ',trim(finidat_rtm) @@ -238,6 +258,9 @@ subroutine iac_init(EClock) allocate(lnd2iac_vars%npp(iac_ctl%nlon,iac_ctl%nlat,iac_ctl%npft)) allocate(lnd2iac_vars%hr(iac_ctl%nlon,iac_ctl%nlat,iac_ctl%npft)) allocate(lnd2iac_vars%pftwgt(iac_ctl%nlon,iac_ctl%nlat,iac_ctl%npft)) + allocate(lnd2iac_vars%forc_hdm(iac_ctl%nlon,iac_ctl%nlat)) + allocate(lnd2iac_vars%hdd(iac_ctl%nlon,iac_ctl%nlat)) + allocate(lnd2iac_vars%cdd(iac_ctl%nlon,iac_ctl%nlat)) allocate(iac2lnd_vars%pct_pft(iac_ctl%nlon,iac_ctl%nlat,iac_ctl%npft)) allocate(iac2lnd_vars%pct_pft_prev(iac_ctl%nlon,iac_ctl%nlat,iac_ctl%npft)) diff --git a/iac/gcam b/iac/gcam index ed35249..4d02db1 160000 --- a/iac/gcam +++ b/iac/gcam @@ -1 +1 @@ -Subproject commit ed35249e6892b33fc3418c9aa4a92641e852f3a1 +Subproject commit 4d02db1ed8cf08b4f23a4408aff96fd3091b4a60