<p><b>mpetersen@lanl.gov</b> 2011-05-24 10:56:23 -0600 (Tue, 24 May 2011)</p><p>BRANCH COMMIT: Contains partially completed split Higdon, and baroclinic Stage 1 with splitting on G rather than provisional u.<br>
</p><hr noshade><pre><font color="gray">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-20 16:44:36 UTC (rev 849)
+++ branches/ocean_projects/timesplitting_mrp/src/core_ocean/module_time_integration.F        2011-05-24 16:56:23 UTC (rev 850)
@@ -609,7 +609,8 @@
! Stage 2: Barotropic velocity (2D) prediction, explicitly subcycled
!
- if (trim(config_time_integration) == 'higdon_unsplit') then
+ if (trim(config_time_integration) == 'higdon_unsplit') then
+
block => domain % blocklist
do while (associated(block))
block % state % time_levs(2) % state % uBtr % array(:) = 0.0
@@ -618,24 +619,21 @@
elseif (trim(config_time_integration) == 'higdon_split') then
+ ! note: eventually, put this in subroutine mpas_core_init
+ ! using an input variable to choose the filter type.
+ allocate(aBtrSumCoef(0:2*config_n_btr_subcycles), &
+ bBtrSumCoef(0:2*config_n_btr_subcycles) )
+ aBtrSumCoef = 1.0/(2.0*config_n_btr_subcycles+1)
+ bBtrSumCoef(2*config_n_btr_subcycles) = aBtrSumCoef(2*config_n_btr_subcycles)
+ do j=2*config_n_btr_subcycles-1,0,-1
+ ! mrp not sure - do I need to divide by 2*config_n_btr_subcycles here?
+ ! or do I divide by sum(bBtrSumCoef) at the end?
+ bBtrSumCoef(j) = bBtrSumCoef(j+1) + aBtrSumCoef(j)
+ end do
-
block => domain % blocklist
do while (associated(block))
- ! note: eventually, put this in subroutine mpas_core_init
- ! using an input variable to choose the filter type.
- allocate(aBtrSumCoef(0:2*config_n_btr_subcycles), &
- bBtrSumCoef(0:2*config_n_btr_subcycles) )
- aBtrSumCoef = 1.0/(2.0*config_n_btr_subcycles+1)
-
- bBtrSumCoef(2*config_n_btr_subcycles) = aBtrSumCoef(2*config_n_btr_subcycles)
- do j=2*config_n_btr_subcycles-1,0,-1
- ! mrp not sure - do I need to divide by 2*config_n_btr_subcycles here?
- ! or do I divide by sum(bBtrSumCoef) at the end?
- bBtrSumCoef(j) = bBtrSumCoef(j+1) + aBtrSumCoef(j)
- end do
-
do iCell=1,block % mesh % nCells
! sshSubcycleOld = sshOld
block % state % time_levs(1) % state % sshSubcycle % array(iCell) &
@@ -697,11 +695,9 @@
tend_ssh(cell1) = tend_ssh(cell1) - flux
tend_ssh(cell2) = tend_ssh(cell2) + flux
end do
- do iCell=1,nCells
- tend_ssh(iCell) = tend_ssh(iCell) / block % mesh % areaCell % array (iCell)
- end do
- do iCell=1,block % mesh % nCells ! couple tracers to h
+ do iCell=1,block % mesh % nCells
+
! \zeta_{n+(j+1)/J} = \zeta_{n+j/J} + \frac{\Delta t}{J} \left( -
! </font>
<font color="gray">abla \cdot {\bf F}_j \right)
! sshSubcycleNew(iCell) = sshSubcycleOld(iCell) + dt/J*tend_ssh(iCell)/areaCell(iCell)
@@ -750,10 +746,26 @@
! - grad_ssh + GBtrForcing(iEdge))
! uBtrNew(iEdge) = uBtrNew(iEdge) + aBtrSumCoef(j)*uBtrSubcycleNew(iEdge)
! end do
+
block => domain % blocklist
do while (associated(block))
allocate(tend_uBtr(block % mesh % nEdges))
+
+ do iCell=1,block % mesh % nCells
+
+ block % state % time_levs(2) % state % sshSubcycle % array(iCell) &
+ = block % state % time_levs(1) % state % sshSubcycle % array(iCell) &
+ + dt/config_n_btr_subcycles &
+ * tend_ssh(iCell) / block % mesh % areaCell % array (iCell)
+
+ block % state % time_levs(2) % state % ssh % array(iCell) &
+ = block % state % time_levs(2) % state % ssh % array(iCell) &
+ + aBtrSumCoef(j) &
+ * block % state % time_levs(2) % state % sshSubcycle % array(iCell)
+
+ end do
+
deallocate(tend_uBtr)
block => block % next
end do ! block
</font>
</pre>