<p><b>pwjones@lanl.gov</b> 2011-08-01 12:18:21 -0600 (Mon, 01 Aug 2011)</p><p><br>
new module structure for ocean horizontal tracer mixing<br>
- improves performance<br>
- adds modularity for future additions<br>
- adds encapsulation for these modules<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/performance/mpas/namelist.input.ocean
===================================================================
--- branches/ocean_projects/performance/mpas/namelist.input.ocean        2011-07-28 22:46:56 UTC (rev 932)
+++ branches/ocean_projects/performance/mpas/namelist.input.ocean        2011-08-01 18:18:21 UTC (rev 933)
@@ -26,8 +26,8 @@
&hmix
config_h_mom_eddy_visc2 = 0.0
config_h_mom_eddy_visc4 = 5.0e8
- config_h_tracer_eddy_diff2 = 10.0
- config_h_tracer_eddy_diff4 = 0.0
+ config_OcnHmixTracerDel2Diff = 10.0
+ config_OcnHmixTracerDel4Diff = 0.0
/
&vmix
config_vert_visc_type = 'const'
Modified: branches/ocean_projects/performance/mpas/namelist.input.sw
===================================================================
--- branches/ocean_projects/performance/mpas/namelist.input.sw        2011-07-28 22:46:56 UTC (rev 932)
+++ branches/ocean_projects/performance/mpas/namelist.input.sw        2011-08-01 18:18:21 UTC (rev 933)
@@ -1,31 +1,54 @@
&sw_model
- config_test_case = 5
+ config_test_case = 0
config_time_integration = 'RK4'
- config_dt = 172.8
- config_ntimesteps = 7500
- config_output_interval = 500
- config_stats_interval = 0
- config_h_ScaleWithMesh = .false.
- config_h_mom_eddy_visc2 = 0.0
- config_h_mom_eddy_visc4 = 0.0
- config_h_tracer_eddy_diff2 = 0.0
- config_h_tracer_eddy_diff4 = 0.0
- config_thickness_adv_order = 2
- config_tracer_adv_order = 2
- config_positive_definite = .false.
- config_monotonic = .false.
- config_wind_stress = .false.
- config_bottom_drag = .false.
+ config_dt = 120.0
+ config_ntimesteps = 20
+ config_output_interval = 7200
+ config_stats_interval = 20
/
&io
- config_input_name = 'grid.nc'
- config_output_name = 'output.nc'
- config_restart_name = 'restart.nc'
+ config_input_name = grid.nc
+ config_output_name = output.nc
+ config_restart_name = restart.nc
/
&restart
- config_restart_interval = 3000
+ config_restart_interval = 10000000
config_do_restart = .false.
config_restart_time = 1036800.0
/
+
+&grid
+ config_vert_grid_type = 'zlevel'
+ config_rho0 = 1000
+/
+&hmix
+ config_h_mom_eddy_visc2 = 1.0e5
+ config_h_mom_eddy_visc4 = 0.0
+ config_OcnHmixTracerDel2Diff = 1.0e4
+ config_OcnHmixTracerDel4Diff = 0.0
+/
+&vmix
+ config_vert_visc_type = 'const'
+ config_vert_diff_type = 'const'
+ config_vmixTanhViscMax = 2.5e-1
+ config_vmixTanhViscMin = 1.0e-4
+ config_vmixTanhDiffMax = 2.5e-2
+ config_vmixTanhDiffMin = 1.0e-5
+ config_vmixTanhZMid = -100
+ config_vmixTanhZWidth = 100
+ config_vert_viscosity = 2.5e-5
+ config_vert_diffusion = 2.5e-5
+/
+&eos
+ config_eos_type = 'jm'
+/
+&advection
+ config_vert_tracer_adv = 'spline'
+ config_vert_tracer_adv_order = 3
+ config_tracer_adv_order = 2
+ config_thickness_adv_order = 2
+ config_positive_definite = .false.
+ config_monotonic = .false.
+/
Modified: branches/ocean_projects/performance/mpas/src/core_ocean/Makefile
===================================================================
--- branches/ocean_projects/performance/mpas/src/core_ocean/Makefile        2011-07-28 22:46:56 UTC (rev 932)
+++ branches/ocean_projects/performance/mpas/src/core_ocean/Makefile        2011-08-01 18:18:21 UTC (rev 933)
@@ -3,6 +3,9 @@
OBJS = module_mpas_core.o \
module_test_cases.o \
module_advection.o \
+ module_OcnHmixTracer.o \
+ module_OcnHmixTracerDel2.o \
+ module_OcnHmixTracerDel4.o \
module_time_integration.o \
module_global_diagnostics.o
@@ -15,11 +18,17 @@
module_advection.o:
+module_OcnHmixTracerDel2.o:
+
+module_OcnHmixTracerDel4.o:
+
+module_OcnHmixTracer.o: module_OcnHmixTracerDel2.o module_OcnHmixTracerDel4.o
+
module_time_integration.o:
module_global_diagnostics.o:
-module_mpas_core.o: module_advection.o module_global_diagnostics.o module_test_cases.o module_time_integration.o
+module_mpas_core.o: module_advection.o module_global_diagnostics.o module_test_cases.o module_OcnHmixTracer.o module_OcnHmixTracerDel2.o module_OcnHmixTracerDel4.o module_time_integration.o
clean:
        $(RM) *.o *.mod *.f90 libdycore.a
