<p><b>ringler@lanl.gov</b> 2010-11-27 09:47:30 -0700 (Sat, 27 Nov 2010)</p><p><br>
T and S restoring occurs after final RK4 stage.<br>
<br>
added option to scale del2 and del4 viscosity with grid spacing<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/workspace_tdr/src/core_ocean/Registry
===================================================================
--- branches/ocean_projects/workspace_tdr/src/core_ocean/Registry 2010-11-27 03:57:39 UTC (rev 632)
+++ branches/ocean_projects/workspace_tdr/src/core_ocean/Registry 2010-11-27 16:47:30 UTC (rev 633)
@@ -22,6 +22,7 @@
namelist real hmix config_h_tracer_eddy_diff2 0.0
namelist real hmix config_h_tracer_eddy_diff4 0.0
namelist real hmix config_apvm_upwinding 0.5
+namelist logical hmix config_scaleVisc false
namelist character vmix config_vert_visc_type const
namelist character vmix config_vert_diff_type const
namelist real vmix config_vert_viscosity 2.5e-4
@@ -137,8 +138,9 @@
var persistent integer boundaryVertex ( nVertLevels nVertices ) 0 iro boundaryVertex mesh - -
var persistent integer boundaryCell ( nVertLevels nCells ) 0 iro boundaryCell mesh - -
var persistent real u_src ( nVertLevels nEdges ) 0 iro u_src mesh - -
-var persistent real temperatureRestore ( nCells ) 0 iro temperatureRestore mesh - -
-var persistent real salinityRestore ( nCells ) 0 iro salinityRestore mesh - -
+var persistent real temperatureRestore ( nCells ) 0 i-o temperatureRestore mesh - -
+var persistent real salinityRestore ( nCells ) 0 i-o salinityRestore mesh - -
+var persistent real meshSpacing ( nCells ) 0 i-o meshSpacing mesh - -
# Prognostic variables: read from input, saved in restart, and written to output
var persistent real u ( nVertLevels nEdges Time ) 2 iro u state - -
@@ -146,8 +148,8 @@
var persistent real rho ( nVertLevels nCells Time ) 2 iro rho state - -
var persistent real temperature ( nVertLevels nCells Time ) 2 iro temperature state tracers dynamics
var persistent real salinity ( nVertLevels nCells Time ) 2 iro salinity state tracers dynamics
-var persistent real tracer1 ( nVertLevels nCells Time ) 2 ir- tracer1 state tracers testing
-var persistent real tracer2 ( nVertLevels nCells Time ) 2 ir- tracer2 state tracers testing
+var persistent real tracer1 ( nVertLevels nCells Time ) 2 iro tracer1 state tracers testing
+var persistent real tracer2 ( nVertLevels nCells Time ) 2 iro tracer2 state tracers testing
# Tendency variables: neither read nor written to any files
var persistent real tend_u ( nVertLevels nEdges Time ) 1 - u tend - -
Modified: branches/ocean_projects/workspace_tdr/src/core_ocean/module_time_integration.F
===================================================================
--- branches/ocean_projects/workspace_tdr/src/core_ocean/module_time_integration.F 2010-11-27 03:57:39 UTC (rev 632)
+++ branches/ocean_projects/workspace_tdr/src/core_ocean/module_time_integration.F 2010-11-27 16:47:30 UTC (rev 633)
@@ -71,7 +71,7 @@
type (block_type), pointer :: block
type (state_type) :: provis
- integer :: rk_step
+ integer :: rk_step, index_temperature, index_salinity
real (kind=RKIND), dimension(4) :: rk_weights, rk_substep_weights
@@ -230,6 +230,20 @@
block % state % time_levs(2) % state % tracers % array(:,k,iCell) &
/ block % state % time_levs(2) % state % h % array(k,iCell)
end do
+
+ if(config_surface_restore) then
+ index_temperature = block % state % time_levs(2) % state % index_temperature
+ index_salinity = block % state % time_levs(2) % state % index_salinity
+ block % state % time_levs(2) % state % tracers % array(index_temperature,1,iCell) = &
+ block % state % time_levs(2) % state % tracers % array(index_temperature,1,iCell) - &
+ (block % state % time_levs(2) % state % tracers % array(index_temperature,1,iCell) - &
+ block % mesh % temperatureRestore % array(iCell)) / (config_restore_timescale * 86400.0) * dt
+ block % state % time_levs(2) % state % tracers % array(index_salinity,1,iCell) = &
+ block % state % time_levs(2) % state % tracers % array(index_salinity,1,iCell) - &
+ (block % state % time_levs(2) % state % tracers % array(index_salinity,1,iCell) - &
+ block % mesh % salinityRestore % array(iCell)) / (config_restore_timescale * 86400.0) * dt
+ endif
+
end do
if (config_test_case == 1) then ! For case 1, wind field should be fixed
@@ -272,7 +286,7 @@
upstream_bias, wTopEdge, rho0Inv, r
real (kind=RKIND), dimension(:), pointer :: &
h_s, fVertex, fEdge, dvEdge, dcEdge, areaCell, areaTriangle, &
- zMidZLevel, zTopZLevel
+ zMidZLevel, zTopZLevel, meshSpacing
real (kind=RKIND), dimension(:,:), pointer :: &
weightsOnEdge, kiteAreasOnVertex, h_edge, h, u, v, pressure, &
tend_h, tend_u, circulation, vorticity, ke, ke_edge, pv_edge, &
@@ -329,8 +343,9 @@
zMidZLevel => grid % zMidZLevel % array
zTopZLevel => grid % zTopZLevel % array
maxLevelCell => grid % maxLevelCell % array
- maxLevelEdgeTop => grid % maxLevelEdgeTop % array
- maxLevelVertexBot => grid % maxLevelVertexBot % array
+ maxLevelEdgeTop => grid % maxLevelEdgeTop % array
+ maxLevelVertexBot => grid % maxLevelVertexBot % array
+ meshSpacing => grid % meshSpacing % array
tend_h => tend % h % array
tend_u => tend % u % array
@@ -466,10 +481,13 @@
u_diffusion = ( divergence(k,cell2) - divergence(k,cell1) ) / dcEdge(iEdge) &
-( vorticity(k,vertex2) - vorticity(k,vertex1) ) / dvEdge(iEdge)
- u_diffusion = config_h_mom_eddy_visc2 * u_diffusion
- tend_u(k,iEdge) = tend_u(k,iEdge) + u_diffusion
+ if(config_scaleVisc) then
+ u_diffusion = u_diffusion*0.5*(meshSpacing(cell1)+meshSpacing(cell2))**2
+ endif
+ tend_u(k,iEdge) = tend_u(k,iEdge) + config_h_mom_eddy_visc2 * u_diffusion
+
end do
end do
end if
@@ -558,6 +576,10 @@
-( delsq_vorticity(k,vertex2) &
- delsq_vorticity(k,vertex1) ) / dvEdge(iEdge)
+ if(config_scaleVisc) then
+ u_diffusion = u_diffusion*0.5*(meshSpacing(cell1)+meshSpacing(cell2))**3
+ endif
+
tend_u(k,iEdge) = tend_u(k,iEdge) - config_h_mom_eddy_visc4 * u_diffusion
end do
end do
@@ -707,7 +729,7 @@
integer, dimension(:), pointer :: nEdgesOnCell, nEdgesOnEdge, &
maxLevelCell, maxLevelEdgeTop, maxLevelVertexBot
integer, dimension(:,:), pointer :: cellsOnEdge, boundaryCell
- real (kind=RKIND), dimension(:), pointer :: zTopZLevel, temperatureRestore, salinityRestore
+ real (kind=RKIND), dimension(:), pointer :: zTopZLevel
real (kind=RKIND), dimension(:,:), allocatable:: fluxVertTop, tracerTop, boundaryMask
real (kind=RKIND), dimension(:,:,:), allocatable::tr_flux, tr_div, delsq_tracer
@@ -744,28 +766,12 @@
deriv_two => grid % deriv_two % array
- temperatureRestore => grid % temperatureRestore % array
- salinityRestore => grid % salinityRestore % array
-
!
! initialize tracer tendency (RHS of tracer equation) to zero.
!
tend_tr(:,:,:) = 0.0
!
- ! include surface restore on temperature and salinity
- !
- if(config_surface_restore) then
- r = 1.0 / (config_restore_timescale * 86400.0)
- do iCell=1,nCells
- tend_tr(s % index_temperature,1,iCell) = tend_tr(s % index_temperature,1,iCell) - &
- h(1,iCell)*(tracers(s % index_temperature,1,iCell) - temperatureRestore(iCell))*r
- tend_tr(s % index_salinity,1,iCell) = tend_tr(s % index_salinity,1,iCell) - &
- h(1,iCell)*(tracers(s % index_salinity,1,iCell) - salinityRestore(iCell))*r
- enddo
- endif
-
- !
! tracer tendency: horizontal advection term -div( h \phi u)
!
! mrp 101115 note: in order to include flux boundary conditions, we will need to
</font>
</pre>