Index: /branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/CONFIG/SHARED/field_def_bgc.xml
===================================================================
--- /branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/CONFIG/SHARED/field_def_bgc.xml	(revision 8130)
+++ /branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/CONFIG/SHARED/field_def_bgc.xml	(revision 8131)
@@ -432,4 +432,5 @@
        <field id= "DMS_ARAN"   long_name="DMS Surface Concentration, AranamiT"       unit="nmol/L"      />
        <field id= "DMS_HALL"   long_name="DMS Surface Concentration, Halloran"       unit="nmol/L"      />
+       <field id= "DMS_ANDM"   long_name="DMS Surface Concentration, Anderson modif" unit="nmol/L"      />
        <field id= "ATM_XCO2"   long_name="Atmospheric xCO2"                          unit="ppm"         />
        <field id= "OCN_FCO2"   long_name="Surface ocean fCO2"                        unit="uatm"        />
@@ -769,4 +770,5 @@
       <field field_ref= "DMS_ARAN"   name="DMS_ARAN"   />
       <field field_ref= "DMS_HALL"   name="DMS_HALL"   />
+      <field field_ref= "DMS_ANDM"   name="DMS_ANDM"   />
     </field_group>
 
Index: /branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90
===================================================================
--- /branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90	(revision 8130)
+++ /branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90	(revision 8131)
@@ -20,4 +20,5 @@
    !!  -   !  2015-10  (J. Palm)              Update for diag outputs through iom_use  
    !!  -   !  2016-11  (A. Yool)              Updated diags for CMIP6
+   !!  -   !  2017-05  (A. Yool)              Added extra DMS calculation
    !!----------------------------------------------------------------------
    !!
@@ -335,5 +336,5 @@
       REAL(wp) ::    dms_surf
       !! AXY (13/03/15): add in other DMS calculations
-      REAL(wp) ::    dms_andr, dms_simo, dms_aran, dms_hall
+      REAL(wp) ::    dms_andr, dms_simo, dms_aran, dms_hall, dms_andm, dms_nlim, dms_wtkn
 
 # endif
@@ -371,5 +372,5 @@
       REAL(wp), POINTER, DIMENSION(:,:  ) :: mdetc2d, gmidc2d, gmedc2d, f_pco2a2d, f_pco2w2d, f_co2flux2d
       REAL(wp), POINTER, DIMENSION(:,:  ) :: f_TDIC2d, f_TALK2d, f_kw6602d, f_pp02d, f_o2flux2d, f_o2sat2d
-      REAL(wp), POINTER, DIMENSION(:,:  ) :: dms_andr2d, dms_simo2d, dms_aran2d, dms_hall2d, dms_surf2d
+      REAL(wp), POINTER, DIMENSION(:,:  ) :: dms_andr2d, dms_simo2d, dms_aran2d, dms_hall2d, dms_andm2d, dms_surf2d
       REAL(wp), POINTER, DIMENSION(:,:  ) :: iben_n2d, iben_fe2d, iben_c2d, iben_si2d, iben_ca2d, oben_n2d
       REAL(wp), POINTER, DIMENSION(:,:  ) :: oben_fe2d, oben_c2d, oben_si2d, oben_ca2d, sfr_ocal2d
@@ -1001,4 +1002,8 @@
                CALL wrk_alloc( jpi, jpj,   dms_hall2d    )
                dms_hall2d(:,:)      = 0.0 !!
+            ENDIF
+            IF( med_diag%DMS_ANDM%dgsave ) THEN
+               CALL wrk_alloc( jpi, jpj,   dms_andm2d    )
+               dms_andm2d(:,:)      = 0.0 !!
             ENDIF
          ENDIF   
@@ -1841,16 +1846,31 @@
                   !!                 in trcdms_medusa.F90)
                   !!
+                  !! AXY (25/05/17): amended to additionally pass DIN limitation as well as [DIN];
+                  !!                 accounts for differences in nutrient half-saturations; changes
+                  !!                 also made in trc_dms_medusa; this permits an additional DMS
+                  !!                 calculation while retaining the existing Anderson one
+                  !!
                   IF (jdms .eq. 1) THEN