Modified: branches/ocean_projects/performance/mpas/src/core_ocean/Registry
===================================================================
--- branches/ocean_projects/performance/mpas/src/core_ocean/Registry        2011-07-28 22:46:56 UTC (rev 932)
+++ branches/ocean_projects/performance/mpas/src/core_ocean/Registry        2011-08-01 18:18:21 UTC (rev 933)
@@ -18,8 +18,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_diff2 0.0
-namelist real hmix config_h_tracer_eddy_diff4 0.0
+namelist real hmix config_OcnHmixTracerDel2Diff 0.0
+namelist real hmix config_OcnHmixTracerDel4Diff 0.0
namelist real hmix config_apvm_upwinding 0.5
namelist character vmix config_vert_visc_type const
namelist character vmix config_vert_diff_type const
Added: branches/ocean_projects/performance/mpas/src/core_ocean/module_OcnHmixTracer.F
===================================================================
--- branches/ocean_projects/performance/mpas/src/core_ocean/module_OcnHmixTracer.F         (rev 0)
+++ branches/ocean_projects/performance/mpas/src/core_ocean/module_OcnHmixTracer.F        2011-08-01 18:18:21 UTC (rev 933)
@@ -0,0 +1,164 @@
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
+!
+! OcnHmixTracer
+!
+!> \brief MPAS ocean horizontal tracer mixing driver
+!> \author Phil Jones
+!> \date 26 July 2011
+!> \version SVN:$Id:$
+!> \details
+!> This module contains the main driver routine for computing
+!> horizontal mixing tendencies. It primarily decides which mixing
+!> parameterizations are being used and calls them individually.
+!> Detailed mixing parameterizations are contained in their own
+!> modules.
+!
+!-----------------------------------------------------------------------
+
+module OcnHmixTracer
+
+ use grid_types
+ use configure
+ use OcnHmixTracerDel2
+ use OcnHmixTracerDel4
+
+ implicit none
+ private
+ save
+
+ !--------------------------------------------------------------------
+ !
+ ! Public parameters
+ !
+ !--------------------------------------------------------------------
+
+ !--------------------------------------------------------------------
+ !
+ ! Public member functions
+ !
+ !--------------------------------------------------------------------
+
+ public :: OcnHmixTracerTend, &
+ OcnHmixTracerInit
+
+ !--------------------------------------------------------------------
+ !
+ ! Private module variables
+ !
+ !--------------------------------------------------------------------
+
+
+!***********************************************************************
+
+contains
+
+!***********************************************************************
+!
+! routine OcnHmixTracerTend
+!
+!> \brief Computes tendency term for horizontal tracer mixing
+!> \author Phil Jones
+!> \date 26 July 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine computes the horizontal mixing tendency for tracers
+!> 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 OcnHmixTracerTend(tend, tracers, hEdge, grid)
+
+ !-----------------------------------------------------------------
+ !
+ ! input variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:,:), intent(in) :: &
+ tracers !< Input: array of tracers
+
+ real (kind=RKIND), dimension(:,:), intent(in) :: &
+ hEdge !< Input: layer thickness at cell edges
+
+ type (mesh_type), intent(in) :: &
+ grid !< Input: grid information
+
+ !-----------------------------------------------------------------
+ !
+ ! input/output variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:,:), intent(inout) :: &
+ tend !< Input/Output: tracer tendency
+
+ !-----------------------------------------------------------------
+ !
+ ! output variables
+ !
+ !-----------------------------------------------------------------
+
+ !-----------------------------------------------------------------
+ !
+ ! local variables
+ !
+ !-----------------------------------------------------------------
+
+ !-----------------------------------------------------------------
+ !
+ ! call relevant routines for computing tendencies
+ ! note that the user can choose multiple options and the
+ ! tendencies will be added together
+ !
+ !-----------------------------------------------------------------
+
+ call OcnHmixTracerDel2Tend(tend, tracers, hEdge, grid)
+ call OcnHmixTracerDel4Tend(tend, tracers, hEdge, grid)
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnHmixTracerTend
+
+!***********************************************************************
+!
+! routine OcnHmixTracerInit
+!
+!> \brief Initializes ocean tracer horizontal mixing quantities
+!> \author Phil Jones
+!> \date 26 July 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine initializes a variety of quantities related to
+!> horizontal tracer mixing in the ocean. Since a variety of
+!> parameterizations are available, this routine primarily calls the
+!> individual init routines for each parameterization.
+!
+!-----------------------------------------------------------------------
+
+
+ subroutine OcnHmixTracerInit
+
+ !--------------------------------------------------------------------
+
+ !-----------------------------------------------------------------
+ !
+ ! call individual init routines for each parameterization
+ !
+ !-----------------------------------------------------------------
+
+ call OcnHmixTracerDel2Init
+ call OcnHmixTracerDel4Init
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnHmixTracerInit
+
+!***********************************************************************
+
+end module OcnHmixTracer
+
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Added: branches/ocean_projects/performance/mpas/src/core_ocean/module_OcnHmixTracerDel2.F
===================================================================
--- branches/ocean_projects/performance/mpas/src/core_ocean/module_OcnHmixTracerDel2.F         (rev 0)
+++ branches/ocean_projects/performance/mpas/src/core_ocean/module_OcnHmixTracerDel2.F        2011-08-01 18:18:21 UTC (rev 933)
@@ -0,0 +1,249 @@
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
+!
+! OcnHMixTracerDel2
+!
+!> \brief Ocean horizontal mixing - Laplacian parameterization
+!> \author Phil Jones
+!> \date 1 August 2011
+!> \version SVN:$Id:$
+!> \details
+!> This module contains routines for computing horizontal mixing
+!> tendencies using a Laplacian formulation.
+!
+!-----------------------------------------------------------------------
+
+module OcnHMixTracerDel2
+
+ use grid_types
+ use configure
+
+ implicit none
+ private
+ save
+
+ !--------------------------------------------------------------------
+ !
+ ! Public parameters
+ !
+ !--------------------------------------------------------------------
+
+ !--------------------------------------------------------------------
+ !
+ ! Public member functions
+ !
+ !--------------------------------------------------------------------
+
+ public :: OcnHMixTracerDel2Tend, &
+ OcnHMixTracerDel2Init
+
+ !-------------------------------------------------------------------
+ !
+ ! Private module variables
+ !
+ !--------------------------------------------------------------------
+
+ logical :: &
+ hmixDel2On !< local flag to determine whether del2 chosen
+
+ real (kind=RKIND) :: &
+ eddyDiff2 !< base eddy diffusivity for Laplacian
+
+!***********************************************************************
+
+contains
+
+!***********************************************************************
+!
+! routine OcnHMixTracerDel2Tend
+!
+!> \brief Computes tendency term for Laplacian horizontal tracer mixing
+!> \author Phil Jones
+!> \date 1 August 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine computes the horizontal mixing tendency for tracers
+!> based on a Laplacian form for the mixing. This tendency takes the
+!> form </font>
<font color="black">abla \cdot (h \kappa_2 </font>
<font color="blue">abla \phi), where \phi is the
+!> tracer, h is the thickness and \kappa_2 is a base horizontal
+!> diffusivity.
+!
+!-----------------------------------------------------------------------
+
+ subroutine OcnHMixTracerDel2Tend(tend, tracers, hEdge, grid)
+
+ !-----------------------------------------------------------------
+ !
+ ! input variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:,:), intent(in) :: &
+ tracers !< Input: array of tracers
+
+ real (kind=RKIND), dimension(:,:), intent(in) :: &
+ hEdge !< Input: layer thickness at cell edges
+
+ type (mesh_type), intent(in) :: &
+ grid !< Input: grid information
+
+ !-----------------------------------------------------------------
+ !
+ ! input/output variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:,:), intent(inout) :: &
+ tend !< Input/Output: tracer tendency
+
+ !-----------------------------------------------------------------
+ !
+ ! local variables
+ !
+ !-----------------------------------------------------------------
+
+ integer :: &
+ i, k, iCell, iEdge, iTracer, &! loop counters
+ cell1, cell2, &! cell indices
+ nEdges, nCells, nVertLevels, nTracers ! array sizes
+
+ real (kind=RKIND) :: &
+ invAreaCell1, invAreaCell2, &! reciprocal of cell areas
+ flux ! flux across edge
+
+ real (kind=RKIND), dimension(:), pointer :: &
+ dvEdge, dcEdge, areaCell
+
+ integer, dimension(:,:), pointer :: boundaryEdge
+
+ integer, dimension(:), pointer :: maxLevelEdgeTop
+
+ integer, dimension(:,:), pointer :: cellsOnEdge
+
+ real (kind=RKIND), dimension(:), allocatable:: boundaryMask
+
+ real (kind=RKIND), dimension(:,:,:), allocatable:: del2tracer
+
+ !-----------------------------------------------------------------
+ !
+ ! exit if this mixing is not selected
+ !
+ !-----------------------------------------------------------------
+
+ if (.not. hmixDel2On) return
+
+ !-----------------------------------------------------------------
+ !
+ ! initialize some arrays and sizes
+ !
+ !-----------------------------------------------------------------
+
+ areaCell => grid % areaCell % array
+ cellsOnEdge => grid % cellsOnEdge % array
+ dvEdge => grid % dvEdge % array
+ dcEdge => grid % dcEdge % array
+ boundaryEdge => grid % boundaryEdge % array
+ maxLevelEdgeTop => grid % maxLevelEdgeTop % array
+
+ nEdges = grid % nEdges
+ nCells = grid % nCells
+ nVertLevels = size(tracers,dim=2)
+ nTracers = size(tracers,dim=1)
+
+ !-----------------------------------------------------------------
+ !
+ ! del2 horizontal tracer diffusion
+ ! tracer tendency: div(h \kappa_2 </font>
<font color="blue">abla \phi)
+ !
+ !-----------------------------------------------------------------
+
+ allocate(boundaryMask(nVertLevels))
+
+ do iEdge=1,nEdges
+
+ !*** get cell pairs that share edge
+
+ cell1 = grid % cellsOnEdge % array(1,iEdge)
+ cell2 = grid % cellsOnEdge % array(2,iEdge)
+ invAreaCell1 = 1.d0/areaCell(cell1)
+ invAreaCell2 = 1.d0/areaCell(cell2)
+
+ !*** compute a boundary mask at each vertical level
+ !*** on this edge
+
+ do k=1,maxLevelEdgeTop(iEdge)
+ if (boundaryEdge(k,iEdge) == 1) then
+ boundaryMask(k) = 0.d0
+ else
+ boundaryMask(k) = 1.d0
+ endif
+ end do
+
+ !*** compute \kappa_2 </font>
<font color="blue">abla \phi on edge
+
+ do k=1,maxLevelEdgeTop(iEdge)
+ do iTracer=1,nTracers
+
+ ! div(h \kappa_2 </font>
<font color="gray">abla \phi) at cell center
+
+ flux = dvEdge(iEdge) * hEdge(k,iEdge) * eddyDiff2 * &
+ ( tracers(iTracer,k,cell2) &
+ - tracers(iTracer,k,cell1))/dcEdge(iEdge)* &
+ boundaryMask(k)
+
+ 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)
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnHMixTracerDel2Tend
+
+!***********************************************************************
+!
+! routine OcnHMixTracerDel2Init
+!
+!> \brief Initializes ocean tracer Laplacian horizontal mixing
+!> \author Phil Jones
+!> \date 26 July 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine initializes a variety of quantities related to
+!> Laplacian horizontal tracer mixing in the ocean.
+!
+!-----------------------------------------------------------------------
+
+ subroutine OcnHMixTracerDel2Init
+
+ !--------------------------------------------------------------------
+ !
+ ! set some local module variables based on input config choices
+ !
+ !--------------------------------------------------------------------
+
+ hmixDel2On = .false.
+
+ if (config_OcnHmixTracerDel2Diff > 0.0) then
+
+ hmixDel2On = .true.
+
+ eddyDiff2 = config_OcnHmixTracerDel2Diff
+
+ endif
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnHMixTracerDel2Init
+
+!***********************************************************************
+
+end module OcnHMixTracerDel2
+
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Added: branches/ocean_projects/performance/mpas/src/core_ocean/module_OcnHmixTracerDel4.F
===================================================================
--- branches/ocean_projects/performance/mpas/src/core_ocean/module_OcnHmixTracerDel4.F         (rev 0)
+++ branches/ocean_projects/performance/mpas/src/core_ocean/module_OcnHmixTracerDel4.F        2011-08-01 18:18:21 UTC (rev 933)
@@ -0,0 +1,305 @@
+!|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
+!
+! OcnHMixTracerDel4
+!
+!> \brief Ocean horizontal mixing - biharmonic parameterization
+!> \author Phil Jones
+!> \date 1 August 2011
+!> \version SVN:$Id:$
+!> \details
+!> This module contains routines and variables for computing
+!> horizontal mixing tendencies using a biharmonic formulation.
+!
+!-----------------------------------------------------------------------
+
+module OcnHMixTracerDel4
+
+ use grid_types
+ use configure
+
+ implicit none
+ private
+ save
+
+ !--------------------------------------------------------------------
+ !
+ ! Public parameters
+ !
+ !--------------------------------------------------------------------
+
+ !--------------------------------------------------------------------
+ !
+ ! Public member functions
+ !
+ !--------------------------------------------------------------------
+
+ public :: OcnHMixTracerDel4Tend, &
+ OcnHMixTracerDel4Init
+
+ !--------------------------------------------------------------------
+ !
+ ! Private module variables
+ !
+ !--------------------------------------------------------------------
+
+ logical :: &
+ hmixDel4On !< local flag to determine whether del4 chosen
+
+ real (kind=RKIND) :: &
+ eddyDiff4 !< base eddy diffusivity for biharmonic
+
+
+!***********************************************************************
+
+contains
+
+!***********************************************************************
+!
+! routine OcnHMixTracerDel4Tend
+!
+!> \brief Computes tendency term for biharmonic horizontal tracer mixing
+!> \author Phil Jones
+!> \date 1 August 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine computes the horizontal mixing tendency for tracers
+!> based on a biharmonic form for the mixing. This mixing tendency
+!> takes the form:
+!> </font>
<font color="black">abla \cdot (h \kappa_4 </font>
<font color="black">abla [</font>
<font color="black">abla \cdot (h </font>
<font color="blue">abla \phi)])
+!> where \phi is a tracer, h is the layer thickness and \kappa_4 is
+!> a base diffusivity.
+!
+!-----------------------------------------------------------------------
+
+ subroutine OcnHMixTracerDel4Tend(tend, tracers, hEdge, grid)
+
+ !-----------------------------------------------------------------
+ !
+ ! input variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:,:), intent(in) :: &
+ tracers !< Input: array of tracers
+
+ real (kind=RKIND), dimension(:,:), intent(in) :: &
+ hEdge !< Input: layer thickness at cell edges
+
+ type (mesh_type), intent(in) :: &
+ grid !< Input: grid information
+
+ !-----------------------------------------------------------------
+ !
+ ! input/output variables
+ !
+ !-----------------------------------------------------------------
+
+ real (kind=RKIND), dimension(:,:,:), intent(inout) :: &
+ tend !< Input/Output: tracer tendency
+
+ !-----------------------------------------------------------------
+ !
+ ! local variables
+ !
+ !-----------------------------------------------------------------
+
+ integer :: &
+ i, k, iCell, iEdge, iTracer, &! loop counters
+ cell1, cell2, &! cell indices
+ nEdges, nCells, nVertLevels, nTracers ! array sizes
+
+ real (kind=RKIND) :: &
+ invAreaCell1, invAreaCell2, &! reciprocal of cell areas
+ flux ! flux across edge
+
+ real (kind=RKIND), dimension(:), pointer :: &
+ dvEdge, dcEdge, areaCell
+
+ integer, dimension(:,:), pointer :: boundaryEdge
+
+ integer, dimension(:), pointer :: maxLevelCell, maxLevelEdgeTop
+
+ integer, dimension(:,:), pointer :: cellsOnEdge
+
+ real (kind=RKIND), dimension(:), allocatable:: boundaryMask
+
+ real (kind=RKIND), dimension(:,:,:), allocatable:: del2tracer
+
+ !-----------------------------------------------------------------
+ !
+ ! exit if this mixing is not selected
+ !
+ !-----------------------------------------------------------------
+
+ if (.not. hmixDel4On) return
+
+ !-----------------------------------------------------------------
+ !
+ ! initialize some arrays and sizes
+ !
+ !-----------------------------------------------------------------
+
+ areaCell => grid % areaCell % array
+ cellsOnEdge => grid % cellsOnEdge % array
+ dvEdge => grid % dvEdge % array
+ dcEdge => grid % dcEdge % array
+ boundaryEdge => grid % boundaryEdge % array
+ maxLevelCell => grid % maxLevelCell % array
+ maxLevelEdgeTop => grid % maxLevelEdgeTop % array
+
+ nEdges = grid % nEdges
+ nCells = grid % nCells
+ nVertLevels = size(tracers,dim=2)
+ nTracers = size(tracers,dim=1)
+
+ !-----------------------------------------------------------------
+ !
+ ! del4 horizontal tracer diffusion
+ ! tracer tendency: div(h \kappa_4 </font>
<font color="black">abla [div(h </font>
<font color="blue">abla \phi)])
+ ! this is computed with two applications of the Laplacian
+ !
+ !-----------------------------------------------------------------
+
+ allocate(boundaryMask(nVertLevels))
+ allocate(del2tracer(nTracers,nVertLevels,nCells+1))
+
+ del2tracer(:,:,:) = 0.d0
+
+ !*** first application of Laplacian
+ !*** del2: div(h </font>
<font color="blue">abla \phi) at cell center
+
+ do iEdge=1,nEdges
+
+ !*** get cell pair that share edge
+
+ cell1 = grid % cellsOnEdge % array(1,iEdge)
+ cell2 = grid % cellsOnEdge % array(2,iEdge)
+
+ !*** compute a boundary mask at each vertical level
+ !*** on this edge
+
+ do k=1,maxLevelEdgeTop(iEdge)
+ if (boundaryEdge(k,iEdge) == 1) then
+ boundaryMask(k) = 0.d0
+ else
+ boundaryMask(k) = 1.d0
+ endif
+ end do
+
+ do k=1,maxLevelEdgeTop(iEdge)
+ do iTracer=1,nTracers
+
+ flux = dvEdge(iEdge)*hEdge(k,iEdge) * &
+ (tracers(iTracer,k,cell2) - &
+ tracers(iTracer,k,cell1))/dcEdge(iEdge)* &
+ boundaryMask(k)
+
+ del2tracer(iTracer,k,cell1) = &
+ del2tracer(iTracer,k,cell1) + flux
+
+ del2tracer(iTracer,k,cell2) = &
+ del2tracer(iTracer,k,cell2) - flux
+
+ end do
+ end do
+
+ end do
+
+ do iCell = 1,nCells
+ invAreaCell1 = 1.d0/areaCell(iCell)
+ do k=1,maxLevelCell(iCell)
+ do iTracer=1,nTracers
+ del2tracer(iTracer,k,iCell) = &
+ del2tracer(iTracer,k,iCell) * invAreaCell1
+ end do
+ end do
+ end do
+
+ !*** second application of Laplacian
+ !*** del2: div(h </font>
<font color="gray">abla [delsq_tracer]) at cell center
+
+ do iEdge=1,nEdges
+
+ !*** get cell pair that share edge
+
+ cell1 = grid % cellsOnEdge % array(1,iEdge)
+ cell2 = grid % cellsOnEdge % array(2,iEdge)
+ invAreaCell1 = 1.0 / areaCell(cell1)
+ invAreaCell2 = 1.0 / areaCell(cell2)
+
+ !*** compute a boundary mask at each vertical level
+ !*** on this edge
+
+ do k=1,maxLevelEdgeTop(iEdge)
+ if (boundaryEdge(k,iEdge) == 1) then
+ boundaryMask(k) = 0.d0
+ else
+ boundaryMask(k) = 1.d0
+ endif
+ end do
+
+ do k=1,maxLevelEdgeTop(iEdge)
+ do iTracer=1,nTracers
+
+ flux = dvEdge(iEdge)* eddyDiff4 * &
+ (del2tracer(iTracer,k,cell2) - &
+ del2tracer(iTracer,k,cell1))/dcEdge(iEdge) * &
+ boundaryMask(k)
+
+ tend(iTracer,k,cell1) = &
+ tend(iTracer,k,cell1) - flux * invAreaCell1
+ tend(iTracer,k,cell2) = &
+ tend(iTracer,k,cell2) + flux * invAreaCell2
+
+ enddo
+ enddo
+ end do
+
+ deallocate(del2tracer, boundaryMask)
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnHMixTracerDel4Tend
+
+!***********************************************************************
+!
+! routine OcnHMixTracerDel4Init
+!
+!> \brief Initializes ocean tracer biharmonic horizontal mixing
+!> \author Phil Jones
+!> \date 1 August 2011
+!> \version SVN:$Id$
+!> \details
+!> This routine initializes a variety of quantities related to
+!> biharmonic horizontal tracer mixing in the ocean.
+!
+!-----------------------------------------------------------------------
+
+
+ subroutine OcnHMixTracerDel4Init
+
+ !--------------------------------------------------------------------
+ !
+ ! set some local module variables based on input config choices
+ !
+ !--------------------------------------------------------------------
+
+ hmixDel4On = .false.
+
+ if (config_OcnHmixTracerDel4Diff /= 0.0) then
+
+ hmixDel4On = .true.
+
+ eddyDiff4 = config_OcnHmixTracerDel4Diff
+
+ endif
+
+ !--------------------------------------------------------------------
+
+ end subroutine OcnHMixTracerDel4Init
+
+!***********************************************************************
+
+end module OcnHMixTracerDel4
+
+!||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Modified: branches/ocean_projects/performance/mpas/src/core_ocean/module_mpas_core.F
===================================================================
--- branches/ocean_projects/performance/mpas/src/core_ocean/module_mpas_core.F        2011-07-28 22:46:56 UTC (rev 932)
+++ branches/ocean_projects/performance/mpas/src/core_ocean/module_mpas_core.F        2011-08-01 18:18:21 UTC (rev 933)
@@ -3,6 +3,7 @@
use mpas_framework
use dmpar
use test_cases
+ use OcnHmixTracer
type (io_output_object) :: restart_obj
integer :: restart_frame
@@ -49,6 +50,12 @@
block => block % next
end do
+ !***
+ !*** initialize other core modules
+ !***
+
+ call OcnHmixTracerInit
+
! 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
! no initial statistics write.
Modified: branches/ocean_projects/performance/mpas/src/core_ocean/module_time_integration.F
===================================================================
--- branches/ocean_projects/performance/mpas/src/core_ocean/module_time_integration.F        2011-07-28 22:46:56 UTC (rev 932)
+++ branches/ocean_projects/performance/mpas/src/core_ocean/module_time_integration.F        2011-08-01 18:18:21 UTC (rev 933)
@@ -1,3 +1,5 @@
+
+
module time_integration
use grid_types
@@ -6,6 +8,8 @@
use dmpar
use vector_reconstruction
use spline_interpolation
+ use timer
+ use OcnHmixTracer
contains
@@ -344,6 +348,7 @@
u_src => grid % u_src % array
+ call timer_start("height")
!
! height tendency: start accumulating tendency terms
!
@@ -401,11 +406,14 @@
end do
endif ! coordinate type
+ call timer_stop("height")
+ call timer_start("velocity")
!
! velocity tendency: start accumulating tendency terms
!
tend_u(:,:) = 0.0
+ call timer_start("vel vert advect")
!
! velocity tendency: vertical advection term -w du/dz
!
@@ -434,6 +442,8 @@
deallocate(w_dudzTopEdge)
endif
+ call timer_stop("vel vert advect")
+ call timer_start("press grad")
!
! velocity tendency: pressure gradient
!
@@ -459,6 +469,8 @@
enddo
endif
+ call timer_stop("press grad")
+ call timer_start("vel dissipation")
!
! velocity tendency: del2 dissipation, </font>
<font color="black">u_2 </font>
<font color="black">abla^2 u
! computed as </font>
<font color="black">u( </font>
<font color="black">abla divergence + k \times </font>
<font color="gray">abla vorticity )
@@ -582,6 +594,8 @@
end if
+ call timer_stop("vel dissipation")
+ call timer_start("coriolis")
!
! velocity tendency: nonlinear Coriolis term and grad of kinetic energy
!
@@ -604,6 +618,8 @@
end do
end do
+ call timer_stop("coriolis")
+ call timer_start("vel forcing")
!
! velocity tendency: forcing and bottom drag
!
@@ -637,6 +653,8 @@
enddo
+ call timer_stop("vel forcing")
+ call timer_start("vel vmix")
!
! velocity tendency: vertical mixing d/dz( nu_v du/dz))
!
@@ -679,6 +697,8 @@
enddo
end do
+ call timer_stop("vel vmix")
+ call timer_stop("velocity")
deallocate(fluxVertTop, vertViscTop)
end subroutine compute_tend
@@ -727,6 +747,16 @@
real (kind=RKIND), dimension(:,:,:), pointer :: deriv_two
real (kind=RKIND) :: coef_3rd_order, flux3Coef, cSignWTop
+ real (kind=RKIND), dimension(:),allocatable :: &
+ boundaryMaskA, &
+ dvEdgeA, dcEdgeA
+
+ real (kind=RKIND), dimension(:,:),allocatable :: &
+ h_edgeA
+
+ real (kind=RKIND), dimension(:,:,:),allocatable :: &
+ tracersA, tend_trA
+
u => s % u % array
h => s % h % array
boundaryCell=> grid % boundaryCell % array
@@ -757,11 +787,13 @@
deriv_two => grid % deriv_two % array
+ call timer_start("scalar")
!
! initialize tracer tendency (RHS of tracer equation) to zero.
!
tend_tr(:,:,:) = 0.0
+ call timer_start("scalar horz adv")
!
! tracer tendency: horizontal advection term -div( h \phi u)
!
@@ -889,6 +921,8 @@
endif ! if (config_tracer_adv_order == 2 )
+ call timer_stop("scalar horz adv")
+ call timer_start("scalar vert adv")
!
! tracer tendency: vertical advection term -d/dz( h \phi w)
@@ -1070,121 +1104,21 @@
endif ! ZLevel
- !
- ! tracer tendency: del2 horizontal tracer diffusion, div(h \kappa_2 </font>
<font color="red">abla \phi)
- !
- if ( config_h_tracer_eddy_diff2 > 0.0 ) then
+ call timer_stop("scalar vert adv")
- !
- ! 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
+ !***
+ !*** tracer tendency: horizontal tracer diffusion
+ !***
- 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)
+ call timer_start("scalar hdiff")
+ call OcnHmixTracerTend(tend_tr, tracers, h_edge, grid)
+ call timer_stop("scalar hdiff")
- do k=1,maxLevelEdgeTop(iEdge)
- do iTracer=1,num_tracers
- ! \kappa_2 </font>
<font color="red">abla \phi on edge
- tracer_turb_flux = 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)
-
- 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
-
- !
- ! 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 = 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)
-
- end if
-
- !
! tracer tendency: vertical diffusion h d/dz( \kappa_v d\phi/dz)
!
+ call timer_start("scalar vdiff")
allocate(vertDiffTop(nVertLevels+1))
if (config_vert_diff_type.eq.'const') then
vertDiffTop = config_vert_diffusion
@@ -1232,6 +1166,8 @@
enddo ! iCell loop
deallocate(fluxVertTop, vertDiffTop)
+ call timer_stop("scalar vdiff")
+ call timer_stop("scalar")
! print some diagnostics - for debugging
! print *, 'after vertical mixing',&
! 'iTracer,k, minval(tend_tr(itracer,k,:)), maxval(tend_tr(itracer,k,:))'
</font>
</pre>