<p><b>dwj07@fsu.edu</b> 2011-09-19 11:14:29 -0600 (Mon, 19 Sep 2011)</p><p><br>
-- BRANCH COMMITT --<br>
<br>
Adding modules for Tracer hmix and Restoring for temperature and sainity.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/performance/src/core_ocean/Makefile
===================================================================
--- branches/ocean_projects/performance/src/core_ocean/Makefile 2011-09-16 21:57:58 UTC (rev 1007)
+++ branches/ocean_projects/performance/src/core_ocean/Makefile 2011-09-19 17:14:29 UTC (rev 1008)
@@ -26,6 +26,10 @@
module_OcnTracerHadv2.o \
module_OcnTracerHadv3.o \
module_OcnTracerHadv4.o \
+ module_OcnTracerHmix.o \
+ module_OcnTracerHmixDel2.o \
+ module_OcnTracerHmixDel4.o \
+ module_OcnRestoring.o \
module_time_integration.o \
module_global_diagnostics.o
@@ -88,6 +92,14 @@
module_OcnTracerVadvStencil4.o:
+module_OcnTracerHmix.o: module_OcnTracerHmixDel2.o module_OcnTracerHmixDel4.o
+
+module_OcnTracerHmixDel2.o:
+
+module_OcnTracerHmixDel4.o:
+
+module_OcnRestoring.o:
+
module_mpas_core.o: module_advection.o \
module_OcnThickHadv.o \
module_OcnThickVadv.o \
@@ -114,6 +126,10 @@
module_OcnTracerVadvStencil2.o \
module_OcnTracerVadvStencil3.o \
module_OcnTracerVadvStencil4.o \
+ module_OcnTracerHmix.o \
+ module_OcnTracerHmixDel2.o \
+ module_OcnTracerHmixDel4.o \
+ module_OcnRestoring.o \
module_time_integration.o
clean:
Added: branches/ocean_projects/performance/src/core_ocean/module_OcnRestoring.F
===================================================================
--- branches/ocean_projects/performance/src/core_ocean/module_OcnRestoring.F (rev 0)
+++ branches/ocean_projects/performance/src/core_ocean/module_OcnRestoring.F 2011-09-19 17:14:29 UTC (rev 1008)
@@ -0,0 +1,199 @@
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
+!
+! OcnRestoring
+!
+!> \brief MPAS ocean restoring
+!> \author Doug Jacobsen
+!> \date 19 September 2011
+!> \version SVN:$Id:$
+!> \details
+!> This module contains the main driver routine for computing
+!> tendencies for restoring.
+!
+!-----------------------------------------------------------------------
+
+module OcnRestoring
+
+ use grid_types
+ use configure
+
+ implicit none
+ private
+ save
+
+ !--------------------------------------------------------------------
+ !
+ ! Public parameters
+ !
+ !--------------------------------------------------------------------
+
+ !--------------------------------------------------------------------
+ !
+ ! Public member functions
+ !
+ !--------------------------------------------------------------------
+
+ public :: OcnRestoringTend, &
+ OcnRestoringInit
+
+ !--------------------------------------------------------------------
+ !
+ ! Private module variables
+ !
+ !--------------------------------------------------------------------
+
+ logical :: restoringOn
+
+ real (kind=RKIND) :: temperatureTimeScale, salinityTimeScale
+
+
+!***********************************************************************
+
+contains
+
+!***********************************************************************
+!
+! routine OcnRestoringTend
+!
+!> \brief Computes tendency term for restoring
+!> \author Doug Jacobsen
+!> \date 19 September 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine computes the restoring tendency for tracers
+!> based on current state.
+!
+!-----------------------------------------------------------------------
+
+ subroutine OcnRestoringTend(grid, h, indexT, indexS, tracers, tend, err)
+
+ !-----------------------------------------------------------------
+ !
+ ! input variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:), intent(in) :: &
+ h !< Input: thickness
+
+ type (mesh_type), intent(in) :: &
+ grid !< Input: grid information
+
+ real (kind=RKIND), dimension(:,:,:), intent(in) :: &
+ tracers !< Input: tracer quantities
+
+ integer, intent(in) :: indexT, indexS
+
+ !-----------------------------------------------------------------
+ !
+ ! input/output variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:,:), intent(inout) :: &
+ tend !< Input/Output: velocity tendency
+
+ !-----------------------------------------------------------------
+ !
+ ! output variables
+ !
+ !-----------------------------------------------------------------
+
+ integer, intent(out) :: err
+
+ !-----------------------------------------------------------------
+ !
+ ! local variables
+ !
+ !-----------------------------------------------------------------
+
+ integer :: iCell, nCellsSolve, k
+
+ real (kind=RKIND), dimension(:), pointer :: temperatureRestore, salinityRestore
+
+ !-----------------------------------------------------------------
+ !
+ ! call relevant routines for computing tendencies
+ ! note that the user can choose multiple options and the
+ ! tendencies will be added together
+ !
+ !-----------------------------------------------------------------
+
+ err = 0
+
+ if(.not.restoringOn) return
+
+ nCellsSolve = grid % nCellsSolve
+
+ temperatureRestore => grid % temperatureRestore % array
+ salinityRestore => grid % salinityRestore % array
+
+ k = 1 ! restoring only in top layer
+ do iCell=1,nCellsSolve
+
+ tend(indexT, k, iCell) = tend(indexT, k, iCell) &
+ - h(k,iCell)*(tracers(indexT, k, iCell) - temperatureRestore(iCell)) &
+ / (temperatureTimeScale * 86400.0)
+
+ tend(indexS, k, iCell) = tend(indexS, k, iCell) &
+ - h(k,iCell)*(tracers(indexS, k, iCell) - salinityRestore(iCell)) &
+ / (salinityTimeScale * 86400.0)
+
+! write(6,10) iCell, tracers(indexT, k, iCell), &
+! temperatureRestore(iCell), tracers(indexT, k, iCell), &
+! (tracers(indexT, k, iCell) - temperatureRestore(iCell)) &
+! / (config_restoreT_timescale * 86400.0)
+
+ enddo
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnRestoringTend
+
+!***********************************************************************
+!
+! routine OcnRestoringInit
+!
+!> \brief Initializes ocean tracer horizontal mixing quantities
+!> \author Doug Jacobsen
+!> \date 19 September 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine initializes a variety of quantities related to
+!> horizontal velocity mixing in the ocean. Since a variety of
+!> parameterizations are available, this routine primarily calls the
+!> individual init routines for each parameterization.
+!
+!-----------------------------------------------------------------------
+
+
+ subroutine OcnRestoringInit(err)
+
+ !--------------------------------------------------------------------
+
+ !-----------------------------------------------------------------
+ !
+ ! call individual init routines for each parameterization
+ !
+ !-----------------------------------------------------------------
+
+ integer, intent(out) :: err
+
+ err = 0
+ restoringOn = .false.
+
+ if(config_restoreTS) then
+ restoringOn = .true.
+ temperatureTimeScale = config_restoreT_timescale
+ salinityTimeScale = config_restoreS_timescale
+ endif
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnRestoringInit
+
+!***********************************************************************
+
+end module OcnRestoring
+
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Added: branches/ocean_projects/performance/src/core_ocean/module_OcnTracerHmix.F
===================================================================
--- branches/ocean_projects/performance/src/core_ocean/module_OcnTracerHmix.F (rev 0)
+++ branches/ocean_projects/performance/src/core_ocean/module_OcnTracerHmix.F 2011-09-19 17:14:29 UTC (rev 1008)
@@ -0,0 +1,175 @@
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
+!
+! OcnTracerHmix
+!
+!> \brief MPAS ocean horizontal tracer mixing driver
+!> \author Doug Jacobsen
+!> \date 19 September 2011
+!> \version SVN:$Id:$
+!> \details
+!> This module contains the main driver routine for computing
+!> horizontal mixing tendencies.
+!>
+!> It provides an init and a tend function. Each are described below.
+!
+!-----------------------------------------------------------------------
+
+module OcnTracerHmix
+
+ use grid_types
+ use configure
+ use OcnTracerHmixDel2
+ use OcnTracerHmixDel4
+
+ implicit none
+ private
+ save
+
+ !--------------------------------------------------------------------
+ !
+ ! Public parameters
+ !
+ !--------------------------------------------------------------------
+
+ !--------------------------------------------------------------------
+ !
+ ! Public member functions
+ !
+ !--------------------------------------------------------------------
+
+ public :: OcnTracerHmixTend, &
+ OcnTracerHmixInit
+
+ !--------------------------------------------------------------------
+ !
+ ! Private module variables
+ !
+ !--------------------------------------------------------------------
+
+
+!***********************************************************************
+
+contains
+
+!***********************************************************************
+!
+! routine OcnTracerHmixTend
+!
+!> \brief Computes tendency term for horizontal tracer mixing
+!> \author Doug Jacobsen
+!> \date 19 September 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine computes the horizontal mixing tendency for tracer
+!> based on current state and user choices of mixing parameterization.
+!> Multiple parameterizations may be chosen and added together. These
+!> tendencies are generally computed by calling the specific routine
+!> for the chosen parameterization, so this routine is primarily a
+!> driver for managing these choices.
+!
+!-----------------------------------------------------------------------
+
+ subroutine OcnTracerHmixTend(grid, h_edge, tracers, tend, err)
+
+ !-----------------------------------------------------------------
+ !
+ ! input variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:), intent(in) :: &
+ h_edge !< Input: thickness at edge
+
+ type (mesh_type), intent(in) :: &
+ grid !< Input: grid information
+
+ real (kind=RKIND), dimension(:,:,:), intent(in) :: &
+ tracers !< Input: tracer quantities
+
+ !-----------------------------------------------------------------
+ !
+ ! input/output variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:,:), intent(inout) :: &
+ tend !< Input/Output: velocity tendency
+
+ !-----------------------------------------------------------------
+ !
+ ! output variables
+ !
+ !-----------------------------------------------------------------
+
+ integer, intent(out) :: err
+
+ !-----------------------------------------------------------------
+ !
+ ! local variables
+ !
+ !-----------------------------------------------------------------
+
+ integer :: err1, err2
+
+ !-----------------------------------------------------------------
+ !
+ ! call relevant routines for computing tendencies
+ ! note that the user can choose multiple options and the
+ ! tendencies will be added together
+ !
+ !-----------------------------------------------------------------
+
+ call OcnTracerHmixDel2Tend(grid, h_edge, tracers, tend, err1)
+ call OcnTracerHmixDel4Tend(grid, h_edge, tracers, tend, err2)
+
+ err = err1 .or. err2
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnTracerHmixTend
+
+!***********************************************************************
+!
+! routine OcnTracerHmixInit
+!
+!> \brief Initializes ocean tracer horizontal mixing quantities
+!> \author Doug Jacobsen
+!> \date 19 September 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine initializes a variety of quantities related to
+!> horizontal velocity mixing in the ocean. Since a variety of
+!> parameterizations are available, this routine primarily calls the
+!> individual init routines for each parameterization.
+!
+!-----------------------------------------------------------------------
+
+
+ subroutine OcnTracerHmixInit(err)
+
+ !--------------------------------------------------------------------
+
+ !-----------------------------------------------------------------
+ !
+ ! call individual init routines for each parameterization
+ !
+ !-----------------------------------------------------------------
+
+ integer, intent(out) :: err
+
+ integer :: err1, err2
+
+ call OcnTracerHmixDel2Init(err1)
+ call OcnTracerHmixDel4Init(err2)
+
+ err = err1 .or. err2
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnTracerHmixInit
+
+!***********************************************************************
+
+end module OcnTracerHmix
+
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Added: branches/ocean_projects/performance/src/core_ocean/module_OcnTracerHmixDel2.F
===================================================================
--- branches/ocean_projects/performance/src/core_ocean/module_OcnTracerHmixDel2.F (rev 0)
+++ branches/ocean_projects/performance/src/core_ocean/module_OcnTracerHmixDel2.F 2011-09-19 17:14:29 UTC (rev 1008)
@@ -0,0 +1,232 @@
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
+!
+! OcnTracerHmixDel2
+!
+!> \brief MPAS ocean horizontal tracer mixing driver
+!> \author Doug Jacobsen
+!> \date 19 September 2011
+!> \version SVN:$Id:$
+!> \details
+!> This module contains the main driver routine for computing
+!> horizontal mixing tendencies.
+!>
+!> It provides an init and a tend function. Each are described below.
+!
+!-----------------------------------------------------------------------
+
+module OcnTracerHmixDel2
+
+ use grid_types
+ use configure
+ use timer
+
+ implicit none
+ private
+ save
+
+ !--------------------------------------------------------------------
+ !
+ ! Public parameters
+ !
+ !--------------------------------------------------------------------
+
+ !--------------------------------------------------------------------
+ !
+ ! Public member functions
+ !
+ !--------------------------------------------------------------------
+
+ public :: OcnTracerHmixDel2Tend, &
+ OcnTracerHmixDel2Init
+
+ !--------------------------------------------------------------------
+ !
+ ! Private module variables
+ !
+ !--------------------------------------------------------------------
+
+ logical :: del2On
+
+ real (kind=RKIND) :: eddyDiff2
+
+
+!***********************************************************************
+
+contains
+
+!***********************************************************************
+!
+! routine OcnTracerHmixDel2Tend
+!
+!> \brief Computes laplacian tendency term for horizontal tracer mixing
+!> \author Doug Jacobsen
+!> \date 19 September 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine computes the horizontal mixing tendency for tracers
+!> based on current state using a laplacian parameterization.
+!
+!-----------------------------------------------------------------------
+
+ subroutine OcnTracerHmixDel2Tend(grid, h_edge, tracers, tend, err)
+
+ !-----------------------------------------------------------------
+ !
+ ! input variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:), intent(in) :: &
+ h_edge !< Input: thickness at edge
+
+ type (mesh_type), intent(in) :: &
+ grid !< Input: grid information
+
+ real (kind=RKIND), dimension(:,:,:), intent(in) :: &
+ tracers !< Input: tracer quantities
+
+ !-----------------------------------------------------------------
+ !
+ ! input/output variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:,:), intent(inout) :: &
+ tend !< Input/Output: velocity tendency
+
+ !-----------------------------------------------------------------
+ !
+ ! output variables
+ !
+ !-----------------------------------------------------------------
+
+ integer, intent(out) :: err
+
+ !-----------------------------------------------------------------
+ !
+ ! local variables
+ !
+ !-----------------------------------------------------------------
+
+ integer :: iEdge, nEdges, nVertLevels, cell1, cell2
+ integer :: k, iTracer, num_tracers
+
+ integer, dimension(:,:), allocatable :: boundaryMask
+
+ integer, dimension(:), pointer :: maxLevelEdgeTop
+ integer, dimension(:,:), pointer :: cellsOnEdge, boundaryEdge
+
+ real (kind=RKIND) :: invAreaCell1, invAreaCell2
+ real (kind=RKIND) :: tracer_turb_flux, flux
+
+ real (kind=RKIND), dimension(:), pointer :: areaCell, dvEdge, dcEdge
+ real (kind=RKIND), dimension(:), pointer :: meshScalingDel2
+
+ !-----------------------------------------------------------------
+ !
+ ! call relevant routines for computing tendencies
+ ! note that the user can choose multiple options and the
+ ! tendencies will be added together
+ !
+ !-----------------------------------------------------------------
+
+ err = 0
+
+ if (.not.del2On) return
+
+ call timer_start("compute_scalar_tend-horiz diff 2")
+
+ nEdges = grid % nEdges
+ nVertLevels = grid % nVertLevels
+ num_tracers = size(tracers, dim=1)
+
+ maxLevelEdgeTop => grid % maxLevelEdgeTop % array
+ cellsOnEdge => grid % cellsOnEdge % array
+ boundaryEdge => grid % boundaryEdge % array
+ areaCell => grid % areaCell % array
+ dvEdge => grid % dvEdge % array
+ dcEdge => grid % dcEdge % array
+ meshScalingDel2 => grid % meshScalingDel2 % array
+
+ !
+ ! compute a boundary mask to enforce insulating boundary conditions in the horizontal
+ !
+ allocate(boundaryMask(nVertLevels, nEdges+1))
+ boundaryMask = 1.0
+ where(boundaryEdge.eq.1) boundaryMask=0.0
+
+ do iEdge=1,nEdges
+ cell1 = cellsOnEdge(1,iEdge)
+ cell2 = cellsOnEdge(2,iEdge)
+ invAreaCell1 = 1.0/areaCell(cell1)
+ invAreaCell2 = 1.0/areaCell(cell2)
+
+ do k=1,maxLevelEdgeTop(iEdge)
+ do iTracer=1,num_tracers
+ ! \kappa_2 </font>
<font color="blue">abla \phi on edge
+ tracer_turb_flux = meshScalingDel2(iEdge) * eddyDiff2 &
+ *( 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(iTracer,k,cell1) = tend(iTracer,k,cell1) + flux * invAreaCell1
+ tend(iTracer,k,cell2) = tend(iTracer,k,cell2) - flux * invAreaCell2
+ end do
+ end do
+
+ end do
+
+ deallocate(boundaryMask)
+ call timer_stop("compute_scalar_tend-horiz diff 2")
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnTracerHmixDel2Tend
+
+!***********************************************************************
+!
+! routine OcnTracerHmixDel2Init
+!
+!> \brief Initializes ocean tracer horizontal mixing quantities
+!> \author Doug Jacobsen
+!> \date 19 September 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine initializes a variety of quantities related to
+!> laplacian horizontal velocity mixing in the ocean.
+!
+!-----------------------------------------------------------------------
+
+
+ subroutine OcnTracerHmixDel2Init(err)
+
+ !--------------------------------------------------------------------
+
+ !-----------------------------------------------------------------
+ !
+ ! call individual init routines for each parameterization
+ !
+ !-----------------------------------------------------------------
+
+ integer, intent(out) :: err
+
+ err = 0
+
+ del2on = .false.
+
+ if ( config_h_tracer_eddy_diff2 > 0.0 ) then
+ del2On = .true.
+ eddyDiff2 = config_h_tracer_eddy_diff2
+ endif
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnTracerHmixDel2Init
+
+!***********************************************************************
+
+end module OcnTracerHmixDel2
+
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Added: branches/ocean_projects/performance/src/core_ocean/module_OcnTracerHmixDel4.F
===================================================================
--- branches/ocean_projects/performance/src/core_ocean/module_OcnTracerHmixDel4.F (rev 0)
+++ branches/ocean_projects/performance/src/core_ocean/module_OcnTracerHmixDel4.F 2011-09-19 17:14:29 UTC (rev 1008)
@@ -0,0 +1,263 @@
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
+!
+! OcnTracerHmixDel4
+!
+!> \brief MPAS ocean horizontal tracer mixing driver
+!> \author Doug Jacobsen
+!> \date 19 September 2011
+!> \version SVN:$Id:$
+!> \details
+!> This module contains the main driver routine for computing
+!> horizontal mixing tendencies.
+!>
+!> It provides an init and a tend function. Each are described below.
+!
+!-----------------------------------------------------------------------
+
+module OcnTracerHmixDel4
+
+ use grid_types
+ use configure
+ use timer
+
+ implicit none
+ private
+ save
+
+ !--------------------------------------------------------------------
+ !
+ ! Public parameters
+ !
+ !--------------------------------------------------------------------
+
+ !--------------------------------------------------------------------
+ !
+ ! Public member functions
+ !
+ !--------------------------------------------------------------------
+
+ public :: OcnTracerHmixDel4Tend, &
+ OcnTracerHmixDel4Init
+
+ !--------------------------------------------------------------------
+ !
+ ! Private module variables
+ !
+ !--------------------------------------------------------------------
+
+ logical :: Del4On
+
+ real (kind=RKIND) :: eddyDiff4
+
+
+!***********************************************************************
+
+contains
+
+!***********************************************************************
+!
+! routine OcnTracerHmixDel4Tend
+!
+!> \brief Computes biharmonic tendency term for horizontal tracer mixing
+!> \author Doug Jacobsen
+!> \date 19 September 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine computes the horizontal mixing tendency for tracers
+!> based on current state using a biharmonic parameterization.
+!
+!-----------------------------------------------------------------------
+
+ subroutine OcnTracerHmixDel4Tend(grid, h_edge, tracers, tend, err)
+
+ !-----------------------------------------------------------------
+ !
+ ! input variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:), intent(in) :: &
+ h_edge !< Input: thickness at edge
+
+ type (mesh_type), intent(in) :: &
+ grid !< Input: grid information
+
+ real (kind=RKIND), dimension(:,:,:), intent(in) :: &
+ tracers !< Input: tracer quantities
+
+ !-----------------------------------------------------------------
+ !
+ ! input/output variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:,:), intent(inout) :: &
+ tend !< Input/Output: velocity tendency
+
+ !-----------------------------------------------------------------
+ !
+ ! output variables
+ !
+ !-----------------------------------------------------------------
+
+ integer, intent(out) :: err
+
+ !-----------------------------------------------------------------
+ !
+ ! local variables
+ !
+ !-----------------------------------------------------------------
+
+ integer :: iEdge, nEdges, num_tracers, nVertLevels, nCells
+ integer :: iTracer, k, iCell, cell1, cell2
+
+ integer, dimension(:,:), allocatable :: boundaryMask
+
+ integer, dimension(:), pointer :: maxLevelEdgeTop, maxLevelCell
+ integer, dimension(:,:), pointer :: boundaryEdge, cellsOnEdge
+
+ real (kind=RKIND) :: invAreaCell1, invAreaCell2, r, tracer_turb_flux, flux
+
+ real (kind=RKIND), dimension(:,:,:), allocatable :: delsq_tracer
+
+ real (kind=RKIND), dimension(:), pointer :: dcEdge, dvEdge, areaCell, meshScalingDel4
+
+
+ !-----------------------------------------------------------------
+ !
+ ! call relevant routines for computing tendencies
+ ! note that the user can choose multiple options and the
+ ! tendencies will be added together
+ !
+ !-----------------------------------------------------------------
+
+ err = 0
+
+ if (.not.Del4On) return
+
+ call timer_start("compute_scalar_tend-horiz diff 4")
+
+ nEdges = grid % nEdges
+ nCells = grid % nCells
+ num_tracers = size(tracers, dim=1)
+ nVertLevels = grid % nVertLevels
+
+ maxLevelEdgeTop => grid % maxLevelEdgeTop % array
+ maxLevelCell => grid % maxLevelCell % array
+ boundaryEdge => grid % boundaryEdge % array
+ cellsOnEdge => grid % cellsOnEdge % array
+
+ dcEdge => grid % dcEdge % array
+ dvEdge => grid % dvEdge % array
+ areaCell => grid % areaCell % array
+ meshScalingDel4 => grid % meshScalingDel4 % array
+
+ allocate(boundaryMask(nVertLevels, nEdges+1))
+ boundaryMask = 1.0
+ where(boundaryEdge.eq.1) boundaryMask=0.0
+
+ allocate(delsq_tracer(num_tracers,nVertLevels, nCells+1))
+
+ delsq_tracer(:,:,:) = 0.
+
+ ! first del2: div(h </font>
<font color="blue">abla \phi) at cell center
+ do iEdge=1,nEdges
+ cell1 = cellsOnEdge(1,iEdge)
+ cell2 = cellsOnEdge(2,iEdge)
+
+ do k=1,maxLevelEdgeTop(iEdge)
+ 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)
+ 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)
+ end do
+ end do
+ end do
+
+ do iCell = 1,nCells
+ r = 1.0 / areaCell(iCell)
+ do k=1,maxLevelCell(iCell)
+ do iTracer=1,num_tracers
+ delsq_tracer(iTracer,k,iCell) = delsq_tracer(iTracer,k,iCell) * r
+ end do
+ 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)
+ invAreaCell1 = 1.0 / areaCell(cell1)
+ invAreaCell2 = 1.0 / areaCell(cell2)
+
+ do k=1,maxLevelEdgeTop(iEdge)
+ do iTracer=1,num_tracers
+ tracer_turb_flux = meshScalingDel4(iEdge) * eddyDiff4 &
+ *( delsq_tracer(iTracer,k,cell2) &
+ - delsq_tracer(iTracer,k,cell1))/dcEdge(iEdge)
+ flux = dvEdge (iEdge) * tracer_turb_flux
+
+ tend(iTracer,k,cell1) = tend(iTracer,k,cell1) &
+ - flux * invAreaCell1 * boundaryMask(k,iEdge)
+ tend(iTracer,k,cell2) = tend(iTracer,k,cell2) &
+ + flux * invAreaCell2 * boundaryMask(k,iEdge)
+
+ enddo
+ enddo
+ end do
+
+ deallocate(delsq_tracer)
+ call timer_stop("compute_scalar_tend-horiz diff 4")
+ !--------------------------------------------------------------------
+
+ end subroutine OcnTracerHmixDel4Tend
+
+!***********************************************************************
+!
+! routine OcnTracerHmixDel4Init
+!
+!> \brief Initializes ocean tracer horizontal mixing quantities
+!> \author Doug Jacobsen
+!> \date 19 September 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine initializes a variety of quantities related to
+!> biharmonic horizontal velocity mixing in the ocean.
+!
+!-----------------------------------------------------------------------
+
+
+ subroutine OcnTracerHmixDel4Init(err)
+
+ !--------------------------------------------------------------------
+
+ !-----------------------------------------------------------------
+ !
+ ! call individual init routines for each parameterization
+ !
+ !-----------------------------------------------------------------
+
+ integer, intent(out) :: err
+
+ err = 0
+ Del4on = .false.
+
+ if ( config_h_tracer_eddy_diff4 > 0.0 ) then
+ Del4On = .true.
+ eddyDiff4 = config_h_tracer_eddy_diff4
+ endif
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnTracerHmixDel4Init
+
+!***********************************************************************
+
+end module OcnTracerHmixDel4
+
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Modified: branches/ocean_projects/performance/src/core_ocean/module_mpas_core.F
===================================================================
--- branches/ocean_projects/performance/src/core_ocean/module_mpas_core.F 2011-09-16 21:57:58 UTC (rev 1007)
+++ branches/ocean_projects/performance/src/core_ocean/module_mpas_core.F 2011-09-19 17:14:29 UTC (rev 1008)
@@ -12,6 +12,8 @@
use OcnTracerHadv
use OcnTracerVadv
+ use OcnTracerHmix
+ use OcnRestoring
type (io_output_object) :: restart_obj
integer :: restart_frame
@@ -77,6 +79,9 @@
call mpas_init_block(block, block % mesh, dt)
block % state % time_levs(1) % state % xtime % scalar = startTimeStamp
block => block % next
+
+ !dwj 110919 This allows the restorings to grab the indices for
+ ! temperature and salinity tracers from state.
end do
call OcnVelPressureGradInit(err)
@@ -86,6 +91,8 @@
call OcnTracerHadvInit(err)
call OcnTracerVadvInit(err)
+ call OcnTracerHmixInit(err)
+ call OcnRestoringInit(err)
! mrp 100316 In order for this to work, we need to pass domain % dminfo as an
! input arguement into mpas_init. Ask about that later. For now, there will be
Modified: branches/ocean_projects/performance/src/core_ocean/module_time_integration.F
===================================================================
--- branches/ocean_projects/performance/src/core_ocean/module_time_integration.F 2011-09-16 21:57:58 UTC (rev 1007)
+++ branches/ocean_projects/performance/src/core_ocean/module_time_integration.F 2011-09-19 17:14:29 UTC (rev 1008)
@@ -19,6 +19,8 @@
use OcnTracerHadv
use OcnTracerVadv
+ use OcnTracerHmix
+ use OcnRestoring
contains
@@ -2326,13 +2328,6 @@
deriv_two => grid % deriv_two % array
- if(config_restoreTS) then
- index_temperature = s % index_temperature
- index_salinity = s % index_salinity
- temperatureRestore => grid % temperatureRestore % array
- salinityRestore => grid % salinityRestore % array
- endif
-
!
! initialize tracer tendency (RHS of tracer equation) to zero.
!
@@ -2367,118 +2362,9 @@
! tracer tendency: del2 horizontal tracer diffusion, div(h \kappa_2 </font>
<font color="red">abla \phi)
!
call timer_start("compute_scalar_tend-horiz diff")
- if ( config_h_tracer_eddy_diff2 > 0.0 ) then
- call timer_start("compute_scalar_tend-horiz diff 2")
- !
- ! compute a boundary mask to enforce insulating boundary conditions in the horizontal
- !
- allocate(boundaryMask(nVertLevels, nEdges+1))
- boundaryMask = 1.0
- where(boundaryEdge.eq.1) boundaryMask=0.0
+ call OcnTracerHmixTend(grid, h_edge, tracers, tend_tr, err)
- do iEdge=1,grid % nEdges
- cell1 = grid % cellsOnEdge % array(1,iEdge)
- cell2 = grid % cellsOnEdge % array(2,iEdge)
- invAreaCell1 = 1.0/areaCell(cell1)
- invAreaCell2 = 1.0/areaCell(cell2)
-
- do k=1,maxLevelEdgeTop(iEdge)
- do iTracer=1,num_tracers
- ! \kappa_2 </font>
<font color="red">abla \phi on edge
- tracer_turb_flux = meshScalingDel2(iEdge) * config_h_tracer_eddy_diff2 &
- *( tracers(iTracer,k,cell2) &
- - tracers(iTracer,k,cell1))/dcEdge(iEdge)
-
- ! div(h \kappa_2 </font>
<font color="red">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
- end do
-
- end do
-
- deallocate(boundaryMask)
- call timer_stop("compute_scalar_tend-horiz diff 2")
-
- end if
-
- !
- ! tracer tendency: del4 horizontal tracer diffusion, &
- ! div(h \kappa_4 </font>
<font color="black">abla [div(h </font>
<font color="red">abla \phi)])
- !
- if ( config_h_tracer_eddy_diff4 > 0.0 ) then
- call timer_start("compute_scalar_tend-horiz diff 4")
-
- !
- ! compute a boundary mask to enforce insulating boundary conditions in the horizontal
- !
- allocate(boundaryMask(nVertLevels, nEdges+1))
- boundaryMask = 1.0
- where(boundaryEdge.eq.1) boundaryMask=0.0
-
- allocate(delsq_tracer(num_tracers,nVertLevels, nCells+1))
-
- delsq_tracer(:,:,:) = 0.
-
- ! first del2: div(h </font>
<font color="red">abla \phi) at cell center
- do iEdge=1,grid % nEdges
- cell1 = grid % cellsOnEdge % array(1,iEdge)
- cell2 = grid % cellsOnEdge % array(2,iEdge)
-
- do k=1,maxLevelEdgeTop(iEdge)
- 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)
- 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)
- end do
- end do
-
- end do
-
- do iCell = 1,nCells
- r = 1.0 / areaCell(iCell)
- do k=1,maxLevelCell(iCell)
- do iTracer=1,num_tracers
- delsq_tracer(iTracer,k,iCell) = delsq_tracer(iTracer,k,iCell) * r
- end do
- 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)
- invAreaCell1 = 1.0 / areaCell(cell1)
- invAreaCell2 = 1.0 / areaCell(cell2)
-
- do k=1,maxLevelEdgeTop(iEdge)
- do iTracer=1,num_tracers
- tracer_turb_flux = meshScalingDel4(iEdge) * config_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)
-
- enddo
- enddo
- end do
-
- deallocate(delsq_tracer)
- call timer_stop("compute_scalar_tend-horiz diff 4")
-
- end if
call timer_stop("compute_scalar_tend-horiz diff")
! mrp 110516 printing
@@ -2533,28 +2419,12 @@
!
! add restoring to T and S in top model layer
!
- if(config_restoreTS) then
- call timer_start("compute_scalar_tend-restoring")
- k = 1 ! restoring only in top layer
- do iCell=1,nCellsSolve
+ call timer_start("compute_scalar_tend-restoring")
- tend_tr(index_temperature, k, iCell) = tend_tr(index_temperature, k, iCell) &
- - h(k,iCell)*(tracers(index_temperature, k, iCell) - temperatureRestore(iCell)) &
- / (config_restoreT_timescale * 86400.0)
+ call OcnRestoringTend(grid, h, s%index_temperature, s%index_salinity, tracers, tend_tr, err)
- tend_tr(index_salinity, k, iCell) = tend_tr(index_salinity, k, iCell) &
- - h(k,iCell)*(tracers(index_salinity, k, iCell) - salinityRestore(iCell)) &
- / (config_restoreS_timescale * 86400.0)
+ call timer_stop("compute_scalar_tend-restoring")
- ! write(6,10) iCell, tracers(index_temperature, k, iCell), &
- ! temperatureRestore(iCell), tracers(index_temperature, k, iCell), &
- ! (tracers(index_temperature, k, iCell) - temperatureRestore(iCell)) &
- ! / (config_restoreT_timescale * 86400.0)
-
- enddo
-
- call timer_stop("compute_scalar_tend-restoring")
- endif
10 format(2i8,10e20.10)
call timer_stop("compute_scalar_tend")
</font>
</pre>