+                     !!
+                     !! calculate weighted half-saturation for DIN uptake
+                     dms_wtkn = ((zphn * xnln) + (zphd * xnld)) / (zphn + zphd)
                      !!
                      !! feed in correct inputs
                      if (jdms_input .eq. 0) then
                         !! use instantaneous inputs
-                        CALL trc_dms_medusa( zchn, zchd, hmld(ji,jj), qsr(ji,jj), zdin, &  ! inputs
-                        dms_andr, dms_simo, dms_aran, dms_hall )                           ! outputs
+                        dms_nlim = zdin / (zdin + dms_wtkn)
+                        !!
+                        CALL trc_dms_medusa( zchn, zchd,                           &  ! inputs
+                        hmld(ji,jj), qsr(ji,jj),                                   &  ! inputs
+                        zdin, dms_nlim,                                            &  ! inputs
+                        dms_andr, dms_simo, dms_aran, dms_hall, dms_andm )            ! outputs
                      else
                         !! use diel-average inputs
+                        dms_nlim = zn_dms_din(ji,jj) / (zn_dms_din(ji,jj) + dms_wtkn)
+                        !!
                         CALL trc_dms_medusa( zn_dms_chn(ji,jj), zn_dms_chd(ji,jj), &  ! inputs
-                        zn_dms_mld(ji,jj), zn_dms_qsr(ji,jj), zn_dms_din(ji,jj),   &  ! inputs
-                        dms_andr, dms_simo, dms_aran, dms_hall )                      ! outputs
+                        zn_dms_mld(ji,jj), zn_dms_qsr(ji,jj),                      &  ! inputs
+                        zn_dms_din(ji,jj), dms_nlim,                               &  ! inputs
+                        dms_andr, dms_simo, dms_aran, dms_hall, dms_andm )            ! outputs
                      endif
                      !!
@@ -1864,4 +1884,6 @@
                      elseif (jdms_model .eq. 4) then
                         dms_surf = dms_hall
+                     elseif (jdms_model .eq. 5) then
+                        dms_surf = dms_andm
                      endif
                      !!
@@ -1882,4 +1904,7 @@
                        IF( med_diag%DMS_HALL%dgsave ) THEN
                          dms_hall2d(ji,jj) = dms_hall
+                       ENDIF
+                       IF( med_diag%DMS_ANDM%dgsave ) THEN
+                         dms_andm2d(ji,jj) = dms_andm
                        ENDIF
 #   if defined key_debug_medusa
@@ -2527,7 +2552,7 @@
                   !! temperature-dependent
                   fdd  = xmd  * fun_Q10 * zdet
-#if defined key_roam
+# if defined key_roam
                   fddc = xmdc * fun_Q10 * zdtc
-#endif
+# endif
                else
                   !! temperature-independent
@@ -4914,4 +4939,8 @@
                          CALL wrk_dealloc( jpi, jpj,   dms_hall2d   )
                        ENDIF
+                       IF( med_diag%DMS_ANDM%dgsave ) THEN
+                         CALL iom_put( "DMS_ANDM"  , dms_andm2d )
+                         CALL wrk_dealloc( jpi, jpj,   dms_andm2d   )
+                       ENDIF
                      ENDIF
                      !! AXY (24/11/16): extra MOCSY diagnostics
Index: /branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcdms_medusa.F90
===================================================================
--- /branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcdms_medusa.F90	(revision 8130)
+++ /branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcdms_medusa.F90	(revision 8131)
@@ -6,4 +6,5 @@
    !! History :
    !!  -   !  2014-08  (J. Palmieri - A. Yool)    added for UKESM1 project
+   !!  -   !  2017-05  (A. Yool)                  add extra Anderson scheme
    !!----------------------------------------------------------------------
 #if defined key_medusa && defined key_roam
@@ -38,6 +39,6 @@
 !=======================================================================
 !
