<p><b>qchen3@fsu.edu</b> 2011-12-20 10:12:07 -0700 (Tue, 20 Dec 2011)</p><p>BRANCH COMMIT<br>
<br>
The following variables are added to Registry:<br>
config_h_kappa<br>
config_h_kappa_q (for specifying a constant GM coefficient)<br>
h_kappa<br>
h_kappa_q (as a state array variable, to accommodate the idea of spatially <br>
           varying GM coefficient)<br>
uBolus<br>
uTransport<br>
<br>
Two subroutines are added in mpas_ocn_tendency.F:<br>
ocn_get_h_kappa and ocn_get_h_kappa_q.<br>
Right now, these two subroutines assign config_h_kappa and config_h_kappa_q <br>
to h_kappa and h_kappa_q, respectively. <br>
<br>
The computation of uBolus and uTransport is also implemented<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/gmvar/Makefile
===================================================================
--- branches/ocean_projects/gmvar/Makefile        2011-12-20 15:52:45 UTC (rev 1264)
+++ branches/ocean_projects/gmvar/Makefile        2011-12-20 17:12:07 UTC (rev 1265)
@@ -111,9 +111,9 @@
         &quot;CC = gcc&quot; \
         &quot;SFC = ifort&quot; \
         &quot;SCC = gcc&quot; \
-        &quot;FFLAGS = -real-size 64 -O3 -convert big_endian -FR&quot; \
-        &quot;CFLAGS = -O3 -m64&quot; \
-        &quot;LDFLAGS = -O3&quot; \
+        &quot;FFLAGS = -real-size 64 -g -convert big_endian -FR&quot; \
+        &quot;CFLAGS = -g -m64&quot; \
+        &quot;LDFLAGS = -g&quot; \
         &quot;CORE = $(CORE)&quot; \
         &quot;CPPFLAGS = -DRKIND=8 $(MODEL_FORMULATION) -DUNDERSCORE -m64 $(FILE_OFFSET) $(ZOLTAN_DEFINE)&quot; )
 

Modified: branches/ocean_projects/gmvar/src/core_ocean/Registry
===================================================================
--- branches/ocean_projects/gmvar/src/core_ocean/Registry        2011-12-20 15:52:45 UTC (rev 1264)
+++ branches/ocean_projects/gmvar/src/core_ocean/Registry        2011-12-20 17:12:07 UTC (rev 1265)
@@ -47,7 +47,9 @@
 namelist logical   hmix     config_include_KE_vertex false
 namelist real      hmix     config_h_tracer_eddy_diff2  0.0
 namelist real      hmix     config_h_tracer_eddy_diff4  0.0
-namelist real      hmix     config_apvm_upwinding       0.5
+namelist real      hmix     config_apvm_upwinding       0.0
+namelist real      hmix     config_h_kappa              0.0
+namelist real      hmix     config_h_kappa_q            0.0
 namelist logical   hmix     config_mom_decay            false 
 namelist real      hmix     config_mom_decay_time       3600.0
 namelist character vmix     config_vert_visc_type       const
@@ -188,7 +190,7 @@
 
 # Prognostic variables: read from input, saved in restart, and written to output
 var persistent real    u ( nVertLevels nEdges Time ) 2 ir u state - -
-var persistent real    h ( nVertLevels nCells Time ) 2 ir h state - -
+var persistent real    h ( nVertLevels nCells Time ) 2 iro h state - -
 var persistent real    rho ( nVertLevels nCells Time ) 2 iro rho state - -
 var persistent real    temperature ( nVertLevels nCells Time ) 2 iro temperature state tracers dynamics
 var persistent real    salinity ( nVertLevels nCells Time ) 2 iro salinity state tracers dynamics
@@ -217,6 +219,16 @@
 
 # Diagnostic fields: only written to output
 var persistent real    v ( nVertLevels nEdges Time ) 2 - v state - -
+var persistent real    uBolus ( nVertLevels nEdges Time ) 2 o uBolus state - -
+var persistent real    uBolusX ( nVertLevels nEdges Time ) 2 o uBolusX state - -
+var persistent real    uBolusY ( nVertLevels nEdges Time ) 2 o uBolusY state - -
+var persistent real    uBolusZ ( nVertLevels nEdges Time ) 2 o uBolusZ state - -
+var persistent real    uBolusZonal ( nVertLevels nEdges Time ) 2 o uBolusZonal state - -
+var persistent real    uBolusMeridional ( nVertLevels nEdges Time ) 2 o uBolusMeridional state - -
+var persistent real    uTransport ( nVertLevels nEdges Time ) 2 o uTransport state - -
+var persistent real    h_kappa ( nVertLevels nEdges Time ) 2 o h_kappa state - -
+var persistent real    h_kappa_q ( nVertLevels nEdges Time ) 2 o h_kappa_q state - -
+var persistent real    hu ( nVertLevels nEdges Time ) 2 0 hu state - -
 var persistent real    divergence ( nVertLevels nCells Time ) 2 o divergence state - -
 var persistent real    vorticity ( nVertLevels nVertices Time ) 2 o vorticity state - -
 var persistent real    pv_edge ( nVertLevels nEdges Time ) 2 - pv_edge state - -

