<p><b>mpetersen@lanl.gov</b> 2010-11-19 10:30:41 -0700 (Fri, 19 Nov 2010)</p><p>BRANCH COMMIT For testing, set all land values to -1e34.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/test_topography_mrp/src/core_ocean/module_time_integration.F
===================================================================
--- branches/ocean_projects/test_topography_mrp/src/core_ocean/module_time_integration.F        2010-11-18 23:08:48 UTC (rev 628)
+++ branches/ocean_projects/test_topography_mrp/src/core_ocean/module_time_integration.F        2010-11-19 17:30:41 UTC (rev 629)
@@ -67,7 +67,7 @@
       type (domain_type), intent(inout) :: domain
       real (kind=RKIND), intent(in) :: dt
 
-      integer :: iCell, k, i
+      integer :: iCell, iEdge, k, i
       type (block_type), pointer :: block
       type (state_type) :: provis
 
@@ -211,6 +211,19 @@
               end do
            end do
 
+ ! mrp 101118 debug test
+           do iEdge=1,block % mesh % nEdges
+              ! u=-1e34 gives NaN
+!               block % state % time_levs(2) % state % u % array(block % mesh % maxLevelEdgeTop % array(iEdge)+1:block % mesh % nVertLevels,iEdge) = -1e34
+           end do
+           do iCell=1,block % mesh % nCells
+              ! h=-1e34 gives slightly different answer
+!              block % state % time_levs(2) % state % h % array(block % mesh % maxLevelCell % array(iCell)+1:block % mesh % nVertLevels,iCell) = -1e34
+              ! tracer=-1e34 gives bit-for-bit match
+              block % state % time_levs(2) % state % tracers % array(:,block % mesh % maxLevelCell % array(iCell)+1:block % mesh % nVertLevels,iCell) =  -1e34
+           end do
+ ! mrp 101118 debug test end
+
            block =&gt; block % next
         end do
 
@@ -358,20 +371,31 @@
 
       if (config_vert_grid_type.eq.'isopycnal') then
         kmax = nVertLevels
+      do iEdge=1,nEdges
+         cell1 = cellsOnEdge(1,iEdge)
+         cell2 = cellsOnEdge(2,iEdge)
+         do k=1,kmax
+            flux = u(k,iEdge) * dvEdge(iEdge) * h_edge(k,iEdge)
+            tend_h(k,cell1) = tend_h(k,cell1) - flux
+            tend_h(k,cell2) = tend_h(k,cell2) + flux
+         end do
+      end do
+
       elseif (config_vert_grid_type.eq.'zlevel') then
+! mrp changed 101118.  Clean up if this stays
         kmax = 1
-      endif
- 
       do iEdge=1,nEdges
          cell1 = cellsOnEdge(1,iEdge)
          cell2 = cellsOnEdge(2,iEdge)
-         do k=1,kmax
+         do k=1,min(1,maxLevelEdgeTop(iEdge))
             flux = u(k,iEdge) * dvEdge(iEdge) * h_edge(k,iEdge)
             tend_h(k,cell1) = tend_h(k,cell1) - flux
             tend_h(k,cell2) = tend_h(k,cell2) + flux
          end do
       end do
 
+      endif
+ 
       do iCell=1,nCells
          do k=1,kmax
             tend_h(k,iCell) = tend_h(k,iCell) / areaCell(iCell)
@@ -1320,7 +1344,7 @@
       do iEdge=1,nEdges
          vertex1 = verticesOnEdge(1,iEdge)
          vertex2 = verticesOnEdge(2,iEdge)
-         do k=1,maxLevelEdgeBot(iEdge)
+         do k=1,maxLevelEdgeTop(iEdge) ! mrp 101118 changed from maxLevelEdgeBot for testing
             circulation(k,vertex1) = circulation(k,vertex1) - dcEdge(iEdge) * u(k,iEdge)
             circulation(k,vertex2) = circulation(k,vertex2) + dcEdge(iEdge) * u(k,iEdge)
          end do
@@ -1338,7 +1362,7 @@
       do iEdge=1,nEdges
          cell1 = cellsOnEdge(1,iEdge)
          cell2 = cellsOnEdge(2,iEdge)
-         do k=1,maxLevelEdgeBot(iEdge)
+         do k=1,maxLevelEdgeTop(iEdge) ! mrp 101118 changed from maxLevelEdgeBot for testing
              divergence(k,cell1) = divergence(k,cell1) + u(k,iEdge)*dvEdge(iEdge)
              divergence(k,cell2) = divergence(k,cell2) - u(k,iEdge)*dvEdge(iEdge)
          enddo
@@ -1373,7 +1397,7 @@
       do iEdge = 1,nEdges
          do i=1,nEdgesOnEdge(iEdge)
             eoe = edgesOnEdge(i,iEdge)
-            do k = 1,maxLevelEdgeBot(iEdge)
+            do k = 1,maxLevelEdgeTop(iEdge) ! mrp 101118 changed from maxLevelEdgeBot for testing
                v(k,iEdge) = v(k,iEdge) + weightsOnEdge(i,iEdge) * u(k, eoe)
             end do
          end do
@@ -1433,7 +1457,7 @@
       do iVertex = 1,nVertices
          do i=1,vertexDegree
             iEdge = edgesOnVertex(i,iVertex)
-            do k=1,maxLevelEdgeBot(iEdge)
+            do k=1,maxLevelEdgeTop(iEdge) ! mrp 101118 changed from maxLevelEdgeBot for testing
                pv_edge(k,iEdge) =  pv_edge(k,iEdge) + 0.5 * pv_vertex(k,iVertex)
             enddo
         end do
@@ -1457,7 +1481,7 @@
       !   ( this computes gradPVt at all edges bounding real cells and distance-1 ghost cells )
       !
       do iEdge = 1,nEdges
-         do k = 1,maxLevelEdgeBot(iEdge)
+         do k = 1,maxLevelEdgeTop(iEdge) ! mrp 101118 changed from maxLevelEdgeBot for testing
            gradPVt(k,iEdge) = (  pv_vertex(k,verticesOnEdge(2,iEdge)) &amp;
                                - pv_vertex(k,verticesOnEdge(1,iEdge))) &amp;
                                  /dvEdge(iEdge)
@@ -1468,7 +1492,7 @@
       ! Modify PV edge with upstream bias.
       !
       do iEdge = 1,nEdges
-         do k = 1,maxLevelEdgeBot(iEdge)
+         do k = 1,maxLevelEdgeTop(iEdge) ! mrp 101118 changed from maxLevelEdgeBot for testing
            pv_edge(k,iEdge) = pv_edge(k,iEdge) &amp;
              - 0.5 * dt* (  u(k,iEdge) * gradPVn(k,iEdge) &amp;
                           + v(k,iEdge) * gradPVt(k,iEdge) )
@@ -1556,7 +1580,7 @@
         do iEdge=1,nEdges
            cell1 = cellsOnEdge(1,iEdge)
            cell2 = cellsOnEdge(2,iEdge)
-           do k=1,maxLevelEdgeBot(iEdge)
+           do k=1,maxLevelEdgeTop(iEdge) ! mrp 101118 changed from maxLevelEdgeBot for testing
               flux = u(k,iEdge) * dvEdge(iEdge) 
               div_u(k,cell1) = div_u(k,cell1) + flux
               div_u(k,cell2) = div_u(k,cell2) - flux

</font>
</pre>