-   SUBROUTINE trc_dms_medusa( chn, chd, mld, xqsr, xdin,   &  !! inputs
-     &  dms_andr, dms_simo, dms_aran, dms_hall )              !! outputs
+   SUBROUTINE trc_dms_medusa( chn, chd, mld, xqsr, xdin, xlim,  &  !! inputs
+     &  dms_andr, dms_simo, dms_aran, dms_hall dms_andm)           !! outputs
 !      
 !=======================================================================
@@ -69,4 +70,10 @@
       !!                    and is included because it is formally
       !!                    published (and different from the above)
+      !!
+      !! AXY (25/05/17): add extra "corrected" Anderson scheme
+      !!
+      !!                 5. As Anderson et al. (2001) but modified to
+      !!                    more accurately reflect nutrient limitation
+      !!                    status of phytoplankton community
       !!
       !! AXY (08/07/15): amend to remove Julien's original calculation
@@ -84,8 +91,10 @@
       REAL(wp), INTENT( in )    :: xqsr                 !! surface irradiance        (W/m2)
       REAL(wp), INTENT( in )    :: xdin                 !! surface DIN               (mmol N/m3)
+      REAL(wp), INTENT( in )    :: xlim                 !! surface DIN limitation    (mmol N/m3)
       REAL(wp), INTENT( inout ) :: dms_andr             !! DMS surface concentration (nmol/L) 
       REAL(wp), INTENT( inout ) :: dms_simo             !! DMS surface concentration (nmol/L) 
       REAL(wp), INTENT( inout ) :: dms_aran             !! DMS surface concentration (nmol/L) 
       REAL(wp), INTENT( inout ) :: dms_hall             !! DMS surface concentration (nmol/L) 
+      REAL(wp), INTENT( inout ) :: dms_andm             !! DMS surface concentration (nmol/L) 
 !
       REAL(wp) :: CHL, cmr, sw_dms
@@ -147,4 +156,15 @@
            dms_hall = fq3
         endif
+!
+! AXY (25/05/17): modified Anderson et al. (2001)
+        Jterm = xqsr + 1.0e-6
+        !! this version fixes the hard-coded assumption above
+        Qterm = xlim
+        fq1 = log10(CHL * Jterm * Qterm)
+        if (fq1 > 1.72) then
+           dms_andm = (8.24 * (fq1 - 1.72)) + 2.29
+        else
+           dms_andm = 2.29
+        endif
 
   END SUBROUTINE trc_dms_medusa
Index: /branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcnam_medusa.F90
===================================================================
--- /branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcnam_medusa.F90	(revision 8130)
+++ /branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcnam_medusa.F90	(revision 8131)
@@ -493,5 +493,5 @@
 !!	 jdms_input  :  use instant (0) or diel-avg (1) inputs
 !!       jdms_model  :  choice of DMS model passed to atmosphere
-!!                      1 = ANDR, 2 = SIMO, 3 = ARAN, 4 = HALL
+!!                      1 = ANDR, 2 = SIMO, 3 = ARAN, 4 = HALL, 5 = ANDM
 !!
       IF(lwp) THEN
@@ -1026,4 +1026,7 @@
             WRITE(numout,*)     &
             &   ' choice of DMS model passed to atmosphere,                  jdms_model  = Halloran et al. (2010)'
+	 elseif (jdms_model .eq. 5) then
+            WRITE(numout,*)     &
+            &   ' choice of DMS model passed to atmosphere,                  jdms_model  = Anderson et al. (2001; modified)'
          endif
 !!
@@ -2095,4 +2098,9 @@
           med_diag%DMS_HALL%dgsave = .FALSE.
       ENDIF
+      IF  (iom_use("DMS_ANDM")) THEN 
+          med_diag%DMS_ANDM%dgsave = .TRUE.
+      ELSE 
+          med_diag%DMS_ANDM%dgsave = .FALSE.
+      ENDIF
       IF  (iom_use("ATM_XCO2")) THEN 
           med_diag%ATM_XCO2%dgsave = .TRUE.
