<p><b>duda</b> 2011-05-19 15:41:14 -0600 (Thu, 19 May 2011)</p><p>Compute meshScalingDel2 and meshScalingDel4 once during core initialization <br>
through a call to a new routine, compute_mesh_scaling(), from mpas_init_block(). <br>
Also, add meshScalingDel2 and meshScalingDel4 to the Registry.<br>
<br>
<br>
M    src/core_sw/module_mpas_core.F<br>
M    src/core_sw/Registry<br>
M    src/core_sw/module_time_integration.F<br>
</p><hr noshade><pre><font color="gray">Modified: trunk/mpas/src/core_sw/Registry
===================================================================
--- trunk/mpas/src/core_sw/Registry        2011-05-19 17:01:39 UTC (rev 846)
+++ trunk/mpas/src/core_sw/Registry        2011-05-19 21:41:14 UTC (rev 847)
@@ -53,7 +53,6 @@
 var persistent real    xCell ( nCells ) 0 iro xCell mesh - -
 var persistent real    yCell ( nCells ) 0 iro yCell mesh - -
 var persistent real    zCell ( nCells ) 0 iro zCell mesh - -
-var persistent real    meshDensity ( nCells ) 0 iro meshDensity mesh - -
 var persistent integer indexToCellID ( nCells ) 0 iro indexToCellID mesh - -
 
 var persistent real    latEdge ( nEdges ) 0 iro latEdge mesh - -
@@ -70,6 +69,10 @@
 var persistent real    zVertex ( nVertices ) 0 iro zVertex mesh - -
 var persistent integer indexToVertexID ( nVertices ) 0 iro indexToVertexID mesh - -
 
+var persistent real    meshDensity ( nCells ) 0 iro meshDensity mesh - -
+var persistent real    meshScalingDel2 ( nEdges ) 0 ro meshScalingDel2 mesh - -
+var persistent real    meshScalingDel4 ( nEdges ) 0 ro meshScalingDel4 mesh - -
+
 var persistent integer cellsOnEdge ( TWO nEdges ) 0 iro cellsOnEdge mesh - -
 var persistent integer nEdgesOnCell ( nCells ) 0 iro nEdgesOnCell mesh - -
 var persistent integer nEdgesOnEdge ( nEdges ) 0 iro nEdgesOnEdge mesh - -

Modified: trunk/mpas/src/core_sw/module_mpas_core.F
===================================================================
--- trunk/mpas/src/core_sw/module_mpas_core.F        2011-05-19 17:01:39 UTC (rev 846)
+++ trunk/mpas/src/core_sw/module_mpas_core.F        2011-05-19 21:41:14 UTC (rev 847)
@@ -55,7 +55,8 @@
    
    
       call compute_solve_diagnostics(dt, block % state % time_levs(1) % state, mesh)
-   
+      call compute_mesh_scaling(mesh) 
+
       call rbfInterp_initialize(mesh)
       call init_reconstruct(mesh)
       call reconstruct(block % state % time_levs(1) % state, mesh)
@@ -214,4 +215,36 @@
    
    end subroutine mpas_core_finalize
 
+
+   subroutine compute_mesh_scaling(mesh)
+
+      use grid_types
+
+      implicit none
+
+      type (mesh_type), intent(inout) :: mesh
+
+      integer :: iEdge, cell1, cell2
+      real (kind=RKIND), dimension(:), pointer :: meshDensity, meshScalingDel2, meshScalingDel4
+
+      meshDensity =&gt; mesh % meshDensity % array
+      meshScalingDel2 =&gt; mesh % meshScalingDel2 % array
+      meshScalingDel4 =&gt; mesh % meshScalingDel4 % array
+
+      !
+      ! Compute the scaling factors to be used in the del2 and del4 dissipation
+      !
+      meshScalingDel2(:) = 1.0
+      meshScalingDel4(:) = 1.0
+      if (config_h_ScaleWithMesh) then
+         do iEdge=1,mesh%nEdges
+            cell1 = mesh % cellsOnEdge % array(1,iEdge)
+            cell2 = mesh % cellsOnEdge % array(2,iEdge)
+            meshScalingDel2(iEdge) = 1.0 / ( (meshDensity(cell1) + meshDensity(cell2) )/2.0)**(5.0/12.0)
+            meshScalingDel4(iEdge) = 1.0 / ( (meshDensity(cell1) + meshDensity(cell2) )/2.0)**(5.0/6.0)
+         end do
+      end if
+
+   end subroutine compute_mesh_scaling
+
 end module mpas_core

Modified: trunk/mpas/src/core_sw/module_time_integration.F
===================================================================
--- trunk/mpas/src/core_sw/module_time_integration.F        2011-05-19 17:01:39 UTC (rev 846)
+++ trunk/mpas/src/core_sw/module_time_integration.F        2011-05-19 21:41:14 UTC (rev 847)
@@ -258,7 +258,8 @@
       real (kind=RKIND) :: flux, vorticity_abs, workpv, q, upstream_bias
 
       integer :: nCells, nEdges, nVertices, nVertLevels
-      real (kind=RKIND), dimension(:), pointer :: h_s, fVertex, fEdge, dvEdge, dcEdge, areaCell, areaTriangle, meshDensity
+      real (kind=RKIND), dimension(:), pointer :: h_s, fVertex, fEdge, dvEdge, dcEdge, areaCell, areaTriangle, &amp;
+                                                  meshScalingDel2, meshScalingDel4
       real (kind=RKIND), dimension(:,:), pointer :: vh, weightsOnEdge, kiteAreasOnVertex, h_edge, h, u, v, tend_h, tend_u, &amp;
                                                     circulation, vorticity, ke, pv_edge, divergence, h_vertex
       integer, dimension(:,:), pointer :: cellsOnEdge, cellsOnVertex, verticesOnEdge, edgesOnCell, edgesOnEdge, edgesOnVertex
@@ -273,7 +274,6 @@
       real (kind=RKIND), parameter :: rho_ref = 1000.0
       real (kind=RKIND) :: ke_edge
 
-      real (kind=RKIND), allocatable, dimension(:) :: meshScalingDel2, meshScalingDel4
 
       h           =&gt; s % h % array
       u           =&gt; s % u % array
@@ -314,23 +314,9 @@
 
       u_src =&gt; grid % u_src % array
 
-      meshDensity =&gt; grid % meshDensity % array
+      meshScalingDel2 =&gt; grid % meshScalingDel2 % array
+      meshScalingDel4 =&gt; grid % meshScalingDel4 % array
 
-      !
-      ! Compute the scaling factors to be used in the del2 and del4 dissipation
-      !
-      allocate(meshScalingDel2(nEdges))
-      allocate(meshScalingDel4(nEdges))
-      meshScalingDel2(:) = 1.0
-      meshScalingDel4(:) = 1.0
-      if(config_h_ScaleWithMesh) then
-      do iEdge=1,nEdges
-         cell1 = cellsOnEdge(1,iEdge)
-         cell2 = cellsOnEdge(2,iEdge)
-         meshScalingDel2(iEdge) = 1.0 / ( (meshDensity(cell1) + meshDensity(cell2) )/2.0)**(5.0/12.0)
-         meshScalingDel4(iEdge) = 1.0 / ( (meshDensity(cell1) + meshDensity(cell2) )/2.0)**(5.0/6.0)
-      enddo
-      endif
 
       !
       ! Compute height tendency for each cell
@@ -542,9 +528,6 @@
          end do
      endif
  
-     deallocate(meshScalingDel2)
-     deallocate(meshScalingDel4)
-
    end subroutine compute_tend
 
 

</font>
</pre>