<p><b>mpetersen@lanl.gov</b> 2010-04-27 16:18:36 -0600 (Tue, 27 Apr 2010)</p><p>Added vertical viscous term, wind forcing just at surface, bottom drag.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/z_level_mrp/mpas/src/core_ocean/Registry
===================================================================
--- branches/ocean_projects/z_level_mrp/mpas/src/core_ocean/Registry 2010-04-27 18:41:17 UTC (rev 215)
+++ branches/ocean_projects/z_level_mrp/mpas/src/core_ocean/Registry 2010-04-27 22:18:36 UTC (rev 216)
@@ -85,8 +85,8 @@
var real coeffs_reconstruct ( R3 maxEdges nCells ) - coeffs_reconstruct - -
# Arrays for z-level version of mpas-ocean
-var integer kmaxCell ( nCells ) iro kmaxCell - -
-var integer kmaxEdge ( nEdges ) iro kmaxEdge - -
+var integer maxLevelsCell ( nCells ) iro kmaxCell - -
+var integer maxLevelsEdge ( nEdges ) iro kmaxEdge - -
var real hZLevel ( nVertLevels ) iro hZLevel - -
var real zmidZLevel ( nVertLevels ) iro zmidZLevel - -
var real zbotZLevel ( nVertLevels ) iro zbotZLevel - -
@@ -120,6 +120,9 @@
var real zmid ( nVertLevels nCells Time ) o zmid - -
var real zbot ( nVertLevels nCells Time ) o zbot - -
var real zSurface ( nCells Time ) o zSurface - -
+var real zmidEdge ( nVertLevels nEdges Time ) o zmidEdge - -
+var real zbotEdge ( nVertLevels nEdges Time ) o zbotEdge - -
+var real zSurfaceEdge ( nEdges Time ) o zSurfaceEdge - -
var real pmid ( nVertLevels nCells Time ) o pmid - -
var real pbot ( nVertLevels nCells Time ) o pbot - -
var real pmidZLevel ( nVertLevels nCells Time ) o pmidZLevel - -
Modified: branches/ocean_projects/z_level_mrp/mpas/src/core_ocean/module_test_cases.F
===================================================================
--- branches/ocean_projects/z_level_mrp/mpas/src/core_ocean/module_test_cases.F 2010-04-27 18:41:17 UTC (rev 215)
+++ branches/ocean_projects/z_level_mrp/mpas/src/core_ocean/module_test_cases.F 2010-04-27 22:18:36 UTC (rev 216)
@@ -111,22 +111,37 @@
u_src = 0.0
elseif (config_test_case == 0 ) then
! for rectangular basin:
- do iEdge=1,nEdges
- u_src(1:nVertLevels, iEdge) = u_src(1, iEdge) / nVertLevels
- end do
+ !do iEdge=1,nEdges
+ ! u_src(1:nVertLevels, iEdge) = u_src(1, iEdge) / nVertLevels
+ !end do
endif
+ ! set rho directly:
+ !delta_rho=0.0
+ !do iLevel = 1,nVertLevels
+ ! rho(iLevel,1) = delta_rho*(iLevel-1)
+ !enddo
+ !rho(:,1) = rho(:,1) - sum(rho(:,1))/nVertLevels + 1000.0
+ !do iLevel = 1,nVertLevels
+ ! rho(iLevel,:) = rho(iLevel,1)
+ !enddo
+
tracers = 0.0
do iCell = 1,nCells
+ ! for three layer test
+ tracers(index_salinity,1,iCell) = 2.0 ! salinity
+ tracers(index_salinity,2,iCell) = 6.0 ! salinity
+ tracers(index_salinity,3,iCell) = 13.0 ! salinity
do iLevel = 1,nVertLevels
- !tracers(1,iLevel,iCell) = 2.0
- !tracers(2,iLevel,iCell) = 1.0
- !tracers(3,iLevel,iCell) = xCell(iCell)
!if (xCell(iCell)<500*1e3) tracers(4,iLevel,iCell) = 1.0
! tracers(index_temperature,iLevel,iCell) = 5.0 ! temperature
! tracers(index_temperature,iLevel,iCell) = 6.0-iLevel ! temperature
- tracers(index_temperature,iLevel,iCell) = 6.0-iLevel+yCell(iCell)/4000.e3 ! temperature
- tracers(index_salinity,iLevel,iCell) = 35.0 ! salinity
+! tracers(index_temperature,iLevel,iCell) = 6.0-iLevel+yCell(iCell)/4000.e3 ! temperature
+! tracers(index_salinity,iLevel,iCell) = 35.0 ! salinity
+
+ ! for three layer test
+ tracers(index_temperature,iLevel,iCell) = 5.0 ! temperature
+
! tracers(index_temperature,iLevel,iCell) = &
! 10.0* xCell(iCell)/2500.e3
! tracers(index_salinity,iLevel,iCell) = &
@@ -134,18 +149,12 @@
tracers(index_tracer1,iLevel,iCell) = 1.0
tracers(index_tracer2,iLevel,iCell) = &
(yCell(iCell)/4000.e3 + xCell(iCell)/2500.e3 )/2.0
+ rho(iLevel,iCell) = 1000.0*( 1.0 &
+ - 2.5e-4*tracers(index_temperature,iLevel,iCell) &
+ + 7.6e-4*tracers(index_salinity,iLevel,iCell))
enddo
enddo
- delta_rho=0.0
- do iLevel = 1,nVertLevels
- rho(iLevel,1) = delta_rho*(iLevel-1)
- enddo
- rho(:,1) = rho(:,1) - sum(rho(:,1))/nVertLevels + 1000.0
- do iLevel = 1,nVertLevels
- rho(iLevel,:) = rho(iLevel,1)
- enddo
-
#ifdef EXPAND_LEVELS
print '(10a)', 'EXPAND_LEVELS compiler flag is on.', &
' Copying h and u from k=1 to other levels.'
Modified: branches/ocean_projects/z_level_mrp/mpas/src/core_ocean/module_time_integration.F
===================================================================
--- branches/ocean_projects/z_level_mrp/mpas/src/core_ocean/module_time_integration.F 2010-04-27 18:41:17 UTC (rev 215)
+++ branches/ocean_projects/z_level_mrp/mpas/src/core_ocean/module_time_integration.F 2010-04-27 22:18:36 UTC (rev 216)
@@ -285,17 +285,18 @@
real (kind=RKIND) :: flux, vorticity_abs, h_vertex, workpv, q, &
upstream_bias, wHat, areaSumInv, rho0Inv
real (kind=RKIND), dimension(:), pointer :: &
- h_s, fVertex, fEdge, dvEdge, dcEdge, areaCell, areaTriangle
+ h_s, fVertex, fEdge, dvEdge, dcEdge, areaCell, areaTriangle, zSurfaceEdge
real (kind=RKIND), dimension(grid % nVertLevels) :: w_dudz
real (kind=RKIND), dimension(:,:), pointer :: &
vh, weightsOnEdge, kiteAreasOnVertex, h_edge, h, u, v, &
tend_h, tend_u, Fv, circulation, vorticity, ke, pv_edge, &
- divergence, MontPot, pmidZLevel
+ divergence, MontPot, pmidZLevel, zMidEdge, zbotEdge
integer, dimension(:), pointer :: nEdgesOnCell, nEdgesOnEdge
integer, dimension(:,:), pointer :: &
cellsOnEdge, cellsOnVertex, verticesOnEdge, edgesOnCell, &
edgesOnEdge, edgesOnVertex
- real (kind=RKIND) :: u_diffusion, visc
+ real (kind=RKIND) :: u_diffusion, visc, dist
+ real (kind=RKIND), dimension(:), allocatable:: fluxVert
!ocean
real (kind=RKIND), dimension(:,:), pointer :: u_src
@@ -318,6 +319,9 @@
!mrp 100112:
MontPot => s % MontPot % array
pmidZLevel => s % pmidZLevel % array
+ zmidEdge => s % zmidEdge % array
+ zbotEdge => s % zbotEdge % array
+ zSurfaceEdge=> s % zSurfaceEdge % array
!mrp 100112 end
weightsOnEdge => grid % weightsOnEdge % array
@@ -441,6 +445,7 @@
!
tend_u(:,:) = 0.0
rho0Inv = 1.0/config_rho0
+ allocate(fluxVert(0:nVertLevels))
do iEdge=1,grid % nEdgesSolve
cell1 = cellsOnEdge(1,iEdge)
cell2 = cellsOnEdge(2,iEdge)
@@ -504,12 +509,36 @@
+ u_diffusion &
- ( ke(k,cell2) - ke(k,cell1) ) / dcEdge(iEdge)
-!ocean
- tend_u(k,iEdge) = tend_u(k,iEdge) + u_src(k,iEdge)/rho_ref/h_edge(k,iEdge)
-!ocean
+ end do
- end do
+ ! forcing in top layer only
+ tend_u(1,iEdge) = tend_u(1,iEdge) + u_src(1,iEdge)/rho_ref/h_edge(1,iEdge)
+
+ ! bottom drag - change from nVertLevels to KMT later.
+ ! du/dt = u/tau where tau=100 days.
+ tend_u(nVertLevels,iEdge) = tend_u(nVertLevels,iEdge) &
+ - u(nVertLevels,iEdge)/(100.0*86400.0)
+
+ ! vertical mixing
+ fluxVert(0) = 0.0
+ fluxVert(nVertLevels) = 0.0
+ do k=nVertLevels-1,1,-1
+ fluxVert(k) = ( u(k,iEdge) - u(k+1,iEdge) ) &
+ / (zMidEdge(k,iEdge) -zMidEdge(k+1,iEdge))
+ enddo
+ ! note vertical viscosity is hardwired as 2.5e-4 here.
+ fluxVert = 2.5e-4 * fluxVert
+ do k=1,nVertLevels
+ if(k.eq.1) then
+ dist = zSurfaceEdge(iEdge) - zBotEdge(k,iEdge)
+ else
+ dist = zBotEdge(k-1,iEdge) - zBotEdge(k,iEdge)
+ endif
+ tend_u(k,iEdge) = tend_u(k,iEdge) - (fluxVert(k-1) - fluxVert(k))/dist
+ enddo
+
end do
+ deallocate(fluxVert)
#endif
#ifdef NCAR_FORMULATION
@@ -682,12 +711,12 @@
integer :: nCells, nEdges, nVertices, nVertLevels
real (kind=RKIND), dimension(:), pointer :: &
h_s, fVertex, fEdge, dvEdge, dcEdge, areaCell, areaTriangle, &
- zSurface, hZLevel
+ zSurface, hZLevel, zSurfaceEdge
real (kind=RKIND), dimension(:,:), pointer :: &
vh, weightsOnEdge, kiteAreasOnVertex, h_edge, h, u, v, w, &
tend_h, tend_u, circulation, vorticity, ke, &
pv_edge, pv_vertex, pv_cell, gradPVn, gradPVt, divergence, &
- zmid, zbot, pmid, pbot, MontPot, rho, temperature, salinity, pmidZLevel
+ zmid, zbot, zmidEdge, zbotEdge, pmid, pbot, MontPot, rho, temperature, salinity, pmidZLevel
real (kind=RKIND), dimension(:,:,:), pointer :: tracers
character :: c1*6
@@ -718,11 +747,14 @@
tracers => s % tracers % array
zmid => s % zmid % array
zbot => s % zbot % array
+ zmidEdge => s % zmidEdge % array
+ zbotEdge => s % zbotEdge % array
pmid => s % pmid % array
pmidZLevel => s % pmidZLevel % array
pbot => s % pbot % array
MontPot => s % MontPot % array
zSurface => s % zSurface % array
+ zSurfaceEdge=> s % zSurfaceEdge % array
!mrp 100112 end
weightsOnEdge => grid % weightsOnEdge % array
@@ -1046,7 +1078,19 @@
end do
!mrp 100112 end
+ do iEdge=1,nEdges
+ cell1 = cellsOnEdge(1,iEdge)
+ cell2 = cellsOnEdge(2,iEdge)
+ if(cell1.gt.0.and.cell2.gt.0) then
+ do k=1,nVertLevels
+ zBotEdge(k,iEdge) = (zBot(k,cell1)+zBot(k,cell2))/2.0
+ zMidEdge(k,iEdge) = (zMid(k,cell1)+zMid(k,cell2))/2.0
+ enddo
+ zSurfaceEdge(iEdge) = (zSurface(cell1) + zSurface(cell2))/2.0
+ endif
+ enddo
+
!mrp 100426:
!
! For z-level model.
</font>
</pre>