<p><b>mpetersen@lanl.gov</b> 2012-10-19 16:24:29 -0600 (Fri, 19 Oct 2012)</p><p>branch commit, leith_mrp: dd computation of viscosity for output.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/leith_mrp/src/core_ocean/mpas_ocn_vel_hmix.F
===================================================================
--- branches/ocean_projects/leith_mrp/src/core_ocean/mpas_ocn_vel_hmix.F        2012-10-19 22:17:06 UTC (rev 2234)
+++ branches/ocean_projects/leith_mrp/src/core_ocean/mpas_ocn_vel_hmix.F        2012-10-19 22:24:29 UTC (rev 2235)
@@ -126,8 +126,10 @@
       !
       !-----------------------------------------------------------------
 
+      viscosity = 0.0
+
       call mpas_timer_start(&quot;del2&quot;, .false., del2Timer)
-      call ocn_vel_hmix_del2_tend(grid, divergence, vorticity, tend, err1)
+      call ocn_vel_hmix_del2_tend(grid, divergence, vorticity, viscosity, tend, err1)
       call mpas_timer_stop(&quot;del2&quot;, del2Timer)
 
       call mpas_timer_start(&quot;leith&quot;, .false., leithTimer)

Modified: branches/ocean_projects/leith_mrp/src/core_ocean/mpas_ocn_vel_hmix_del2.F
===================================================================
--- branches/ocean_projects/leith_mrp/src/core_ocean/mpas_ocn_vel_hmix_del2.F        2012-10-19 22:17:06 UTC (rev 2234)
+++ branches/ocean_projects/leith_mrp/src/core_ocean/mpas_ocn_vel_hmix_del2.F        2012-10-19 22:24:29 UTC (rev 2235)
@@ -70,7 +70,7 @@
 !
 !-----------------------------------------------------------------------
 
-   subroutine ocn_vel_hmix_del2_tend(grid, divergence, vorticity, tend, err)!{{{
+   subroutine ocn_vel_hmix_del2_tend(grid, divergence, vorticity, viscosity, tend, err)!{{{
 
       !-----------------------------------------------------------------
       !
@@ -96,6 +96,8 @@
       real (kind=RKIND), dimension(:,:), intent(inout) :: &amp;
          tend             !&lt; Input/Output: velocity tendency
 
+      real (kind=RKIND), dimension(:,:), intent(inout) :: &amp;
+         viscosity       !&lt; Input: viscosity
 
       !-----------------------------------------------------------------
       !
@@ -113,10 +115,10 @@
 
       integer :: iEdge, nEdgesSolve, cell1, cell2, vertex1, vertex2
       integer :: k
-      integer, dimension(:), pointer :: maxLevelEdgeTop
+      integer, dimension(:), pointer :: maxLevelEdgeTop, nEdgesOnCell
       integer, dimension(:,:), pointer :: cellsOnEdge, verticesOnEdge, edgeMask
 
-      real (kind=RKIND) :: u_diffusion, invLength1, invLength2
+      real (kind=RKIND) :: u_diffusion, invLength1, invLength2, visc2
       real (kind=RKIND), dimension(:), pointer :: meshScalingDel2, &amp;
               dcEdge, dvEdge
 
@@ -132,6 +134,7 @@
 
       nEdgesSolve = grid % nEdgesSolve
       maxLevelEdgeTop =&gt; grid % maxLevelEdgeTop % array
+      nEdgesOnCell =&gt; grid % nEdgesOnCell % array
       cellsOnEdge =&gt; grid % cellsOnEdge % array
       verticesOnEdge =&gt; grid % verticesOnEdge % array
       meshScalingDel2 =&gt; grid % meshScalingDel2 % array
@@ -158,10 +161,13 @@
                           -viscVortCoef &amp;
                           *( vorticity(k,vertex2) - vorticity(k,vertex1) ) * invLength2
 
-            u_diffusion = meshScalingDel2(iEdge) * eddyVisc2 * u_diffusion
+            visc2 = meshScalingDel2(iEdge) * eddyVisc2
 
-            tend(k,iEdge) = tend(k,iEdge) + edgeMask(k, iEdge) * u_diffusion
+            tend(k,iEdge) = tend(k,iEdge) + edgeMask(k, iEdge) * visc2 * u_diffusion
 
+            viscosity(k,cell1) = viscosity(k,cell1) + visc2/nEdgesOnCell(cell1)
+            viscosity(k,cell2) = viscosity(k,cell2) + visc2/nEdgesOnCell(cell2)
+
          end do
       end do
 

Modified: branches/ocean_projects/leith_mrp/src/core_ocean/mpas_ocn_vel_hmix_leith.F
===================================================================
--- branches/ocean_projects/leith_mrp/src/core_ocean/mpas_ocn_vel_hmix_leith.F        2012-10-19 22:17:06 UTC (rev 2234)
+++ branches/ocean_projects/leith_mrp/src/core_ocean/mpas_ocn_vel_hmix_leith.F        2012-10-19 22:24:29 UTC (rev 2235)
@@ -115,7 +115,7 @@
 
       integer :: iEdge, nEdgesSolve, cell1, cell2, vertex1, vertex2
       integer :: k
-      integer, dimension(:), pointer :: maxLevelEdgeTop
+      integer, dimension(:), pointer :: maxLevelEdgeTop, nEdgesOnCell
       integer, dimension(:,:), pointer :: cellsOnEdge, verticesOnEdge, edgeMask
 
       real (kind=RKIND) :: u_diffusion, invLength1, invLength2, visc2
@@ -134,6 +134,7 @@
 
       nEdgesSolve = grid % nEdgesSolve
       maxLevelEdgeTop =&gt; grid % maxLevelEdgeTop % array
+      nEdgesOnCell =&gt; grid % nEdgesOnCell % array
       cellsOnEdge =&gt; grid % cellsOnEdge % array
       verticesOnEdge =&gt; grid % verticesOnEdge % array
       meshScaling =&gt; grid % meshScaling % array
@@ -141,8 +142,6 @@
       dcEdge =&gt; grid % dcEdge % array
       dvEdge =&gt; grid % dvEdge % array
 
-      viscosity = 1.0
-
       do iEdge=1,nEdgesSolve
          cell1 = cellsOnEdge(1,iEdge)
          cell2 = cellsOnEdge(2,iEdge)
@@ -165,12 +164,11 @@
             visc2 = ( config_leith_parameter * config_leith_dx * meshScaling(iEdge) / 3.14)**3 &amp;
                      * abs( vorticity(k,vertex2) - vorticity(k,vertex1) ) * invLength1 * sqrt(3.0)
             visc2 = min(visc2, config_leith_visc2_max)
-            u_diffusion = visc2 * u_diffusion
 
-            viscosity(k,cell1) = visc2
-            viscosity(k,cell2) = visc2
+            tend(k,iEdge) = tend(k,iEdge) + edgeMask(k, iEdge) * visc2 * u_diffusion
 
-            tend(k,iEdge) = tend(k,iEdge) + edgeMask(k, iEdge) * u_diffusion
+            viscosity(k,cell1) = viscosity(k,cell1) + visc2/nEdgesOnCell(cell1)
+            viscosity(k,cell2) = viscosity(k,cell2) + visc2/nEdgesOnCell(cell2)
 
          end do
       end do

</font>
</pre>