<p><b>dwj07@fsu.edu</b> 2011-11-04 12:53:55 -0600 (Fri, 04 Nov 2011)</p><p><br>
-- BRANCH COMMIT --<br>
<br>
Loop masking in ocn_diagnostic_solve.<br>
<br>
Also, added a routine to mpas_ocn_tendency to intialize flags used for loop masking, and ke cell/vertex computations.<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-11-04 18:11:18 UTC (rev 1170)
+++ branches/ocean_projects/performance/src/core_ocean/mpas_ocn_mpas_core.F 2011-11-04 18:53:55 UTC (rev 1171)
@@ -83,6 +83,9 @@
call ocn_equation_of_state_init(err_tmp)
err = err .or. err_tmp
+ call ocn_tendency_init(err_tmp)
+ err = err .or. err_tmp
+
if(err.eq.1) then
call mpas_dmpar_abort(dminfo)
endif
Modified: branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tendency.F
===================================================================
--- branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tendency.F 2011-11-04 18:11:18 UTC (rev 1170)
+++ branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tendency.F 2011-11-04 18:53:55 UTC (rev 1171)
@@ -64,7 +64,8 @@
ocn_tend_scalar, &
ocn_diagnostic_solve, &
ocn_wtop, &
- ocn_fuperp
+ ocn_fuperp, &
+ ocn_tendency_init
!--------------------------------------------------------------------
!
@@ -72,7 +73,11 @@
!
!--------------------------------------------------------------------
+ integer :: hadv2nd, hadv3rd, hadv4th
+ integer :: ke_cell_flag, ke_vertex_flag
+ real (kind=RKIND) :: coef_3rd_order, fCoef
+
!***********************************************************************
contains
@@ -416,7 +421,7 @@
integer :: iEdge, iCell, iVertex, k, cell1, cell2, vertex1, vertex2, eoe, i, j
- integer :: boundaryMask, velMask, nCells, nEdges, nVertices, nVertLevels, vertexDegree, fCoef, err
+ integer :: boundaryMask, velMask, nCells, nEdges, nVertices, nVertLevels, vertexDegree, err
integer, dimension(:), pointer :: nEdgesOnCell, nEdgesOnEdge, &
maxLevelCell, maxLevelEdgeTop, maxLevelEdgeBot, &
@@ -426,6 +431,8 @@
real (kind=RKIND) :: d2fdx2_cell1, d2fdx2_cell2, coef_3rd_order, r_tmp, invAreaCell1, invAreaCell2, invAreaTri1, invAreaTri2, invLength, h_vertex
+ real (kind=RKIND), dimension(:), allocatable:: pTop
+
real (kind=RKIND), dimension(:), pointer :: &
h_s, fVertex, dvEdge, dcEdge, areaCell, areaTriangle, &
hZLevel
@@ -435,7 +442,6 @@
pv_edge, pv_vertex, pv_cell, gradPVn, gradPVt, divergence, &
rho, temperature, salinity, kev, kevc
real (kind=RKIND), dimension(:,:,:), pointer :: tracers, deriv_two
- real (kind=RKIND), dimension(:), allocatable:: pTop
real (kind=RKIND), dimension(:,:), allocatable:: div_u
character :: c1*6
@@ -499,12 +505,12 @@
! mrp 110516 efficiency note: For z-level, only do this on level 1. h_edge for all
! lower levels is defined by hZlevel.
- coef_3rd_order = 0.
- if (config_thickness_adv_order == 3) coef_3rd_order = 1.0
- if (config_thickness_adv_order == 3 .and. config_monotonic) coef_3rd_order = 0.25
+! coef_3rd_order = 0.
+! if (config_thickness_adv_order == 3) coef_3rd_order = 1.0
+! if (config_thickness_adv_order == 3 .and. config_monotonic) coef_3rd_order = 0.25
- if (config_thickness_adv_order == 2) then
- do iEdge=1,nEdges
+! if (config_thickness_adv_order == 2) then
+ do iEdge=1,nEdges*hadv2nd
cell1 = cellsOnEdge(1,iEdge)
cell2 = cellsOnEdge(2,iEdge)
do k=1,maxLevelEdgeTop(iEdge)
@@ -512,9 +518,9 @@
end do
end do
- else if (config_thickness_adv_order == 3) then
+! else if (config_thickness_adv_order == 3) then
- do iEdge=1,nEdges
+ do iEdge=1,nEdges*hadv3rd
cell1 = cellsOnEdge(1,iEdge)
cell2 = cellsOnEdge(2,iEdge)
@@ -548,9 +554,9 @@
end do ! do k
end do ! do iEdge
- else if (config_thickness_adv_order == 4) then
+! else if (config_thickness_adv_order == 4) then
- do iEdge=1,nEdges
+ do iEdge=1,nEdges*hadv4th
cell1 = cellsOnEdge(1,iEdge)
cell2 = cellsOnEdge(2,iEdge)
@@ -583,7 +589,7 @@
end do ! do k
end do ! do iEdge
- endif ! if(config_thickness_adv_order == 2)
+! endif ! if(config_thickness_adv_order == 2)
!
! set the velocity and height at dummy address
@@ -649,22 +655,24 @@
! Compute kinetic energy in each vertex
!
kev(:,:) = 0.0; kevc(:,:) = 0.0
- do iEdge=1,nEdges
+ do iEdge=1,nEdges*ke_vertex_flag
do k=1,nVertLevels
- kev(k,verticesOnEdge(1,iEdge)) = kev(k,verticesOnEdge(1,iEdge)) + dcEdge(iEdge) * dvEdge(iEdge) * u(k,iEdge)**2
- kev(k,verticesOnEdge(2,iEdge)) = kev(k,verticesOnEdge(2,iEdge)) + dcEdge(iEdge) * dvEdge(iEdge) * u(k,iEdge)**2
+ r_tmp = dcEdge(iEdge) * dvEdge(iEdge) * u(k, iEdge)**2
+ kev(k,verticesOnEdge(1,iEdge)) = kev(k,verticesOnEdge(1,iEdge)) + r_tmp
+ kev(k,verticesOnEdge(2,iEdge)) = kev(k,verticesOnEdge(2,iEdge)) + r_tmp
end do
end do
- do iVertex = 1,nVertices
+ do iVertex = 1,nVertices*ke_vertex_flag
do k=1,nVertLevels
- kev(k,iVertex) = kev(k,iVertex) / areaTriangle(iVertex) / 4.0
+ kev(k,iVertex) = kev(k,iVertex) / areaTriangle(iVertex) * 0.25
enddo
enddo
- do iVertex = 1, nVertices
+ do iVertex = 1, nVertices*ke_vertex_flag
do i=1,grid % vertexDegree
iCell = cellsOnVertex(i,iVertex)
+ invAreaCell1 = 1.0 / areaCell(iCell)
do k=1,nVertLevels
- kevc(k,iCell) = kevc(k,iCell) + kiteAreasOnVertex(i, iVertex) * kev(k, iVertex) / areaCell(iCell)
+ kevc(k,iCell) = kevc(k,iCell) + kiteAreasOnVertex(i, iVertex) * kev(k, iVertex) * invAreaCell1
enddo
enddo
enddo
@@ -672,20 +680,17 @@
!
! Compute kinetic energy in each cell by blending ke and kevc
!
- if(config_include_KE_vertex) then
- do iCell=1,nCells
+ do iCell=1,nCells*ke_vertex_flag
do k=1,nVertLevels
ke(k,iCell) = 5.0/8.0*ke(k,iCell) + 3.0/8.0*kevc(k,iCell)
end do
end do
- endif
!
! Compute ke on cell edges at velocity locations for quadratic bottom drag.
!
! mrp 101025 efficiency note: we could get rid of ke_edge completely by
! using sqrt(u(k,iEdge)**2 + v(k,iEdge)**2) in its place elsewhere.
- ke_edge = 0.0 !mrp remove 0 for efficiency
do iEdge=1,nEdges
cell1 = cellsOnEdge(1,iEdge)
cell2 = cellsOnEdge(2,iEdge)
@@ -698,25 +703,14 @@
! Compute height at vertices, pv at vertices, and average pv to edge locations
! ( this computes pv_vertex at all vertices bounding real cells and distance-1 ghost cells )
!
- if (trim(config_time_integration) == 'RK4') then
- ! for RK4, PV is really PV = (eta+f)/h
- fCoef = 1
- elseif (trim(config_time_integration) == 'split_explicit' &
- .or.trim(config_time_integration) == 'unsplit_explicit') then
- ! for split explicit, PV is eta/h because f is added separately to the momentum forcing.
-! mrp temp, new should be:
- fCoef = 0
-! old, for testing:
-! fCoef = 1
- end if
-
do iVertex = 1,nVertices
+ invAreaTri1 = 1.0 / areaTriangle(iVertex)
do k=1,maxLevelVertexBot(iVertex)
h_vertex = 0.0
do i=1,vertexDegree
h_vertex = h_vertex + h(k,cellsOnVertex(i,iVertex)) * kiteAreasOnVertex(i,iVertex)
end do
- h_vertex = h_vertex / areaTriangle(iVertex)
+ h_vertex = h_vertex * invAreaTri1
pv_vertex(k,iVertex) = (fCoef*fVertex(iVertex) + vorticity(k,iVertex)) / h_vertex
end do
@@ -745,8 +739,8 @@
enddo
enddo
- gradPVn(:,:) = 0.0
- gradPVt(:,:) = 0.0
+! gradPVn(:,:) = 0.0
+! gradPVt(:,:) = 0.0
do iEdge = 1,nEdges
cell1 = cellsOnEdge(1, iEdge)
cell2 = cellsOnEdge(2, iEdge)
@@ -1036,7 +1030,71 @@
end subroutine ocn_fuperp!}}}
!***********************************************************************
+!
+! routine ocn_tendency_init
+!
+!> \brief Initializes flags used within tendency routines.
+!> \author Doug Jacobsen
+!> \date 4 November 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine initializes flags related to quantities computed within
+!> other tendency routines.
+!
+!-----------------------------------------------------------------------
+ subroutine ocn_tendency_init(err)!{{{
+ integer, intent(out) :: err
+
+ err = 0
+
+ coef_3rd_order = 0.
+
+ if (config_thickness_adv_order == 2) then
+ hadv2nd = 1
+ hadv3rd = 0
+ hadv4th = 0
+ else if (config_thickness_adv_order == 3) then
+ hadv2nd = 0
+ hadv3rd = 1
+ hadv4th = 0
+
+ if(config_monotonic) then
+ coef_3rd_order = 0.25
+ else
+ coef_3rd_order = 1.0
+ endif
+ else if (config_thickness_adv_order == 4) then
+ hadv2nd = 0
+ hadv3rd = 0
+ hadv4th = 1
+ end if
+
+
+ if(config_include_KE_vertex) then
+ ke_vertex_flag = 1
+ ke_cell_flag = 0
+ else
+ ke_vertex_flag = 0
+ ke_cell_flag = 1
+ endif
+
+ if (trim(config_time_integration) == 'RK4') then
+ ! for RK4, PV is really PV = (eta+f)/h
+ fCoef = 1
+ elseif (trim(config_time_integration) == 'split_explicit' &
+ .or.trim(config_time_integration) == 'unsplit_explicit') then
+ ! for split explicit, PV is eta/h because f is added separately to the momentum forcing.
+ ! mrp temp, new should be:
+ fCoef = 0
+ ! old, for testing:
+ ! fCoef = 1
+ end if
+
+ end subroutine ocn_tendency_init!}}}
+
+!***********************************************************************
+
end module ocn_tendency
!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
</font>
</pre>