Modified: branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_tendency.F
===================================================================
--- branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_tendency.F        2011-12-20 15:52:45 UTC (rev 1264)
+++ branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_tendency.F        2011-12-20 17:12:07 UTC (rev 1265)
@@ -620,7 +620,8 @@
         weightsOnEdge, kiteAreasOnVertex, h_edge, h, u, v, w, pressure,&amp;
         circulation, vorticity, ke, ke_edge, MontPot, wTop, &amp;
         pv_edge, pv_vertex, pv_cell, gradPVn, gradPVt, divergence, &amp;
-        rho, temperature, salinity, kev, kevc
+        rho, temperature, salinity, kev, kevc, uBolus, uTransport, &amp;
+        h_kappa, h_kappa_q
       real (kind=RKIND), dimension(:,:,:), pointer :: tracers
       real (kind=RKIND), dimension(:), allocatable:: pTop
       real (kind=RKIND), dimension(:,:), allocatable:: div_u
@@ -642,6 +643,8 @@
       h           =&gt; s % h % array
       u           =&gt; s % u % array
       v           =&gt; s % v % array
+      uBolus      =&gt; s % uBolus % array
+      uTransport  =&gt; s % uTransport % array
       wTop        =&gt; s % wTop % array
       h_edge      =&gt; s % h_edge % array
       circulation =&gt; s % circulation % array
@@ -660,6 +663,8 @@
       tracers     =&gt; s % tracers % array
       MontPot     =&gt; s % MontPot % array
       pressure    =&gt; s % pressure % array
+      h_kappa     =&gt; s % h_kappa % array
+      h_kappa_q   =&gt; s % h_kappa_q % array
 
       weightsOnEdge     =&gt; grid % weightsOnEdge % array
       kiteAreasOnVertex =&gt; grid % kiteAreasOnVertex % array
@@ -695,6 +700,7 @@
       boundaryEdge =&gt; grid % boundaryEdge % array
       boundaryCell =&gt; grid % boundaryCell % array
 
+
       !
       ! Compute height on cell edges at velocity locations
       !   Namelist options control the order of accuracy of the reconstructed h_edge value
@@ -906,6 +912,24 @@
       enddo
 
       !
+      ! Compute Bolus velocity as bolus == - /kappa (/grad h) / h
+      !
+      uBolus(:,:) = 0.0
+      call ocn_get_h_kappa(s, grid)
+      do iEdge = 1,nEdges
+         cell1 = cellsOnEdge(1,iEdge)
+         cell2 = cellsOnEdge(2,iEdge)
+         do k=1,maxLevelEdgeTop(iEdge)
+            uBolus(k,iEdge) = -h_kappa(k,iEdge) * (h(k,cell2) - h(k,cell1))/dcEdge(iEdge)/h_edge(k,iEdge)
+         end do
+      end do
+
+      !
+      ! Compute the transport velocity (this velocity transports thickness and PV)
+      !
+      uTransport(:,:) = u(:,:) + uBolus(:,:)
+
+      !
       ! Compute kinetic energy in each cell by blending ke and kevc
       !
       if(config_include_KE_vertex) then
@@ -1034,7 +1058,7 @@
       do iEdge = 1,nEdges
          do k = 1,maxLevelEdgeBot(iEdge)
            pv_edge(k,iEdge) = pv_edge(k,iEdge) &amp;
-             - 0.5 * dt* (  u(k,iEdge) * gradPVn(k,iEdge) &amp;
+             - config_apvm_upwinding * dt* (  u(k,iEdge) * gradPVn(k,iEdge) &amp;
                           + v(k,iEdge) * gradPVt(k,iEdge) )
          enddo
       enddo
@@ -1346,6 +1370,38 @@
 
    end subroutine ocn_fuperp!}}}
 
+
+   subroutine ocn_get_h_kappa(s, grid)
+
+      type (state_type), intent(inout) :: s
+      type (mesh_type), intent(in) :: grid
+
+      real(kind=RKIND), dimension(:,:), pointer    :: h_kappa
+
+
+      h_kappa  =&gt; s % h_kappa % array
+
+      h_kappa(:,:) = config_h_kappa
+
+
+   end subroutine ocn_get_h_kappa
+
+
+   subroutine ocn_get_h_kappa_q(s, grid)
+
+      type (state_type), intent(inout) :: s
+      type (mesh_type), intent(in) :: grid
+
+      real(kind=RKIND), dimension(:,:), pointer    :: h_kappa_q
+
+
+      h_kappa_q  =&gt; s % h_kappa_q % array
+
+      h_kappa_q(:,:) = config_h_kappa_q
+
+
+   end subroutine ocn_get_h_kappa_q
+
 !***********************************************************************
 
 end module ocn_tendency

</font>
</pre>