<p><b>mpetersen@lanl.gov</b> 2012-10-09 13:27:28 -0600 (Tue, 09 Oct 2012)</p><p>-- BRANCH COMMIT --<br>
<br>
partial_bottom_cells branch<br>
<br>
1. Change h_s variable to bottomDepth.<br>
2. Change usage of refBottomDepth to bottomDepth.<br>
3. Add config_set_bottomDepth_to_ref to initialize bottomDepth for non-pbc runs.<br>
4. Merged trunk to branch.<br>
<br>
Tested with 120km global and 3x80x10 overflow domain. Bit-for-bit match if<br>
config_set_bottomDepth_to_ref=true for branch runs.<br>
</p><hr noshade><pre><font color="gray">Index: branches/ocean_projects/partial_bottom_cells
===================================================================
--- branches/ocean_projects/partial_bottom_cells 2012-10-05 15:40:59 UTC (rev 2198)
+++ branches/ocean_projects/partial_bottom_cells 2012-10-09 19:27:28 UTC (rev 2199)
Property changes on: branches/ocean_projects/partial_bottom_cells
___________________________________________________________________
Modified: svn:mergeinfo
## -20,3 +20,4 ##
/branches/omp_blocks/multiple_blocks:1803-2084
/branches/source_renaming:1082-1113
/branches/time_manager:924-962
+/trunk/mpas:2172-2198
\ No newline at end of property
Modified: branches/ocean_projects/partial_bottom_cells/src/core_ocean/Registry
===================================================================
--- branches/ocean_projects/partial_bottom_cells/src/core_ocean/Registry 2012-10-05 15:40:59 UTC (rev 2198)
+++ branches/ocean_projects/partial_bottom_cells/src/core_ocean/Registry 2012-10-09 19:27:28 UTC (rev 2199)
@@ -31,6 +31,7 @@
namelist character grid config_pressure_type pressure
namelist real grid config_rho0 1028
namelist logical grid config_enforce_zstar_at_restart false
+namelist logical grid config_set_bottomDepth_to_ref false
namelist integer split_explicit_ts config_n_ts_iter 2
namelist integer split_explicit_ts config_n_bcl_iter_beg 2
namelist integer split_explicit_ts config_n_bcl_iter_mid 2
@@ -162,7 +163,7 @@
var persistent real kiteAreasOnVertex ( vertexDegree nVertices ) 0 iro kiteAreasOnVertex mesh - -
var persistent real fEdge ( nEdges ) 0 iro fEdge mesh - -
var persistent real fVertex ( nVertices ) 0 iro fVertex mesh - -
-var persistent real h_s ( nCells ) 0 iro h_s mesh - -
+var persistent real bottomDepth ( nCells ) 0 iro bottomDepth mesh - -
% Space needed for advection
var persistent real deriv_two ( maxEdges2 TWO nEdges ) 0 - deriv_two mesh - -
@@ -192,8 +193,8 @@
var persistent integer maxLevelEdgeBot ( nEdges ) 0 - maxLevelEdgeBot mesh - -
var persistent integer maxLevelVertexTop ( nVertices ) 0 - maxLevelVertexTop mesh - -
var persistent integer maxLevelVertexBot ( nVertices ) 0 - maxLevelVertexBot mesh - -
-var persistent real referenceBottomDepth ( nVertLevels ) 0 iro referenceBottomDepth mesh - -
-var persistent real referenceBottomDepthTopOfCell ( nVertLevelsP1 ) 0 - referenceBottomDepthTopOfCell mesh - -
+var persistent real refBottomDepth ( nVertLevels ) 0 iro refBottomDepth mesh - -
+var persistent real refBottomDepthTopOfCell ( nVertLevelsP1 ) 0 - refBottomDepthTopOfCell mesh - -
var persistent real hZLevel ( nVertLevels ) 0 iro hZLevel mesh - -
var persistent real zstarWeight ( nVertLevels ) 0 - zstarWeight mesh - -
@@ -208,10 +209,10 @@
var persistent real temperatureRestore ( nCells ) 0 iro temperatureRestore mesh - -
var persistent real salinityRestore ( nCells ) 0 iro salinityRestore mesh - -
-% mrp trying to figure out why these do not appear
-var persistent real windStressMonthly ( nMonths nEdges ) 0 iro windStressMonthly mesh - -
-var persistent real temperatureRestoreMonthly ( nMonths nCells ) 0 iro temperatureRestoreMonthly mesh - -
-var persistent real salinityRestoreMonthly ( nMonths nCells ) 0 iro salinityRestoreMonthly mesh - -
+% Monthly forcing variables. Change to iro when monthly forcing is used.
+var persistent real windStressMonthly ( nMonths nEdges ) 0 - windStressMonthly mesh - -
+var persistent real temperatureRestoreMonthly ( nMonths nCells ) 0 - temperatureRestoreMonthly mesh - -
+var persistent real salinityRestoreMonthly ( nMonths nCells ) 0 - salinityRestoreMonthly mesh - -
% Prognostic variables: read from input, saved in restart, and written to output
var persistent real u ( nVertLevels nEdges Time ) 2 ir u state - -
Modified: branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_equation_of_state_jm.F
===================================================================
--- branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_equation_of_state_jm.F 2012-10-05 15:40:59 UTC (rev 2198)
+++ branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_equation_of_state_jm.F 2012-10-09 19:27:28 UTC (rev 2199)
@@ -96,7 +96,7 @@
real (kind=RKIND), dimension(:), pointer :: &
- referenceBottomDepth, pRefEOS
+ refBottomDepth, pRefEOS
real (kind=RKIND), dimension(:,:), intent(inout) :: &
rho
real (kind=RKIND), dimension(:,:,:), intent(in) :: tracers
@@ -197,7 +197,7 @@
nCells = grid % nCells
maxLevelCell => grid % maxLevelCell % array
nVertLevels = grid % nVertLevels
- referenceBottomDepth => grid % referenceBottomDepth % array
+ refBottomDepth => grid % refBottomDepth % array
! Jackett and McDougall
@@ -214,14 +214,14 @@
allocate(pRefEOS(nVertLevels),p(nVertLevels),p2(nVertLevels))
! This could be put in the init routine.
- ! Note I am using referenceBottomDepth, so pressure on top level does
+ ! Note I am using refBottomDepth, so pressure on top level does
! not include SSH contribution. I am not sure if that matters, but
! POP does it the same way.
- depth = 0.5*referenceBottomDepth(1)
+ depth = 0.5*refBottomDepth(1)
pRefEOS(1) = 0.059808*(exp(-0.025*depth) - 1.0) &
+ 0.100766*depth + 2.28405e-7*depth**2
do k = 2,nVertLevels
- depth = 0.5*(referenceBottomDepth(k)+referenceBottomDepth(k-1))
+ depth = 0.5*(refBottomDepth(k)+refBottomDepth(k-1))
pRefEOS(k) = 0.059808*(exp(-0.025*depth) - 1.0) &
+ 0.100766*depth + 2.28405e-7*depth**2
enddo
Modified: branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_monthly_forcing.F
===================================================================
--- branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_monthly_forcing.F 2012-10-05 15:40:59 UTC (rev 2198)
+++ branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_monthly_forcing.F 2012-10-09 19:27:28 UTC (rev 2199)
@@ -176,6 +176,8 @@
if(config_use_monthly_forcing) then
monthlyForcingOn = .true.
+
+ write (0,'(a)') " Monthly forcing is on. Make sure monthly forcing variables include iro in Registry, and are in your initial condition or restart file."
end if
!--------------------------------------------------------------------
Modified: branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_mpas_core.F
===================================================================
--- branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_mpas_core.F 2012-10-05 15:40:59 UTC (rev 2198)
+++ branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_mpas_core.F 2012-10-09 19:27:28 UTC (rev 2199)
@@ -541,13 +541,15 @@
type (domain_type), intent(inout) :: domain
- integer :: i, iCell, iEdge, iVertex, k
+ integer :: i, iCell, iEdge, iVertex, k, nCells
type (block_type), pointer :: block
integer :: iTracer, cell, cell1, cell2
real (kind=RKIND) :: uhSum, hSum, hEdge1
- real (kind=RKIND), dimension(:), pointer :: referenceBottomDepth, &
- referenceBottomDepthTopOfCell, zstarWeight, hZLevel
+
+ integer, dimension(:), pointer :: maxLevelCell
+ real (kind=RKIND), dimension(:), pointer :: refBottomDepth, &
+ refBottomDepthTopOfCell, zstarWeight, hZLevel, bottomDepth
real (kind=RKIND), dimension(:,:), pointer :: h
integer :: nVertLevels
@@ -557,28 +559,41 @@
do while (associated(block))
h => block % state % time_levs(1) % state % h % array
- referenceBottomDepth => block % mesh % referenceBottomDepth % array
- referenceBottomDepthTopOfCell => block % mesh % referenceBottomDepthTopOfCell % array
+ refBottomDepth => block % mesh % refBottomDepth % array
+ refBottomDepthTopOfCell => block % mesh % refBottomDepthTopOfCell % array
+ bottomDepth => block % mesh % bottomDepth % array
zstarWeight => block % mesh % zstarWeight % array
hZLevel => block % mesh % hZLevel % array
+ maxLevelCell => block % mesh % maxLevelCell % array
+
+ nCells = block % mesh % nCells
nVertLevels = block % mesh % nVertLevels
! mrp 120208 right now hZLevel is in the grid.nc file.
- ! We would like to transition to using referenceBottomDepth
+ ! We would like to transition to using refBottomDepth
! as the defining variable instead, and will transition soon.
! When the transition is done, hZLevel can be removed from
! registry and the following four lines deleted.
- referenceBottomDepth(1) = hZLevel(1)
+ refBottomDepth(1) = hZLevel(1)
do k = 2,nVertLevels
- referenceBottomDepth(k) = referenceBottomDepth(k-1) + hZLevel(k)
+ refBottomDepth(k) = refBottomDepth(k-1) + hZLevel(k)
end do
! TopOfCell needed where zero depth for the very top may be referenced.
- referenceBottomDepthTopOfCell(1) = 0.0
+ refBottomDepthTopOfCell(1) = 0.0
do k = 1,nVertLevels
- referenceBottomDepthTopOfCell(k+1) = referenceBottomDepth(k)
+ refBottomDepthTopOfCell(k+1) = refBottomDepth(k)
end do
+ ! Some older grids do not have the bottomDepth variable.
+ ! If so, set config_set_bottomDepth_to_ref=.true. and bottomDepth will
+ ! simply be initialized with refBottomDepth. This does not use pbcs.
+ if (config_set_bottomDepth_to_ref) then
+ do iCell = 1,nCells
+ bottomDepth(iCell) = refBottomDepth(maxLevelCell(iCell))
+ enddo
+ endif
+
! Initialization of zstarWeights. This determines how SSH perturbations
! are distributed throughout the column.
if (config_vert_grid_type.eq.'zlevel') then
@@ -605,6 +620,29 @@
endif
+
+ ! Check that bottomDepth and maxLevelCell match forSome older grids do not have the bottomDepth variable.
+ if (config_vert_grid_type.eq.'zlevel'.or. &
+ config_vert_grid_type.eq.'zstar'.or. &
+ config_vert_grid_type.eq.'zstarweights') then
+ do iCell = 1,nCells
+ if (bottomDepth(iCell) > refBottomDepth(maxLevelCell(iCell)).or. &
+ bottomDepth(iCell) < refBottomDepthTopOfCell(maxLevelCell(iCell))) then
+ write (0,*) ' fatal error: bottomDepth and maxLevelCell do not match:'
+ write (0,*) ' iCell, maxLevelCell(iCell), bottomDepth(iCell): ', &
+ iCell, maxLevelCell(iCell), bottomDepth(iCell)
+ write (0,*) ' refBottomDepth(maxLevelCell(iCell)), refBottomDepthTopOfCell(maxLevelCell(iCell)): ', &
+ refBottomDepth(maxLevelCell(iCell)), refBottomDepthTopOfCell(maxLevelCell(iCell))
+ endif
+ ! mrp, for debugging. Delete later:
+ write (0,'(a,2i5,3f10.3)') ' iCell, maxLevelCell(iCell), bottomDepth(iCell), refb, reft: ', &
+ iCell, maxLevelCell(iCell), bottomDepth(iCell), &
+ refBottomDepth(maxLevelCell(iCell)), refBottomDepthTopOfCell(maxLevelCell(iCell))
+ enddo
+ endif
+
+
+
block => block % next
end do
@@ -625,7 +663,7 @@
integer :: iTracer, cell, cell1, cell2
real (kind=RKIND) :: uhSum, hSum, hEdge1
- real (kind=RKIND), dimension(:), pointer :: referenceBottomDepth
+ real (kind=RKIND), dimension(:), pointer :: refBottomDepth
real (kind=RKIND), dimension(:,:), pointer :: h
integer :: nVertLevels
@@ -635,7 +673,7 @@
do while (associated(block))
h => block % state % time_levs(1) % state % h % array
- referenceBottomDepth => block % mesh % referenceBottomDepth % array
+ refBottomDepth => block % mesh % refBottomDepth % array
nVertLevels = block % mesh % nVertLevels
! Compute barotropic velocity at first timestep
@@ -651,7 +689,7 @@
if (config_filter_btr_mode) then
do iCell=1,block % mesh % nCells
block % state % time_levs(1) % state % h % array(1,iCell) &
- = block % mesh % referenceBottomDepth % array(1)
+ = block % mesh % refBottomDepth % array(1)
enddo
endif
@@ -735,7 +773,7 @@
real (kind=RKIND) :: hSum, sumZstarWeights
real (kind=RKIND), dimension(:), pointer :: hZLevel, zstarWeight, &
- referenceBottomDepth
+ refBottomDepth
real (kind=RKIND), dimension(:,:), pointer :: h
! Initialize z-level grid variables from h, read in from input file.
@@ -747,7 +785,7 @@
hZLevel => block % mesh % hZLevel % array
maxLevelCell => block % mesh % maxLevelCell % array
zstarWeight => block % mesh % zstarWeight % array
- referenceBottomDepth => block % mesh % referenceBottomDepth % array
+ refBottomDepth => block % mesh % refBottomDepth % array
do iCell=1,block % mesh % nCells
! Compute the total column thickness, hSum, and the sum of zstar weights.
@@ -762,7 +800,7 @@
! where zeta is SSH and W_k are weights
do k = 1,maxLevelCell(iCell)
h(k,iCell) = hZLevel(k) &
- + (hSum - referenceBottomDepth(maxLevelCell(iCell))) &
+ + (hSum - refBottomDepth(maxLevelCell(iCell))) &
* zstarWeight(k)/sumZstarWeights
enddo
Modified: branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_tendency.F
===================================================================
--- branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_tendency.F 2012-10-05 15:40:59 UTC (rev 2198)
+++ branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_tendency.F 2012-10-09 19:27:28 UTC (rev 2199)
@@ -411,8 +411,7 @@
real (kind=RKIND), dimension(:), allocatable:: pTop
real (kind=RKIND), dimension(:), pointer :: &
- h_s, fVertex, dvEdge, dcEdge, areaCell, areaTriangle, &
- referenceBottomDepth, ssh
+ bottomDepth, fVertex, dvEdge, dcEdge, areaCell, areaTriangle, ssh
real (kind=RKIND), dimension(:,:), pointer :: &
weightsOnEdge, kiteAreasOnVertex, h_edge, h, u, v, pressure,&
circulation, vorticity, ke, ke_edge, MontPot, wTop, zMid, &
@@ -435,11 +434,11 @@
kev => s % kev % array
kevc => s % kevc % array
ke_edge => s % ke_edge % array
- Vor_edge => s % Vor_edge % array
- Vor_vertex => s % Vor_vertex % array
- Vor_cell => s % Vor_cell % array
- gradVor_n => s % gradVor_n % array
- gradVor_t => s % gradVor_t % array
+ Vor_edge => s % Vor_edge % array
+ Vor_vertex => s % Vor_vertex % array
+ Vor_cell => s % Vor_cell % array
+ gradVor_n => s % gradVor_n % array
+ gradVor_t => s % gradVor_t % array
rho => s % rho % array
MontPot => s % MontPot % array
pressure => s % pressure % array
@@ -460,9 +459,8 @@
dvEdge => grid % dvEdge % array
areaCell => grid % areaCell % array
areaTriangle => grid % areaTriangle % array
- h_s => grid % h_s % array
+ bottomDepth => grid % bottomDepth % array
fVertex => grid % fVertex % array
- referenceBottomDepth => grid % referenceBottomDepth % array
deriv_two => grid % deriv_two % array
maxLevelCell => grid % maxLevelCell % array
maxLevelEdgeTop => grid % maxLevelEdgeTop % array
@@ -779,9 +777,9 @@
pTop(1) = 0.0
! For isopycnal mode, p is the Montgomery Potential.
! At top layer it is g*SSH, where SSH may be off by a
- ! constant (ie, h_s can be relative to top or bottom)
+ ! constant (ie, bottomDepth can be relative to top or bottom)
MontPot(1,iCell) = gravity &
- * (h_s(iCell) + sum(h(1:nVertLevels,iCell)))
+ * (bottomDepth(iCell) + sum(h(1:nVertLevels,iCell)))
do k=2,nVertLevels
pTop(k) = pTop(k-1) + rho(k-1,iCell)*gravity* h(k-1,iCell)
@@ -810,10 +808,10 @@
! Compute zMid, the z-coordinate of the middle of the layer.
! This is used for the rho g grad z momentum term.
- ! Note the negative sign, since referenceBottomDepth is positive
+ ! Note the negative sign, since bottomDepth is positive
! and z-coordinates are negative below the surface.
k = maxLevelCell(iCell)
- zMid(k:nVertLevels,iCell) = -referenceBottomDepth(k) + 0.5*h(k,iCell)
+ zMid(k:nVertLevels,iCell) = -bottomDepth(iCell) + 0.5*h(k,iCell)
do k=maxLevelCell(iCell)-1, 1, -1
zMid(k,iCell) = zMid(k+1,iCell) &
@@ -830,13 +828,11 @@
!
do iCell=1,nCells
! Start at the bottom where we know the depth, and go up.
- ! The bottom depth for this cell is
- ! referenceBottomDepth(maxLevelCell(iCell)).
- ! Note the negative sign, since referenceBottomDepth is positive
+ ! The bottom depth for this cell is bottomDepth(iCell).
+ ! Note the negative sign, since bottomDepth is positive
! and z-coordinates are negative below the surface.
- ssh(iCell) = -referenceBottomDepth(maxLevelCell(iCell)) &
- + sum(h(1:maxLevelCell(iCell),iCell))
+ ssh(iCell) = - bottomDepth(iCell) + sum(h(1:maxLevelCell(iCell),iCell))
end do
@@ -878,7 +874,7 @@
real (kind=RKIND), dimension(:), pointer :: &
- h_s, fVertex, fEdge, dvEdge, dcEdge, areaCell, areaTriangle, zstarWeight
+ dvEdge, areaCell, zstarWeight
real (kind=RKIND), dimension(:,:), pointer :: uTransport,h,wTop, h_edge
real (kind=RKIND), dimension(:,:), allocatable:: div_hu
real (kind=RKIND), dimension(:), allocatable:: div_hu_btr, h_tend_col
Modified: branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_test_cases.F
===================================================================
--- branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_test_cases.F 2012-10-05 15:40:59 UTC (rev 2198)
+++ branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_test_cases.F 2012-10-09 19:27:28 UTC (rev 2199)
@@ -342,10 +342,10 @@
do iCell=1,grid % nCells
if (grid % lonCell % array(iCell) < 0.0) grid % lonCell % array(iCell) = grid % lonCell % array(iCell) + 2.0 * pii
r = sqrt(min(rr**2.0, (grid % lonCell % array(iCell) - lambda_c)**2.0 + (grid % latCell % array(iCell) - theta_c)**2.0))
- grid % h_s % array(iCell) = hs0 * (1.0 - r/rr)
+ grid % bottomDepth % array(iCell) = hs0 * (1.0 - r/rr)
end do
! output about mountain
-print *, 'h_s',minval(grid % h_s % array),sum(grid % h_s % array)/grid % nCells, maxval(grid % h_s % array)
+print *, 'bottomDepth',minval(grid % bottomDepth % array),sum(grid % bottomDepth % array)/grid % nCells, maxval(grid % bottomDepth % array)
!
! Initialize tracer fields
@@ -372,7 +372,7 @@
)**2.0 &
) / &
gravity
- state % h % array(1,iCell) = state % h % array(1,iCell) - grid % h_s % array(iCell)
+ state % h % array(1,iCell) = state % h % array(1,iCell) - grid % bottomDepth % array(iCell)
end do
end subroutine sw_test_case_5
Modified: branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_time_integration_rk4.F
===================================================================
--- branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_time_integration_rk4.F 2012-10-05 15:40:59 UTC (rev 2198)
+++ branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_time_integration_rk4.F 2012-10-09 19:27:28 UTC (rev 2199)
@@ -236,9 +236,8 @@
do iCell=1,block % mesh % nCells
do k=1,block % mesh % maxLevelCell % array(iCell)
block % state % time_levs(2) % state % tracers % array(:,k,iCell) = &
- ( block % state % time_levs(2) % state % tracers % array(:,k,iCell) &
- + rk_weights(rk_step) * block % tend % tracers % array(:,k,iCell) )&
- / block % state % time_levs(2) % state % h % array(k, iCell)
+ block % state % time_levs(2) % state % tracers % array(:,k,iCell) &
+ + rk_weights(rk_step) * block % tend % tracers % array(:,k,iCell)
end do
end do
@@ -256,6 +255,20 @@
! A little clean up at the end: decouple new scalar fields and compute diagnostics for new state
!
call mpas_timer_start("RK4-cleaup phase")
+
+ ! Rescale tracers
+ block => domain % blocklist
+ do while(associated(block))
+ do iCell = 1, block % mesh % nCells
+ do k = 1, block % mesh % maxLevelCell % array(iCell)
+ block % state % time_levs(2) % state % tracers % array(:, k, iCell) = block % state % time_levs(2) % state % tracers % array(:, k, iCell) &
+ / block % state % time_levs(2) % state % h % array(k, iCell)
+ end do
+ end do
+ block => block % next
+ end do
+
+
if (config_implicit_vertical_mix) then
call mpas_timer_start("RK4-implicit vert mix")
block => domain % blocklist
Modified: branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_time_integration_split.F
===================================================================
--- branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_time_integration_split.F 2012-10-05 15:40:59 UTC (rev 2198)
+++ branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_time_integration_split.F 2012-10-09 19:27:28 UTC (rev 2199)
@@ -411,10 +411,24 @@
do iEdge=1,block % mesh % nEdges
cell1 = block % mesh % cellsOnEdge % array(1,iEdge)
cell2 = block % mesh % cellsOnEdge % array(2,iEdge)
-
+
+ ! method 0: orig, works only without pbc:
sshEdge = 0.5 * (block % state % time_levs(oldBtrSubcycleTime) % state % sshSubcycle % array(cell1) &
+ block % state % time_levs(oldBtrSubcycleTime) % state % sshSubcycle % array(cell2) )
- hSum = sshEdge + block % mesh % referenceBottomDepthTopOfCell % array (block % mesh % maxLevelEdgeTop % array(iEdge)+1)
+ hSum = sshEdge + block % mesh % refBottomDepthTopOfCell % array (block % mesh % maxLevelEdgeTop % array(iEdge)+1)
+
+ ! method 1, should match before:
+! sshEdge = 0.5 *( block % state % time_levs(oldBtrSubcycleTime) % state % sshSubcycle % array(cell1) &
+! + block % state % time_levs(oldBtrSubcycleTime) % state % sshSubcycle % array(cell2) )
+! hSum = sshEdge + min(block % mesh % bottomDepth % array(cell1), &
+! block % mesh % bottomDepth % array(cell2))
+
+ ! method 2, better, I think.
+ ! take average of full thickness at two neighboring cells
+! hsum = 0.5 *( block % state % time_levs(oldBtrSubcycleTime) % state % sshSubcycle % array(cell1) &
+! + block % state % time_levs(oldBtrSubcycleTime) % state % sshSubcycle % array(cell2) &
+! + block % mesh % bottomDepth % array(cell1) &
+! + block % mesh % bottomDepth % array(cell2) )
flux = ((1.0-config_btr_gam1_uWt1) * block % state % time_levs(oldBtrSubcycleTime) % state % uBtrSubcycle % array(iEdge) &
+ config_btr_gam1_uWt1 * block % state % time_levs(newBtrSubcycleTime) % state % uBtrSubcycle % array(iEdge)) &
@@ -513,7 +527,7 @@
+ config_btr_gam2_SSHWt1 *block % state % time_levs(newBtrSubcycleTime) % state % sshSubcycle % array(cell2)
sshEdge = 0.5 * (sshCell1 + sshCell2)
- hSum = sshEdge + block % mesh % referenceBottomDepthTopOfCell % array (block % mesh % maxLevelEdgeTop % array(iEdge)+1)
+ hSum = sshEdge + block % mesh % refBottomDepthTopOfCell % array (block % mesh % maxLevelEdgeTop % array(iEdge)+1)
flux = ((1.0-config_btr_gam3_uWt2) * block % state % time_levs(oldBtrSubcycleTime) % state % uBtrSubcycle % array(iEdge) &
+ config_btr_gam3_uWt2 * block % state % time_levs(newBtrSubcycleTime) % state % uBtrSubcycle % array(iEdge)) &
Modified: branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F
===================================================================
--- branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F 2012-10-05 15:40:59 UTC (rev 2198)
+++ branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F 2012-10-09 19:27:28 UTC (rev 2199)
@@ -177,22 +177,22 @@
integer :: k, nVertLevels
- real (kind=RKIND), dimension(:), pointer :: referenceBottomDepth
+ real (kind=RKIND), dimension(:), pointer :: refBottomDepth
err = 0
if(.not.tanhViscOn) return
nVertLevels = grid % nVertLevels
- referenceBottomDepth => grid % referenceBottomDepth % array
+ refBottomDepth => grid % refBottomDepth % array
- ! referenceBottomDepth is used here for simplicity. Using zMid and h, which
+ ! refBottomDepth is used here for simplicity. Using zMid and h, which
! vary in time, would give the exact location of the top, but it
! would only change the diffusion value very slightly.
vertViscTopOfEdge = 0.0
do k=2,nVertLevels
vertViscTopOfEdge(k,:) = -(config_max_visc_tanh-config_min_visc_tanh)/2.0 &
- *tanh((referenceBottomDepth(k-1)+config_ZMid_tanh) &
+ *tanh((refBottomDepth(k-1)+config_ZMid_tanh) &
/config_zWidth_tanh) &
+ (config_max_visc_tanh+config_min_visc_tanh)/2
end do
@@ -250,22 +250,22 @@
integer :: k, nVertLevels
- real (kind=RKIND), dimension(:), pointer :: referenceBottomDepth
+ real (kind=RKIND), dimension(:), pointer :: refBottomDepth
err = 0
if(.not.tanhDiffOn) return
nVertLevels = grid % nVertLevels
- referenceBottomDepth => grid % referenceBottomDepth % array
+ refBottomDepth => grid % refBottomDepth % array
- ! referenceBottomDepth is used here for simplicity. Using zMid and h, which
+ ! refBottomDepth is used here for simplicity. Using zMid and h, which
! vary in time, would give the exact location of the top, but it
! would only change the diffusion value very slightly.
vertDiffTopOfCell = 0.0
do k=2,nVertLevels
vertDiffTopOfCell(k,:) = -(config_max_diff_tanh-config_min_diff_tanh)/2.0 &
- *tanh((referenceBottomDepth(k-1)+config_ZMid_tanh) &
+ *tanh((refBottomDepth(k-1)+config_ZMid_tanh) &
/config_zWidth_tanh) &
+ (config_max_diff_tanh+config_min_diff_tanh)/2
end do
Modified: branches/ocean_projects/partial_bottom_cells/src/framework/mpas_io_input.F
===================================================================
--- branches/ocean_projects/partial_bottom_cells/src/framework/mpas_io_input.F 2012-10-05 15:40:59 UTC (rev 2198)
+++ branches/ocean_projects/partial_bottom_cells/src/framework/mpas_io_input.F 2012-10-09 19:27:28 UTC (rev 2199)
@@ -262,7 +262,7 @@
call mpas_dmpar_abort(domain % dminfo)
end if
!write(0,*) 'MGD DEBUGGING time = ', input_obj % time
- write(0,*) 'Restarting model from time ', timeStamp
+ write(0,*) 'Restarting model from time ', trim(timeStamp)
end if
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
</font>
</pre>