<p><b>duda</b> 2011-04-06 18:18:11 -0600 (Wed, 06 Apr 2011)</p><p>BRANCH COMMIT<br>
<br>
Add correction to the computation of ke in non-hydrostatic <br>
model to avoid Hollingsworth instability.<br>
<br>
Note: Once we have subTriangleAreasOnVertex in our regular<br>
mesh files, the corrections can be made more exact.<br>
<br>
<br>
M    module_time_integration.F<br>
</p><hr noshade><pre><font color="gray">Modified: branches/atmos_physics/src/core_nhyd_atmos/module_time_integration.F
===================================================================
--- branches/atmos_physics/src/core_nhyd_atmos/module_time_integration.F        2011-04-06 00:07:21 UTC (rev 785)
+++ branches/atmos_physics/src/core_nhyd_atmos/module_time_integration.F        2011-04-07 00:18:11 UTC (rev 786)
@@ -3028,7 +3028,12 @@
       integer, dimension(:,:), pointer :: cellsOnEdge, cellsOnVertex, verticesOnEdge, edgesOnCell, edgesOnEdge, edgesOnVertex
       integer, dimension(:), pointer :: nEdgesOnCell, nEdgesOnEdge
 
+      !WCS-instability
+      logical, parameter :: hollingsworth=.true.
+      real (kind=RKIND), allocatable, dimension(:,:) :: ke_vertex
+      real (kind=RKIND)  :: ke_fact
 
+
       h           =&gt; s % rho % array
       u           =&gt; s % u % array
       v           =&gt; diag % v % array
@@ -3147,7 +3152,51 @@
          end do
       end do
 
+      !WCS-instability
+      ! Compute ke at cell vertices - AG's new KE construction, part 1
+      ! *** approximation here because we don't have inner triangle areas
       !
+      if (hollingsworth) then
+      allocate (ke_vertex(nVertLevels,nVertices))
+      do iVertex=1,nVertices
+         do k=1,nVertLevels
+!            ke_vertex(k,iVertex) = (  subTriangleAreasOnVertex(1,iVertex)*u(k,EdgesOnVertex(1,iVertex))**2.0  &amp;
+!                                    + subTriangleAreasOnVertex(2,iVertex)*u(k,EdgesOnVertex(2,iVertex))**2.0  &amp;
+!                                    + subTriangleAreasOnVertex(3,iVertex)*u(k,EdgesOnVertex(3,iVertex))**2.0  &amp;
+!                                                 ) / AreaTriangle(iVertex)
+
+            ke_vertex(k,iVertex) = (  dcEdge(EdgesOnVertex(1,iVertex))*dvEdge(EdgesOnVertex(1,iVertex))*u(k,EdgesOnVertex(1,iVertex))**2.0  &amp;
+                                     +dcEdge(EdgesOnVertex(2,iVertex))*dvEdge(EdgesOnVertex(2,iVertex))*u(k,EdgesOnVertex(2,iVertex))**2.0  &amp;
+                                     +dcEdge(EdgesOnVertex(3,iVertex))*dvEdge(EdgesOnVertex(3,iVertex))*u(k,EdgesOnVertex(3,iVertex))**2.0  &amp;
+                                                 ) * 0.25 / AreaTriangle(iVertex)
+
+         end do
+      end do
+
+      ! adjust ke at cell vertices - AG's new KE construction, part 2
+      !
+
+      ke_fact = 1.0 - .375
+
+      do iCell=1,nCells
+         do k=1,nVertLevels
+            ke(k,iCell) = ke_fact*ke(k,iCell)
+         end do
+      end do
+
+      do iVertex = 1, nVertices
+       do i=1,grid % vertexDegree
+          iCell = cellsOnVertex(i,iVertex)
+          do k = 1,nVertLevels
+             ke(k,iCell) = ke(k,iCell) + (1.-ke_fact)*kiteAreasOnVertex(i, iVertex) * ke_vertex(k, iVertex) / areaCell(iCell)
+          end do
+       end do
+      end do
+      deallocate (ke_vertex)
+      end if
+      !END of WCS-instability
+
+      !
       ! Compute v (tangential) velocities
       !
       v(:,:) = 0.0

</font>
</pre>