<p><b>mpetersen@lanl.gov</b> 2012-02-02 08:52:46 -0700 (Thu, 02 Feb 2012)</p><p>Remove old lines that compute h from SSH. Those are specific to the old z-level implementation. Correct a bug where h_edge is indexed by iCell rather than iEdge.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/ale_split_exp/src/core_ocean/mpas_ocn_mpas_core.F
===================================================================
--- branches/ocean_projects/ale_split_exp/src/core_ocean/mpas_ocn_mpas_core.F 2012-02-01 22:51:17 UTC (rev 1449)
+++ branches/ocean_projects/ale_split_exp/src/core_ocean/mpas_ocn_mpas_core.F 2012-02-02 15:52:46 UTC (rev 1450)
@@ -113,13 +113,6 @@
call ocn_init_z_level(domain)
endif
- if (trim(config_new_btr_variables_from) == 'btr_avg' &
- .and.trim(config_time_integration) == 'unsplit_explicit') then
- print *, ' unsplit_explicit option must use',&
- ' config_new_btr_variables_from==last_subcycle'
- call mpas_dmpar_abort(dminfo)
- endif
-
!
! Initialize core
!
Modified: branches/ocean_projects/ale_split_exp/src/core_ocean/mpas_ocn_time_integration_split.F
===================================================================
--- branches/ocean_projects/ale_split_exp/src/core_ocean/mpas_ocn_time_integration_split.F 2012-02-01 22:51:17 UTC (rev 1449)
+++ branches/ocean_projects/ale_split_exp/src/core_ocean/mpas_ocn_time_integration_split.F 2012-02-02 15:52:46 UTC (rev 1450)
@@ -88,15 +88,9 @@
eoe, oldBtrSubcycleTime, newBtrSubcycleTime, uPerpTime, BtrCorIter, &
n_bcl_iter(config_n_ts_iter), &
vertex1, vertex2, iVertex
-! mrp 120131 remove h computed from SSH in current step
-!, sshSwapFlag
-
type (block_type), pointer :: block
real (kind=RKIND) :: uhSum, hSum, sshEdge, flux, &
- uPerp, uCorr, temp, coef, FBtr_coeff, sshCell1, sshCell2
-! mrp 120131 remove h computed from SSH in current step
-! real (kind=RKIND), dimension(:), pointer :: sshNew
-
+ uPerp, uCorr, temp, temp_h, coef, FBtr_coeff, sshCell1, sshCell2
integer :: num_tracers, ucorr_coef, err
real (kind=RKIND), dimension(:,:), pointer :: &
u, h, h_edge, ke_edge, vertViscTopOfEdge, vertDiffTopOfCell
@@ -106,9 +100,6 @@
real (kind=RKIND), dimension(:), allocatable:: A,C,uTemp
real (kind=RKIND), dimension(:,:), allocatable:: tracersTemp
-! mrp 120131 remove h computed from SSH in current step
-! sshSwapFlag = 1
-
call mpas_timer_start("se timestep", .false., timer_main)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
@@ -246,24 +237,12 @@
/block % mesh % dcEdge % array(iEdge) )
enddo
-! mrp 120131 remove h computed from SSH in current step
-! ! Compute GBtrForcing, the vertically averaged forcing
-! sshEdge = 0.5*( block % state % time_levs(1) % state % ssh % array(cell1) &
-! + block % state % time_levs(1) % state % ssh % array(cell2) )
-! uhSum = (sshEdge + block % mesh % hZLevel % array(1)) * uTemp(1)
-! hSum = sshEdge + block % mesh % hZLevel % array(1)
+ uhSum = block % state % time_levs(2) % state % h_edge % array(1,iEdge) * uTemp(1)
+ hSum = block % state % time_levs(2) % state % h_edge % array(1,iEdge)
- uhSum = block % state % time_levs(2) % state % h_edge % array(1,iCell) * uTemp(1)
- hSum = block % state % time_levs(2) % state % h_edge % array(1,iCell)
-
do k=2,block % mesh % maxLevelEdgeTop % array(iEdge)
-! mrp 120131 remove h computed from SSH in current step
-! uhSum = uhSum + block % mesh % hZLevel % array(k) * uTemp(k)
-! hSum = hSum + block % mesh % hZLevel % array(k)
-
- uhSum = uhSum + block % state % time_levs(2) % state % h_edge % array(k,iCell) * uTemp(k)
- hSum = hSum + block % state % time_levs(2) % state % h_edge % array(k,iCell)
-
+ uhSum = uhSum + block % state % time_levs(2) % state % h_edge % array(k,iEdge) * uTemp(k)
+ hSum = hSum + block % state % time_levs(2) % state % h_edge % array(k,iEdge)
enddo
block % state % time_levs(1) % state % GBtrForcing % array(iEdge) = split*uhSum/hSum/dt
@@ -323,7 +302,8 @@
! For Split_Explicit unsplit, simply set uBtrNew=0, uBtrSubcycle=0, and uNew=uBclNew
block % state % time_levs(2) % state % uBtr % array(:) = 0.0
- block % state % time_levs(newBtrSubcycleTime) % state % uBtrSubcycle % array(:) = 0.0
+! mrp 120131 remove h computed from SSH in current step
+! block % state % time_levs(newBtrSubcycleTime) % state % uBtrSubcycle % array(:) = 0.0
block % state % time_levs(2) % state % u % array(:,:) = block % state % time_levs(2) % state % uBcl % array(:,:)
@@ -343,10 +323,6 @@
do iCell=1,block % mesh % nCells
! sshSubcycleOld = sshOld
block % state % time_levs(oldBtrSubcycleTime) % state % sshSubcycle % array(iCell) = block % state % time_levs(1) % state % ssh % array(iCell)
-
-! mrp 120131 remove h computed from SSH in current step
-! ! sshNew = sshOld This is the first for the summation
-! block % state % time_levs(2) % state % ssh % array(iCell) = block % state % time_levs(1) % state % ssh % array(iCell)
enddo
do iEdge=1,block % mesh % nEdges
@@ -441,7 +417,7 @@
! config_btr_gam1_uWt1= 1 flux = uBtrNew*H
! config_btr_gam1_uWt1=0.5 flux = 1/2*(uBtrNew+uBtrOld)*H
! config_btr_gam1_uWt1= 0 flux = uBtrOld*H
- ! mrp 120201 efficiency: could we combint the following edge and cell loops?
+ ! mrp 120201 efficiency: could we combine the following edge and cell loops?
do iEdge=1,block % mesh % nEdges
cell1 = block % mesh % cellsOnEdge % array(1,iEdge)
cell2 = block % mesh % cellsOnEdge % array(2,iEdge)
@@ -550,7 +526,7 @@
! config_btr_gam3_uWt2= 1 flux = uBtrNew*H
! config_btr_gam3_uWt2=0.5 flux = 1/2*(uBtrNew+uBtrOld)*H
! config_btr_gam3_uWt2= 0 flux = uBtrOld*H
- ! mrp 120201 efficiency: could we combint the following edge and cell loops?
+ ! mrp 120201 efficiency: could we combine the following edge and cell loops?
do iEdge=1,block % mesh % nEdges
cell1 = block % mesh % cellsOnEdge % array(1,iEdge)
cell2 = block % mesh % cellsOnEdge % array(2,iEdge)
@@ -599,13 +575,6 @@
block => domain % blocklist
do while (associated(block))
-! mrp 120131 remove h computed from SSH in current step
-! ! Accumulate SSH in running sum over the subcycles.
-! do iCell=1,block % mesh % nCells
-! block % state % time_levs(2) % state % ssh % array(iCell) = block % state % time_levs(2) % state % ssh % array(iCell) &
-! + block % state % time_levs(newBtrSubcycleTime) % state % sshSubcycle % array(iCell)
-! end do
-
! uBtrNew = uBtrNew + uBtrSubcycleNEW
! This accumulates the sum.
! If the Barotropic Coriolis iteration is limited to one, this could
@@ -641,19 +610,6 @@
/ (config_n_btr_subcycles*config_btr_subcycle_loop_factor + 1)
end do
-! mrp 120131 remove h computed from SSH in current step
-! if (config_SSH_from=='avg_of_SSH_subcycles') then
-! do iCell=1,block % mesh % nCells
-! block % state % time_levs(2) % state % ssh % array(iCell) = block % state % time_levs(2) % state % ssh % array(iCell) &
-! / (config_n_btr_subcycles*config_btr_subcycle_loop_factor + 1)
-! end do
-! elseif (config_SSH_from=='avg_flux') then
-! ! see below
-! else
-! write(0,*) 'Abort: Unknown config_SSH_from option: ' //trim(config_SSH_from)
-! call mpas_dmpar_abort(dminfo)
-! endif
-
block => block % next
end do ! block
@@ -703,8 +659,6 @@
! end do
!!!!!!!!!!!!!!!! END do not delete this completely yet.
-! mrp 120131 remove h computed from SSH in current step
-! endif
! Correction velocity uCorr = (Flux - Sum(h u*))/H
! or, for the full latex version:
@@ -719,29 +673,17 @@
endif
do iEdge=1,block % mesh % nEdges
-! mrp 120131 remove h computed from SSH in current step
-! cell1 = block % mesh % cellsOnEdge % array(1,iEdge)
-! cell2 = block % mesh % cellsOnEdge % array(2,iEdge)
-! mrp 120131 remove h computed from SSH in current step
-! sshEdge = 0.5 * (block % state % time_levs(2) % state % ssh % array(cell1) &
-! + block % state % time_levs(2) % state % ssh % array(cell2))
-
! This is u^{avg}
uTemp(:) = block % state % time_levs(2) % state % uBtr % array(iEdge) &
+ block % state % time_levs(2) % state % uBcl % array(:,iEdge)
- uhSum = block % state % time_levs(2) % state % h_edge % array(1,iCell) * uTemp(1)
- hSum = block % state % time_levs(2) % state % h_edge % array(1,iCell)
+ uhSum = block % state % time_levs(2) % state % h_edge % array(1,iEdge) * uTemp(1)
+ hSum = block % state % time_levs(2) % state % h_edge % array(1,iEdge)
do k=2,block % mesh % maxLevelEdgeTop % array(iEdge)
-! mrp 120131 remove h computed from SSH in current step
-! uhSum = uhSum + block % mesh % hZLevel % array(k) * uTemp(k)
-! hSum = hSum + block % mesh % hZLevel % array(k)
-
- uhSum = uhSum + block % state % time_levs(2) % state % h_edge % array(k,iCell) * uTemp(k)
- hSum = hSum + block % state % time_levs(2) % state % h_edge % array(k,iCell)
-
+ uhSum = uhSum + block % state % time_levs(2) % state % h_edge % array(k,iEdge) * uTemp(k)
+ hSum = hSum + block % state % time_levs(2) % state % h_edge % array(k,iEdge)
enddo
uCorr = ucorr_coef*(( block % state % time_levs(1) % state % FBtr % array(iEdge) - uhSum)/hSum)
@@ -759,34 +701,16 @@
enddo
endif
-! mrp 120131 remove h computed from SSH in current step
-! ! Put new sshEdge values in h_edge array, for the OcnTendScalar call below.
-! block % state % time_levs(2) % state % h_edge % array(1,iEdge) = sshEdge + block % mesh % hZLevel % array(1)
-! do k=2,block % mesh % nVertLevels
-! block % state % time_levs(2) % state % h_edge % array(k,iEdge) = block % mesh % hZLevel % array(k)
-! enddo
-
end do ! iEdge
-! mrp 120131 remove h computed from SSH in current step
-! ! Put new SSH values in h array, for the OcnTendScalar call below.
-! do iCell=1,block % mesh % nCells
-! block % state % time_levs(2) % state % h % array(1,iCell) = block % state % time_levs(2) % state % ssh % array(iCell) &
-! + block % mesh % hZLevel % array(1)
-! ! mrp 110601 efficiency note: Since h just moves back and forth between pointers,
-! ! this is not necessary once initialized.
-! do k=2,block % mesh % nVertLevels * sshSwapFlag
-! block % state % time_levs(2) % state % h % array(k,iCell) = block % mesh % hZLevel % array(k)
-! enddo
-! enddo ! iCell
-
deallocate(uTemp)
block => block % next
end do ! block
-! mrp 120131 remove h computed from SSH in current step
-! sshSwapFlag = 0
- endif ! split_explicit
+
+
+ endif ! split_explicit
+
call mpas_timer_stop("se btr vel", timer_btr_vel)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
@@ -799,30 +723,13 @@
do while (associated(block))
call ocn_wtop(block % state % time_levs(1) % state,block % state % time_levs(2) % state, block % mesh)
-! mrp 120131 remove h computed from SSH in current step
-! if (trim(config_time_integration) == 'unsplit_explicit') then
-! call ocn_tend_h(block % tend, block % state % time_levs(2) % state , block % diagnostics, block % mesh)
-! endif
call ocn_tend_h (block % tend, block % state % time_levs(2) % state , block % diagnostics, block % mesh)
call ocn_tend_scalar(block % tend, block % state % time_levs(2) % state , block % diagnostics, block % mesh)
block => block % next
end do
-! mrp 120131 remove h computed from SSH in current step
-! ! update halo for thickness for unsplit only
-! if (trim(config_time_integration) == 'unsplit_explicit') then
-! call mpas_timer_start("se halo h", .false., timer_halo_h)
-! block => domain % blocklist
-! do while (associated(block))
-! call mpas_dmpar_exch_halo_field2d_real(domain % dminfo, block % tend % h % array(:,:), &
-! block % mesh % nVertLevels, block % mesh % nCells, block % parinfo % cellsToSend, block % parinfo % cellsToRecv)
-! block => block % next
-! end do
-! call mpas_timer_stop("se halo h", timer_halo_h)
-! endif ! unsplit_explicit
-
- ! update halo for thickness and tracers
+ ! update halo for thickness and tracer tendencies
call mpas_timer_start("se halo h", .false., timer_halo_h)
block => domain % blocklist
do while (associated(block))
@@ -840,74 +747,34 @@
block => domain % blocklist
do while (associated(block))
-! mrp 120131 remove h computed from SSH in current step
-! allocate(hNew(block % mesh % nVertLevels))
-! mrp 120131 remove h computed from SSH in current step
-! if (trim(config_new_btr_variables_from) == 'last_subcycle') then
-! ! This points to the last barotropic SSH subcycle
-! sshNew => block % state % time_levs(oldBtrSubcycleTime) % state % sshSubcycle % array
-! elseif (trim(config_new_btr_variables_from) == 'btr_avg') then
-! ! This points to the tendency variable SSH*
-! sshNew => block % state % time_levs(2) % state % ssh % array
-! endif
-
-! mrp 120131 remove h computed from SSH in current step
-! if (trim(config_time_integration) == 'unsplit_explicit') then
-! do iCell=1,block % mesh % nCells
-! ! this is h_{n+1}
-! do k = 1,maxLevelCell(iCell)
-! block % state % time_levs(2) % state % h % array(k,iCell) &
-! = block % state % time_levs(1) % state % h % array(k,iCell) &
-! + dt* block % tend % h % array(k,iCell)
-! enddo
-
-! mrp 120131 remove h computed from SSH in current step
-! ! this is only for the hNew computation below, so there is the correct
-! ! value in the ssh variable for unsplit_explicit case.
-! block % state % time_levs(oldBtrSubcycleTime) % state % sshSubcycle % array(iCell) &
-! = block % state % time_levs(2) % state % h % array(1,iCell) &
-! - block % mesh % hZLevel % array(1)
-! end do ! iCell
-
-! mrp 120131 remove h computed from SSH in current step
-! endif ! unsplit_explicit
-
-
!!!!!!!!!!!!!!! not on last iteration: reset variables for next iteration
if (split_explicit_step < config_n_ts_iter) then
! Only need T & S for earlier iterations,
! then all the tracers needed the last time through.
-! mrp 120131 remove h computed from SSH in current step
-! hNew(:) = block % mesh % hZLevel % array(:)
do iCell=1,block % mesh % nCells
! sshNew is a pointer, defined above.
-! mrp 120131 remove h computed from SSH in current step
-! hNew(1) = sshNew(iCell) + block % mesh % hZLevel % array(1)
do k=1,block % mesh % maxLevelCell % array(iCell)
! this is h_{n+1}
- temp &
+ temp_h &
= block % state % time_levs(1) % state % h % array(k,iCell) &
+ dt* block % tend % h % array(k,iCell)
! this is h_{n+1/2}
block % state % time_levs(2) % state % h % array(k,iCell) &
= 0.5*( block % state % time_levs(1) % state % h % array(k,iCell) &
- + temp)
+ + temp_h)
do i=1,2
! This is Phi at n+1
temp = ( block % state % time_levs(1) % state % tracers % array(i,k,iCell) &
* block % state % time_levs(1) % state % h % array(k,iCell) &
+ dt * block % tend % tracers % array(i,k,iCell)) &
- / block % state % time_levs(2) % state % h % array(k,iCell)
+ / temp_h
-! mrp 120131 remove h computed from SSH in current step
-! + dt * block % tend % tracers % array(i,k,iCell)) / hNew(k)
-
! This is Phi at n+1/2
block % state % time_levs(2) % state % tracers % array(i,k,iCell) = 0.5*( &
block % state % time_levs(1) % state % tracers % array(i,k,iCell) + temp )
@@ -917,18 +784,6 @@
! mrp 120131 remove h computed from SSH in current step
-! if (trim(config_time_integration) == 'unsplit_explicit') then
-! ! compute h*, which is h at n+1/2 and put into array hNew
-! ! on last iteration, hNew remains at n+1
-! do iCell=1,block % mesh % nCells
-! block % state % time_levs(2) % state % h % array(1,iCell) = 0.5*( &
-! block % state % time_levs(2) % state % h % array(1,iCell) &
-! + block % state % time_levs(1) % state % h % array(1,iCell) )
-! end do ! iCell
-! mrp 120131 remove h computed from SSH in current step
-! endif ! unsplit_explicit
-
-! mrp 120131 remove h computed from SSH in current step
! this next section is probably not needed - just u^{tr} should be fine.
! compute u*, the velocity for tendency terms. Put in uNew.
! uBclNew is at time n+1/2 here.
@@ -952,12 +807,7 @@
!!!!!!!!!!!!!!! On last iteration: Compute all variables at time n+1
elseif (split_explicit_step == config_n_ts_iter) then
-! mrp 120131 remove h computed from SSH in current step
-! hNew(:) = block % mesh % hZLevel % array(:)
do iCell=1,block % mesh % nCells
-! mrp 120131 remove h computed from SSH in current step
-! ! sshNew is a pointer, defined above.
-! hNew(1) = sshNew(iCell) + block % mesh % hZLevel % array(1)
do k=1,block % mesh % maxLevelCell % array(iCell)
! this is h_{n+1}
@@ -973,9 +823,6 @@
+ dt * block % tend % tracers % array(i,k,iCell)) &
/ block % state % time_levs(2) % state % h % array(k,iCell)
-! mrp 120131 remove h computed from SSH in current step
-! + dt * block % tend % tracers % array(i,k,iCell)) / hNew(k)
-
enddo
end do
end do
@@ -999,8 +846,6 @@
enddo ! iEdges
endif ! split_explicit_step
-! mrp 120131 remove h computed from SSH in current step
- ! deallocate(hNew)
block => block % next
end do
@@ -1014,6 +859,9 @@
do while (associated(block))
call mpas_dmpar_exch_halo_field3d_real(domain % dminfo, block % state % time_levs(2) % state % tracers % array(:,:,:), &
block % tend % num_tracers, block % mesh % nVertLevels, block % mesh % nCells, block % parinfo % cellsToSend, block % parinfo % cellsToRecv)
+
+
+
block => block % next
end do
call mpas_timer_stop("se halo tracers", timer_halo_tracers)
@@ -1036,54 +884,11 @@
! END large iteration loop
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-! mrp 120131 remove h computed from SSH in current step
-! !
-! ! A little clean up at the end: decouple new scalar fields and compute diagnostics for new state
-! !
block => domain % blocklist
do while (associated(block))
-! mrp 120131 remove h computed from SSH in current step
-! if (trim(config_new_btr_variables_from) == 'last_subcycle') then
-! do iEdge=1,block % mesh % nEdges
-! ! uBtrNew = uBtrSubcycleNew (old here is because counter already flipped)
-! ! This line is not needed if u is resplit at the beginning of the timestep.
-! block % state % time_levs(2) % state % uBtr % array(iEdge) &
-! = block % state % time_levs(oldBtrSubcycleTime) % state % uBtrSubcycle % array(iEdge)
-! enddo ! iEdges
-! elseif (trim(config_new_btr_variables_from) == 'btr_avg') then
-! ! uBtrNew from u*. this is done above, so u* is already in
-! ! block % state % time_levs(2) % state % uBtr % array(iEdge)
-! else
-! write(0,*) 'Abort: Unknown config_new_btr_variables_from: '//trim(config_time_integration)
-! call mpas_dmpar_abort(dminfo)
-! endif
-
-! mrp 120131 remove h computed from SSH in current step
-! if (trim(config_time_integration) == 'split_explicit') then
-! if (trim(config_new_btr_variables_from) == 'last_subcycle') then
-! do iCell=1,block % mesh % nCells
-! ! SSH for the next step is from the end of the barotropic subcycle.
-! block % state % time_levs(2) % state % ssh % array(iCell) = block % state % time_levs(oldBtrSubcycleTime) % state % sshSubcycle % array(iCell)
-! end do ! iCell
-! elseif (trim(config_new_btr_variables_from) == 'btr_avg') then
-! ! sshNew from ssh*. This is done above, so ssh* is already in
-! ! block % state % time_levs(2) % state % ssh % array(iCell)
-! endif
-! do iCell=1,block % mesh % nCells
-! ! Put new SSH values in h array, for the OcnTendScalar call below.
-! block % state % time_levs(2) % state % h % array(1,iCell) = block % state % time_levs(2) % state % ssh % array(iCell) + block % mesh % hZLevel % array(1)
-! ! mrp 110601 efficiency note: Since h just moves back and forth between pointers,
-! ! this is not necessary once initialized.
-! do k=2,block % mesh % nVertLevels * sshSwapFlag
-! block % state % time_levs(2) % state % h % array(k,iCell) = block % mesh % hZLevel % array(k)
-! end do
-! end do ! iCell
-! sshSwapFlag = 0
-! end if ! split_explicit
-
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
! Implicit vertical mixing, done after timestep is complete
@@ -1133,6 +938,8 @@
block => block % next
end do
call mpas_timer_stop("se timestep", timer_main)
+
+
end subroutine ocn_time_integrator_split!}}}
subroutine filter_btr_mode_tend_u(tend, s, d, grid)!{{{
</font>
</pre>