<p><b>mpetersen@lanl.gov</b> 2010-11-02 16:11:00 -0600 (Tue, 02 Nov 2010)</p><p>BRANCH COMMIT: Remove unneeded 3D variables zMid, zTop, zMidEdge, zTopEdge, ssh. There are redundant with h and h_edge variables.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/topography2_mrp/src/core_ocean/Registry
===================================================================
--- branches/ocean_projects/topography2_mrp/src/core_ocean/Registry 2010-11-02 12:52:48 UTC (rev 589)
+++ branches/ocean_projects/topography2_mrp/src/core_ocean/Registry 2010-11-02 22:11:00 UTC (rev 590)
@@ -167,16 +167,10 @@
var persistent real uReconstructZ ( nVertLevels nCells Time ) 2 o uReconstructZ state - -
var persistent real uReconstructZonal ( nVertLevels nCells Time ) 2 o uReconstructZonal state - -
var persistent real uReconstructMeridional ( nVertLevels nCells Time ) 2 o uReconstructMeridional state - -
-var persistent real zMid ( nVertLevels nCells Time ) 2 - zMid state - -
-var persistent real zTop ( nVertLevelsP1 nCells Time ) 2 - zTop state - -
-var persistent real zMidEdge ( nVertLevels nEdges Time ) 2 - zMidEdge state - -
-var persistent real zTopEdge ( nVertLevelsP1 nEdges Time ) 2 - zTopEdge state - -
var persistent real p ( nVertLevels nCells Time ) 2 o p state - -
-var persistent real pTop ( nVertLevelsP1 nCells Time ) 2 - pTop state - -
var persistent real pZLevel ( nVertLevels nCells Time ) 2 - pZLevel state - -
var persistent real MontPot ( nVertLevels nCells Time ) 2 o MontPot state - -
var persistent real wTop ( nVertLevelsP1 nCells Time ) 2 o wTop state - -
-var persistent real ssh ( nCells Time ) 2 o ssh state - -
var persistent real pgrad ( nVertLevels nCells Time ) 2 o pgrad state - -
var persistent real pgrad_edge ( nVertLevels nEdges Time ) 2 o pgrad_edge state - -
Modified: branches/ocean_projects/topography2_mrp/src/core_ocean/module_global_diagnostics.F
===================================================================
--- branches/ocean_projects/topography2_mrp/src/core_ocean/module_global_diagnostics.F 2010-11-02 12:52:48 UTC (rev 589)
+++ branches/ocean_projects/topography2_mrp/src/core_ocean/module_global_diagnostics.F 2010-11-02 22:11:00 UTC (rev 590)
@@ -37,7 +37,7 @@
real (kind=RKIND) :: areaCellGlobal, areaEdgeGlobal, areaTriangleGlobal
real (kind=RKIND), dimension(:), pointer :: areaCell, dcEdge, dvEdge, areaTriangle, areaEdge
real (kind=RKIND), dimension(:,:), pointer :: h, u, v, h_edge, circulation, vorticity, ke, pv_edge, pv_vertex, &
- pv_cell, gradPVn, gradPVt, zMid, zTop, p, pTop, MontPot, wTop, rho, tracerTemp
+ pv_cell, gradPVn, gradPVt, p, MontPot, wTop, rho, tracerTemp
real (kind=RKIND), dimension(:,:,:), pointer :: tracers
real (kind=RKIND) :: volumeCellGlobal, volumeEdgeGlobal, CFLNumberGlobal
@@ -80,10 +80,7 @@
pv_cell => state % pv_cell % array
gradPVn => state % gradPVn % array
gradPVt => state % gradPVt % array
- zMid => state % zMid % array
- zTop => state % zTop % array
p => state % p % array
- pTop => state % pTop % array
MontPot => state % MontPot % array
variableIndex = 0
@@ -156,30 +153,12 @@
gradPVt(:,1:nEdgesSolve), sums(variableIndex), mins(variableIndex), maxes(variableIndex), verticalSumMins(variableIndex), &
verticalSumMaxes(variableIndex))
- ! zMid
- variableIndex = variableIndex + 1
- call computeFieldAreaWeightedLocalStats(dminfo, nVertLevels, nCellsSolve, areaCell(1:nCellsSolve), &
- zMid(:,1:nCellsSolve), sums(variableIndex), mins(variableIndex), maxes(variableIndex), verticalSumMins(variableIndex), &
- verticalSumMaxes(variableIndex))
-
- ! zTop
- variableIndex = variableIndex + 1
- call computeFieldAreaWeightedLocalStats(dminfo, nVertLevels, nCellsSolve, areaCell(1:nCellsSolve), &
- zTop(:,1:nCellsSolve), sums(variableIndex), mins(variableIndex), maxes(variableIndex), verticalSumMins(variableIndex), &
- verticalSumMaxes(variableIndex))
-
! p
variableIndex = variableIndex + 1
call computeFieldVolumeWeightedLocalStats(dminfo, nVertLevels, nCellsSolve, areaCell(1:nCellsSolve), h(:,1:nCellsSolve), &
p(:,1:nCellsSolve), sums(variableIndex), mins(variableIndex), maxes(variableIndex), verticalSumMins(variableIndex), &
verticalSumMaxes(variableIndex))
- ! pTop
- variableIndex = variableIndex + 1
- call computeFieldVolumeWeightedLocalStats(dminfo, nVertLevels, nCellsSolve, areaCell(1:nCellsSolve), h(:,1:nCellsSolve), &
- pTop(:,1:nCellsSolve), sums(variableIndex), mins(variableIndex), maxes(variableIndex), verticalSumMins(variableIndex), &
- verticalSumMaxes(variableIndex))
-
! MontPot
variableIndex = variableIndex + 1
call computeFieldVolumeWeightedLocalStats(dminfo, nVertLevels, nCellsSolve, areaCell(1:nCellsSolve), h(:,1:nCellsSolve), &
@@ -309,22 +288,10 @@
variableIndex = variableIndex + 1
averages(variableIndex) = sums(variableIndex)/volumeEdgeGlobal
- ! zMid
- variableIndex = variableIndex + 1
- averages(variableIndex) = sums(variableIndex)/(areaCellGlobal*nVertLevels)
-
- ! zTop
- variableIndex = variableIndex + 1
- averages(variableIndex) = sums(variableIndex)/(areaCellGlobal*nVertLevels)
-
! p
variableIndex = variableIndex + 1
averages(variableIndex) = sums(variableIndex)/volumeCellGlobal
- ! pTop
- variableIndex = variableIndex + 1
- averages(variableIndex) = sums(variableIndex)/volumeCellGlobal
-
! MontPot
variableIndex = variableIndex + 1
averages(variableIndex) = sums(variableIndex)/volumeCellGlobal
Modified: branches/ocean_projects/topography2_mrp/src/core_ocean/module_time_integration.F
===================================================================
--- branches/ocean_projects/topography2_mrp/src/core_ocean/module_time_integration.F 2010-11-02 12:52:48 UTC (rev 589)
+++ branches/ocean_projects/topography2_mrp/src/core_ocean/module_time_integration.F 2010-11-02 22:11:00 UTC (rev 590)
@@ -277,7 +277,7 @@
real (kind=RKIND), dimension(:,:), pointer :: &
weightsOnEdge, kiteAreasOnVertex, h_edge, h, u, v, wTop, &
tend_h, tend_u, circulation, vorticity, ke, ke_edge, pv_edge, &
- divergence, MontPot, pZLevel, zMidEdge, zTopEdge
+ divergence, MontPot, pZLevel
type (dm_info) :: dminfo
integer, dimension(:), pointer :: nEdgesOnCell, nEdgesOnEdge, &
@@ -309,8 +309,6 @@
pv_edge => s % pv_edge % array
MontPot => s % MontPot % array
pZLevel => s % pZLevel % array
- zTopEdge => s % zTopEdge % array
- zMidEdge => s % zMidEdge % array
weightsOnEdge => grid % weightsOnEdge % array
kiteAreasOnVertex => grid % kiteAreasOnVertex % array
@@ -660,18 +658,22 @@
allocate(fluxVertTop(nVertLevels+1))
fluxVertTop(1) = 0.0
do iEdge=1,grid % nEdgesSolve
-
+
do k=2,maxLevelEdgeTop(iEdge)
fluxVertTop(k) = vertViscTop(k) &
* ( u(k-1,iEdge) - u(k,iEdge) ) &
- / (zMidEdge(k-1,iEdge) - zMidEdge(k,iEdge))
+ * 2 / (h_edge(k-1,iEdge) + h_edge(k,iEdge))
+! mrp previous line was this, but z variables are redundant with h.
+! / (zMidEdge(k-1,iEdge) - zMidEdge(k,iEdge))
enddo
fluxVertTop(maxLevelEdgeTop(iEdge)+1) = 0.0
do k=1,maxLevelEdgeTop(iEdge)
tend_u(k,iEdge) = tend_u(k,iEdge) &
+ (fluxVertTop(k) - fluxVertTop(k+1)) &
- /(zTopEdge(k,iEdge) - zTopEdge(k+1,iEdge))
+ / h_edge(k,iEdge)
+! mrp previous line was this, but z variables are redundant with h.
+! /(zTopEdge(k,iEdge) - zTopEdge(k+1,iEdge))
enddo
end do
@@ -698,11 +700,11 @@
integer :: i, k, iCell, iEdge, iTracer, cell1, cell2, upwindCell,&
nEdges, nCells, nVertLevels, num_tracers
real (kind=RKIND) :: invAreaCell1, invAreaCell2, tracer_turb_flux
- real (kind=RKIND) :: flux, tracer_edge, r, dist
+ real (kind=RKIND) :: flux, tracer_edge, r
real (kind=RKIND), dimension(:), pointer :: &
h_s, fVertex, fEdge, dvEdge, dcEdge, areaCell, areaTriangle
real (kind=RKIND), dimension(:,:), pointer :: &
- u,h,wTop, h_edge, zMid, zTop
+ u,h,wTop, h_edge
real (kind=RKIND), dimension(:,:,:), pointer :: &
tracers, tend_tr
integer, dimension(:,:), pointer :: boundaryEdge
@@ -728,8 +730,6 @@
wTop => s % wTop % array
tracers => s % tracers % array
h_edge => s % h_edge % array
- zMid => s % zMid % array
- zTop => s % zTop % array
tend_tr => tend % tracers % array
@@ -1065,16 +1065,17 @@
! compute \kappa_v d\phi/dz
fluxVertTop(iTracer,k) = vertDiffTop(k) &
* (tracers(iTracer,k-1,iCell) - tracers(iTracer,k,iCell) )&
- / (zMid(k-1,iCell) -zMid(k,iCell))
+ * 2 / (h(k-1,iCell) + h(k,iCell))
enddo
enddo
fluxVertTop(:,maxLevelCell(iCell)+1) = 0.0
do k=1,maxLevelCell(iCell)
- dist = zTop(k,iCell) - zTop(k+1,iCell)
do iTracer=1,num_tracers
+ ! This is h d/dz( fluxVertTop) but h and dz cancel, so
+ ! reduces to delta( fluxVertTop)
tend_tr(iTracer,k,iCell) = tend_tr(iTracer,k,iCell) &
- + h(k,iCell)*(fluxVertTop(iTracer,k) - fluxVertTop(iTracer,k+1))/dist
+ + fluxVertTop(iTracer,k) - fluxVertTop(iTracer,k+1)
enddo
enddo
@@ -1119,14 +1120,15 @@
real (kind=RKIND), dimension(:), pointer :: &
h_s, fVertex, fEdge, dvEdge, dcEdge, areaCell, areaTriangle, &
- hZLevel, ssh, zMidZLevel, zTopZLevel
+ hZLevel
real (kind=RKIND), dimension(:,:), pointer :: &
weightsOnEdge, kiteAreasOnVertex, h_edge, h, u, v, w, wTop, &
circulation, vorticity, ke, ke_edge, pgrad, pgrad_edge, &
pv_edge, pv_vertex, pv_cell, gradPVn, gradPVt, divergence, &
- zMid, zTop, zMidEdge, zTopEdge, p, pTop, MontPot, rho, temperature, salinity, pZLevel
+ MontPot, rho, temperature, salinity, pZLevel
real (kind=RKIND), dimension(:,:,:), pointer :: tracers
real (kind=RKIND), dimension(:,:), allocatable:: div_u
+ real (kind=RKIND), dimension(:), allocatable:: pTop
character :: c1*6
integer, dimension(:,:), pointer :: cellsOnEdge, cellsOnVertex, &
@@ -1157,15 +1159,8 @@
gradPVt => s % gradPVt % array
rho => s % rho % array
tracers => s % tracers % array
- zMid => s % zMid % array
- zTop => s % zTop % array
- zMidEdge => s % zMidEdge % array
- zTopEdge => s % zTopEdge % array
- p => s % p % array
pZLevel => s % pZLevel % array
- pTop => s % pTop % array
MontPot => s % MontPot % array
- ssh => s % ssh % array
pgrad => s % pgrad % array
pgrad_edge => s % pgrad_edge % array
@@ -1192,8 +1187,6 @@
maxLevelEdgeTop => grid % maxLevelEdgeTop % array
maxLevelEdgeBot => grid % maxLevelEdgeBot % array
maxLevelVertexBot => grid % maxLevelVertexBot % array
- zMidZLevel => grid % zMidZLevel % array
- zTopZLevel => grid % zTopZLevel % array
nCells = grid % nCells
nEdges = grid % nEdges
@@ -1487,13 +1480,6 @@
enddo
!
- ! Compute sea surface height
- !
- do iCell=1,nCells
- ssh(iCell) = h(1,iCell) - hZLevel(1)
- enddo
-
- !
! equation of state
!
! For a isopycnal model, density should remain constant.
@@ -1515,94 +1501,36 @@
endif
!
- ! Z locations
+ ! Pressure
+ ! This section must be after computing rho
!
if (config_vert_grid_type.eq.'isopycnal') then
+ allocate(pTop(nVertLevels+1))
- ! Compute mid- and top-depths of each layer, from bottom up.
- do iCell=1,nCells
- ! h_s is ocean topography: elevation above lowest point,
- ! and is positive. z coordinates are pos upward.
- ! h_s is only used for isopycnal coordinates.
- zTop(nVertLevels+1,iCell) = h_s(iCell)
- do k=nVertLevels,1,-1
- zMid(k,iCell) = zTop(k+1,iCell) + 0.5*h(k,iCell)
- zTop(k,iCell) = zTop(k+1,iCell) + h(k,iCell)
- end do
- end do
-
- do iEdge=1,nEdges
- cell1 = cellsOnEdge(1,iEdge)
- cell2 = cellsOnEdge(2,iEdge)
- do k=1,maxLevelEdgeBot(iEdge)
- zTopEdge(k,iEdge) = (zTop(k,cell1)+zTop(k,cell2))/2.0
- zMidEdge(k,iEdge) = (zMid(k,cell1)+zMid(k,cell2))/2.0
- enddo
- zTopEdge(nVertLevels+1,iEdge) = ( zTop(nVertLevels+1,cell1) &
- + zTop(nVertLevels+1,cell2))/2.0
- enddo
-
- elseif (config_vert_grid_type.eq.'zlevel') then
-
- ! For z-level coordinates, z variables only change for the top
- ! layer, so z variables for k=2:nVertLevels+1 can be computed
- ! from 1D ZLevel variables
- do iCell=1,nCells
- do k=2,nVertLevels
- zMid(k,iCell) = zMidZLevel(k)
- zTop(k,iCell) = zTopZLevel(k)
- end do
- zTop(nVertLevels+1,iCell) = zTopZLevel(nVertLevels+1)
-
- zMid(1,iCell) = zTopZLevel(2) + 0.5*h(1,iCell)
- zTop(1,iCell) = zTopZLevel(2) + h(1,iCell)
- enddo
- zMid(1,nCells+1)=0
- zTop(1,nCells+1)=0
-
- do iEdge=1,nEdges
- cell1 = cellsOnEdge(1,iEdge)
- cell2 = cellsOnEdge(2,iEdge)
- zTopEdge(1,iEdge) = (zTop(1,cell1)+zTop(1,cell2))/2.0
- zMidEdge(1,iEdge) = (zMid(1,cell1)+zMid(1,cell2))/2.0
-
- do k=2,nVertLevels
- zMidEdge(k,iEdge) = zMidZLevel(k)
- zTopEdge(k,iEdge) = zTopZLevel(k)
- end do
- zTopEdge(nVertLevels+1,iEdge) = zTopZLevel(nVertLevels+1)
- enddo
-
- endif
-
- !
- ! Pressure.
- ! This section must be after computing rho and zTop.
- !
- if (config_vert_grid_type.eq.'isopycnal') then
-
! For Isopycnal model.
- ! Compute mid- and bottom-depth of each layer, from bottom up.
- ! Then compute mid- and bottom-pressure of each layer, and
- ! Montgomery Potential, from top down
+ ! Compute pressure at top of each layer, and then
+ ! Montgomery Potential.
do iCell=1,nCells
! assume atmospheric pressure at the surface is zero for now.
- pTop(1,iCell) = 0.0
+ pTop(1) = 0.0
do k=1,nVertLevels
delta_p = rho(k,iCell)*gravity* h(k,iCell)
- p(k ,iCell) = pTop(k,iCell) + 0.5*delta_p
- pTop(k+1,iCell) = pTop(k,iCell) + delta_p
+ pTop(k+1) = pTop(k) + delta_p
end do
- MontPot(1,iCell) = gravity * zTop(1,iCell)
+ ! MontPot at top layer is g*SSH, where SSH may be off by a
+ ! constant i the horizontal.
+ MontPot(1,iCell) = gravity &
+ * (h_s(iCell) + sum(h(1:nVertLevels,iCell)))
do k=2,nVertLevels
! from delta M = p delta / rho
MontPot(k,iCell) = MontPot(k-1,iCell) &
- + pTop(k,iCell)*(1.0/rho(k,iCell) - 1.0/rho(k-1,iCell))
+ + pTop(k)*(1.0/rho(k,iCell) - 1.0/rho(k-1,iCell))
end do
end do
+ deallocate(pTop)
elseif (config_vert_grid_type.eq.'zlevel') then
@@ -1629,7 +1557,9 @@
endif
- ! compute vertical velocity
+ !
+ ! vertical velocity
+ !
if (config_vert_grid_type.eq.'isopycnal') then
! set vertical velocity to zero in isopycnal case
wTop=0.0
Modified: branches/ocean_projects/topography2_mrp/src/core_ocean/mpas_interface.F
===================================================================
--- branches/ocean_projects/topography2_mrp/src/core_ocean/mpas_interface.F 2010-11-02 12:52:48 UTC (rev 589)
+++ branches/ocean_projects/topography2_mrp/src/core_ocean/mpas_interface.F 2010-11-02 22:11:00 UTC (rev 590)
@@ -186,42 +186,12 @@
end do
maxLevelVertexTop(nVertices+1) = 0
- ! mrp temp printing, can remove later
- print *, 'min/max MaxLevelCell, min/max MaxLevelEdge, ',&
- 'min/max MaxLevelVertex'
- print '(20i5)', &
- minval(maxLevelCell(1:nCells)), maxval(maxLevelCell(1:nCells)), &
- minval(maxLevelEdgeTop(1:nEdges)), maxval(maxLevelEdgeTop(1:nEdges)), &
- minval(maxLevelVertexBot(1:nVertices)), maxval(maxLevelVertexBot(1:nVertices))
+ block => block % next
+ end do
- block => block % next
- end do
+ ! Note: We do not update halos on maxLevel* variables. I want the
+ ! outside edge of a halo to be zero on each processor.
- ! update halos for maxLevel variables
-! mrp 101028 actually, dont update halos. I want maxLevel* on the
-! outside edge of a halo to be zero, so we do not loop over those.
-! block => domain % blocklist
-! do while (associated(block))
-
-! call dmpar_exch_halo_field1dInteger(domain % dminfo, &
-! block % mesh % maxLevelCell % array, block % mesh % nCells, &
-! block % parinfo % cellsToSend, block % parinfo % cellsToRecv)
-! call dmpar_exch_halo_field1dInteger(domain % dminfo, &
-! block % mesh % maxLevelEdgeTop % array, block % mesh % nEdges, &
-! block % parinfo % edgesToSend, block % parinfo % edgesToRecv)
-! call dmpar_exch_halo_field1dInteger(domain % dminfo, &
-! block % mesh % maxLevelEdgeBot % array, block % mesh % nEdges, &
-! block % parinfo % edgesToSend, block % parinfo % edgesToRecv)
-! call dmpar_exch_halo_field1dInteger(domain % dminfo, &
-! block % mesh % maxLevelVertexTop % array, block % mesh % nVertices, &
-! block % parinfo % verticesToSend, block % parinfo % verticesToRecv)
-! call dmpar_exch_halo_field1dInteger(domain % dminfo, &
-! block % mesh % maxLevelVertexBot % array, block % mesh % nVertices, &
-! block % parinfo % verticesToSend, block % parinfo % verticesToRecv)
-
-! block => block % next
-! end do
-
end subroutine mpas_init_domain
</font>
</pre>