<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 => mesh % meshDensity % array
+ meshScalingDel2 => mesh % meshScalingDel2 % array
+ meshScalingDel4 => 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, &
+ meshScalingDel2, meshScalingDel4
real (kind=RKIND), dimension(:,:), pointer :: vh, weightsOnEdge, kiteAreasOnVertex, h_edge, h, u, v, tend_h, tend_u, &
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 => s % h % array
u => s % u % array
@@ -314,23 +314,9 @@
u_src => grid % u_src % array
- meshDensity => grid % meshDensity % array
+ meshScalingDel2 => grid % meshScalingDel2 % array
+ meshScalingDel4 => 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>