<p><b>mpetersen@lanl.gov</b> 2011-05-31 11:14:20 -0600 (Tue, 31 May 2011)</p><p>BRANCH COMMIT: Use uBcl_{n+1/2} for Coriolis term in Baroclinic iteration. Still use uBcl_{n+1} to proceed to next time step. This compiles and works, but does not exactly match previous commit after 10 timesteps, as expected.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/timesplitting_mrp/src/core_ocean/Registry
===================================================================
--- branches/ocean_projects/timesplitting_mrp/src/core_ocean/Registry        2011-05-31 14:42:34 UTC (rev 862)
+++ branches/ocean_projects/timesplitting_mrp/src/core_ocean/Registry        2011-05-31 17:14:20 UTC (rev 863)
@@ -187,7 +187,7 @@
var persistent real sshSubcycle ( nCells Time ) 2 o sshSubcycle state - -
var persistent real FBtr ( nEdges Time ) 1 o FBtr state - -
var persistent real GBtrForcing ( nEdges Time ) 1 o GBtrForcing state - -
-var persistent real uBcl ( nVertLevels nEdges Time ) 2 o uBcl state - -
+var persistent real uBcl ( nVertLevels nEdges Time ) 3 o uBcl state - -
# Diagnostic fields: only written to output
var persistent real v ( nVertLevels nEdges Time ) 2 - v state - -
Modified: branches/ocean_projects/timesplitting_mrp/src/core_ocean/module_time_integration.F
===================================================================
--- branches/ocean_projects/timesplitting_mrp/src/core_ocean/module_time_integration.F        2011-05-31 14:42:34 UTC (rev 862)
+++ branches/ocean_projects/timesplitting_mrp/src/core_ocean/module_time_integration.F        2011-05-31 17:14:20 UTC (rev 863)
@@ -571,10 +571,15 @@
block % state % time_levs(1) % state % GBtrForcing % array(iEdge) = split*GSum/hSum/dt
do k=1,block % mesh % maxLevelEdgeTop % array(iEdge)
- !{\bf u}'_{k,n+1} = {\bf u}'_{k,n} + \Delta t \left({\bf G}_k - {\overline {\bf G}}\right)
+ ! These two steps are together here:
+ !{\bf u}'_{k,n+1} = {\bf u}'_{k,n} - \Delta t {\overline {\bf G}}
+ !{\bf u}'_{k,n+1/2} = \frac{1}{2}\left({\bf u}^{'}_{k,n} +{\bf u}'_{k,n+1}\right)
+ ! so that uBclNew is at time n+1/2
block % state % time_levs(2) % state % uBcl % array(k,iEdge) &
- = uTemp(k) &
- - dt * block % state % time_levs(1) % state % GBtrForcing % array(iEdge)
+ = 0.5*( &
+ block % state % time_levs(1) % state % uBcl % array(k,iEdge) &
+ + uTemp(k) - dt * block % state % time_levs(1) % state % GBtrForcing % array(iEdge))
+
enddo
enddo ! iEdge
@@ -613,6 +618,9 @@
block => domain % blocklist
do while (associated(block))
block % state % time_levs(2) % state % uBtr % array(:) = 0.0
+ ! mrp 110531 change following lines once I have pointers set up for subcycling
+ block % state % time_levs(2) % state % uBtrSubcycle % array(:) = 0.0
+ block % state % time_levs(1) % state % uBtrSubcycle % array(:) = 0.0
block => block % next
end do ! block
@@ -824,11 +832,11 @@
do while (associated(block))
do iEdge=1,block % mesh % nEdges
+ ! compute u*, the velocity for tendency terms. Put in uNew.
+ ! uBclNew is at time n+1/2 here.
block % state % time_levs(2) % state % u % array(:,iEdge) &
= block % state % time_levs(2) % state % uBtr % array(iEdge) &
- + 0.5*( &
- block % state % time_levs(1) % state % uBcl % array(:,iEdge) &
- + block % state % time_levs(2) % state % uBcl % array(:,iEdge) )
+ + block % state % time_levs(2) % state % uBcl % array(:,iEdge)
end do ! iEdge
@@ -947,6 +955,21 @@
! block % state % time_levs(2) % state % tracers % array(3,1,1:nCells) = domain % dminfo % my_proc_id
! printing end
+ ! Recompute final u to go on to next step.
+ ! u_{n+1} = uBtr_{n+1} + uBcl_{n+1}
+ ! Right now uBclNew is at time n+1/2, so back compute to get uBcl at time n+1
+ ! using uBcl_{n+1/2} = 1/2*(uBcl_n + u_Bcl_{n+1})
+ ! so the following lines are
+ ! u_{n+1} = uBtr_{n+1} + 2*uBcl_{n+1/2} - uBcl_n
+
+ ! mrp 110531 change once I have pointers set up for subcycling
+ do iEdge=1,block % mesh % nEdges
+ block % state % time_levs(2) % state % u % array(:,iEdge) &
+ = block % state % time_levs(2) % state % uBtrSubcycle % array(iEdge) &
+ +2*block % state % time_levs(2) % state % uBcl % array(:,iEdge) &
+ - block % state % time_levs(1) % state % uBcl % array(:,iEdge)
+ enddo
+
u => block % state % time_levs(2) % state % u % array
tracers => block % state % time_levs(2) % state % tracers % array
h => block % state % time_levs(2) % state % h % array
@@ -957,10 +980,6 @@
maxLevelCell => block % mesh % maxLevelCell % array
maxLevelEdgeTop => block % mesh % maxLevelEdgeTop % array
- nCells = block % mesh % nCells
- nEdges = block % mesh % nEdges
- nVertLevels = block % mesh % nVertLevels
-
! dividing by h above for higdon.
! do iCell=1,nCells
! do k=1,maxLevelCell(iCell)
</font>
</pre>