Index: /branches/UKMO/dev_r5518_GO6_package_asm_3d_bgc/NEMOGCM/NEMO/OPA_SRC/ASM/asmbgc.F90
===================================================================
--- /branches/UKMO/dev_r5518_GO6_package_asm_3d_bgc/NEMOGCM/NEMO/OPA_SRC/ASM/asmbgc.F90	(revision 9381)
+++ /branches/UKMO/dev_r5518_GO6_package_asm_3d_bgc/NEMOGCM/NEMO/OPA_SRC/ASM/asmbgc.F90	(revision 9382)
@@ -896,6 +896,6 @@
             ! reset to zero at the start of trc_stp, called after this routine
 #if defined key_medusa && defined key_foam_medusa
-            WHERE( phyto2d_balinc(:,:,:,:) > 0.0_wp .OR. &
-               &   trn(:,:,:,jp_msa0:jp_msa1) + phyto2d_balinc(:,:,:,:) * zincwgt > 0.0_wp )
+            WHERE( phyto2d_balinc(:,:,:,jp_msa0:jp_msa1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_msa0:jp_msa1) + phyto2d_balinc(:,:,:,jp_msa0:jp_msa1) * zincwgt > 0.0_wp )
                trn(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1) + &
                   &                         phyto2d_balinc(:,:,:,jp_msa0:jp_msa1) * zincwgt
@@ -904,6 +904,6 @@
             END WHERE
 #elif defined key_hadocc
-            WHERE( phyto2d_balinc(:,:,:,:) > 0.0_wp .OR. &
-               &   trn(:,:,:,jp_had0:jp_had1) + phyto2d_balinc(:,:,:,:) * zincwgt > 0.0_wp )
+            WHERE( phyto2d_balinc(:,:,:,jp_had0:jp_had1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_had0:jp_had1) + phyto2d_balinc(:,:,:,jp_had0:jp_had1) * zincwgt > 0.0_wp )
                trn(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1) + &
                   &                         phyto2d_balinc(:,:,:,jp_had0:jp_had1) * zincwgt
@@ -932,6 +932,6 @@
                &           ' Background state is taken from model rather than background file' )
 #if defined key_medusa && defined key_foam_medusa
-            WHERE( phyto2d_balinc(:,:,:,:) > 0.0_wp .OR. &
-               &   trn(:,:,:,jp_msa0:jp_msa1) + phyto2d_balinc(:,:,:,:) * zincwgt > 0.0_wp )
+            WHERE( phyto2d_balinc(:,:,:,jp_msa0:jp_msa1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_msa0:jp_msa1) + phyto2d_balinc(:,:,:,jp_msa0:jp_msa1) > 0.0_wp )
                trn(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1) + &
                   &                         phyto2d_balinc(:,:,:,jp_msa0:jp_msa1)
@@ -939,6 +939,6 @@
             END WHERE
 #elif defined key_hadocc
-            WHERE( phyto2d_balinc(:,:,:,:) > 0.0_wp .OR. &
-               &   trn(:,:,:,jp_had0:jp_had1) + phyto2d_balinc(:,:,:,:) * zincwgt > 0.0_wp )
+            WHERE( phyto2d_balinc(:,:,:,jp_had0:jp_had1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_had0:jp_had1) + phyto2d_balinc(:,:,:,jp_had0:jp_had1) > 0.0_wp )
                trn(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1) + &
                   &                         phyto2d_balinc(:,:,:,jp_had0:jp_had1)
@@ -1083,6 +1083,6 @@
             ! reset to zero at the start of trc_stp, called after this routine
 #if defined key_medusa && defined key_foam_medusa
-            WHERE( phyto3d_balinc(:,:,:,:) > 0.0_wp .OR. &
-               &   trn(:,:,:,jp_msa0:jp_msa1) + phyto3d_balinc(:,:,:,:) * zincwgt > 0.0_wp )
+            WHERE( phyto3d_balinc(:,:,:,jp_msa0:jp_msa1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_msa0:jp_msa1) + phyto3d_balinc(:,:,:,jp_msa0:jp_msa1) * zincwgt > 0.0_wp )
                trn(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1) + &
                   &                         phyto3d_balinc(:,:,:,jp_msa0:jp_msa1) * zincwgt
@@ -1091,6 +1091,6 @@
             END WHERE
 #elif defined key_hadocc
