From 69f30fef51c9f43145a77b8d4c95c0b3d0a87999 Mon Sep 17 00:00:00 2001
From: Tomas Torsvik <tomas.torsvik.work@gmail.com>
Date: Tue, 24 Jan 2023 13:13:58 +0100
Subject: [PATCH] Read parameters for alkalinization scenario from namelist
 file.

---
 cime_config/buildnml    | 19 +++++++++++++++++
 hamocc/mo_read_oafx.F90 | 45 ++++++++++++++++++++++-------------------
 2 files changed, 43 insertions(+), 21 deletions(-)

diff --git a/cime_config/buildnml b/cime_config/buildnml
index 57cc45eb..5c4c4f09 100755
--- a/cime_config/buildnml
+++ b/cime_config/buildnml
@@ -1495,6 +1495,25 @@ cat >>! $RUNDIR/ocn_in$inststr << EOF
   SEDPORFILE   = $SEDPORFILE
 /
 
+! NAMELIST FOR ALKALINIZATION SCENARIO
+!
+! CONTENTS:
+!
+! ADDALK_TOT    :  Pmol alkalinity/yr added in the scenarios.
+! CDRMIP_LATMAX : Max latitude where alkalinity is added according to the
+!                   CDRMIP protocol
+! CDRMIP_LATMIN : Min latitude where alkalinity is added according to the
+!                   CDRMIP protocol
+! RAMP_START    : Start year for ramp up in 'ramp' scenario
+! RAMP_END      : End year for 'ramp' scenario
+&BGCNML
+  ADDALK_TOT    = 0.135
+  CDRMIP_LATMAX = 70.0
+  CDRMIP_LATMIN = -60.0
+  RAMP_START    = 2025
+  RAMP_END      = 2035
+/
+
 ! IO-NAMELIST FOR DIAGNOSTIC iHAMOCC OUTPUT
 !
 ! Namelist acronyms:
diff --git a/hamocc/mo_read_oafx.F90 b/hamocc/mo_read_oafx.F90
index d6f2cb1b..be978f99 100644
--- a/hamocc/mo_read_oafx.F90
+++ b/hamocc/mo_read_oafx.F90
@@ -86,18 +86,21 @@ module mo_read_oafx
   !                ALK/yr-1 from 2025 to 2035 over the ice-free surface ocean  
   !                (assumed to be between 60S and 70N)
   !
-  real, parameter :: addalk_0p14   = 0.14  ! Pmol alkalinity/yr added in the
-  real, parameter :: addalk_0p56   = 0.56  ! 'const_0p14' and 'const_0p56' 
-                                           ! scenarios
-  real, parameter :: cdrmip_latmax =  70.0 ! Min and max latitude where
-  real, parameter :: cdrmip_latmin = -60.0 ! alkalinity is added according
+  real, protected :: addalk_tot    = 0.56  ! Pmol alkalinity/yr added in the
+                                           ! scenarios. Read from namelist file
+                                           ! to overwrite default value.
+  real, protected :: cdrmip_latmax =  70.0 ! Min and max latitude where
+  real, protected :: cdrmip_latmin = -60.0 ! alkalinity is added according
                                            ! to the CDRMIP protocol
-  real, parameter :: addalk_ramp   = 0.135 ! Max Pmol alkalinity/yr added
-  integer, parameter :: ramp_start = 2025  ! in 2035 in the 'ramp' scenario,
-  integer, parameter :: ramp_end   = 2035  ! starting at 0 Pmol/yr in 2025.
+  integer, protected :: ramp_start = 2025  ! In 'ramp' scenario, start at
+  integer, protected :: ramp_end   = 2035  ! 0 Pmol/yr in ramp_start, and max
+                                           ! addalk_const Pmol/yr in ramp_end.
 
   logical,   save :: lini = .false.
 
+  namelist /bgcoafx/ addalk_tot,cdrmip_latmax,cdrmip_latmin,                  &
+       &   ramp_start,ramp_end
+
 !******************************************************************************
 contains
 
@@ -126,17 +129,18 @@ subroutine ini_read_oafx(kpie,kpje,pdlxp,pdlyp,pglat,omask)
 !
 !******************************************************************************
   use mod_xc,         only: xcsum,xchalt,mnproc,nbdy,ips
-  use mo_control_bgc, only: io_stdo_bgc,do_oalk
+  use mo_control_bgc, only: io_stdo_bgc,do_oalk,bgc_namelist,get_bgc_namelist
 
-  implicit none 
+  implicit none
 
   integer, intent(in) :: kpie,kpje
   real,    intent(in) :: pdlxp(kpie,kpje), pdlyp(kpie,kpje)
-  real,    intent(in) :: pglat(1-nbdy:kpie+nbdy,1-nbdy:kpje+nbdy)  
+  real,    intent(in) :: pglat(1-nbdy:kpie+nbdy,1-nbdy:kpje+nbdy)
   real,    intent(in) :: omask(kpie,kpje)
 
   integer :: i,j,errstat
-  real    :: avflx,ztotarea,addalk_tot
+  integer :: iounit
+  real    :: avflx,ztotarea
   real    :: ztmp1(1-nbdy:kpie+nbdy,1-nbdy:kpje+nbdy)
 
   ! Return if alkalinization is turned off
@@ -166,6 +170,13 @@ subroutine ini_read_oafx(kpie,kpje,pdlxp,pdlyp,pglat,omask)
         write(io_stdo_bgc,*)' '
       endif
 
+      ! Read parameters for alkalinization fluxes from namelist file
+      if(.not. allocated(bgc_namelist)) call get_bgc_namelist
+      open (newunit=iounit, file=bgc_namelist, status='old'                   &
+           &   ,action='read')
+      read (unit=iounit, nml=BGCOAFX)
+      close (unit=iounit)
+
       ! Allocate field to hold constant alkalinization fluxes
       if(mnproc.eq.1) then
         write(io_stdo_bgc,*)'Memory allocation for variable oalkflx ...'
@@ -189,15 +200,7 @@ subroutine ini_read_oafx(kpie,kpje,pdlxp,pdlyp,pglat,omask)
       enddo
 
       call xcsum(ztotarea,ztmp1,ips)
-      
-      if( trim(oalkscen)=='const_0p14') then
-        addalk_tot = addalk_0p14
-      else if( trim(oalkscen)=='const_0p56') then
-        addalk_tot = addalk_0p56
-      else
-        addalk_tot = addalk_ramp
-      endif
-    
+
       ! Calculate alkalinity flux (kmol m^2 yr-1) to be applied
       avflx = addalk_tot/ztotarea*Pmol2kmol
       if(mnproc.eq.1) then