<p><b>jgraham@lanl.gov</b> 2012-11-26 15:18:07 -0700 (Mon, 26 Nov 2012)</p><p>ready for first attempt, everything is coded<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/turbulent_transfer/src/core_ocean/mpas_ocn_tendency.F
===================================================================
--- branches/ocean_projects/turbulent_transfer/src/core_ocean/mpas_ocn_tendency.F        2012-11-26 21:02:18 UTC (rev 2319)
+++ branches/ocean_projects/turbulent_transfer/src/core_ocean/mpas_ocn_tendency.F        2012-11-26 22:18:07 UTC (rev 2320)
@@ -271,16 +271,30 @@
type (diagnostics_type), intent(in) :: d !< Input: Diagnostic information
type (mesh_type), intent(in) :: grid !< Input: Grid information
+
real (kind=RKIND), dimension(:,:), pointer :: &
h_edge, h, u, rho, zMid, pressure, &
circulation, vorticity, viscosity, ke, ke_edge, Vor_edge, &
MontPot, wTop, divergence, vertViscTopOfEdge, tend_u_cor, &
- tend_u_press, tend_u_hmix, tend_u_force, tend_u_vadv, tend_u_vmix
+ tend_u_press, tend_u_hmix, tend_u_force, tend_u_vadv, tend_u_vmix, &
+ u_weighted
real (kind=RKIND), dimension(:,:), pointer :: u_src
+ real (kind=RKIND) :: invAreaCell1
+
+ real (kind=RKIND), dimension(:), pointer :: &
+ dvEdge, dcEdge, areaCell
+
integer :: err
+ integer :: iEdge, iCell, k, i
+ integer :: boundaryMask, velMask, nCells, nEdges, nVertices, nVertLevels, vertexDegree, err
+
+ integer, dimension(:), pointer :: nEdgesOnCell, maxLevelCell
+ integer, dimension(:,:), pointer :: edgesOnCell
+
+
call mpas_timer_start("ocn_tend_u_diagnostic")
u => s % u % array
@@ -298,12 +312,23 @@
pressure => s % pressure % array
vertViscTopOfEdge => d % vertViscTopOfEdge % array
+
+ u_weighted => s % uWeighted % array
tend_u_cor => s % uCoriolisTend % array
tend_u_press => s % uPressureGradTend % array
tend_u_hmix => s % uHMixTend % array
tend_u_force => s % uForcingTend % array
tend_u_vadv => s % uVAdvTend % array
tend_u_vmix => s % uVMixTend % array
+
+ nEdgesOnCell => grid % nEdgesOnCell % array
+ edgesOnCell => grid % edgesOnCell % array
+ dcEdge => grid % dcEdge % array
+ dvEdge => grid % dvEdge % array
+ areaCell => grid % areaCell % array
+ maxLevelCell => grid % maxLevelCell % array
+ nCells = grid % nCells
+
u_src => grid % u_src % array
@@ -311,6 +336,7 @@
! velocity tendency terms
!
! jpg 121101 efficiency: I'm guessing these memory addresses need to be zeroed each call, but I'm not sure how time state 2 is setup
+ u_weighted(:,:) = 0.0
tend_u_cor(:,:) = 0.0
tend_u_press(:,:) = 0.0
tend_u_hmix(:,:) = 0.0
@@ -318,6 +344,24 @@
tend_u_vadv(:,:) = 0.0
tend_u_vmix(:,:) = 0.0
+
+ ! jpg 121126 attempt to follow ke weighting given in Ringler, Thuburn
+ ! et al to give 2 * "balanced" velocity
+ ! at edges so that spectral conservation will follow actual energy
+ ! conservation in MPAS
+
+ do iCell = 1, nCells
+ invAreaCell1 = 1.0 / areaCell(iCell)
+ do i = 1, nEdgesOnCell(iCell)
+ iEdge = edgesOnCell(i, iCell)
+ do k = 1, maxLevelCell(iCell)
+ u_weighted(k, iEdge) = u_weighted(k, iEdge) + 2.0 * 0.25 * &
+ dvEdge(iEdge) * dcEdge(iEdge) * u(k, iEdge) * invAreaCell1
+ end do
+ end do
+ end do
+
+
!
! velocity tendency: nonlinear Coriolis term and grad of kinetic energy
!
</font>
</pre>