<p><b>mpetersen@lanl.gov</b> 2012-03-23 11:17:18 -0600 (Fri, 23 Mar 2012)</p><p>Move computation of uBolus from ocn_diagnostic_solve to rk4, in preparation of split explicit changes.  Using 120km zlevel rk4, get bit-for-bit match before and after this commit.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_mpas_core.F
===================================================================
--- branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_mpas_core.F        2012-03-23 00:24:20 UTC (rev 1698)
+++ branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_mpas_core.F        2012-03-23 17:17:18 UTC (rev 1699)
@@ -19,6 +19,7 @@
    use ocn_tracer_hadv
    use ocn_tracer_vadv
    use ocn_tracer_hmix
+   use ocn_gm
    use ocn_restoring
 
    use ocn_equation_of_state
@@ -247,6 +248,20 @@
       call ocn_diagnostic_solve(dt, block % state % time_levs(1) % state, mesh)
       call mpas_timer_stop(&quot;diagnostic solve&quot;, initDiagSolveTimer)
 
+! mrp 120323 add:
+         ! Apply the GM closure as a bolus velocity
+         if (config_h_kappa .GE. epsilon(0D0)) then
+             call ocn_gm_compute_uBolus(block % state % time_levs(1) % state, block % mesh)
+         else
+            ! mrp efficiency note: if uBolusGM is guaranteed to be zero, this can be removed.
+            block % state % time_levs(1) % state % uBolusGM % array(:,:) = 0.0
+         end if
+
+            block % state % time_levs(1) % state % uTransport % array(:,:) &amp;
+          = block % state % time_levs(1) % state % u % array(:,:) &amp;
+          + block % state % time_levs(1) % state % uBolusGM % array(:,:)
+! mrp 120323 end
+
       call ocn_wtop(block % state % time_levs(1) % state,block % state % time_levs(1) % state, mesh)
 
       call ocn_compute_mesh_scaling(mesh)

Modified: branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_tendency.F
===================================================================
--- branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_tendency.F        2012-03-23 00:24:20 UTC (rev 1698)
+++ branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_tendency.F        2012-03-23 17:17:18 UTC (rev 1699)
@@ -841,14 +841,16 @@
 
       end do
 
-      ! Applications of any bolus-like parametrizations should be made after the next line.
-      uTransport(:,:) = u(:,:)
+! mrp 120323 remove:
+!      ! Applications of any bolus-like parametrizations should be made after the next line.
+!      uTransport(:,:) = u(:,:)
+!      ! Apply the GM closure as a bolus velocity
+!      if(config_h_kappa .GE. epsilon(0D0)) then
+!         call ocn_gm_compute_uBolus(s,grid)
+!         uTransport(:,:) = uTransport(:,:) + uBolusGM(:,:)
+!      end if
+! mrp 120323 end
 
-      ! Apply the GM closure as a bolus velocity
-      if(config_h_kappa .GE. epsilon(0D0)) then
-         call ocn_gm_compute_uBolus(s,grid)
-         uTransport(:,:) = uTransport(:,:) + uBolusGM(:,:)
-      end if
 
 
    end subroutine ocn_diagnostic_solve!}}}

Modified: branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_time_integration_rk4.F
===================================================================
--- branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_time_integration_rk4.F        2012-03-23 00:24:20 UTC (rev 1698)
+++ branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_time_integration_rk4.F        2012-03-23 17:17:18 UTC (rev 1699)
@@ -25,6 +25,7 @@
 
    use ocn_equation_of_state
    use ocn_vmix
+   use ocn_gm
    use ocn_time_average
 
    implicit none
@@ -241,6 +242,18 @@
 
               call ocn_diagnostic_solve(dt, provis, block % mesh)
 
+! mrp 120323 add:
+              ! Apply the GM closure as a bolus velocity
+              if(config_h_kappa .GE. epsilon(0D0)) then
+                 call ocn_gm_compute_uBolus(provis, block % mesh)
+              else
+                 ! mrp efficiency note: if uBolusGM is guaranteed to be zero, this can be removed.
+                 provis % uBolusGM % array(:,:) = 0.0
+              end if
+
+              provis % uTransport % array(:,:) = provis % u % array(:,:) + provis % uBolusGM % array(:,:)
+! mrp 120323 end
+
               block =&gt; block % next
            end do
         end if
@@ -343,6 +356,20 @@
 
          call ocn_diagnostic_solve(dt, block % state % time_levs(2) % state, block % mesh)
 
+! mrp 120323 add:
+         ! Apply the GM closure as a bolus velocity
+         if (config_h_kappa .GE. epsilon(0D0)) then
+             call ocn_gm_compute_uBolus(block % state % time_levs(2) % state, block % mesh)
+         else
+            ! mrp efficiency note: if uBolusGM is guaranteed to be zero, this can be removed.
+            block % state % time_levs(2) % state % uBolusGM % array(:,:) = 0.0
+         end if
+
+            block % state % time_levs(2) % state % uTransport % array(:,:) &amp;
+          = block % state % time_levs(2) % state % u % array(:,:) &amp;
+          + block % state % time_levs(2) % state % uBolusGM % array(:,:)
+! mrp 120323 end
+
          call mpas_reconstruct(block % mesh, block % state % time_levs(2) % state % u % array,          &amp;
                           block % state % time_levs(2) % state % uReconstructX % array,            &amp;
                           block % state % time_levs(2) % state % uReconstructY % array,            &amp;

</font>
</pre>