<p><b>mpetersen@lanl.gov</b> 2012-10-09 14:46:25 -0600 (Tue, 09 Oct 2012)</p><p>partial_bottom_cells branch. I have three methods to compute the barotropic thickness at the cell edge in the code right now. Method 0 is the original, and does not work for pbcs. Method 1: use the min thickness of the neighboring cells. Method 1 works for pbcs, and matches bit-for-bit with method 0. Method 2: use average thickness of two neighboring cells. This was tested, but is not a bit-for-bit match with other methods, as expected.<br>
</p><hr noshade><pre><font color="gray">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-09 19:27:28 UTC (rev 2199)
+++ branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_mpas_core.F 2012-10-09 20:46:25 UTC (rev 2200)
@@ -620,12 +620,26 @@
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
+
+ ! mrp, for debugging. Delete later:
+ write (0,'(a,2i5,10f10.3)') ' iCell, maxLevelCell, bdepth, refb, reft, sumh-bottomDepth: ', &
+ iCell, maxLevelCell(iCell), bottomDepth(iCell), &
+ refBottomDepth(maxLevelCell(iCell)), refBottomDepthTopOfCell(maxLevelCell(iCell)), &
+ sum(h(1:maxLevelCell(iCell),iCell))-bottomDepth(iCell)
+
+ ! Check if abs(ssh)>2m. If so, print warning.
+ if (abs(sum(h(1:maxLevelCell(iCell),iCell))-bottomDepth(iCell))>2.0) then
+ write (0,*) ' Warning: abs(sum(h)-bottomDepth)>2m. Most likely, initial h does not match bottomDepth.'
+ write (0,*) ' iCell, K=maxLevelCell(iCell), bottomDepth(iCell),sum(h),bottomDepth,hZLevel(K),h(K): ', &
+ iCell, maxLevelCell(iCell), bottomDepth(iCell),sum(h(1:maxLevelCell(iCell),iCell)),bottomDepth(iCell), &
+ hZLevel(maxLevelCell(iCell)), h(maxLevelCell(iCell),iCell)
+ endif
+
+ ! Check that bottomDepth and maxLevelCell match forSome older grids do not have the bottomDepth variable.
if (bottomDepth(iCell) > refBottomDepth(maxLevelCell(iCell)).or. &
bottomDepth(iCell) < refBottomDepthTopOfCell(maxLevelCell(iCell))) then
write (0,*) ' fatal error: bottomDepth and maxLevelCell do not match:'
@@ -633,16 +647,12 @@
iCell, maxLevelCell(iCell), bottomDepth(iCell)
write (0,*) ' refBottomDepth(maxLevelCell(iCell)), refBottomDepthTopOfCell(maxLevelCell(iCell)): ', &
refBottomDepth(maxLevelCell(iCell)), refBottomDepthTopOfCell(maxLevelCell(iCell))
+ call mpas_dmpar_finalize(domain % dminfo)
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
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-09 19:27:28 UTC (rev 2199)
+++ branches/ocean_projects/partial_bottom_cells/src/core_ocean/mpas_ocn_time_integration_split.F 2012-10-09 20:46:25 UTC (rev 2200)
@@ -412,24 +412,21 @@
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 % refBottomDepthTopOfCell % array (block % mesh % maxLevelEdgeTop % array(iEdge)+1)
+ ! method 0: orig, works only without pbc:
+ !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))
+ 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) )
-
+ !hSum = sshEdge + 0.5 *( 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)) &
* hSum
@@ -525,9 +522,19 @@
+ config_btr_gam2_SSHWt1 *block % state % time_levs(newBtrSubcycleTime) % state % sshSubcycle % array(cell1)
sshCell2 = (1-config_btr_gam2_SSHWt1)*block % state % time_levs(oldBtrSubcycleTime) % state % sshSubcycle % array(cell2) &
+ config_btr_gam2_SSHWt1 *block % state % time_levs(newBtrSubcycleTime) % state % sshSubcycle % array(cell2)
+ sshEdge = 0.5 * (sshCell1 + sshCell2)
- sshEdge = 0.5 * (sshCell1 + sshCell2)
- hSum = sshEdge + block % mesh % refBottomDepthTopOfCell % array (block % mesh % maxLevelEdgeTop % array(iEdge)+1)
+ ! method 0: orig, works only without pbc:
+ !hSum = sshEdge + block % mesh % refBottomDepthTopOfCell % array (block % mesh % maxLevelEdgeTop % array(iEdge)+1)
+
+ ! method 1, should match before:
+ 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 = sshEdge + 0.5 *( block % mesh % bottomDepth % array(cell1) &
+ ! + block % mesh % bottomDepth % array(cell2) )
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)) &
</font>
</pre>