Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Output sea flux species dry deposition velocity with DryDep collection #2606

Open
wants to merge 5 commits into
base: dev/14.6.0
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -11,9 +11,11 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),

### Changed
- Moved PINO3H to be in alphabetical order in `species_database.yml`
- Moved dry dep velocity diagnostic outputs for sea flux and satellite diagnostic species into the `DryDep` collection

### Fixed
- Fixed CEDS HEMCO_Config.rc entries to emit TMB into the TMB species (and not HCOOH)
- In `hco_interface_gc_mod.F90`, update `SatDiagnColEmis` and `SatDiagnSurfFlux` arrays with `(I,J,S)` instead of `(:,:,S)`

## [14.5.0] - 2024-11-07
### Added
Expand Down
72 changes: 44 additions & 28 deletions GeosCore/drydep_mod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -36,7 +36,6 @@ MODULE DRYDEP_MOD
#if defined( MODEL_CESM )
PUBLIC :: UPDATE_DRYDEPFREQ
#else

!
! !PRIVATE MEMBER FUNCTIONS:
!
Expand Down Expand Up @@ -610,10 +609,10 @@ SUBROUTINE UPDATE_DRYDEPFREQ( Input_Opt, State_Chm, State_Diag, State_Grid, &
ENDDO
!$OMP END PARALLEL DO

! Set diagnostics - cnsider moving?
IF ( State_Diag%Archive_DryDepVel .or. &
State_Diag%Archive_DryDepVelForALT1 .or. &
State_Diag%Archive_SatDiagnDryDepVel ) THEN
! Set diagnostics - consider moving?
IF ( State_Diag%Archive_DryDepVel .or. &
State_Diag%Archive_DryDepVelForALT1 .or. &
State_Diag%Archive_SatDiagnDryDepVel ) THEN

!$OMP PARALLEL DO &
!$OMP DEFAULT( SHARED )&
Expand All @@ -623,22 +622,38 @@ SUBROUTINE UPDATE_DRYDEPFREQ( Input_Opt, State_Chm, State_Diag, State_Grid, &
! Point to State_Chm%DryDepVel [m/s]
NDVZ = NDVZIND(D)

! Dry dep velocity [cm/s]
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Removing these diagnostics from drydep_mod and only setting in compute_sflx_for_vdiff means they will always be zero if using full mixing.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This should now be fixed in commits d9ee2e2 and c7cfcfa.

IF ( State_Diag%Archive_DryDepVel ) THEN
S = State_Diag%Map_DryDepVel%id2slot(D)
IF ( S > 0 ) THEN
State_Diag%DryDepVel(:,:,S) = &
State_Chm%DryDepVel(:,:,NDVZ) * 100.0_f4
! When using the full PBL mixing option (aka TURBDAY),
! update DryDepVel and SatDiagnDryDepVel diagnostics here.
!
! When using the non-local PBL mixing option (VDIFF), then
! update these diagnostics in Compute_Sflx_for_Vdiff (in
! GeosCore/hco_interface_gc_mod.F90). This is necessary in
! order to capture the air-sea deposition velocities computed
! by the HEMCO "SeaFlux" extension. For more information,
! see https://github.com/geoschem/geos-chem/pull/2606 and
! https://github.com/geoschem/geos-chem/issues/2564.
!
! -- Bob Yantosca (03 Dec 2024)
IF ( .not. Input_Opt%LNLPBL ) THEN

! Dry dep velocity [cm/s]
IF ( State_Diag%Archive_DryDepVel ) THEN
S = State_Diag%Map_DryDepVel%id2slot(D)
IF ( S > 0 ) THEN
State_Diag%DryDepVel(:,:,S) = &
State_Chm%DryDepVel(:,:,NDVZ) * 100.0_f4
ENDIF
ENDIF
ENDIF

! Satellite diagnostic: Dry dep velocity [cm/s]
IF ( State_Diag%Archive_SatDiagnDryDepVel ) THEN
S = State_Diag%Map_SatDiagnDryDepVel%id2slot(D)
IF ( S > 0 ) THEN
State_Diag%SatDiagnDryDepVel(:,:,S) = &
State_Chm%DryDepVel(:,:,NDVZ) * 100.0_f4
! Satellite diagnostic: Dry dep velocity [cm/s]
IF ( State_Diag%Archive_SatDiagnDryDepVel ) THEN
S = State_Diag%Map_SatDiagnDryDepVel%id2slot(D)
IF ( S > 0 ) THEN
State_Diag%SatDiagnDryDepVel(:,:,S) = &
State_Chm%DryDepVel(:,:,NDVZ) * 100.0_f4
ENDIF
ENDIF

ENDIF

! Dry dep velocity [cm/s] for species at altitude (e.g. 10m)
Expand Down Expand Up @@ -1640,9 +1655,9 @@ SUBROUTINE DEPVEL( Input_Opt, State_Chm, State_Diag, State_Grid, &
! - Bob Yantosca (17 May 2023)
IF ( N_SPC == ID_Hg0 ) THEN

! Assume lower reactivity
F0_K = 3.0e-05_f8
! Assume lower reactivity
F0_K = 3.0e-05_f8

! But if this is the rainforest land type and we fall
! within the bounding box of the Amazon rainforest,
! then increase reactivity as inferred from observations.
Expand Down Expand Up @@ -1738,10 +1753,10 @@ SUBROUTINE DEPVEL( Input_Opt, State_Chm, State_Diag, State_Grid, &
RHB(I,J), &
W10(I,J), &
VTSoutput, &
Input_Opt, &
Input_Opt, &
State_Chm )

VTSoutput_(K,LDT) = VTSoutput
VTSoutput_(K,LDT) = VTSoutput

ELSE IF ( SpcInfo%DD_DustDryDep ) THEN

Expand Down Expand Up @@ -2254,6 +2269,7 @@ SUBROUTINE DEPVEL( Input_Opt, State_Chm, State_Diag, State_Grid, &
!** Load array State_Chm%DryDepVel
DO 550 K=1,NUMDEP
IF (.NOT.LDEP(K)) GOTO 550

State_Chm%DryDepVel(I,J,K) = VD(K)

! Now check for negative deposition velocity before returning to
Expand Down Expand Up @@ -3217,7 +3233,7 @@ FUNCTION AERO_SFCRSII( K, II, PRESS, TEMP, USTAR, RHB, W10, VTSout, Input_Opt, S
TYPE(ChmState), INTENT(IN) :: State_Chm ! Chemistry State object

!
! !OUTPUT PARAMETERS:
! !OUTPUT PARAMETERS:
!
REAL(f8), INTENT(OUT) :: VTSout ! Settling velocity [m/s]

Expand Down Expand Up @@ -3984,7 +4000,7 @@ END FUNCTION DUST_SFCRSI
FUNCTION ADUST_SFCRSII( K, II, PRESS, TEMP, USTAR, &
VTSout, RHB, State_Chm ) RESULT( RS )
!
! !USES:
! !USES:
!
USE Species_Mod, ONLY : Species
USE State_Chm_Mod, ONLY : ChmState
Expand All @@ -4001,7 +4017,7 @@ FUNCTION ADUST_SFCRSII( K, II, PRESS, TEMP, USTAR, &
TYPE(ChmState), INTENT(IN) :: State_Chm ! Chemistry State object

!
! !OUTPUT PARAMETERS:
! !OUTPUT PARAMETERS:
!
REAL(f8), INTENT(OUT) :: VTSout ! Settling velocity [m/s]

Expand Down Expand Up @@ -4217,7 +4233,7 @@ FUNCTION ADUST_SFCRSII( K, II, PRESS, TEMP, USTAR, &
!BC
ELSE IF ( K == idd_BCPI .OR. K == idd_BCPO ) THEN
! DIAM is not changed

!OA
ELSE
DIAM = DIAM * ((1.0_fp + 0.1_fp * RHBL / (1.0_fp - RHBL)) &
Expand Down Expand Up @@ -4339,7 +4355,7 @@ FUNCTION DUST_SFCRSII( K, II, PRESS, TEMP, USTAR, DIAM, &
REAL(f8), INTENT(IN) :: DEN ! Particle density [kg/m3]

!
! !OUTPUT PARAMETERS:
! !OUTPUT PARAMETERS:
!
REAL(f8), INTENT(OUT) :: VTSout ! Settling velocity [m/s]

Expand Down
54 changes: 49 additions & 5 deletions GeosCore/hco_interface_gc_mod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -4563,6 +4563,9 @@ SUBROUTINE Compute_Sflx_for_Vdiff( Input_Opt, State_Chm, State_Diag, &
REAL(fp), TARGET :: dflx(State_Grid%NX, &
State_Grid%NY, &
State_Chm%nAdvect )
REAL(fp), TARGET :: dvel(State_Grid%NX, &
State_Grid%NY, &
State_Chm%nAdvect )

! Pointers and Objects
REAL(f4), POINTER :: Ptr2D(:,:) => NULL()
Expand All @@ -4581,6 +4584,7 @@ SUBROUTINE Compute_Sflx_for_Vdiff( Input_Opt, State_Chm, State_Diag, &
RC = GC_SUCCESS
dflx = 0.0_fp
eflx = 0.0_fp
dvel = 0.0_fp
colEflx = 0.0_fp
ThisSpc => NULL()
errMsg = ''
Expand Down Expand Up @@ -4777,6 +4781,9 @@ SUBROUTINE Compute_Sflx_for_Vdiff( Input_Opt, State_Chm, State_Diag, &
dflx(I,J,NA) = dflx(I,J,NA) + ( dep &
* State_Chm%Species(NA)%Conc(I,J,1) &
/ (AIRMW / ThisSpc%MW_g) )
! Get dry deposition velocity in cm/s:
dvel(I,J,NA) = dvel(I,J,NA) + ( dep * State_Met%BXHEIGHT(I,J,1) * &
1.0e+2_fp )
ENDIF
ENDIF
ENDDO ! I
Expand Down Expand Up @@ -4838,6 +4845,9 @@ SUBROUTINE Compute_Sflx_for_Vdiff( Input_Opt, State_Chm, State_Diag, &
* State_Chm%Species(N)%Conc(I,J,1) &
/ ( AIRMW / ThisSpc%MW_g )

! Get dry deposition velocity in cm/s:
dvel(I,J,N) = dvel(I,J,N) + ( State_Chm%DryDepFreq(I,J,ND) * &
State_Met%BXHEIGHT(I,J,1) * 1.0e+2_fp )

IF ( Input_Opt%ITS_A_MERCURY_SIM .and. ThisSpc%Is_Hg0 ) THEN

Expand All @@ -4858,8 +4868,10 @@ SUBROUTINE Compute_Sflx_for_Vdiff( Input_Opt, State_Chm, State_Diag, &

! Free species database pointer
ThisSpc => NULL()

ENDDO


!=====================================================================
! Convert DFLX from 1/s to kg/m2/s
!
Expand Down Expand Up @@ -4891,12 +4903,11 @@ SUBROUTINE Compute_Sflx_for_Vdiff( Input_Opt, State_Chm, State_Diag, &
!=====================================================================
! Defining Satellite Diagnostics
!=====================================================================

! Define emission satellite diagnostics
IF ( State_Diag%Archive_SatDiagnColEmis ) THEN
DO S = 1, State_Diag%Map_SatDiagnColEmis%nSlots
N = State_Diag%Map_SatDiagnColEmis%slot2id(S)
State_Diag%SatDiagnColEmis(:,:,S) = colEflx(:,:,N)
State_Diag%SatDiagnColEmis(I,J,S) = colEflx(I,J,N)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is a bug fix we should note in the changelog. These diagnostics must have been many times larger than they should be before this fix.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@lizziel: We can test this to make sure we get zero-diff output for this diagnostic

ENDDO
ENDIF

Expand All @@ -4907,7 +4918,7 @@ SUBROUTINE Compute_Sflx_for_Vdiff( Input_Opt, State_Chm, State_Diag, &
IF ( State_Diag%Archive_SatDiagnSurfFlux ) THEN
DO S = 1, State_Diag%Map_SatDiagnSurfFlux%nSlots
N = State_Diag%Map_SatDiagnSurfFlux%slot2id(S)
State_Diag%SatDiagnSurfFlux(:,:,S) = State_Chm%SurfaceFlux(:,:,N)
State_Diag%SatDiagnSurfFlux(I,J,S) = State_Chm%SurfaceFlux(I,J,N)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Same comment as above. Note that both of these diagnostics have bugs that are now fixed.

Copy link
Contributor Author

@yantosca yantosca Dec 3, 2024

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Changelog fixed in commit d9ee2e2.

ENDDO
ENDIF

Expand Down Expand Up @@ -4951,7 +4962,7 @@ SUBROUTINE Compute_Sflx_for_Vdiff( Input_Opt, State_Chm, State_Diag, &
ENDIF

ENDDO
ENDDO
ENDDO
!$OMP END PARALLEL DO

!### Uncomment for debug output
Expand All @@ -4970,7 +4981,9 @@ SUBROUTINE Compute_Sflx_for_Vdiff( Input_Opt, State_Chm, State_Diag, &
! so we don't need to do any further special handling.
!=======================================================================
IF ( Input_Opt%LGTMM .or. Input_Opt%LSOILNOX .or. &
State_Diag%Archive_DryDepMix .or. State_Diag%Archive_DryDep ) THEN
State_Diag%Archive_DryDepMix .or. State_Diag%Archive_DryDep .or. &
State_Diag%Archive_DryDepVel .or. &
State_Diag%Archive_SatDiagnDryDepVel ) THEN

! Loop over only the drydep species
! If drydep is turned off, nDryDep=0 and the loop won't execute
Expand Down Expand Up @@ -5016,6 +5029,36 @@ SUBROUTINE Compute_Sflx_for_Vdiff( Input_Opt, State_Chm, State_Diag, &
ENDIF
ENDIF

!-----------------------------------------------------------------
! HISTORY: Update dry deposition velocity [cm/s]
!
! Here we save the dry deposition velocities (stored in the
! DVEL array) into the DryDepVel and SatDiagnDryDepVel History
! diagnostics. This is necessary in order to capture the
! air-sea deposition velocity computed by the HEMCO "SeaFlux"
! extension for certain species.
!
! When using the full PBL mixing option (aka TURBDAY), the
! DryDepVel and SatDiagnDryDepVel diagnostics will be archived
! in drydep_mod.F90 instead.
!-----------------------------------------------------------------

! Dry deposition velocity [cm/s]
IF ( State_Diag%Archive_DryDepVel ) THEN
S = State_Diag%Map_DryDepVel%id2slot(ND)
IF ( S > 0 ) THEN
State_Diag%DryDepVel(:,:,S) = dvel(:,:,N)
ENDIF
ENDIF

! Satellite diagnostic dry deposition velocity (cm/s):
IF ( State_Diag%Archive_SatDiagnDryDepVel ) THEN
S = State_Diag%Map_SatDiagnDryDepVel%id2slot(ND)
IF ( S > 0 ) THEN
State_Diag%SatDiagnDryDepVel(:,:,S) = dvel(:,:,N)
ENDIF
ENDIF

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

As I noted above, putting these here would only work for non-local mixing. Full mixing does not call this routine.

Copy link
Contributor Author

@yantosca yantosca Dec 3, 2024

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This should now be fixed in commits d9ee2e2 and c7cfcfa.

!-----------------------------------------------------------------
! If Soil NOx is turned on, then call SOIL_DRYDEP to
! archive dry deposition fluxes for nitrogen species
Expand All @@ -5039,6 +5082,7 @@ SUBROUTINE Compute_Sflx_for_Vdiff( Input_Opt, State_Chm, State_Diag, &
ThisSpc => NULL()
ENDDO
!$OMP END PARALLEL DO

ENDIF

!=======================================================================
Expand Down