<p><b>qchen3@fsu.edu</b> 2011-12-20 10:12:07 -0700 (Tue, 20 Dec 2011)</p><p>BRANCH COMMIT<br>
<br>
The following variables are added to Registry:<br>
config_h_kappa<br>
config_h_kappa_q (for specifying a constant GM coefficient)<br>
h_kappa<br>
h_kappa_q (as a state array variable, to accommodate the idea of spatially <br>
varying GM coefficient)<br>
uBolus<br>
uTransport<br>
<br>
Two subroutines are added in mpas_ocn_tendency.F:<br>
ocn_get_h_kappa and ocn_get_h_kappa_q.<br>
Right now, these two subroutines assign config_h_kappa and config_h_kappa_q <br>
to h_kappa and h_kappa_q, respectively. <br>
<br>
The computation of uBolus and uTransport is also implemented<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/gmvar/Makefile
===================================================================
--- branches/ocean_projects/gmvar/Makefile        2011-12-20 15:52:45 UTC (rev 1264)
+++ branches/ocean_projects/gmvar/Makefile        2011-12-20 17:12:07 UTC (rev 1265)
@@ -111,9 +111,9 @@
        "CC = gcc" \
        "SFC = ifort" \
        "SCC = gcc" \
-        "FFLAGS = -real-size 64 -O3 -convert big_endian -FR" \
-        "CFLAGS = -O3 -m64" \
-        "LDFLAGS = -O3" \
+        "FFLAGS = -real-size 64 -g -convert big_endian -FR" \
+        "CFLAGS = -g -m64" \
+        "LDFLAGS = -g" \
        "CORE = $(CORE)" \
        "CPPFLAGS = -DRKIND=8 $(MODEL_FORMULATION) -DUNDERSCORE -m64 $(FILE_OFFSET) $(ZOLTAN_DEFINE)" )
