<p><b>mpetersen@lanl.gov</b> 2010-08-31 16:04:18 -0600 (Tue, 31 Aug 2010)</p><p>Adding some clarifying comments and changing a few namelist items on dissipation branch<br>
</p><hr noshade><pre><font color="gray">Modified: branches/dissipation/src/core_ocean/Registry
===================================================================
--- branches/dissipation/src/core_ocean/Registry 2010-08-31 01:54:38 UTC (rev 484)
+++ branches/dissipation/src/core_ocean/Registry 2010-08-31 22:04:18 UTC (rev 485)
@@ -7,7 +7,6 @@
namelist integer sw_model config_ntimesteps 7500
namelist integer sw_model config_output_interval 500
namelist integer sw_model config_stats_interval 100
-namelist logical sw_model config_eden false
namelist character io config_input_name grid.nc
namelist character io config_output_name output.nc
namelist character io config_restart_name restart.nc
@@ -18,8 +17,8 @@
namelist real grid config_rho0 1028
namelist real hmix config_h_mom_eddy_visc2 0.0
namelist real hmix config_h_mom_eddy_visc4 0.0
-namelist real hmix config_h_tracer_eddy_visc2 0.0
-namelist real hmix config_h_tracer_eddy_visc4 0.0
+namelist real hmix config_h_tracer_eddy_diff2 0.0
+namelist real hmix config_h_tracer_eddy_diff4 0.0
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
Modified: branches/dissipation/src/core_ocean/module_time_integration.F
===================================================================
--- branches/dissipation/src/core_ocean/module_time_integration.F 2010-08-31 01:54:38 UTC (rev 484)
+++ branches/dissipation/src/core_ocean/module_time_integration.F 2010-08-31 22:04:18 UTC (rev 485)
@@ -388,7 +388,6 @@
end do
endif ! coordinate type
-
!
! velocity tendency: start accumulating tendency terms
!
@@ -448,7 +447,9 @@
endif
!
- ! velocity tendency: del2 and/or del4 dissipation
+ ! velocity tendency: del2 dissipation, </font>
<font color="black">u_2 </font>
<font color="blue">abla^2 u
+ ! computed as </font>
<font color="black">u( </font>
<font color="black">abla divergence + k \times </font>
<font color="gray">abla vorticity )
+ ! strictly only valid for h_mom_eddy_visc2 == constant
!
if ( h_mom_eddy_visc2 > 0.0 ) then
do iEdge=1,grid % nEdgesSolve
@@ -459,12 +460,12 @@
do k=1,nVertLevels
- !
- ! Compute diffusion, computed as </font>
<font color="black">abla divergence - k \times </font>
<font color="red">abla vorticity
- ! strictly only valid for h_mom_eddy_visc2 == constant
- !
- u_diffusion = ( divergence(k,cell2) - divergence(k,cell1) ) / dcEdge(iEdge) &
- -( vorticity(k,vertex2) - vorticity(k,vertex1) ) / dvEdge(iEdge)
+ ! Here -( vorticity(k,vertex2) - vorticity(k,vertex1) ) / dvEdge(iEdge)
+ ! is - </font>
<font color="blue">abla vorticity pointing from vertex 2 to vertex 1, or equivalently
+ ! + k \times </font>
<font color="gray">abla vorticity pointing from cell1 to cell2.
+
+ u_diffusion = ( divergence(k,cell2) - divergence(k,cell1) ) / dcEdge(iEdge) &
+ -( vorticity(k,vertex2) - vorticity(k,vertex1) ) / dvEdge(iEdge)
u_diffusion = h_mom_eddy_visc2 * u_diffusion
tend_u(k,iEdge) = tend_u(k,iEdge) + u_diffusion
@@ -473,7 +474,12 @@
end do
end if
-
+ !
+ ! velocity tendency: del4 dissipation, -</font>
<font color="black">u_4 </font>
<font color="blue">abla^4 u
+ ! computed as </font>
<font color="black">abla^2 u = </font>
<font color="black">abla divergence + k \times </font>
<font color="gray">abla vorticity
+ ! applied recursively.
+ ! strictly only valid for h_mom_eddy_visc4 == constant
+ !
if ( h_mom_eddy_visc4 > 0.0 ) then
allocate(delsq_divergence(nVertLevels, nCells+1))
@@ -483,6 +489,7 @@
delsq_u(:,:) = 0.0
+ ! Compute </font>
<font color="black">abla^2 u = </font>
<font color="black">abla divergence + k \times </font>
<font color="gray">abla vorticity
do iEdge=1,grid % nEdges
cell1 = cellsOnEdge(1,iEdge)
cell2 = cellsOnEdge(2,iEdge)
@@ -490,26 +497,25 @@
vertex2 = verticesOnEdge(2,iEdge)
do k=1,nVertLevels
-
- !
- ! Compute diffusion, computed as </font>
<font color="black">abla divergence - k \times </font>
<font color="red">abla vorticity
- ! strictly only valid for h_mom_eddy_visc4 == constant
- !
- u_diffusion = ( divergence(k,cell2) - divergence(k,cell1) ) / dcEdge(iEdge) &
- -( vorticity(k,vertex2) - vorticity(k,vertex1) ) / dvEdge(iEdge)
+
+ u_diffusion = ( divergence(k,cell2) - divergence(k,cell1) ) / dcEdge(iEdge) &
+ -( vorticity(k,vertex2) - vorticity(k,vertex1)) / dvEdge(iEdge)
delsq_u(k,iEdge) = delsq_u(k,iEdge) + u_diffusion
end do
end do
+ ! vorticity using </font>
<font color="gray">abla^2 u
delsq_circulation(:,:) = 0.0
do iEdge=1,nEdges
vertex1 = verticesOnEdge(1,iEdge)
vertex2 = verticesOnEdge(2,iEdge)
do k=1,nVertLevels
- delsq_circulation(k,vertex1) = delsq_circulation(k,vertex1) - dcEdge(iEdge) * delsq_u(k,iEdge)
- delsq_circulation(k,vertex2) = delsq_circulation(k,vertex2) + dcEdge(iEdge) * delsq_u(k,iEdge)
+ delsq_circulation(k,vertex1) = delsq_circulation(k,vertex1) &
+ - dcEdge(iEdge) * delsq_u(k,iEdge)
+ delsq_circulation(k,vertex2) = delsq_circulation(k,vertex2) &
+ + dcEdge(iEdge) * delsq_u(k,iEdge)
end do
end do
do iVertex=1,nVertices
@@ -519,13 +525,16 @@
end do
end do
+ ! Divergence using </font>
<font color="gray">abla^2 u
delsq_divergence(:,:) = 0.0
do iEdge=1,nEdges
cell1 = cellsOnEdge(1,iEdge)
cell2 = cellsOnEdge(2,iEdge)
do k=1,nVertLevels
- delsq_divergence(k,cell1) = delsq_divergence(k,cell1) + delsq_u(k,iEdge)*dvEdge(iEdge)
- delsq_divergence(k,cell2) = delsq_divergence(k,cell2) - delsq_u(k,iEdge)*dvEdge(iEdge)
+ delsq_divergence(k,cell1) = delsq_divergence(k,cell1) &
+ + delsq_u(k,iEdge)*dvEdge(iEdge)
+ delsq_divergence(k,cell2) = delsq_divergence(k,cell2) &
+ - delsq_u(k,iEdge)*dvEdge(iEdge)
end do
end do
do iCell = 1,nCells
@@ -535,6 +544,8 @@
end do
end do
+ ! Compute - \kappa </font>
<font color="blue">abla^4 u
+ ! as </font>
<font color="black">abla div(</font>
<font color="black">abla^2 u) + k \times </font>
<font color="black">abla ( k \cross curl(</font>
<font color="gray">abla^2 u) )
do iEdge=1,grid % nEdgesSolve
cell1 = cellsOnEdge(1,iEdge)
cell2 = cellsOnEdge(2,iEdge)
@@ -543,12 +554,10 @@
do k=1,nVertLevels
- !
- ! Compute diffusion, computed as </font>
<font color="black">abla divergence - k \times </font>
<font color="gray">abla vorticity
- ! only valid for h_mom_eddy_visc4 == constant
- !
- u_diffusion = ( delsq_divergence(k,cell2) - delsq_divergence(k,cell1) ) / dcEdge(iEdge) &
- -( delsq_vorticity(k,vertex2) - delsq_vorticity(k,vertex1) ) / dvEdge(iEdge)
+ u_diffusion = ( delsq_divergence(k,cell2) &
+ - delsq_divergence(k,cell1) ) / dcEdge(iEdge) &
+ -( delsq_vorticity(k,vertex2) &
+ - delsq_vorticity(k,vertex1) ) / dvEdge(iEdge)
tend_u(k,iEdge) = tend_u(k,iEdge) - h_mom_eddy_visc4 * u_diffusion
end do
@@ -667,7 +676,7 @@
integer :: iCell, iEdge, k, iTracer, cell1, cell2, upwindCell,&
nEdges, nCells, nVertLevels
- real (kind=RKIND) :: h_tracer_eddy_visc2, h_tracer_eddy_visc4, invAreaCell1, invAreaCell2, tracer_turb_flux
+ real (kind=RKIND) :: h_tracer_eddy_diff2, h_tracer_eddy_diff4, invAreaCell1, invAreaCell2, tracer_turb_flux
real (kind=RKIND) :: flux, tracer_edge, r
real (kind=RKIND) :: dist
real (kind=RKIND), dimension(:), pointer :: &
@@ -710,8 +719,8 @@
nVertLevels = grid % nVertLevels
- h_tracer_eddy_visc2 = config_h_tracer_eddy_visc2
- h_tracer_eddy_visc4 = config_h_tracer_eddy_visc4
+ h_tracer_eddy_diff2 = config_h_tracer_eddy_diff2
+ h_tracer_eddy_diff4 = config_h_tracer_eddy_diff4
!
! tracer tendency: horizontal advection term -div( h \phi u)
@@ -809,12 +818,10 @@
deallocate(tracerTop)
!
- ! tracer tendency: horizontal tracer diffusion, del2 and/or del4
+ ! tracer tendency: del2 horizontal tracer diffusion, div(h \kappa_2 </font>
<font color="gray">abla \phi)
!
+ if ( h_tracer_eddy_diff2 > 0.0 ) then
- if ( h_tracer_eddy_visc2 > 0.0 ) then
-
-
!
! compute a boundary mask to enforce insulating boundary conditions in the horizontal
!
@@ -822,7 +829,6 @@
boundaryMask = 1.0
where(boundaryEdge.eq.1) boundaryMask=0.0
-
do iEdge=1,grid % nEdges
cell1 = grid % cellsOnEdge % array(1,iEdge)
cell2 = grid % cellsOnEdge % array(2,iEdge)
@@ -831,8 +837,14 @@
do k=1,grid % nVertLevels
do iTracer=1,num_tracers
- tracer_turb_flux = h_tracer_eddy_visc2*(tracers(iTracer,k,cell2) - tracers(iTracer,k,cell1))/dcEdge(iEdge)
- flux = dvEdge (iEdge) * h_edge(k,iEdge) * tracer_turb_flux * boundaryMask(k, iEdge)
+ ! \kappa_2 </font>
<font color="blue">abla \phi on edge
+ tracer_turb_flux = h_tracer_eddy_diff2 &
+ *( tracers(iTracer,k,cell2) &
+ - tracers(iTracer,k,cell1))/dcEdge(iEdge)
+
+ ! div(h \kappa_2 </font>
<font color="gray">abla \phi) at cell center
+ flux = dvEdge (iEdge) * h_edge(k,iEdge) &
+ * tracer_turb_flux * boundaryMask(k, iEdge)
tend_tr(iTracer,k,cell1) = tend_tr(iTracer,k,cell1) + flux * invAreaCell1
tend_tr(iTracer,k,cell2) = tend_tr(iTracer,k,cell2) - flux * invAreaCell2
end do
@@ -844,7 +856,11 @@
end if
- if ( h_tracer_eddy_visc4 > 0.0 ) then
+ !
+ ! tracer tendency: del4 horizontal tracer diffusion, &
+ ! div(h \kappa_4 </font>
<font color="black">abla [div(h </font>
<font color="gray">abla \phi)])
+ !
+ if ( h_tracer_eddy_diff4 > 0.0 ) then
!
! compute a boundary mask to enforce insulating boundary conditions in the horizontal
@@ -857,6 +873,7 @@
delsq_tracer(:,:,:) = 0.
+ ! first del2: div(h </font>
<font color="gray">abla \phi) at cell center
do iEdge=1,grid % nEdges
cell1 = grid % cellsOnEdge % array(1,iEdge)
cell2 = grid % cellsOnEdge % array(2,iEdge)
@@ -864,9 +881,13 @@
do k=1,grid % nVertLevels
do iTracer=1,num_tracers
delsq_tracer(iTracer,k,cell1) = delsq_tracer(iTracer,k,cell1) &
- + dvEdge(iEdge)*h_edge(k,iEdge)*(tracers(iTracer,k,cell2) - tracers(iTracer,k,cell1))/dcEdge(iEdge) * boundaryMask(k,iEdge)
+ + dvEdge(iEdge)*h_edge(k,iEdge) &
+ *(tracers(iTracer,k,cell2) - tracers(iTracer,k,cell1)) &
+ /dcEdge(iEdge) * boundaryMask(k,iEdge)
delsq_tracer(iTracer,k,cell2) = delsq_tracer(iTracer,k,cell2) &
- - dvEdge(iEdge)*h_edge(k,iEdge)*(tracers(iTracer,k,cell2) - tracers(iTracer,k,cell1))/dcEdge(iEdge) * boundaryMask(k,iEdge)
+ - dvEdge(iEdge)*h_edge(k,iEdge) &
+ *(tracers(iTracer,k,cell2) - tracers(iTracer,k,cell1)) &
+ /dcEdge(iEdge) * boundaryMask(k,iEdge)
end do
end do
@@ -881,6 +902,7 @@
end do
end do
+ ! second del2: div(h </font>
<font color="gray">abla [delsq_tracer]) at cell center
do iEdge=1,grid % nEdges
cell1 = grid % cellsOnEdge % array(1,iEdge)
cell2 = grid % cellsOnEdge % array(2,iEdge)
@@ -889,11 +911,15 @@
do k=1,grid % nVertLevels
do iTracer=1,num_tracers
- tracer_turb_flux = h_tracer_eddy_visc4*(delsq_tracer(iTracer,k,cell2) - delsq_tracer(iTracer,k,cell1))/dcEdge(iEdge)
+ tracer_turb_flux = h_tracer_eddy_diff4 &
+ *( delsq_tracer(iTracer,k,cell2) &
+ - delsq_tracer(iTracer,k,cell1))/dcEdge(iEdge)
flux = dvEdge (iEdge) * tracer_turb_flux
- tend_tr(iTracer,k,cell1) = tend_tr(iTracer,k,cell1) - flux * invAreaCell1 * boundaryMask(k,iEdge)
- tend_tr(iTracer,k,cell2) = tend_tr(iTracer,k,cell2) + flux * invAreaCell2 * boundaryMask(k,iEdge)
+ tend_tr(iTracer,k,cell1) = tend_tr(iTracer,k,cell1) &
+ - flux * invAreaCell1 * boundaryMask(k,iEdge)
+ tend_tr(iTracer,k,cell2) = tend_tr(iTracer,k,cell2) &
+ + flux * invAreaCell2 * boundaryMask(k,iEdge)
end do
enddo
</font>
</pre>