<p><b>duda</b> 2009-10-22 11:05:44 -0600 (Thu, 22 Oct 2009)</p><p>Uncomment code for pv upwinding in the normal direction: the<br>
computation of gradPVn is correct for edges bounding real cells, so we<br>
can do upwinding for these same edges. Since we use pv_edge for real<br>
edges and the edgesOnEdge of real edges in compute_tend(), we need to<br>
keep the halo exchange of pv_edge, though.<br>
<br>
M swmodel/src/module_time_integration.F<br>
</p><hr noshade><pre><font color="gray">Modified: trunk/swmodel/src/module_time_integration.F
===================================================================
--- trunk/swmodel/src/module_time_integration.F        2009-10-20 21:44:00 UTC (rev 60)
+++ trunk/swmodel/src/module_time_integration.F        2009-10-22 17:05:44 UTC (rev 61)
@@ -551,7 +551,7 @@
! tdr
!
! Compute height at vertices, pv at vertices, and average pv to edge locations
- ! ( this computes pv_vertex at all vertices bounding real cells )
+ ! ( this computes pv_vertex at all vertices bounding real cells and distance-1 ghost cells )
!
do iVertex = 1,nVertices
if (cellsOnVertex(1,iVertex) > 0 .and. &
@@ -575,7 +575,7 @@
! tdr
!
! Compute gradient of PV in the tangent direction
- ! ( this computes gradPVt at all edges bounding real cells )
+ ! ( this computes gradPVt at all edges bounding real cells and distance-1 ghost cells )
!
do iEdge = 1,nEdges
do k = 1,nVertLevels
@@ -587,7 +587,7 @@
! tdr
!
! Compute pv at the edges
- ! ( this computes pv_edge at all edges bounding real cells )
+ ! ( this computes pv_edge at all edges bounding real cells and distance-1 ghost cells )
!
pv_edge(:,:) = 0.0
do iVertex = 1,nVertices
@@ -616,7 +616,7 @@
! tdr
!
! Compute pv at cell centers
- ! ( this computes pv_cell for all real cells )
+ ! ( this computes pv_cell for all real cells and distance-1 ghost cells )
!
pv_cell(:,:) = 0.0
do iVertex = 1, nVertices
@@ -634,7 +634,7 @@
! tdr
!
! Compute gradient of PV in normal direction
- ! (tdr: 2009-10-02: this is not correct because the pv_cell in the halo is not correct)
+ ! ( this computes gradPVn for all edges bounding real cells )
!
gradPVn(:,:) = 0.0
do iEdge = 1,nEdges
@@ -647,15 +647,14 @@
enddo
! tdr
- ! Modify PV edge with upstream bias. (this can not be done because gradPVn is not correct.)
+ ! Modify PV edge with upstream bias.
!
- ! do iEdge = 1,nEdges
- ! do k = 1,nVertLevels
- ! pv_edge(k,iEdge) = pv_edge(k,iEdge) - 0.5 * u(k,iEdge) * dt * gradPVn(k,iEdge)
- ! enddo
- ! enddo
+ do iEdge = 1,nEdges
+ do k = 1,nVertLevels
+ pv_edge(k,iEdge) = pv_edge(k,iEdge) - 0.5 * u(k,iEdge) * dt * gradPVn(k,iEdge)
+ enddo
+ enddo
-
end subroutine compute_solve_diagnostics
end module time_integration
</font>
</pre>