Modified: branches/ocean_projects/gmvar/src/core_ocean/Registry
===================================================================
--- branches/ocean_projects/gmvar/src/core_ocean/Registry        2011-12-20 15:52:45 UTC (rev 1264)
+++ branches/ocean_projects/gmvar/src/core_ocean/Registry        2011-12-20 17:12:07 UTC (rev 1265)
@@ -47,7 +47,9 @@
namelist logical hmix config_include_KE_vertex false
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 real hmix config_apvm_upwinding 0.0
+namelist real hmix config_h_kappa 0.0
+namelist real hmix config_h_kappa_q 0.0
namelist logical hmix config_mom_decay false
namelist real hmix config_mom_decay_time 3600.0
namelist character vmix config_vert_visc_type const
@@ -188,7 +190,7 @@
# Prognostic variables: read from input, saved in restart, and written to output
var persistent real u ( nVertLevels nEdges Time ) 2 ir u state - -
-var persistent real h ( nVertLevels nCells Time ) 2 ir h state - -
+var persistent real h ( nVertLevels nCells Time ) 2 iro h state - -
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
@@ -217,6 +219,16 @@
# Diagnostic fields: only written to output
var persistent real v ( nVertLevels nEdges Time ) 2 - v state - -
+var persistent real uBolus ( nVertLevels nEdges Time ) 2 o uBolus state - -
+var persistent real uBolusX ( nVertLevels nEdges Time ) 2 o uBolusX state - -
+var persistent real uBolusY ( nVertLevels nEdges Time ) 2 o uBolusY state - -
+var persistent real uBolusZ ( nVertLevels nEdges Time ) 2 o uBolusZ state - -
+var persistent real uBolusZonal ( nVertLevels nEdges Time ) 2 o uBolusZonal state - -
+var persistent real uBolusMeridional ( nVertLevels nEdges Time ) 2 o uBolusMeridional state - -
+var persistent real uTransport ( nVertLevels nEdges Time ) 2 o uTransport state - -
+var persistent real h_kappa ( nVertLevels nEdges Time ) 2 o h_kappa state - -
+var persistent real h_kappa_q ( nVertLevels nEdges Time ) 2 o h_kappa_q state - -
+var persistent real hu ( nVertLevels nEdges Time ) 2 0 hu state - -
var persistent real divergence ( nVertLevels nCells Time ) 2 o divergence state - -
var persistent real vorticity ( nVertLevels nVertices Time ) 2 o vorticity state - -
var persistent real pv_edge ( nVertLevels nEdges Time ) 2 - pv_edge state - -
Modified: branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_tendency.F
===================================================================
--- branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_tendency.F        2011-12-20 15:52:45 UTC (rev 1264)
+++ branches/ocean_projects/gmvar/src/core_ocean/mpas_ocn_tendency.F        2011-12-20 17:12:07 UTC (rev 1265)
@@ -620,7 +620,8 @@
weightsOnEdge, kiteAreasOnVertex, h_edge, h, u, v, w, pressure,&
circulation, vorticity, ke, ke_edge, MontPot, wTop, &
pv_edge, pv_vertex, pv_cell, gradPVn, gradPVt, divergence, &
- rho, temperature, salinity, kev, kevc
+ rho, temperature, salinity, kev, kevc, uBolus, uTransport, &
+ h_kappa, h_kappa_q
real (kind=RKIND), dimension(:,:,:), pointer :: tracers
real (kind=RKIND), dimension(:), allocatable:: pTop
real (kind=RKIND), dimension(:,:), allocatable:: div_u
@@ -642,6 +643,8 @@
h => s % h % array
u => s % u % array
v => s % v % array
+ uBolus => s % uBolus % array
+ uTransport => s % uTransport % array
wTop => s % wTop % array
h_edge => s % h_edge % array
circulation => s % circulation % array
@@ -660,6 +663,8 @@
tracers => s % tracers % array
MontPot => s % MontPot % array
pressure => s % pressure % array
+ h_kappa => s % h_kappa % array
+ h_kappa_q => s % h_kappa_q % array
weightsOnEdge => grid % weightsOnEdge % array
kiteAreasOnVertex => grid % kiteAreasOnVertex % array
@@ -695,6 +700,7 @@
boundaryEdge => grid % boundaryEdge % array
boundaryCell => grid % boundaryCell % array
+
!
! Compute height on cell edges at velocity locations
! Namelist options control the order of accuracy of the reconstructed h_edge value
@@ -906,6 +912,24 @@
enddo
!
+ ! Compute Bolus velocity as bolus == - /kappa (/grad h) / h
+ !
+ uBolus(:,:) = 0.0
+ call ocn_get_h_kappa(s, grid)
+ do iEdge = 1,nEdges
+ cell1 = cellsOnEdge(1,iEdge)
+ cell2 = cellsOnEdge(2,iEdge)
+ do k=1,maxLevelEdgeTop(iEdge)
+ uBolus(k,iEdge) = -h_kappa(k,iEdge) * (h(k,cell2) - h(k,cell1))/dcEdge(iEdge)/h_edge(k,iEdge)
+ end do
+ end do
+
+ !
+ ! Compute the transport velocity (this velocity transports thickness and PV)
+ !
+ uTransport(:,:) = u(:,:) + uBolus(:,:)
+
+ !
! Compute kinetic energy in each cell by blending ke and kevc
!
if(config_include_KE_vertex) then
@@ -1034,7 +1058,7 @@
do iEdge = 1,nEdges
do k = 1,maxLevelEdgeBot(iEdge)
pv_edge(k,iEdge) = pv_edge(k,iEdge) &
- - 0.5 * dt* ( u(k,iEdge) * gradPVn(k,iEdge) &
+ - config_apvm_upwinding * dt* ( u(k,iEdge) * gradPVn(k,iEdge) &
+ v(k,iEdge) * gradPVt(k,iEdge) )
enddo
enddo
@@ -1346,6 +1370,38 @@
end subroutine ocn_fuperp!}}}
+
+ subroutine ocn_get_h_kappa(s, grid)
+
+ type (state_type), intent(inout) :: s
+ type (mesh_type), intent(in) :: grid
+
+ real(kind=RKIND), dimension(:,:), pointer :: h_kappa
+
+
+ h_kappa => s % h_kappa % array
+
+ h_kappa(:,:) = config_h_kappa
+
+
+ end subroutine ocn_get_h_kappa
+
+
+ subroutine ocn_get_h_kappa_q(s, grid)
+
+ type (state_type), intent(inout) :: s
+ type (mesh_type), intent(in) :: grid
+
+ real(kind=RKIND), dimension(:,:), pointer :: h_kappa_q
+
+
+ h_kappa_q => s % h_kappa_q % array
+
+ h_kappa_q(:,:) = config_h_kappa_q
+
+
+ end subroutine ocn_get_h_kappa_q
+
!***********************************************************************
end module ocn_tendency
</font>
</pre>