--- MITgcm/pkg/dic/dic_surfforcing.F	2011/05/05 22:23:27	1.27
+++ MITgcm/pkg/dic/dic_surfforcing.F	2015/01/11 20:02:41	1.30
@@ -1,4 +1,4 @@
-C $Header: /home/ubuntu/mnt/e9_copy/MITgcm/pkg/dic/dic_surfforcing.F,v 1.27 2011/05/05 22:23:27 stephd Exp $
+C $Header: /home/ubuntu/mnt/e9_copy/MITgcm/pkg/dic/dic_surfforcing.F,v 1.30 2015/01/11 20:02:41 jmc Exp $
 C $Name:  $
 
 #include "DIC_OPTIONS.h"
@@ -9,7 +9,7 @@
 
 C !INTERFACE: ==========================================================
       SUBROUTINE DIC_SURFFORCING( PTR_CO2 , PTR_ALK, PTR_PO4, GDC,
-     I           bi,bj,imin,imax,jmin,jmax,
+     I           bi,bj,iMin,iMax,jMin,jMax,
      I           myIter,myTime,myThid)
 
 C !DESCRIPTION:
@@ -46,6 +46,7 @@
 
 C !LOCAL VARIABLES: ====================================================
        INTEGER i,j, kLev
+       _RL co3dummy
 C Number of iterations for pCO2 solvers...
 C Solubility relation coefficients
       _RL SchmidtNoDIC(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
@@ -79,20 +80,19 @@
 c       ENDDO
 c#endif
 
-
 C =================================================================
 C determine inorganic carbon chem coefficients
-        DO j=jmin,jmax
-         DO i=imin,imax
+        DO j=jMin,jMax
+         DO i=iMin,iMax
 
 #ifdef DIC_BIOTIC
 cQQQQ check ptracer numbers
 #ifdef DIC_BOUNDS
-             surfalk(i,j) = max(0.4 _d 0, 
+             surfalk(i,j) = max(0.4 _d 0,
      &                          min(10. _d 0,PTR_ALK(i,j,klev)))
      &                          * maskC(i,j,kLev,bi,bj)
-             surfphos(i,j)  = max(1.0 _d -11, 
-     &                          min(1._d -1, PTR_PO4(i,j,klev)))
+             surfphos(i,j)  = max(1.0 _d -11,
+     &                          min(1. _d -1,PTR_PO4(i,j,klev)))
      &                          * maskC(i,j,kLev,bi,bj)
 #else
              surfalk(i,j) = PTR_ALK(i,j,klev)
@@ -108,11 +108,11 @@
 C FOR NON-INTERACTIVE Si
              surfsi(i,j)   = SILICA(i,j,bi,bj) * maskC(i,j,kLev,bi,bj)
 #ifdef DIC_BOUNDS
-            surftemp(i,j) = max(-4. _d 0, 
+            surftemp(i,j) = max(-4. _d 0,
      &                          min(50. _d 0, theta(i,j,kLev,bi,bj)))
-            surfsalt(i,j) = max(4. _d 0, 
+            surfsalt(i,j) = max(4. _d 0,
      &                          min(50. _d 0, salt(i,j,kLev,bi,bj)))
-            surfdic(i,j)  = max(0.4 _d 0, 
+            surfdic(i,j)  = max(0.4 _d 0,
      &                          min(10. _d 0, PTR_CO2(i,j,kLev)))
 #else
             surftemp(i,j) = theta(i,j,kLev,bi,bj)
@@ -127,8 +127,8 @@
      I                       bi,bj,iMin,iMax,jMin,jMax,myThid)
 C====================================================================
 
-       DO j=jmin,jmax
-        DO i=imin,imax
+       DO j=jMin,jMax
+        DO i=iMin,iMax
 C Compute AtmosP and Kwexch_Pre which are re-used for flux of O2
 
 #ifdef USE_PLOAD
@@ -151,9 +151,9 @@
 
 c pCO2 solver...
 C$TAF LOOP = parallel
-       DO j=jmin,jmax
+       DO j=jMin,jMax
 C$TAF LOOP = parallel
-        DO i=imin,imax
+        DO i=iMin,iMax
 
           IF ( maskC(i,j,kLev,bi,bj).NE.0. _d 0 ) THEN
             CALL CALC_PCO2_APPROX(
@@ -167,7 +167,7 @@
      I        ak0(i,j,bi,bj), fugf(i,j,bi,bj),
      I        ff(i,j,bi,bj),
      I        bt(i,j,bi,bj),st(i,j,bi,bj),ft(i,j,bi,bj),
-     U        pH(i,j,bi,bj),pCO2(i,j,bi,bj),
+     U        pH(i,j,bi,bj),pCO2(i,j,bi,bj),co3dummy,
      I        i,j,kLev,bi,bj,myIter,myThid )
           ELSE
             pCO2(i,j,bi,bj)=0. _d 0
@@ -175,8 +175,8 @@
         ENDDO
        ENDDO
 
-       DO j=jmin,jmax
-        DO i=imin,imax
+       DO j=jMin,jMax
+        DO i=iMin,iMax
 
           IF ( maskC(i,j,kLev,bi,bj).NE.0. _d 0 ) THEN
 C calculate SCHMIDT NO. for CO2
@@ -186,7 +186,7 @@
      &          + sca3 * theta(i,j,kLev,bi,bj)*theta(i,j,kLev,bi,bj)
      &          + sca4 * theta(i,j,kLev,bi,bj)*theta(i,j,kLev,bi,bj)
      &                *theta(i,j,kLev,bi,bj)
-c make sure Schmidt number isn't negative (will happen if temp>39C)
+c make sure Schmidt number is not negative (will happen if temp>39C)
              SchmidtNoDIC(i,j)=max(1.0 _d -2, SchmidtNoDIC(i,j))
 
 C Determine surface flux (FDIC)
@@ -209,7 +209,7 @@
      &         ff(i,j,bi,bj)*pCO2(i,j,bi,bj)
      &         )
 #else
-C Corrected by Val Bennington Nov 2010 per G.A. McKinley's finding
+C Corrected by Val Bennington Nov 2010 per G.A. McKinley s finding
 C of error in application of water vapor correction
 c Flux = kw*rho*(ff*pCO2atm-k0*FugFac*pCO2ocean)
                FluxCO2(i,j,bi,bj) =
@@ -247,8 +247,8 @@
          ENDDO
 
 C update tendency
-         DO j=jmin,jmax
-          DO i=imin,imax
+         DO j=jMin,jMax
+          DO i=iMin,iMax
            GDC(i,j)= recip_drF(kLev)*recip_hFacC(i,j,kLev,bi,bj)
      &              *(FluxCO2(i,j,bi,bj)
 #ifdef ALLOW_OLD_VIRTUALFLUX