From c21bb8c57f20c2fe884ef4f0b9c49a60896f2f3f Mon Sep 17 00:00:00 2001 From: Peter Colarco Date: Mon, 8 Jan 2024 15:28:37 -0500 Subject: [PATCH] prc: updates to include separate explosive and degassing volcanic emissions for SO2 --- CHANGELOG.md | 13 +++ .../AMIP.20C/SU2G_instance_SU.rc | 3 +- .../SU2G_GridComp/AMIP/SU2G_instance_SU.rc | 3 +- .../SU2G_GridComp/SU2G_GridCompMod.F90 | 85 +++++++++++++++---- .../SU2G_GridComp/SU2G_instance_SU.rc | 3 +- Process_Library/GOCART2G_Process.F90 | 23 ++--- 6 files changed, 94 insertions(+), 36 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 25023cc6..4cf021c7 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -13,6 +13,19 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Changed +## [Unreleased - Peter.R.Colarco@nasa.gov] - 2024-01-08 + +### Changed + +- Changed SU2G_instance_SU.rc to now have separate filename inputs for explosive and degassing volcanoes +- Moved present volcanic emission inventories to one or the other line for these new entries; set other + line /dev/null; this is stop gap until next time we update volcanic emission inventories, at which + point will provide (for AMIP and AMIP.20C) separate explosive and degassing emissions +- Made accommodating changes for above in SU2G_GridCompMod.F90 and in the Process Library +- Verified zero diff in current configuration (this is true of tracers and restarts, but not diagnostics: + until an actual split is made in the input emissions then the volcanic emissions are being assigned to + one or the other emission diagnostics (explosive or degassing). + ## [v2.2.1] - 2023-05-30 ### Fixed diff --git a/ESMF/GOCART2G_GridComp/SU2G_GridComp/AMIP.20C/SU2G_instance_SU.rc b/ESMF/GOCART2G_GridComp/SU2G_GridComp/AMIP.20C/SU2G_instance_SU.rc index 1901e354..cdf30ab4 100644 --- a/ESMF/GOCART2G_GridComp/SU2G_GridComp/AMIP.20C/SU2G_instance_SU.rc +++ b/ESMF/GOCART2G_GridComp/SU2G_GridComp/AMIP.20C/SU2G_instance_SU.rc @@ -8,7 +8,8 @@ aerosol_monochromatic_optics_file: ExtData/chemistry/AerosolOptics/v0.0.0/x/opti nbins: 4 # Volcanic pointwise sources -volcano_srcfilen: ExtData/chemistry/CARN/v202106/sfc/so2_volcanic_emissions_Carns.%y4%m2%d2.rc +volcano_srcfilen_explosive: /dev/null +volcano_srcfilen_degassing: ExtData/chemistry/CARN/v202106/sfc/so2_volcanic_emissions_Carns.%y4%m2%d2.rc # Heights [m] of LTO, CDS and CRS aviation emissions layers aviation_vertical_layers: 0.0 100.0 9.0e3 10.0e3 diff --git a/ESMF/GOCART2G_GridComp/SU2G_GridComp/AMIP/SU2G_instance_SU.rc b/ESMF/GOCART2G_GridComp/SU2G_GridComp/AMIP/SU2G_instance_SU.rc index 1901e354..cdf30ab4 100644 --- a/ESMF/GOCART2G_GridComp/SU2G_GridComp/AMIP/SU2G_instance_SU.rc +++ b/ESMF/GOCART2G_GridComp/SU2G_GridComp/AMIP/SU2G_instance_SU.rc @@ -8,7 +8,8 @@ aerosol_monochromatic_optics_file: ExtData/chemistry/AerosolOptics/v0.0.0/x/opti nbins: 4 # Volcanic pointwise sources -volcano_srcfilen: ExtData/chemistry/CARN/v202106/sfc/so2_volcanic_emissions_Carns.%y4%m2%d2.rc +volcano_srcfilen_explosive: /dev/null +volcano_srcfilen_degassing: ExtData/chemistry/CARN/v202106/sfc/so2_volcanic_emissions_Carns.%y4%m2%d2.rc # Heights [m] of LTO, CDS and CRS aviation emissions layers aviation_vertical_layers: 0.0 100.0 9.0e3 10.0e3 diff --git a/ESMF/GOCART2G_GridComp/SU2G_GridComp/SU2G_GridCompMod.F90 b/ESMF/GOCART2G_GridComp/SU2G_GridComp/SU2G_GridCompMod.F90 index 6c4abc37..22366489 100644 --- a/ESMF/GOCART2G_GridComp/SU2G_GridComp/SU2G_GridCompMod.F90 +++ b/ESMF/GOCART2G_GridComp/SU2G_GridComp/SU2G_GridCompMod.F90 @@ -54,6 +54,9 @@ module SU2G_GridCompMod !=========================================================================== ! !Sulfer state type :: ThreadWorkspace + integer :: nymd_last = -1 ! Previous nymd. Updated daily + +! Degassing volcanoes integer :: nVolc = 0 real, allocatable, dimension(:) :: vLat, & vLon, & @@ -62,7 +65,17 @@ module SU2G_GridCompMod vCloud integer, allocatable, dimension(:) :: vStart, & vEnd - integer :: nymd_last = -1 ! Previous nymd. Updated daily +! Explosive volcanoes + integer :: nVolcE = 0 + real, allocatable, dimension(:) :: vLatE, & + vLonE, & + vSO2E, & + vElevE, & + vCloudE + integer, allocatable, dimension(:) :: vStartE, & + vEndE + +! Other point emissions of Sulfate (SO4) integer :: nPts = -1 integer, allocatable, dimension(:) :: pstart, pend real, allocatable, dimension(:) :: pLat, & @@ -86,8 +99,10 @@ module SU2G_GridCompMod !real, pointer :: h2o2_init(:,:,:) ! Special handling for volcanic emissions - character(len=255) :: volcano_srcfilen -! !Workspae for point emissions + character(len=255) :: volcano_srcfilen_degassing + character(len=255) :: volcano_srcfilen_explosive + +! Workspace for point emissions logical :: doing_point_emissions = .false. character(len=255) :: point_emissions_srcfilen ! filename for pointwise emissions type(ThreadWorkspace), allocatable :: workspaces(:) @@ -174,7 +189,8 @@ subroutine SetServices ( GC, RC ) allocate(self%sigma(self%nbins), __STAT__) ! process SU-specific items - call ESMF_ConfigGetAttribute(cfg, self%volcano_srcfilen, label='volcano_srcfilen:', __RC__) + call ESMF_ConfigGetAttribute(cfg, self%volcano_srcfilen_degassing, label='volcano_srcfilen_degassing:', __RC__) + call ESMF_ConfigGetAttribute(cfg, self%volcano_srcfilen_explosive, label='volcano_srcfilen_explosive:', __RC__) call ESMF_ConfigGetAttribute(cfg, self%eAircraftFuel, label='aircraft_fuel_emission_factor:', __RC__) call ESMF_ConfigGetAttribute(cfg, self%fSO4anth, label='so4_anthropogenic_fraction:', __RC__) call ESMF_ConfigGetAttribute(cfg, self%sigma, label='sigma:', __RC__) @@ -718,7 +734,7 @@ subroutine Run1 (GC, import, export, clock, RC) real, dimension(:,:), allocatable :: so2biomass_src, so2biomass_src_, so2anthro_l1_src, & so2anthro_l2_src, so2ship_src, so4ship_src, dmso_conc, & aviation_lto_src, aviation_cds_src, aviation_crs_src - integer, dimension(:), allocatable :: iPointVolc, jPointVolc, iPoint, jPoint + integer, dimension(:), allocatable :: iPoint, jPoint real, dimension(:,:,:), allocatable :: emissions_point character (len=ESMF_MAXSTR) :: fname ! file name for point source emissions logical :: fileExists @@ -812,34 +828,48 @@ subroutine Run1 (GC, import, export, clock, RC) where(1.01*aviation_cds_src > MAPL_UNDEF ) aviation_cds_src = 0. where(1.01*aviation_crs_src > MAPL_UNDEF ) aviation_crs_src = 0. +! Start with a clean emission diagnostic + if(associated(SUEM)) SUEM = 0.0 + + ! Update emissions/production if necessary (daily) ! ----------------------------------------------- thread = MAPL_get_current_thread() workspace => self%workspaces(thread) +! Update Volcanic SO2 Emissions Daily if(workspace%nymd_last /= nymd) then workspace%nymd_last = nymd -! Get pointwise SO2 and altitude of volcanoes from a daily file data base - if(index(self%volcano_srcfilen,'volcanic_') /= 0) then - call StrTemplate(fname, self%volcano_srcfilen, xid='unknown', & +! DEGASSING: Get pointwise SO2 and altitude of volcanoes from a daily file data base + workspace%nVolc = 0 ! case of /dev/null (no volcanoes) or ill-formed filename + if(index(self%volcano_srcfilen_degassing,'volcanic_') /= 0) then + call StrTemplate(fname, self%volcano_srcfilen_degassing, xid='unknown', & nymd=nymd, nhms=120000 ) call ReadPointEmissions (nymd, fname, workspace%nVolc, workspace%vLat, workspace%vLon, & workspace%vElev, workspace%vCloud, workspace%vSO2, workspace%vStart, & workspace%vEnd, label='volcano', __RC__) workspace%vSO2 = workspace%vSO2 * fMassSO2 / fMassSulfur -! Special possible case - if(self%volcano_srcfilen(1:9) == '/dev/null') workspace%nVolc = 0 + end if + +! EXPLOSIVE: Get pointwise SO2 and altitude of volcanoes from a daily file data base + workspace%nVolcE = 0 ! case of /dev/null (no volcanoes) or ill-formed filename + if(index(self%volcano_srcfilen_explosive,'volcanic_') /= 0) then + call StrTemplate(fname, self%volcano_srcfilen_explosive, xid='unknown', & + nymd=nymd, nhms=120000 ) + call ReadPointEmissions (nymd, fname, workspace%nVolcE, workspace%vLatE, workspace%vLonE, & + workspace%vElevE, workspace%vCloudE, workspace%vSO2E, workspace%vStartE, & + workspace%vEndE, label='volcano', __RC__) + workspace%vSO2 = workspace%vSO2 * fMassSO2 / fMassSulfur end if end if -! Apply volcanic emissions -! ------------------------ +! DEGASSING: Apply volcanic emissions +! ----------------------------------- if (workspace%nVolc > 0) then - if (associated(SO2EMVE)) SO2EMVE=0.0 if (associated(SO2EMVN)) SO2EMVN=0.0 - allocate(iPointVolc(workspace%nVolc), jPointVolc(workspace%nVolc), __STAT__) - call MAPL_GetHorzIJIndex(workspace%nVolc, iPointVolc, jPointVolc, & + allocate(iPoint(workspace%nVolc), jPoint(workspace%nVolc), __STAT__) + call MAPL_GetHorzIJIndex(workspace%nVolc, iPoint, jPoint, & grid = grid, & lon = workspace%vLon/real(MAPL_RADIANS_TO_DEGREES), & lat = workspace%vLat/real(MAPL_RADIANS_TO_DEGREES), & @@ -850,8 +880,30 @@ subroutine Run1 (GC, import, export, clock, RC) end if call SUvolcanicEmissions (workspace%nVolc, workspace%vStart, workspace%vEnd, workspace%vSO2, workspace%vElev, & - workspace%vCloud, iPointVolc, jPointVolc, nhms, SO2EMVN, SO2EMVE, SO2, nSO2, SUEM, & + workspace%vCloud, iPoint, jPoint, nhms, SO2EMVN, SO2, nSO2, SUEM, & self%km, self%cdt, MAPL_GRAV, zle, delp, area, workspace%vLat, workspace%vLon, __RC__) + deallocate(iPoint, jPoint, __STAT__) + end if + +! EXPLOSIVE: Apply volcanic emissions +! ----------------------------------- + if (workspace%nVolcE > 0) then + if (associated(SO2EMVE)) SO2EMVE=0.0 + allocate(iPoint(workspace%nVolcE), jPoint(workspace%nVolcE), __STAT__) + call MAPL_GetHorzIJIndex(workspace%nVolcE, iPoint, jPoint, & + grid = grid, & + lon = workspace%vLon/real(MAPL_RADIANS_TO_DEGREES), & + lat = workspace%vLat/real(MAPL_RADIANS_TO_DEGREES), & + rc = status) + if ( status /= 0 ) then + if (mapl_am_i_root()) print*, trim(Iam), ' - cannot get indices for point emissions' + VERIFY_(status) + end if + + call SUvolcanicEmissions (workspace%nVolcE, workspace%vStartE, workspace%vEndE, workspace%vSO2E, workspace%vElevE, & + workspace%vCloudE, iPoint, jPoint, nhms, SO2EMVE, SO2, nSO2, SUEM, & + self%km, self%cdt, MAPL_GRAV, zle, delp, area, workspace%vLatE, workspace%vLonE, __RC__) + deallocate(iPoint, jPoint, __STAT__) end if ! Apply diurnal cycle if so desired @@ -926,6 +978,7 @@ subroutine Run1 (GC, import, export, clock, RC) workspace%pStart, workspace%pEnd, zle, & area, iPoint, jPoint, nhms, emissions_point, __RC__) + deallocate(iPoint, jPoint, __STAT__) SO4 = SO4 + self%cdt * MAPL_GRAV / delp * emissions_point end if diff --git a/ESMF/GOCART2G_GridComp/SU2G_GridComp/SU2G_instance_SU.rc b/ESMF/GOCART2G_GridComp/SU2G_GridComp/SU2G_instance_SU.rc index a842807a..4fd84a32 100644 --- a/ESMF/GOCART2G_GridComp/SU2G_GridComp/SU2G_instance_SU.rc +++ b/ESMF/GOCART2G_GridComp/SU2G_GridComp/SU2G_instance_SU.rc @@ -8,7 +8,8 @@ aerosol_monochromatic_optics_file: ExtData/chemistry/AerosolOptics/v0.0.0/x/opti nbins: 4 # Volcanic pointwise sources -volcano_srcfilen: ExtData/chemistry/CARN/v202106/sfc/so2_volcanic_emissions_CARN_v202106.degassing_only.rc +volcano_srcfilen_explosive: /dev/null +volcano_srcfilen_degassing: ExtData/chemistry/CARN/v202106/sfc/so2_volcanic_emissions_CARN_v202106.degassing_only.rc # Heights [m] of LTO, CDS and CRS aviation emissions layers aviation_vertical_layers: 0.0 100.0 9.0e3 10.0e3 diff --git a/Process_Library/GOCART2G_Process.F90 b/Process_Library/GOCART2G_Process.F90 index 6892082c..d8ff3672 100644 --- a/Process_Library/GOCART2G_Process.F90 +++ b/Process_Library/GOCART2G_Process.F90 @@ -5673,7 +5673,7 @@ end subroutine DMSemission ! !IROUTINE: SUvolcanicEmissions subroutine SUvolcanicEmissions (nVolc, vStart, vEnd, vSO2, vElev, vCloud, iPoint, & - jPoint, nhms, SO2EMVN, SO2EMVE, SO2, nSO2, SU_emis, km, cdt, grav,& + jPoint, nhms, SO2EMVol, SO2, nSO2, SU_emis, km, cdt, grav,& hghte, delp, area, vLat, vLon, rc) ! !USES: implicit NONE @@ -5697,8 +5697,7 @@ subroutine SUvolcanicEmissions (nVolc, vStart, vEnd, vSO2, vElev, vCloud, iPoint real, dimension(:), intent(in) :: vLat ! latitude specified in file [degree] real, dimension(:), intent(in) :: vLon ! longitude specified in file [degree] ! !INOUT PARAMETERS: - real, pointer, dimension(:,:), intent(inout) :: SO2EMVN ! non-explosive volcanic emissions [kg m-2 s-1] - real, pointer, dimension(:,:), intent(inout) :: SO2EMVE ! explosive volcanic emissions [kg m-2 s-1] + real, pointer, dimension(:,:), intent(inout) :: SO2EMVol ! volcanic emissions [kg m-2 s-1] real, pointer, dimension(:,:,:), intent(inout) :: SO2 ! SO2 [kg kg-1] real, pointer, dimension(:,:,:), intent(inout) :: SU_emis ! SU emissions, kg/m2/s real, dimension(:), intent(inout) :: vElev ! bottom elevation of emissions [m] @@ -5722,7 +5721,6 @@ subroutine SUvolcanicEmissions (nVolc, vStart, vEnd, vSO2, vElev, vCloud, iPoint real :: deltaSO2v real, dimension(:,:), allocatable :: z0 real, allocatable, dimension(:,:) :: srcSO2volc - real, allocatable, dimension(:,:) :: srcSO2volce !EOP !------------------------------------------------------------------------- @@ -5731,13 +5729,9 @@ subroutine SUvolcanicEmissions (nVolc, vStart, vEnd, vSO2, vElev, vCloud, iPoint if (nVolc > 0) then allocate(srcSO2volc, mold=area) - allocate(srcSO2volce, mold=area) srcSO2volc = 0. - srcSO2volce = 0. - if (associated(SU_emis)) SU_emis = 0.0 - if (associated(SO2EMVN)) SO2EMVN = 0. - if (associated(SO2EMVE)) SO2EMVE = 0. + if (associated(SO2EMVol)) SO2EMVol = 0. allocate(z0, mold=area) z0 = hghte(:,:,km) @@ -5776,11 +5770,7 @@ subroutine SUvolcanicEmissions (nVolc, vStart, vEnd, vSO2, vElev, vCloud, iPoint ! Diagnostic - sum of volcanos ! ---------------------------- - if (hup .eq. hlow) then - srcSO2volc(i,j) = srcSO2volc(i,j) + so2volcano - else - srcSO2volce(i,j) = srcSO2volce(i,j) + so2volcano - endif + srcSO2volc(i,j) = srcSO2volc(i,j) + so2volcano dzvolc = hup-hlow do k = km, 1, -1 @@ -5827,9 +5817,8 @@ subroutine SUvolcanicEmissions (nVolc, vStart, vEnd, vSO2, vElev, vCloud, iPoint enddo ! it end if ! nVolc > 0 - if (associated(SO2EMVN)) SO2EMVN = SO2EMVN + srcSO2volc - if (associated(SO2EMVE)) SO2EMVE = SO2EMVE + srcSO2volce - if (associated(SU_emis)) SU_emis(:,:,nSO2) = SU_emis(:,:,nSO2) + srcSO2volc + srcSO2volce + if (associated(SO2EMVol)) SO2EMVol = SO2EMVol + srcSO2volc + if (associated(SU_emis)) SU_emis(:,:,nSO2) = SU_emis(:,:,nSO2) + srcSO2volc __RETURN__(__SUCCESS__) end subroutine SUvolcanicEmissions