<p><b>dwj07@fsu.edu</b> 2012-10-19 15:02:11 -0600 (Fri, 19 Oct 2012)</p><p><br>
        -- BRANCH COMMIT --<br>
<br>
        Cleaner version of ocn_wtop at Mark's request.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/option3_b4b_test/src/core_ocean/mpas_ocn_tendency.F
===================================================================
--- branches/ocean_projects/option3_b4b_test/src/core_ocean/mpas_ocn_tendency.F        2012-10-19 20:39:30 UTC (rev 2230)
+++ branches/ocean_projects/option3_b4b_test/src/core_ocean/mpas_ocn_tendency.F        2012-10-19 21:02:11 UTC (rev 2231)
@@ -869,8 +869,8 @@
       real (kind=RKIND), dimension(:), pointer :: &amp;
         dvEdge, areaCell, zstarWeight
       real (kind=RKIND), dimension(:,:), pointer :: uTransport,h,wTop, h_edge
-      real (kind=RKIND), dimension(:,:), allocatable:: div_hu
-      real (kind=RKIND), dimension(:), allocatable:: div_hu_btr, h_tend_col
+      real (kind=RKIND), dimension(:), allocatable:: div_hu, h_tend_col
+      real (kind=RKIND) :: div_hu_btr
 
       integer, dimension(:,:), pointer :: cellsOnEdge, cellsOnVertex, &amp;
         verticesOnEdge, edgesOnCell, edgesOnEdge, edgesOnVertex, &amp;
@@ -898,69 +898,57 @@
       nEdges      = grid % nEdges
       nVertLevels = grid % nVertLevels
 
-      allocate(div_hu(nVertLevels,nCells+1), div_hu_btr(nCells+1), &amp;
-          h_tend_col(nVertLevels))
 
+      if (config_vert_grid_type.eq.'isopycnal') then
+        ! set vertical velocity to zero in isopycnal case
+        wTop=0.0_RKIND
+        return
+      end if
+
+      allocate(div_hu(nVertLevels), h_tend_col(nVertLevels))
+
       !
       ! Compute div(h^{edge} u) for each cell
       ! See Ringler et al. (2010) jcp paper, eqn 19, 21, and fig. 3.
       !
-      div_hu(:,:) = 0.0_RKIND
 
       do iCell=1,nCells
+        div_hu(:) = 0.0_RKIND
+        div_hu_btr = 0.0_RKIND
+        hSum = 0.0_RKIND
         invAreaCell = 1.0_RKIND / areaCell(iCell)
+
         do i = 1, nEdgesOnCell(iCell)
           iEdge = edgesOnCell(i, iCell)
 
           do k = 1, maxLevelEdgeBot(iEdge)
             flux = uTransport(k, iEdge) * dvEdge(iEdge) * h_edge(k, iEdge)
-            div_hu(k, iCell) = div_hu(k, iCell) - edgeSignOnCell(i, iCell) * flux * invAreaCell
+            flux = edgeSignOnCell(i, iCell) * flux * invAreaCell
+            div_hu(k) = div_hu(k) - flux
+            div_hu_btr = div_hu_btr - flux
           end do
         end do
-      end do
 
-      do iCell=1,nCells
-         div_hu_btr(iCell) = 0.0_RKIND
-         do k=1,maxLevelCell(iCell)
-            div_hu_btr(iCell) = div_hu_btr(iCell) + div_hu(k,iCell)
-         end do
-      end do
+        do k = 1, maxLevelCell(iCell)
+           h_tend_col(k) = - zstarWeight(k) * h(k, iCell) * div_hu_btr
+           hSum = hSum + zstarWeight(k) * h(k, iCell)
+        end do
 
-      !
-      ! vertical velocity through layer interface
-      !
-      !dwj: 10/25/2011 - Need to explore isopycnal vs zlevel flags
-      if (config_vert_grid_type.eq.'isopycnal') then
-        ! set vertical velocity to zero in isopycnal case
-        wTop=0.0_RKIND
+        if(hSum &gt; 0.0) then
+           h_tend_col = h_tend_col / hSum
+        end if
 
-      else  ! zlevel or zstar type vertical grid
-
-        do iCell=1,nCells
-
-           hSum = 0.0_RKIND
-           do k=1,maxLevelCell(iCell)
-              h_tend_col(k) = - zstarWeight(k)*h(k,iCell)*div_hu_btr(iCell)
-              hSum = hSum + zstarWeight(k)*h(k,iCell)
-           end do
-           if(hSum &gt; 0.0) then
-             h_tend_col = h_tend_col / hSum
-           else
-           end if
-
-           ! Vertical velocity through layer interface at top and 
-           ! bottom is zero.
-           wTop(1,iCell) = 0.0_RKIND
-           wTop(maxLevelCell(iCell)+1,iCell) = 0.0_RKIND
-           do k=maxLevelCell(iCell),2,-1
-              wTop(k,iCell) = wTop(k+1,iCell) - div_hu(k,iCell) - h_tend_col(k)
-           end do
+        ! Vertical velocity through layer interface at top and 
+        ! bottom is zero.
+        wTop(1,iCell) = 0.0_RKIND
+        wTop(maxLevelCell(iCell)+1,iCell) = 0.0_RKIND
+        do k=maxLevelCell(iCell),2,-1
+           wTop(k,iCell) = wTop(k+1,iCell) - div_hu(k) - h_tend_col(k)
         end do
+      end do
 
-      endif
+      deallocate(div_hu, h_tend_col)
 
-      deallocate(div_hu, div_hu_btr, h_tend_col)
-
    end subroutine ocn_wtop!}}}
 
 !***********************************************************************

</font>
</pre>