Skip to content

added capability to save Hplus concentration #144

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.

Sign up for GitHub

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

Open
yutiPML wants to merge 2 commits into
base: master
Choose a base branch
from
Open

added capability to save Hplus concentration #144

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
676f109
added capability to save Hplus concentration
yutiPML Mar 27, 2025
cfd3666
Merge branch 'master' into 143-hydrogen-ions-concentration-is-not-a-d…
wathen Nov 3, 2025
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
8 changes: 5 additions & 3 deletions src/benthic_carbonate.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,7 @@ module ersem_benthic_carbonate
type (type_dependency_id) :: id_ETW, id_X1X, id_dens, id_pres
type (type_horizontal_dependency_id) :: id_Carb_in,id_pco2_in

type (type_horizontal_diagnostic_variable_id) :: id_ph,id_pco2,id_CarbA, id_BiCarb, id_Carb
type (type_horizontal_diagnostic_variable_id) :: id_ph,id_pco2,id_CarbA, id_BiCarb, id_Carb, id_Hplus
type (type_horizontal_diagnostic_variable_id) :: id_Om_cal,id_Om_arg
integer :: phscale

Expand Down Expand Up @@ -45,6 +45,7 @@ subroutine initialize(self,configunit)
call self%register_diagnostic_variable(self%id_CarbA, 'CarbA', 'mmol/m^3','carbonic acid concentration',source=source_do_bottom)
call self%register_diagnostic_variable(self%id_BiCarb,'BiCarb','mmol/m^3','bicarbonate concentration',source=source_do_bottom)
call self%register_diagnostic_variable(self%id_Carb, 'Carb', 'mmol/m^3','carbonate concentration',source=source_do_bottom)
call self%register_diagnostic_variable(self%id_Carb, 'Hplus', 'mmol/m^3','hydrogen ions concentration',source=source_do_bottom)
call self%register_diagnostic_variable(self%id_Om_cal,'Om_cal','-','calcite saturation',source=source_do_bottom)
call self%register_diagnostic_variable(self%id_Om_arg,'Om_arg','-','aragonite saturation',source=source_do_bottom)
call self%register_dependency(self%id_ETW, standard_variables%temperature)
Expand All @@ -62,7 +63,7 @@ subroutine do_bottom(self,_ARGUMENTS_DO_BOTTOM_)

real(rk) :: O3c,ETW,X1X,density,pres
real(rk) :: TA,Ctot
real(rk) :: pH,PCO2,H2CO3,HCO3,CO3,k0co2
real(rk) :: pH,PCO2,H2CO3,HCO3,CO3,hplus,k0co2
real(rk) :: Om_cal,Om_arg
logical :: success

Expand All @@ -77,7 +78,7 @@ subroutine do_bottom(self,_ARGUMENTS_DO_BOTTOM_)

TA = TA/1.0e6_rk ! from umol kg-1 to mol kg-1
Ctot = O3C / 1.e3_rk / density ! from mmol m-3 to mol kg-1
call co2dyn (ETW,X1X,pres*0.1_rk,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,k0co2,success,self%phscale) ! NB pressure from dbar to bar
call co2dyn (ETW,X1X,pres*0.1_rk,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,Hplus,k0co2,success,self%phscale) ! NB pressure from dbar to bar
if (success) then
! Carbonate system iterative scheme converged - save associated diagnostics.
! Convert outputs from fraction to ppm (pCO2) and from mol kg-1 to mmol m-3 (concentrations).
Expand All @@ -86,6 +87,7 @@ subroutine do_bottom(self,_ARGUMENTS_DO_BOTTOM_)
_SET_HORIZONTAL_DIAGNOSTIC_(self%id_CarbA, H2CO3*1.e3_rk*density)
_SET_HORIZONTAL_DIAGNOSTIC_(self%id_Bicarb,HCO3*1.e3_rk*density)
_SET_HORIZONTAL_DIAGNOSTIC_(self%id_Carb, CO3*1.e3_rk*density)
_SET_HORIZONTAL_DIAGNOSTIC_(self%id_Hplus, Hplus) !*1.e3_rk*density)
else
! Carbonate system iterative scheme did not converge.
! All diagnostics retain their previous value.
Expand Down
31 changes: 18 additions & 13 deletions src/carbonate.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -13,10 +13,10 @@ module ersem_carbonate
! Variable identifiers
type (type_state_variable_id) :: id_O3c,id_TA,id_bioalk
type (type_dependency_id) :: id_ETW, id_X1X, id_dens, id_pres
type (type_dependency_id) :: id_Carb_in,id_BiCarb_in,id_CarbA_in,id_pH_in,id_pco2_in
type (type_dependency_id) :: id_Carb_in,id_BiCarb_in,id_CarbA_in,id_pH_in,id_pco2_in, id_Hplus_in
type (type_horizontal_dependency_id) :: id_wnd,id_PCO2A

