<p><b>dwj07@fsu.edu</b> 2011-10-27 09:26:27 -0600 (Thu, 27 Oct 2011)</p><p><br>
        -- BRANCH COMMIT --<br>
<br>
        Making edgeMask, cellMask, and vertexMask be initialized when boundaryEdge/Cell/Vertex are. In compute_max_z_level.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/performance/src/core_ocean/mpas_ocn_mpas_core.F
===================================================================
--- branches/ocean_projects/performance/src/core_ocean/mpas_ocn_mpas_core.F        2011-10-27 15:02:51 UTC (rev 1149)
+++ branches/ocean_projects/performance/src/core_ocean/mpas_ocn_mpas_core.F        2011-10-27 15:26:27 UTC (rev 1150)
@@ -218,7 +218,6 @@
call mpas_timer_stop("diagnostic solve", initDiagSolveTimer)
call ocn_compute_mesh_scaling(mesh)
- call build_boundary_masks(mesh) !dwj: need to fix this routine
call mpas_rbf_interp_initialize(mesh)
call mpas_init_reconstruct(mesh)
@@ -640,7 +639,7 @@
maxLevelVertexTop, maxLevelVertexBot
integer, dimension(:,:), pointer :: &
cellsOnEdge, cellsOnVertex, boundaryEdge, boundaryCell, &
- boundaryVertex, verticesOnEdge
+ boundaryVertex, verticesOnEdge, edgeMask, cellMask, vertexMask
! Initialize z-level grid variables from h, read in from input file.
block => domain % blocklist
@@ -657,6 +656,9 @@
boundaryEdge => block % mesh % boundaryEdge % array
boundaryCell => block % mesh % boundaryCell % array
boundaryVertex => block % mesh % boundaryVertex % array
+ edgeMask => block % mesh % edgeMask % array
+ cellMask => block % mesh % cellMask % array
+ vertexMask => block % mesh % vertexMask % array
nCells = block % mesh % nCells
nEdges = block % mesh % nEdges
@@ -711,21 +713,30 @@
! set boundary edge
boundaryEdge=1
+ edgeMask=0
do iEdge=1,nEdges
boundaryEdge(1:maxLevelEdgeTop(iEdge),iEdge)=0
+ edgeMask(1:maxLevelEdgeTop(iEdge),iEdge)=1
end do
!
! Find cells and vertices that have an edge on the boundary
!
- boundaryCell(:,:) = 0
+ boundaryCell = 0
+ boundaryVertex = 0
+ cellMask = 1
+ vertexMask = 1
do iEdge=1,nEdges
do k=1,nVertLevels
if (boundaryEdge(k,iEdge).eq.1) then
boundaryCell(k,cellsOnEdge(1,iEdge)) = 1
boundaryCell(k,cellsOnEdge(2,iEdge)) = 1
+ cellMask(k,cellsOnEdge(1,iEdge)) = 0
+ cellMask(k,cellsOnEdge(2,iEdge)) = 0
boundaryVertex(k,verticesOnEdge(1,iEdge)) = 1
boundaryVertex(k,verticesOnEdge(2,iEdge)) = 1
+ vertexMask(k,verticesOnEdge(1,iEdge)) = 0
+ vertexMask(k,verticesOnEdge(2,iEdge)) = 0
endif
end do
end do
@@ -786,50 +797,6 @@
end subroutine ocn_compute_mesh_scaling!}}}
- subroutine build_boundary_masks(grid)!{{{
-
- use mpas_grid_types
-
- type (mesh_type), intent(inout) :: grid
-
- integer :: i, k
- integer :: nCells, nEdges, nVertices, nVertLevels
- integer, dimension(:,:), pointer :: edgeMask, cellMask, vertexMask
- integer, dimension(:,:), pointer :: boundaryEdge, boundaryCell, boundaryVertex
-
- nCells = grid % nCells
- nEdges = grid % nEdges
- nVertices = grid % nVertices
- nVertLevels = grid % nVertLevels
-
- boundaryCell => grid % boundaryCell % array
- boundaryEdge => grid % boundaryEdge % array
- boundaryVertex => grid % boundaryVertex % array
-
- edgeMask => grid % edgeMask % array
- cellMask => grid % cellMask % array
- vertexMask => grid % vertexMask % array
-
- do i = 1, nCells
- do k = 1, nVertLevels
- cellMask(k, i) = transfer(abs(.not.boundaryCell(k, i)), cellMask(k, i))
- end do
- end do
-
- do i = 1, nEdges
- do k = 1, nVertLevels
- edgeMask(k, i) = transfer(abs(.not.boundaryEdge(k, i)), edgeMask(k, i))
- end do
- end do
-
- do i = 1, nVertices
- do k = 1, nVertLevels
- vertexMask(k, i) = transfer(abs(.not.boundaryVertex(k, i)), vertexMask(k, i))
- end do
- end do
-
- end subroutine build_boundary_masks!}}}
-
end module mpas_core
! vim: foldmethod=marker
</font>
</pre>