<p><b>mpetersen@lanl.gov</b> 2010-04-09 11:15:43 -0600 (Fri, 09 Apr 2010)</p><p>Change compute_scalar_tend to use pointers.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/z_level_mrp/mpas/src/core_ocean/module_time_integration.F
===================================================================
--- branches/ocean_projects/z_level_mrp/mpas/src/core_ocean/module_time_integration.F        2010-04-09 17:14:34 UTC (rev 187)
+++ branches/ocean_projects/z_level_mrp/mpas/src/core_ocean/module_time_integration.F        2010-04-09 17:15:43 UTC (rev 188)
@@ -483,29 +483,53 @@
type (grid_state), intent(in) :: s
type (grid_meta), intent(in) :: grid
- integer :: iCell, iEdge, k, iTracer, cell1, cell2
+ integer :: iCell, iEdge, k, iTracer, cell1, cell2, &
+ nEdges, nCells, nVertLevels
real (kind=RKIND) :: flux, tracer_edge
- tend % tracers % array(:,:,:) = 0.0
- do iEdge=1,grid % nEdges
- cell1 = grid % cellsOnEdge % array(1,iEdge)
- cell2 = grid % cellsOnEdge % array(2,iEdge)
+ real (kind=RKIND), dimension(:), pointer :: &
+ h_s, fVertex, fEdge, dvEdge, dcEdge, areaCell, areaTriangle
+ real (kind=RKIND), dimension(:,:), pointer :: &
+ u,h,Fv, h_edge
+ real (kind=RKIND), dimension(:,:,:), pointer :: &
+ tracers, tend_tr
+ integer, dimension(:), pointer :: nEdgesOnCell, nEdgesOnEdge
+ integer, dimension(:,:), pointer :: cellsOnEdge
+
+ u => s % u % array
+ tracers => s % tracers % array
+ h_edge => s % h_edge % array
+
+ tend_tr => tend % tracers % array
+ Fv => tend % w % array
+
+ areaCell => grid % areaCell % array
+ cellsOnEdge => grid % cellsOnEdge % array
+ dvEdge => grid % dvEdge % array
+ nEdges = grid % nEdges
+ nCells = grid % nCells
+ nVertLevels = grid % nVertLevels
+
+ tend_tr(:,:,:) = 0.0
+ do iEdge=1,nEdges
+ cell1 = cellsOnEdge(1,iEdge)
+ cell2 = cellsOnEdge(2,iEdge)
if (cell1 > 0 .and. cell2 > 0) then
- do k=1,grid % nVertLevels
+ do k=1,nVertLevels
do iTracer=1,num_tracers
- tracer_edge = 0.5 * (s % tracers % array(iTracer,k,cell1) + s % tracers % array(iTracer,k,cell2))
- flux = s % u % array(k,iEdge) * grid % dvEdge % array(iEdge) * s % h_edge % array(k,iEdge) * tracer_edge
- tend % tracers % array(iTracer,k,cell1) = tend % tracers % array(iTracer,k,cell1) - flux
- tend % tracers % array(iTracer,k,cell2) = tend % tracers % array(iTracer,k,cell2) + flux
- end do
- end do
+ tracer_edge = 0.5 * (tracers(iTracer,k,cell1) + tracers(iTracer,k,cell2))
+ flux = u(k,iEdge) * dvEdge(iEdge) * h_edge(k,iEdge) * tracer_edge
+ tend_tr(iTracer,k,cell1) = tend_tr(iTracer,k,cell1) - flux
+ tend_tr(iTracer,k,cell2) = tend_tr(iTracer,k,cell2) + flux
+ end do
+ end do
end if
- end do
+ end do
do iCell=1,grid % nCellsSolve
do k=1,grid % nVertLevelsSolve
do iTracer=1,num_tracers
- tend % tracers % array(iTracer,k,iCell) = tend % tracers % array(iTracer,k,iCell) / grid % areaCell % array(iCell)
+ tend_tr(iTracer,k,iCell) = tend_tr(iTracer,k,iCell) / areaCell(iCell)
end do
end do
end do
</font>
</pre>