type (type_diagnostic_variable_id) :: id_ph,id_pco2,id_CarbA, id_BiCarb, id_Carb, id_TA_diag
type (type_diagnostic_variable_id) :: id_ph,id_pco2,id_CarbA, id_BiCarb, id_Carb, id_TA_diag, id_Hplus
type (type_diagnostic_variable_id) :: id_Om_cal,id_Om_arg

type (type_horizontal_diagnostic_variable_id) :: id_fair,id_wnd_diag
Expand Down Expand Up @@ -87,7 +87,7 @@ subroutine initialize(self,configunit)
call self%register_diagnostic_variable(self%id_CarbA, 'CarbA', 'mmol/m^3','carbonic acid concentration',missing_value=0._rk)
call self%register_diagnostic_variable(self%id_BiCarb,'BiCarb','mmol/m^3','bicarbonate concentration',missing_value=0._rk)
call self%register_diagnostic_variable(self%id_Carb, 'Carb', 'mmol/m^3','carbonate concentration',standard_variable=standard_variables%mole_concentration_of_carbonate_expressed_as_carbon,missing_value=0._rk)

call self%register_diagnostic_variable(self%id_Hplus, 'Hplus', 'mmol/m^3','hydrogen ions concentration',missing_value=0._rk)
call self%register_diagnostic_variable(self%id_Om_cal,'Om_cal','-','calcite saturation',missing_value=4._rk)
call self%register_diagnostic_variable(self%id_Om_arg,'Om_arg','-','aragonite saturation',missing_value=3._rk)
end if
Expand All @@ -108,6 +108,7 @@ subroutine initialize(self,configunit)
call self%register_dependency(self%id_CarbA_in,'CarbA','mmol/m^3','previous carbonic acid concentration')
call self%register_dependency(self%id_BiCarb_in,'BiCarb','mmol/m^3','previous bicarbonate concentration')
call self%register_dependency(self%id_pH_in,'pH','-','previous pH')
call self%register_dependency(self%id_Hplus_in,'Hplus','mmol/m^3','previous hydrogen ions concentration')
end if

if (self%iswASFLUX>=1) then
Expand Down Expand Up @@ -155,7 +156,7 @@ subroutine do(self,_ARGUMENTS_DO_)

real(rk) :: O3c,ETW,X1X,density,pres
real(rk) :: TA,bioalk,Ctot
real(rk) :: pH,PCO2,H2CO3,HCO3,CO3,k0co2
real(rk) :: pH,PCO2,H2CO3,HCO3,CO3,k0co2,Hplus
real(rk) :: Om_cal,Om_arg
logical :: success

Expand Down Expand Up @@ -189,7 +190,7 @@ subroutine do(self,_ARGUMENTS_DO_)
TA = TA / 1.0e3_rk / density ! from mmol m-3 to mol kg-1
Ctot = O3C / 1.e3_rk / density ! from mmol m-3 to mol kg-1

CALL CO2DYN (ETW,X1X,pres*0.1_rk,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,k0co2,success,self%phscale) ! NB pressure from dbar to bar
CALL CO2DYN (ETW,X1X,pres*0.1_rk,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,Hplus,k0co2,success,self%phscale) ! NB pressure from dbar to bar

if (.not.success) then
! Carbonate system iterative scheme did not converge.
Expand All @@ -200,16 +201,19 @@ subroutine do(self,_ARGUMENTS_DO_)
_GET_(self%id_Carb_in,CO3)
_GET_(self%id_pCO2_in,pCO2)
_GET_(self%id_pH_in,pH)
_GET_(self%id_Hplus_in,Hplus)
CO3 = CO3/1.e3_rk/density ! from mmol/m3 to mol/kg
HCO3 = HCO3/1.e3_rk/density ! from mmol/m3 to mol/kg
H2CO3 = H2CO3/1.e3_rk/density ! from mmol/m3 to mol/kg
pCO2 = pCO2/1.e6_rk ! from uatm to atm
Hplus= Hplus/1.e3_rk/density ! from mmol/m3 to mol/kg
endif
_SET_DIAGNOSTIC_(self%id_ph,pH)
_SET_DIAGNOSTIC_(self%id_pco2,PCO2*1.e6_rk)
_SET_DIAGNOSTIC_(self%id_CarbA, H2CO3*1.e3_rk*density)
_SET_DIAGNOSTIC_(self%id_Bicarb,HCO3*1.e3_rk*density)
_SET_DIAGNOSTIC_(self%id_Carb, CO3*1.e3_rk*density)
_SET_DIAGNOSTIC_(self%id_Hplus, Hplus*1.e3_rk*density)

