<p><b>dwj07@fsu.edu</b> 2011-10-11 13:45:50 -0600 (Tue, 11 Oct 2011)</p><p><br>
-- BRANCH COMMIT --<br>
<br>
Lots of loop fusion in ocn_diagnostic_solve.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tendency.F
===================================================================
--- branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tendency.F 2011-10-11 17:36:58 UTC (rev 1058)
+++ branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tendency.F 2011-10-11 19:45:50 UTC (rev 1059)
@@ -572,7 +572,7 @@
integer :: iEdge, iCell, iVertex, k, cell1, cell2, vertex1, vertex2, eoe, i, j, cov
- real (kind=RKIND) :: flux, vorticity_abs, h_vertex, workpv, rho0Inv
+ real (kind=RKIND) :: flux, vorticity_abs, h_vertex, workpv, rho0Inv, invAreaTri1, invAreaTri2, invAreaCell1, invAreaCell2
integer :: nCells, nEdges, nVertices, nVertLevels, vertexDegree, fCoef, err
@@ -770,66 +770,45 @@
tracers(s % index_temperature,:,nCells+1) = -1e34
tracers(s % index_salinity,:,nCells+1) = -1e34
- !
- ! Compute circulation and relative vorticity at each vertex
- !
circulation(:,:) = 0.0
+ vorticity(:,:) = 0.0
+ divergence(:,:) = 0.0
+ ke(:,:) = 0.0
+ v(:,:) = 0.0
do iEdge=1,nEdges
vertex1 = verticesOnEdge(1,iEdge)
vertex2 = verticesOnEdge(2,iEdge)
- do k=1,maxLevelEdgeBot(iEdge)
- circulation(k,vertex1) = circulation(k,vertex1) - dcEdge(iEdge) * u(k,iEdge)
- circulation(k,vertex2) = circulation(k,vertex2) + dcEdge(iEdge) * u(k,iEdge)
- end do
- end do
- do iVertex=1,nVertices
- do k=1,maxLevelVertexBot(iVertex)
- vorticity(k,iVertex) = circulation(k,iVertex) / areaTriangle(iVertex)
- end do
- end do
- !
- ! Compute the divergence at each cell center
- !
- divergence(:,:) = 0.0
- do iEdge=1,nEdges
cell1 = cellsOnEdge(1,iEdge)
cell2 = cellsOnEdge(2,iEdge)
- do k=1,maxLevelEdgeBot(iEdge)
- divergence(k,cell1) = divergence(k,cell1) + u(k,iEdge)*dvEdge(iEdge)
- divergence(k,cell2) = divergence(k,cell2) - u(k,iEdge)*dvEdge(iEdge)
- enddo
- end do
- do iCell = 1,nCells
- r = 1.0 / areaCell(iCell)
- do k = 1,maxLevelCell(iCell)
- divergence(k,iCell) = divergence(k,iCell) * r
- enddo
- enddo
- !
- ! Compute kinetic energy in each cell
- !
- ke(:,:) = 0.0
- do iEdge=1,nEdges
- cell1 = cellsOnEdge(1,iEdge)
- cell2 = cellsOnEdge(2,iEdge)
+ invAreaTri1 = 1.0 / areaTriangle(vertex1)
+ invAreaTri2 = 1.0 / areaTriangle(vertex2)
+
+ invAreaCell1 = 1.0 / areaCell(cell1)
+ invAreaCell2 = 1.0 / areaCell(cell2)
do k=1,maxLevelEdgeBot(iEdge)
- ke(k,cell1) = ke(k,cell1) + 0.25 * dcEdge(iEdge) * dvEdge(iEdge) * u(k,iEdge)**2.0
- ke(k,cell2) = ke(k,cell2) + 0.25 * dcEdge(iEdge) * dvEdge(iEdge) * u(k,iEdge)**2.0
- enddo
- end do
- do iCell = 1,nCells
- do k = 1,maxLevelCell(iCell)
- ke(k,iCell) = ke(k,iCell) / areaCell(iCell)
- enddo
- enddo
- !
- ! Compute v (tangential) velocities
- !
- v(:,:) = 0.0
- do iEdge = 1,nEdges
+ ! Compute circulation and relative vorticity at each vertex
+ r = dcEdge(iEdge) * u(k,iEdge)
+ circulation(k,vertex1) = circulation(k,vertex1) - r
+ circulation(k,vertex2) = circulation(k,vertex2) + r
+
+ vorticity(k,vertex1) = vorticity(k,vertex1) - r * invAreaTri1
+ vorticity(k,vertex2) = vorticity(k,vertex2) + r * invAreaTri2
+
+ ! Compute the divergence at each cell center
+ r = dvEdge(iEdge) * u(k,iEdge)
+ divergence(k,cell1) = divergence(k,cell1) + r * invAreaCell1
+ divergence(k,cell2) = divergence(k,cell2) - r * invAreaCell2
+
+ ! Compute kinetic energy in each cell
+ r = r * u(k, iEdge) * dcEdge(iEdge)
+ ke(k,cell1) = ke(k,cell1) + 0.25 * r * invAreaCell1
+ ke(k,cell2) = ke(k,cell2) + 0.25 * r * invAreaCell2
+ end do
+
+ ! Compute v (tangential) velocities
do i=1,nEdgesOnEdge(iEdge)
eoe = edgesOnEdge(i,iEdge)
! mrp 101115 note: in order to include flux boundary conditions,
@@ -838,6 +817,7 @@
v(k,iEdge) = v(k,iEdge) + weightsOnEdge(i,iEdge) * u(k, eoe)
end do
end do
+
end do
!
@@ -882,12 +862,11 @@
end do
end do
- !
- ! Compute pv at cell centers
- ! ( this computes pv_cell for all real cells and distance-1 ghost cells )
- !
pv_cell(:,:) = 0.0
+ pv_edge(:,:) = 0.0
do iVertex = 1,nVertices
+ ! Compute pv at cell centers
+ ! ( this computes pv_cell for all real cells and distance-1 ghost cells )
do i=1,vertexDegree
iCell = cellsOnVertex(i,iVertex)
do k = 1,maxLevelCell(iCell)
@@ -896,21 +875,16 @@
/ areaCell(iCell)
enddo
enddo
- enddo
- !
- ! Compute pv at the edges
- ! ( this computes pv_edge at all edges bounding real cells )
- !
- pv_edge(:,:) = 0.0
- do iVertex = 1,nVertices
+ ! Compute pv at the edges
+ ! ( this computes pv_edge at all edges bounding real cells )
do i=1,vertexDegree
iEdge = edgesOnVertex(i,iVertex)
do k=1,maxLevelEdgeBot(iEdge)
pv_edge(k,iEdge) = pv_edge(k,iEdge) + 0.5 * pv_vertex(k,iVertex)
enddo
end do
- end do
+ enddo
!
! Compute gradient of PV in normal direction
@@ -923,25 +897,18 @@
- pv_cell(k,cellsOnEdge(1,iEdge))) &
/ dcEdge(iEdge)
enddo
- enddo
!
! Compute gradient of PV in the tangent direction
! ( this computes gradPVt at all edges bounding real cells and distance-1 ghost cells )
!
- do iEdge = 1,nEdges
do k = 1,maxLevelEdgeBot(iEdge)
gradPVt(k,iEdge) = ( pv_vertex(k,verticesOnEdge(2,iEdge)) &
- pv_vertex(k,verticesOnEdge(1,iEdge))) &
/dvEdge(iEdge)
- enddo
- enddo
-
!
! Modify PV edge with upstream bias.
!
- do iEdge = 1,nEdges
- do k = 1,maxLevelEdgeBot(iEdge)
pv_edge(k,iEdge) = pv_edge(k,iEdge) &
- 0.5 * dt* ( u(k,iEdge) * gradPVn(k,iEdge) &
+ v(k,iEdge) * gradPVt(k,iEdge) )
</font>
</pre>