<p><b>xylar@lanl.gov</b> 2012-08-26 05:33:11 -0600 (Sun, 26 Aug 2012)</p><p>* BRANCH COMMIT *<br>
<br>
Added the surfacePressure variable to be read in from the input file; added to the total pressure<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/surface_pressure_ice_shelves/src/core_ocean/Registry
===================================================================
--- branches/ocean_projects/surface_pressure_ice_shelves/src/core_ocean/Registry        2012-08-25 12:32:42 UTC (rev 2124)
+++ branches/ocean_projects/surface_pressure_ice_shelves/src/core_ocean/Registry        2012-08-26 11:33:11 UTC (rev 2125)
@@ -274,6 +274,7 @@
var persistent real uSrcReconstructMeridional ( nVertLevels nCells Time ) 2 o uSrcReconstructMeridional state - -
var persistent real MontPot ( nVertLevels nCells Time ) 2 - MontPot state - -
var persistent real pressure ( nVertLevels nCells Time ) 2 - pressure state - -
+var persistent real surfacePressure ( nCells ) 1 iro surfacePressure mesh - -
var persistent real wTop ( nVertLevelsP1 nCells Time ) 2 - wTop state - -
var persistent real rhoDisplaced ( nVertLevels nCells Time ) 2 - rhoDisplaced state - -
Modified: branches/ocean_projects/surface_pressure_ice_shelves/src/core_ocean/mpas_ocn_tendency.F
===================================================================
--- branches/ocean_projects/surface_pressure_ice_shelves/src/core_ocean/mpas_ocn_tendency.F        2012-08-25 12:32:42 UTC (rev 2124)
+++ branches/ocean_projects/surface_pressure_ice_shelves/src/core_ocean/mpas_ocn_tendency.F        2012-08-26 11:33:11 UTC (rev 2125)
@@ -412,7 +412,7 @@
real (kind=RKIND), dimension(:), pointer :: &
h_s, fVertex, dvEdge, dcEdge, areaCell, areaTriangle, &
- referenceBottomDepth, ssh
+ referenceBottomDepth, ssh, surfacePressure
real (kind=RKIND), dimension(:,:), pointer :: &
weightsOnEdge, kiteAreasOnVertex, h_edge, h, u, v, pressure,&
circulation, vorticity, ke, ke_edge, MontPot, wTop, zMid, &
@@ -468,7 +468,8 @@
maxLevelEdgeTop => grid % maxLevelEdgeTop % array
maxLevelEdgeBot => grid % maxLevelEdgeBot % array
maxLevelVertexBot => grid % maxLevelVertexBot % array
-
+ surfacePressure => grid % surfacePressure % array
+
nCells = grid % nCells
nEdges = grid % nEdges
nVertices = grid % nVertices
@@ -799,9 +800,11 @@
do iCell=1,nCells
! pressure for generalized coordinates
! assume atmospheric pressure at the surface is zero for now.
- pressure(1,iCell) = rho(1,iCell)*gravity &
- * 0.5*h(1,iCell)
+ ! mrp 120724 use same rho as top cell right now.
+ ! Later this would be a rho that represents the land ice.
+ pressure(1,iCell) = surfacePressure(iCell) + rho(1,iCell)*gravity*0.5*h(1,iCell)
+
do k=2,maxLevelCell(iCell)
pressure(k,iCell) = pressure(k-1,iCell) &
+ 0.5*gravity*( rho(k-1,iCell)*h(k-1,iCell) &
</font>
</pre>