-            WHERE( phyto3d_balinc(:,:,:,:) > 0.0_wp .OR. &
-               &   trn(:,:,:,jp_had0:jp_had1) + phyto3d_balinc(:,:,:,:) * zincwgt > 0.0_wp )
+            WHERE( phyto3d_balinc(:,:,:,jp_had0:jp_had1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_had0:jp_had1) + phyto3d_balinc(:,:,:,jp_had0:jp_had1) * zincwgt > 0.0_wp )
                trn(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1) + &
                   &                         phyto3d_balinc(:,:,:,jp_had0:jp_had1) * zincwgt
@@ -1119,6 +1119,6 @@
                &           ' Background state is taken from model rather than background file' )
 #if defined key_medusa && defined key_foam_medusa
-            WHERE( phyto3d_balinc(:,:,:,:) > 0.0_wp .OR. &
-               &   trn(:,:,:,jp_msa0:jp_msa1) + phyto3d_balinc(:,:,:,:) * zincwgt > 0.0_wp )
+            WHERE( phyto3d_balinc(:,:,:,jp_msa0:jp_msa1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_msa0:jp_msa1) + phyto3d_balinc(:,:,:,jp_msa0:jp_msa1) > 0.0_wp )
                trn(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1) + &
                   &                         phyto3d_balinc(:,:,:,jp_msa0:jp_msa1)
@@ -1126,6 +1126,6 @@
             END WHERE
 #elif defined key_hadocc
-            WHERE( phyto3d_balinc(:,:,:,:) > 0.0_wp .OR. &
-               &   trn(:,:,:,jp_had0:jp_had1) + phyto3d_balinc(:,:,:,:) * zincwgt > 0.0_wp )
+            WHERE( phyto3d_balinc(:,:,:,jp_had0:jp_had1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_had0:jp_had1) + phyto3d_balinc(:,:,:,jp_had0:jp_had1) > 0.0_wp )
                trn(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1) + &
                   &                         phyto3d_balinc(:,:,:,jp_had0:jp_had1)
@@ -1353,17 +1353,19 @@
             ! reset to zero at the start of trc_stp, called after this routine
 #if defined key_medusa && defined key_foam_medusa
-            DO jk = 1, jpkm1
-               trn(:,:,jk,jp_msa0:jp_msa1) = trn(:,:,jk,jp_msa0:jp_msa1) + &
-                  &                          pco2_balinc(:,:,jk,jp_msa0:jp_msa1) * zincwgt
-               trb(:,:,jk,jp_msa0:jp_msa1) = trb(:,:,jk,jp_msa0:jp_msa1) + &
-                  &                          pco2_balinc(:,:,jk,jp_msa0:jp_msa1) * zincwgt
-            END DO
+            WHERE( pco2_balinc(:,:,:,jp_msa0:jp_msa1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_msa0:jp_msa1) + pco2_balinc(:,:,:,jp_msa0:jp_msa1) * zincwgt > 0.0_wp )
+               trn(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1) + &
+                  &                         pco2_balinc(:,:,:,jp_msa0:jp_msa1) * zincwgt
+               trb(:,:,:,jp_msa0:jp_msa1) = trb(:,:,:,jp_msa0:jp_msa1) + &
+                  &                         pco2_balinc(:,:,:,jp_msa0:jp_msa1) * zincwgt
+            END WHERE
 #elif defined key_hadocc
-            DO jk = 1, jpkm1
-               trn(:,:,jk,jp_had0:jp_had1) = trn(:,:,jk,jp_had0:jp_had1) + &
-                  &                          pco2_balinc(:,:,jk,jp_had0:jp_had1) * zincwgt
-               trb(:,:,jk,jp_had0:jp_had1) = trb(:,:,jk,jp_had0:jp_had1) + &
-                  &                          pco2_balinc(:,:,jk,jp_had0:jp_had1) * zincwgt
-            END DO
+            WHERE( pco2_balinc(:,:,:,jp_had0:jp_had1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_had0:jp_had1) + pco2_balinc(:,:,:,jp_had0:jp_had1) * zincwgt > 0.0_wp )
+               trn(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1) + &
+                  &                         pco2_balinc(:,:,:,jp_had0:jp_had1) * zincwgt
+               trb(:,:,:,jp_had0:jp_had1) = trb(:,:,:,jp_had0:jp_had1) + &
+                  &                         pco2_balinc(:,:,:,jp_had0:jp_had1) * zincwgt
+            END WHERE
 #endif
 
@@ -1383,20 +1385,22 @@
             neuler = 0                    ! Force Euler forward step
 
-#if defined key_medusa && defined key_foam_medusa
             ! Initialize the now fields with the background + increment
             ! Background currently is what the model is initialised with
