<p><b>mpetersen@lanl.gov</b> 2011-02-23 14:50:59 -0700 (Wed, 23 Feb 2011)</p><p>BRANCH COMMIT: Enforce wTop(1,:)=0 for boundary condition. Changed some spacing.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/vert_adv_mrp/src/core_ocean/module_time_integration.F
===================================================================
--- branches/ocean_projects/vert_adv_mrp/src/core_ocean/module_time_integration.F        2011-02-23 16:55:38 UTC (rev 744)
+++ branches/ocean_projects/vert_adv_mrp/src/core_ocean/module_time_integration.F        2011-02-23 21:50:59 UTC (rev 745)
@@ -376,18 +376,18 @@
elseif (config_vert_grid_type.eq.'zlevel') then
- do iEdge=1,nEdges
- cell1 = cellsOnEdge(1,iEdge)
- cell2 = cellsOnEdge(2,iEdge)
- do k=1,min(1,maxLevelEdgeTop(iEdge))
- flux = u(k,iEdge) * dvEdge(iEdge) * h_edge(k,iEdge)
- tend_h(k,cell1) = tend_h(k,cell1) - flux
- tend_h(k,cell2) = tend_h(k,cell2) + flux
+ do iEdge=1,nEdges
+ cell1 = cellsOnEdge(1,iEdge)
+ cell2 = cellsOnEdge(2,iEdge)
+ do k=1,min(1,maxLevelEdgeTop(iEdge))
+ flux = u(k,iEdge) * dvEdge(iEdge) * h_edge(k,iEdge)
+ tend_h(k,cell1) = tend_h(k,cell1) - flux
+ tend_h(k,cell2) = tend_h(k,cell2) + flux
+ end do
end do
- end do
- do iCell=1,nCells
- tend_h(1,iCell) = tend_h(1,iCell) / areaCell(iCell)
- end do
+ do iCell=1,nCells
+ tend_h(1,iCell) = tend_h(1,iCell) / areaCell(iCell)
+ end do
endif ! config_vert_grid_type
@@ -940,7 +940,7 @@
hRatioZLevelK(k) *tracers(iTracer,k-1,iCell) &
+ hRatioZLevelKm1(k)*tracers(iTracer,k ,iCell)
end do
- do k=3,maxLevelCell(iCell)-1
+ do k=3,maxLevelCell(iCell)-1
cSignWTop = sign(flux3Coef,wTop(k,iCell))
do iTracer=1,num_tracers
tracerTop(iTracer,k,iCell) = &
@@ -1057,15 +1057,9 @@
endif ! vertical tracer advection method
do iCell=1,nCellsSolve
- ! For top layer, contribution is from bottom only.
- k=1
+ do k=1,maxLevelCell(iCell)
do iTracer=1,num_tracers
tend_tr(iTracer,k,iCell) = tend_tr(iTracer,k,iCell) &
- + wTop(k+1,iCell)*tracerTop(iTracer,k+1,iCell)
- end do
- do k=2,maxLevelCell(iCell)
- do iTracer=1,num_tracers
- tend_tr(iTracer,k,iCell) = tend_tr(iTracer,k,iCell) &
- ( wTop(k ,iCell)*tracerTop(iTracer,k ,iCell) &
- wTop(k+1,iCell)*tracerTop(iTracer,k+1,iCell))
end do
@@ -1701,7 +1695,7 @@
endif
!
- ! vertical velocity
+ ! vertical velocity through layer interface
!
if (config_vert_grid_type.eq.'isopycnal') then
! set vertical velocity to zero in isopycnal case
@@ -1718,7 +1712,7 @@
do iEdge=1,nEdges
cell1 = cellsOnEdge(1,iEdge)
cell2 = cellsOnEdge(2,iEdge)
- do k=1,maxLevelEdgeBot(iEdge)
+ do k=2,maxLevelEdgeBot(iEdge)
flux = u(k,iEdge) * dvEdge(iEdge)
div_u(k,cell1) = div_u(k,cell1) + flux
div_u(k,cell2) = div_u(k,cell2) - flux
@@ -1726,16 +1720,14 @@
end do
do iCell=1,nCells
- do k=1,maxLevelCell(iCell)
- div_u(k,iCell) = div_u(k,iCell) / areaCell(iCell)
- end do
-
- ! Vertical velocity at bottom is zero.
+ ! Vertical velocity through layer interface at top and
+ ! bottom is zero.
+ wTop(1,iCell) = 0.0
wTop(maxLevelCell(iCell)+1,iCell) = 0.0
- do k=maxLevelCell(iCell),1,-1
- wTop(k,iCell) = wTop(k+1,iCell) - div_u(k,iCell)*h(k,iCell)
+ do k=maxLevelCell(iCell),2,-1
+ wTop(k,iCell) = wTop(k+1,iCell) &
+ - div_u(k,iCell)/areaCell(iCell)*h(k,iCell)
end do
-
end do
deallocate(div_u)
</font>
</pre>