! Call carbonate saturation state subroutine
CALL CaCO3_Saturation (ETW, X1X, pres*1.e4_rk, CO3, Om_cal, Om_arg) ! NB pressure from dbar to Pa
Expand All @@ -227,7 +231,7 @@ subroutine do_surface(self,_ARGUMENTS_DO_SURFACE_)
real(rk) :: wnd,PCO2A
real(rk) :: sc,fwind,UPTAKE,FAIRCO2

real(rk) :: ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,k0co2
real(rk) :: ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,Hplus,k0co2
logical :: success

if (self%iswASFLUX<=0) return
Expand Down Expand Up @@ -265,7 +269,7 @@ subroutine do_surface(self,_ARGUMENTS_DO_SURFACE_)
! for surface box only calculate air-sea flux
!..Only call after 2 days, because the derivation of instability in the
!..
CALL CO2dyn(T, S, PRSS*0.1_rk,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,k0co2,success,self%phscale)
CALL CO2dyn(T, S, PRSS*0.1_rk,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,Hplus,k0co2,success,self%phscale)
if (.not.success) then
_GET_(self%id_pco2_in,PCO2)
PCO2 = PCO2*1.e-6_rk
Expand Down Expand Up @@ -324,13 +328,13 @@ subroutine do_surface(self,_ARGUMENTS_DO_SURFACE_)
!\\
!\\
! !INTERFACE:
SUBROUTINE CO2dyn ( T, S, PRSS,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,k0co2,success,hscale)
SUBROUTINE CO2dyn ( T, S, PRSS,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,Hplus,k0co2,success,hscale)
!
! !LOCAL VARIABLES:
! ! TODO - SORT THESE!
real(rk),intent(in) :: T, S, PRSS ! NB PRSS is pressure in bar
real(rk),intent(inout) :: ctot,TA
real(rk),intent(out) :: pH,PCO2,H2CO3,HCO3,CO3,k0co2
real(rk),intent(out) :: pH,PCO2,H2CO3,HCO3,CO3,Hplus,k0co2
logical, intent(out) :: success
integer, intent(in) :: hscale

Expand Down Expand Up @@ -363,9 +367,10 @@ SUBROUTINE CO2dyn ( T, S, PRSS,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,k0co2,success,hsca
H2CO3 = 0._rk
HCO3 = 0._rk
CO3 = 0._rk
Hplus = 0._rk

CALL CO2SET(PRSS,Tmax,S,k0co2,k1co2,k2co2,kb,hscale)
CALL CO2CLC(k0co2,k1co2,k2co2,kb,ICALC,BORON,BTOT,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,success)
CALL CO2CLC(k0co2,k1co2,k2co2,kb,ICALC,BORON,BTOT,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,Hplus,success)

END SUBROUTINE CO2DYN
!
Expand Down Expand Up @@ -564,19 +569,19 @@ END SUBROUTINE CO2SET
!\\
!\\
! !INTERFACE:
SUBROUTINE CO2CLC(k0co2,k1co2,k2co2,kb,ICALC,BORON,BTOT,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,success)
SUBROUTINE CO2CLC(k0co2,k1co2,k2co2,kb,ICALC,BORON,BTOT,ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,AHplus,success)
!
! !USES:
!
! !LOCAL VARIABLES:
real(rk),intent(in) :: k0co2,k1co2,k2co2,kb
real(rk),intent(inout) :: ctot,TA,pH,PCO2,H2CO3,HCO3,CO3
real(rk),intent(inout) :: ctot,TA,pH,PCO2,H2CO3,HCO3,CO3,AHPLUS
logical,intent(out) :: success
INTEGER ICALC, KARL, LQ
!put counter in to check duration in convergence loop
INTEGER :: COUNTER,C_CHECK,C_SW
real(rk) :: ALKC, ALKB,BTOT
real(rk) :: AKR,AHPLUS
real(rk) :: AKR
real(rk) :: PROD,tol1,tol2,tol3,tol4,steg,fak
real(rk) :: STEGBY,Y,X,W,X1,Y1,X2,Y2,FACTOR,TERM,Z
LOGICAL :: BORON,DONE
Expand Down