-            CALL ctl_warn( ' Doing direct initialisation of MEDUSA with pCO2 assimilation', &
+            CALL ctl_warn( ' Doing direct initialisation with pCO2 assimilation', &
                &           ' Background state is taken from model rather than background file' )
-            trn(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1) + &
-               &                         pco2_balinc(:,:,:,jp_msa0:jp_msa1)
-            trb(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1)
+#if defined key_medusa && defined key_foam_medusa
+            WHERE( pco2_balinc(:,:,:,jp_msa0:jp_msa1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_msa0:jp_msa1) + pco2_balinc(:,:,:,jp_msa0:jp_msa1) > 0.0_wp )
+               trn(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1) + &
+                  &                         pco2_balinc(:,:,:,jp_msa0:jp_msa1)
+               trb(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1)
+            END WHERE
 #elif defined key_hadocc
-            ! Initialize the now fields with the background + increment
-            ! Background currently is what the model is initialised with
-            CALL ctl_warn( ' Doing direct initialisation of HadOCC with pCO2 assimilation', &
-               &           ' Background state is taken from model rather than background file' )
-            trn(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1) + &
-               &                         pco2_balinc(:,:,:,jp_had0:jp_had1)
-            trb(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1)
+            WHERE( pco2_balinc(:,:,:,jp_had0:jp_had1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_had0:jp_had1) + pco2_balinc(:,:,:,jp_had0:jp_had1) > 0.0_wp )
+               trn(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1) + &
+                  &                         pco2_balinc(:,:,:,jp_had0:jp_had1)
+               trb(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1)
+            END WHERE
 #endif
  
@@ -1589,10 +1593,11 @@
             ! Add directly to trn and trb, rather than to tra, because tra gets
             ! reset to zero at the start of trc_stp, called after this routine
-            DO jk = 1, jpkm1
-               trn(:,:,jk,jp_msa0:jp_msa1) = trn(:,:,jk,jp_msa0:jp_msa1) + &
-                  &                          ph_balinc(:,:,jk,jp_msa0:jp_msa1) * zincwgt
-               trb(:,:,jk,jp_msa0:jp_msa1) = trb(:,:,jk,jp_msa0:jp_msa1) + &
-                  &                          ph_balinc(:,:,jk,jp_msa0:jp_msa1) * zincwgt
-            END DO
+            WHERE( ph_balinc(:,:,:,jp_msa0:jp_msa1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_msa0:jp_msa1) + ph_balinc(:,:,:,jp_msa0:jp_msa1) * zincwgt > 0.0_wp )
+               trn(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1) + &
+                  &                         ph_balinc(:,:,:,jp_msa0:jp_msa1) * zincwgt
+               trb(:,:,:,jp_msa0:jp_msa1) = trb(:,:,:,jp_msa0:jp_msa1) + &
+                  &                         ph_balinc(:,:,:,jp_msa0:jp_msa1) * zincwgt
+            END WHERE
 
             ! Do not deallocate arrays - needed by asm_bgc_bal_wri
@@ -1615,7 +1620,10 @@
             CALL ctl_warn( ' Doing direct initialisation of MEDUSA with pH assimilation', &
                &           ' Background state is taken from model rather than background file' )
-            trn(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1) + &
-               &                         ph_balinc(:,:,:,jp_msa0:jp_msa1)
-            trb(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1)
+            WHERE( ph_balinc(:,:,:,jp_msa0:jp_msa1) > 0.0_wp .OR. &
+               &   trn(:,:,:,jp_msa0:jp_msa1) + ph_balinc(:,:,:,jp_msa0:jp_msa1) > 0.0_wp )
+               trn(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1) + &
+                  &                         ph_balinc(:,:,:,jp_msa0:jp_msa1)
+               trb(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1)
+            END WHERE
  
             ! Do not deallocate arrays - needed by asm_bgc_bal_wri
@@ -1654,4 +1662,118 @@
       !!----------------------------------------------------------------------
 
+      IF ( kt <= nit000 ) THEN
+
+         !----------------------------------------------------------------------
+         ! Remove any other balancing increments
+         !----------------------------------------------------------------------
+
+         IF ( ln_pno3inc ) THEN
+#if defined key_hadocc
+            it = jp_had_nut
+#elif defined key_medusa && defined key_foam_medusa
+            it = jpdin
+#else
+            CALL ctl_stop ( ' bgc3d_asm_inc: no compatible BGC model defined' )
+#endif
+            IF ( ln_phytobal ) THEN
+               pno3_bkginc(:,:,:) = pno3_bkginc(:,:,:) - phyto2d_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_plchltotinc .OR. ln_pchltotinc ) THEN
+               pno3_bkginc(:,:,:) = pno3_bkginc(:,:,:) - phyto3d_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_sfco2inc .OR. ln_spco2inc ) THEN
+               pno3_bkginc(:,:,:) = pno3_bkginc(:,:,:) - pco2_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_pphinc ) THEN
+               pno3_bkginc(:,:,:) = pno3_bkginc(:,:,:) - ph_balinc(:,:,:,it)
+            ENDIF
+         ENDIF
+
+         IF ( ln_psi4inc ) THEN
+#if defined key_medusa && defined key_foam_medusa
+            it = jpsil
+#else
+            CALL ctl_stop ( ' bgc3d_asm_inc: no compatible BGC model defined' )
+#endif
+            IF ( ln_phytobal ) THEN
+               psi4_bkginc(:,:,:) = psi4_bkginc(:,:,:) - phyto2d_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_plchltotinc .OR. ln_pchltotinc ) THEN
+               psi4_bkginc(:,:,:) = psi4_bkginc(:,:,:) - phyto3d_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_sfco2inc .OR. ln_spco2inc ) THEN
+               psi4_bkginc(:,:,:) = psi4_bkginc(:,:,:) - pco2_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_pphinc ) THEN
+               psi4_bkginc(:,:,:) = psi4_bkginc(:,:,:) - ph_balinc(:,:,:,it)
+            ENDIF
+         ENDIF
+
+         IF ( ln_pdicinc ) THEN
+#if defined key_hadocc
+            it = jp_had_dic
+#elif defined key_medusa && defined key_foam_medusa
+            it = jpdic
+#else
+            CALL ctl_stop ( ' bgc3d_asm_inc: no compatible BGC model defined' )
+#endif
+            IF ( ln_phytobal ) THEN
+               pdic_bkginc(:,:,:) = pdic_bkginc(:,:,:) - phyto2d_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_plchltotinc .OR. ln_pchltotinc ) THEN
+               pdic_bkginc(:,:,:) = pdic_bkginc(:,:,:) - phyto3d_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_sfco2inc .OR. ln_spco2inc ) THEN
+               pdic_bkginc(:,:,:) = pdic_bkginc(:,:,:) - pco2_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_pphinc ) THEN
+               pdic_bkginc(:,:,:) = pdic_bkginc(:,:,:) - ph_balinc(:,:,:,it)
+            ENDIF
+         ENDIF
+
+         IF ( ln_palkinc ) THEN
+#if defined key_hadocc
+            it = jp_had_alk
+#elif defined key_medusa && defined key_foam_medusa
+            it = jpalk
+#else
+            CALL ctl_stop ( ' bgc3d_asm_inc: no compatible BGC model defined' )
+#endif
+            IF ( ln_phytobal ) THEN
+               palk_bkginc(:,:,:) = palk_bkginc(:,:,:) - phyto2d_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_plchltotinc .OR. ln_pchltotinc ) THEN
+               palk_bkginc(:,:,:) = palk_bkginc(:,:,:) - phyto3d_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_sfco2inc .OR. ln_spco2inc ) THEN
+               palk_bkginc(:,:,:) = palk_bkginc(:,:,:) - pco2_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_pphinc ) THEN
+               palk_bkginc(:,:,:) = palk_bkginc(:,:,:) - ph_balinc(:,:,:,it)
+            ENDIF
+         ENDIF
+
+         IF ( ln_po2inc ) THEN
+#if defined key_medusa && defined key_foam_medusa
+            it = jpoxy
+#else
+            CALL ctl_stop ( ' bgc3d_asm_inc: no compatible BGC model defined' )
+#endif
+            IF ( ln_phytobal ) THEN
+               po2_bkginc(:,:,:) = po2_bkginc(:,:,:) - phyto2d_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_plchltotinc .OR. ln_pchltotinc ) THEN
+               po2_bkginc(:,:,:) = po2_bkginc(:,:,:) - phyto3d_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_sfco2inc .OR. ln_spco2inc ) THEN
+               po2_bkginc(:,:,:) = po2_bkginc(:,:,:) - pco2_balinc(:,:,:,it)
+            ENDIF
+            IF ( ln_pphinc ) THEN
+               po2_bkginc(:,:,:) = po2_bkginc(:,:,:) - ph_balinc(:,:,:,it)
+            ENDIF
+         ENDIF
+
+      ENDIF
+
       IF ( ll_asmiau ) THEN
 
