<p><b>dwj07@fsu.edu</b> 2010-09-02 11:14:50 -0600 (Thu, 02 Sep 2010)</p><p><br>
        namelist.input.sw: Added config_stats_interval<br>
<br>
        module_global_diagnostics.F: modified variables for more readable names.<br>
                                                                 modified and checked for accurate computations.<br>
<br>
        mpas_interface.F: config_stats_interval = 0 causes no global diagnostics.<br>
</p><hr noshade><pre><font color="gray">Modified: trunk/mpas/namelist.input.sw
===================================================================
--- trunk/mpas/namelist.input.sw        2010-09-01 23:21:01 UTC (rev 492)
+++ trunk/mpas/namelist.input.sw        2010-09-02 17:14:50 UTC (rev 493)
@@ -4,6 +4,7 @@
config_dt = 172.8
config_ntimesteps = 7500
config_output_interval = 500
+ config_stats_interval = 0
config_visc = 0.0
/
Modified: trunk/mpas/src/core_sw/module_global_diagnostics.F
===================================================================
--- trunk/mpas/src/core_sw/module_global_diagnostics.F        2010-09-01 23:21:01 UTC (rev 492)
+++ trunk/mpas/src/core_sw/module_global_diagnostics.F        2010-09-02 17:14:50 UTC (rev 493)
@@ -28,7 +28,7 @@
! 1. Define the array to integrate, and the variable for the value above.
! 2. Allocate the array with the correct dimensions.
! 3. Fill the array with the data to be integrated.
- ! eg. G_h = Sum(h dA)/Sum(dA), See below for array filling
+ ! eg. GlobalFluidThickness = Sum(h dA)/Sum(dA), See below for array filling
! 4. Call Function to compute Global Stat that you want.
! 5. Finish computing the global stat/integral
! 6. Write out your global stat to the file
@@ -46,26 +46,36 @@
integer :: nCells
! Step 1
- ! 1. Define the array to integrate, and the variable for the value above.
- real (kind=RKIND) :: areaCellGlobal, areaEdgeGlobal, areaTriangleGlobal
- real (kind=RKIND), dimension(:), pointer :: areaCell, dcEdge, dvEdge, areaTriangle, areaEdge, h_s
- real (kind=RKIND), dimension(:,:), pointer :: h, u, v, h_edge, pv_edge, pv_vertex, pv_cell, h_vertex
+ ! 1. Define the array to integrate, and the variable for the value to be stored in after the integration
+ real (kind=RKIND), dimension(:), pointer :: areaCell, dcEdge, dvEdge, areaTriangle, h_s, fCell, fEdge
+ real (kind=RKIND), dimension(:,:), pointer :: h, u, v, h_edge, pv_edge, pv_vertex, pv_cell, h_vertex, weightsOnEdge
- real (kind=RKIND), dimension(:), allocatable :: eres, h_ref
- real (kind=RKIND), dimension(:,:), allocatable :: pe_vertex, u2, h2, hb, pv_temp, pve_bot, h_top, h_bot, cor_energy, sshr, ke_sink, pe_sink
real (kind=RKIND), dimension(:,:,:), pointer :: tracers
- real (kind=RKIND) :: global_temp, gh, gpv, gpe, ge, ger, e_cor, sshr_val, h_ref_val, ke_sink_val, pe_sink_val, ke
+ real (kind=RKIND), dimension(:), allocatable :: volumeWeightedPotentialEnergyReservoir, averageThickness
+ real (kind=RKIND), dimension(:), allocatable :: potentialEnstrophyReservior, areaEdge, h_s_edge
+
+ real (kind=RKIND), dimension(:,:), allocatable :: cellVolume, cellArea, volumeWeightedPotentialVorticity
+ real (kind=RKIND), dimension(:,:), allocatable :: volumeWeightedPotentialEnstrophy, vertexVolume, volumeWeightedKineticEnergy
+ real (kind=RKIND), dimension(:,:), allocatable :: volumeWeightedPotentialEnergy, volumeWeightedPotentialEnergyTopography
+ real (kind=RKIND), dimension(:,:), allocatable :: keTend_CoriolisForce, keTend_PressureGradient
+ real (kind=RKIND), dimension(:,:), allocatable ::peTend_DivThickness, refAreaWeightedSurfaceHeight, refAreaWeightedSurfaceHeight_edge
+
+ real (kind=RKIND) :: sumCellVolume, sumCellArea, sumVertexVolume, sumrefAreaWeightedSurfaceHeight
+
+ real (kind=RKIND) :: globalFluidThickness, globalPotentialVorticity, globalPotentialEnstrophy, globalEnergy
+ real (kind=RKIND) :: globalCoriolisEnergyTendency, globalKEPETendency, globalPotentialEnstrophyReservoir
+ real (kind=RKIND) :: globalKineticEnergy, globalPotentialEnergy, globalPotentialEnergyReservoir
+ real (kind=RKIND) :: globalKineticEnergyTendency, globalPotentialEnergyTendency
+ real (kind=RKIND) :: global_temp, workpv, q
real (kind=RKIND) :: volumeCellGlobal, volumeEdgeGlobal, CFLNumberGlobal
+
integer :: elementIndex, variableIndex, nVariables, nSums, nMaxes, nMins
- integer :: timeLevel,k,i
+ integer :: timeLevel, eoe, iLevel, iCell, iEdge, iVertex
+ integer :: fileID, iCell1, iCell2, j
- integer, dimension(:,:), pointer :: cellsOnEdge, edgesOnCell
-
- integer, parameter :: kMaxVariables = 1024 ! this must be a little more than double the number of variables to be reduced
- real (kind=RKIND), dimension(kMaxVariables) :: sums, mins, maxes, averages, verticalSumMins, verticalSumMaxes, reductions
-
- integer :: fileID, iCell1, iCell2, iEdge
+ integer, dimension(:,:), pointer :: cellsOnEdge, edgesOnCell, edgesOnEdge
+ integer, dimension(:), pointer :: nEdgesOnEdge
cellsOnEdge => grid % cellsOnEdge % array
edgesOnCell => grid % edgesOnCell % array
@@ -81,8 +91,14 @@
dcEdge => grid % dcEdge % array
dvEdge => grid % dvEdge % array
areaTriangle => grid % areaTriangle % array
+ fCell => grid % fCell % array
+ fEdge => grid % fEdge % array
+ edgesOnEdge => grid % edgesOnEdge % array
+ nEdgesOnEdge => grid % nEdgesOnEdge % array
+
allocate(areaEdge(1:nEdgesSolve))
areaEdge = dcEdge(1:nEdgesSolve)*dvEdge(1:nEdgesSolve)
+ weightsOnEdge => grid % weightsOnEdge % array
h => state % h % array
u => state % u % array
@@ -96,114 +112,171 @@
! Step 2
! 2. Allocate the array with the correct dimensions.
- allocate(h_top(nVertLevels,nCellsSolve))
- allocate(h_bot(nVertLevels,nCellsSolve))
- allocate(h2(nVertLevels,nCellsSolve))
- allocate(hb(nVertLevels,nCellsSolve))
- allocate(sshr(nVertLevels,nCellsSolve))
+ allocate(cellVolume(nVertLevels,nCellsSolve))
+ allocate(cellArea(nVertLevels,nCellsSolve))
+ allocate(refAreaWeightedSurfaceHeight(nVertLevels,nCellsSolve))
+ allocate(refAreaWeightedSurfaceHeight_edge(nVertLevels,nEdgesSolve))
+ allocate(volumeWeightedPotentialVorticity(nVertLevels,nVerticesSolve))
+ allocate(volumeWeightedPotentialEnstrophy(nVertLevels,nVerticesSolve))
+ allocate(potentialEnstrophyReservior(nCellsSolve))
+ allocate(vertexVolume(nVertLevels,nVerticesSolve))
+ allocate(volumeWeightedKineticEnergy(nVertLevels,nEdgesSolve))
+ allocate(volumeWeightedPotentialEnergy(nVertLevels,nCellsSolve))
+ allocate(volumeWeightedPotentialEnergyTopography(nVertLevels,nCellsSolve))
+ allocate(volumeWeightedPotentialEnergyReservoir(nCellsSolve))
+ allocate(keTend_CoriolisForce(nVertLevels,nEdgesSolve))
+ allocate(keTend_PressureGradient(nVertLevels,nEdgesSolve))
+ allocate(peTend_DivThickness(nVertLevels,nCells))
- allocate(h_ref(nCellsSolve))
+ allocate(averageThickness(nCellsSolve))
- allocate(pv_temp(nVertLevels,nVerticesSolve))
- allocate(pe_vertex(nVertLevels,nVerticesSolve))
- allocate(pve_bot(nVertLevels,nVerticesSolve))
+ allocate(h_s_edge(nEdgesSOlve))
- allocate(u2(nVertLevels,nEdgesSolve))
- allocate(cor_energy(nVertLevels,nEdgesSolve))
- allocate(ke_sink(nVertLevels,nEdgesSolve))
- allocate(pe_sink(nVertLevels,nCells))
+ cellVolume = 0
+ refAreaWeightedSurfaceHeight = 0
+ refAreaWeightedSurfaceHeight_edge = 0
+ vertexVolume = 0
+ cellArea = 0
+ averageThickness = 0
+ volumeWeightedPotentialVorticity = 0
+ volumeWeightedPotentialEnstrophy = 0
+ volumeWeightedKineticEnergy = 0
+ volumeWeightedPotentialEnergy = 0
+ volumeWeightedPotentialEnergyTopography = 0
+ volumeWeightedPotentialEnergyReservoir = 0
+ keTend_PressureGradient = 0
+ peTend_DivThickness = 0
+ keTend_CoriolisForce = 0
+ h_s_edge = 0
- allocate(eres(nCellsSolve))
-
- pe_sink = 0
-
! Build Arrays for Global Integrals
! Step 3
! 3. Fill the array with the data to be integrated.
- ! eg. G_h = Sum(h dA)/Sum(dA), See below for array filling
- do i = 1,nVertLevels
- ! G_h Top = Sum(h dA)
- h_top(i,:) = h(i,1:nCellsSolve)*areaCell(1:nCellsSolve)
- ! G_h Bot = Sum(dA)
- h_bot(i,:) = areaCell(1:nCellsSolve)
- pv_temp(i,:) = pv_vertex(i,1:nVerticesSolve)*h_vertex(i,1:nVerticesSolve)*areaTriangle(1:nVerticesSolve)
- pe_vertex(i,:) = pv_vertex(i,1:nVerticesSolve)*pv_vertex(i,1:nVerticesSolve)*h_vertex(i,1:nVerticesSolve)*areaTriangle(1:nVerticesSolve)
- pve_bot(i,:) = h_vertex(i,1:nVerticesSolve)*areaTriangle(1:nVerticesSolve)
- u2(i,:) = u(i,1:nEdgesSolve)*u(i,1:nEdgesSolve)*h_edge(i,1:nEdgesSolve)*areaEdge(1:nEdgesSolve)*0.5
- h2(i,:) = gravity*h(i,1:nCellsSolve)*h(i,1:nCellsSolve)*areaCell(1:nCellsSolve)*0.5
- hb(i,:) = gravity*h(i,1:nCellsSolve)*h_s(1:nCellsSolve)*areaCell(1:nCellsSolve)
- cor_energy(i,:) = h_edge(i,1:nEdgesSolve)*u(i,1:nEdgesSolve)*pv_edge(i,1:nEdgesSolve)*areaEdge(1:nEdgesSolve)
- sshr(i,:) = areaCell(1:nCellsSolve)*(h(i,1:nCellsSolve)+h_s(1:nCellsSolve))
+ ! eg. GlobalFluidThickness = Sum(h dA)/Sum(dA), See below for array filling
+ do iLevel = 1,nVertLevels
+ ! eg. GlobalFluidThickness top (Sum( h dA)) = Sum(cellVolume)
+ cellVolume(iLevel,:) = h(iLevel,1:nCellsSolve)*areaCell(1:nCellsSolve)
+ ! eg. GlobalFluidThickness bot (Sum(dA)) = Sum(cellArea)
+ cellArea(iLevel,:) = areaCell(1:nCellsSolve)
+ volumeWeightedPotentialVorticity(iLevel,:) = pv_vertex(iLevel,1:nVerticesSolve) &
+ *h_vertex(iLevel,1:nVerticesSolve)*areaTriangle(1:nVerticesSolve)
+ volumeWeightedPotentialEnstrophy(iLevel,:) = pv_vertex(iLevel,1:nVerticesSolve) &
+ *pv_vertex(iLevel,1:nVerticesSolve)*h_vertex(iLevel,1:nVerticesSolve)*areaTriangle(1:nVerticesSolve)
+ vertexVolume(iLevel,:) = h_vertex(iLevel,1:nVerticesSolve)*areaTriangle(1:nVerticesSolve)
+ volumeWeightedKineticEnergy(iLevel,:) = u(iLevel,1:nEdgesSolve)*u(iLevel,1:nEdgesSolve) &
+ *h_edge(iLevel,1:nEdgesSolve)*areaEdge(1:nEdgesSolve)*0.5
+ volumeWeightedPotentialEnergy(iLevel,:) = gravity*h(iLevel,1:nCellsSolve)*h(iLevel,1:nCellsSolve)*areaCell(1:nCellsSolve)*0.5
+ volumeWeightedPotentialEnergyTopography(iLevel,:) = gravity*h(iLevel,1:nCellsSolve)*h_s(1:nCellsSolve)*areaCell(1:nCellsSolve)
+ refAreaWeightedSurfaceHeight(iLevel,:) = areaCell(1:nCellsSolve)*(h(iLevel,1:nCellsSolve)+h_s(1:nCellsSolve))
- do k = 1,nEdgesSolve
- iCell1 = cellsOnEdge(k,1)
- iCell2 = cellsOnEdge(k,2)
- ke_sink(i,k) = areaEdge(k)*h_edge(i,k)*u(i,k)*gravity*(h(i,iCell2)+h_s(iCell2) - h(i,iCell1)-h_s(iCell1))/dcEdge(k)
- pe_sink(i,iCell1) = pe_sink(i,iCell1) + h_edge(i,k)*u(i,k)*dvEdge(k)
- pe_sink(i,iCell2) = pe_sink(i,iCell2) - h_edge(i,k)*u(i,k)*dvEdge(k)
+ do iEdge = 1,nEdgesSolve
+ q = 0.0
+ do j = 1,nEdgesOnEdge(iEdge)
+ eoe = edgesOnEdge(j,iEdge)
+ workpv = 0.5 * (pv_edge(iLevel,iEdge) + pv_edge(iLevel,eoe))
+ q = q + weightsOnEdge(j,iEdge) * u(iLevel,eoe) * workpv * h_edge(iLevel,eoe)
+ end do
+ keTend_CoriolisForce(iLevel,iEdge) = h_edge(iLevel,iEdge) * u(iLevel,iEdge) * q * areaEdge(iEdge)
+
+ iCell1 = cellsOnEdge(iEdge,1)
+ iCell2 = cellsOnEdge(iEdge,2)
+
+ refAreaWeightedSurfaceHeight_edge(iLevel,iEdge) = areaEdge(iEdge)*(h_edge(iLevel,iEdge) + 0.5*(h_s(iCell1) + h_s(iCell2)))
+
+ keTend_PressureGradient(iLevel,iEdge) = areaEdge(iEdge)*h_edge(iLevel,iEdge)*u(iLevel,iEdge) &
+ *gravity*(h(iLevel,iCell2)+h_s(iCell2) - h(iLevel,iCell1)-h_s(iCell1))/dcEdge(iEdge)
+ peTend_DivThickness(iLevel,iCell1) = peTend_DivThickness(iLevel,iCell1) &
+ + h_edge(iLevel,iEdge)*u(iLevel,iEdge)*dvEdge(iEdge)
+ peTend_DivThickness(iLevel,iCell2) = peTend_DivThickness(iLevel,iCell2) &
+ - h_edge(iLevel,iEdge)*u(iLevel,iEdge)*dvEdge(iEdge)
end do
- pe_sink(i,:) = pe_sink(i,1:nCellsSolve)*gravity*(h(i,1:nCellsSolve)+h_s(1:nCellsSolve))
+
+ peTend_DivThickness(iLevel,:) = peTend_DivThickness(iLevel,1:nCells)*gravity &
+ *(h(iLevel,1:nCells)+h_s(1:nCells))
end do
+ do iEdge = 1,nEdgesSolve
+ iCell1 = cellsOnEdge(iEdge,1)
+ iCell2 = cellsOnEdge(iEdge,2)
+
+ h_s_edge(iEdge) = 0.5*(h_s(iCell1) + h_s(iCell2))
+ end do
+
! Step 4
! 4. Call Function to compute Global Stat that you want.
- ! KE and PE Sink Terms
- call computeGlobalSum(dminfo, nVertLevels, nEdgesSolve, ke_sink, ke_sink_val)
- call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, pe_sink, pe_sink_val)
+ ! Computing Kinetic and Potential Energy Tendency Terms
+ call computeGlobalSum(dminfo, nVertLevels, nEdgesSolve, keTend_PressureGradient, globalKineticEnergyTendency)
+ call computeGlobalSum(dminfo, nVertLevels, nCells, peTend_DivThickness, globalPotentialEnergyTendency)
- call computeGlobalSum(dminfo, 1, nCellsSolve, h_ref, h_ref_val)
+ ! Computing top and bottom of global mass integral
+ call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, cellVolume, sumCellVolume)
+ call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, cellArea, sumCellArea)
- ! Global Area Average of h
- call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, h_top, gh)
- call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, h_bot, global_temp)
+ globalKineticEnergyTendency = globalKineticEnergyTendency / sumCellVolume
+ globalPotentialEnergyTendency = globalPotentialEnergyTendency / sumCellVolume
! Step 5
! 5. Finish computing the global stat/integral
- gh = gh/global_temp
+ globalFluidThickness = sumCellVolume/sumCellArea
- ! Mean SSH for PE Res
- call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, sshr, sshr_val)
+ ! Compute Average Sea Surface Height for Potential Energy and Enstrophy
+ ! Reservoir computations
+ call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, refAreaWeightedSurfaceHeight, sumrefAreaWeightedSurfaceHeight)
- h_ref(:) = (sshr_val/global_temp)-h_s(1:nCellsSolve)
+ averageThickness(:) = (sumrefAreaWeightedSurfaceHeight/sumCellArea)-h_s(1:nCellsSolve)
+ ! Compute Volume Weighted Averages of Potential Vorticity and Potential Enstrophy
+ call computeGlobalSum(dminfo, nVertLevels, nVerticesSolve, volumeWeightedPotentialVorticity, globalPotentialVorticity)
+ call computeGlobalSum(dminfo, nVertLevels, nVerticesSolve, volumeWeightedPotentialEnstrophy, globalPotentialEnstrophy)
+ call computeGlobalSum(dminfo, nVertLevels, nVerticesSolve, vertexVolume, sumVertexVolume)
- ! Gloabl Area-Thickness Average of pv and pe
- call computeGlobalSum(dminfo, nVertLevels, nVerticesSolve, pv_temp, gpv)
- call computeGlobalSum(dminfo, nVertLevels, nVerticesSolve, pe_vertex, gpe)
- call computeGlobalSum(dminfo, nVertLevels, nVerticesSolve, pve_bot, global_temp)
+ globalPotentialVorticity = globalPotentialVorticity/sumVertexVolume
+ globalPotentialEnstrophy = globalPotentialEnstrophy/sumVertexVolume
- gpv = gpv/global_temp
- gpe = gpe/global_temp
+ ! Compte Potential Enstrophy Reservior
+ potentialEnstrophyReservior(:) = areaCell(:)*fCell(:)*fCell(:)/averageThickness
+ call computeGlobalSum(dminfo, 1, nCellsSolve, potentialEnstrophyReservior, globalPotentialEnstrophyReservoir)
+ globalPotentialEnstrophyReservoir = globalPotentialEnstrophyReservoir/sumCellVolume
- ! Global Integral of Total Energy with subtracting reference res
- call computeGlobalSum(dminfo, nVertLevels, nEdgesSolve, u2, ke)
- call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, h2, global_temp)
- ge = ke + global_temp
- call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, hb, global_temp)
- ge = ge + global_temp
+ globalPotentialEnstrophy = globalPotentialEnstrophy - globalPotentialEnstrophyReservoir
- eres(1:nCellsSolve) = areaCell(1:nCellsSolve)*h_ref*h_ref*gravity*0.5
- call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, eres, global_temp)
- ger = global_temp
- eres(1:nCellsSolve) = areaCell(1:nCellsSolve)*h_ref*h_s(1:nCellsSolve)*gravity
- call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, eres, global_temp)
- ger = ger + global_temp
- ge = ge - ger
+ ! Compute Kinetic and Potential Energy terms to be combined into total energy
+ call computeGlobalSum(dminfo, nVertLevels, nEdgesSolve, volumeWeightedKineticEnergy, globalKineticEnergy)
+ call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, volumeWeightedPotentialEnergy, globalPotentialEnergy)
+ call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, volumeWeightedPotentialEnergyTopography, global_temp)
- ! Global Integral of spurious Coriolis energy for Time to KE doubling
- call computeGlobalSum(dminfo, nVertLevels, nEdgesSolve, cor_energy, e_cor)
+ globalKineticEnergy = globalKineticEnergy/sumCellVolume
+ globalPotentialEnergy = (globalPotentialEnergy + global_temp)/sumCellVolume
+ ! Compute Potential energy reservoir to be subtracted from potential energy term
+ volumeWeightedPotentialEnergyReservoir(1:nCellsSolve) = areaCell(1:nCellsSolve)*averageThickness*averageThickness*gravity*0.5
+ call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, volumeWeightedPotentialEnergyReservoir, globalPotentialEnergyReservoir)
+ volumeWeightedPotentialEnergyReservoir(1:nCellsSolve) = areaCell(1:nCellsSolve)*averageThickness*h_s(1:nCellsSolve)*gravity
+ call computeGlobalSum(dminfo, nVertLevels, nCellsSolve, volumeWeightedPotentialEnergyReservoir, global_temp)
+
+ globalPotentialEnergyReservoir = (globalPotentialEnergyReservoir + global_temp)/sumCellVolume
+
+ globalPotentialEnergy = globalPotentialEnergy - globalPotentialEnergyReservoir
+ globalEnergy = globalKineticEnergy + globalPotentialEnergy
+
+ ! Compute Coriolis energy tendency term
+ call computeGlobalSum(dminfo, nVertLevels, nEdgesSolve, keTend_CoriolisForce, globalCoriolisEnergyTendency)
+ globalCoriolisEnergyTendency = globalCoriolisEnergyTendency/sumCellVolume
+
! Step 6
! 6. Write out your global stat to the file
if (dminfo % my_proc_id == IO_NODE) then
fileID = getFreeUnit()
+
if (timeIndex/config_stats_interval == 1) then
open(fileID, file='GlobalIntegrals.txt',STATUS='unknown')
else
- open(fileID, file='GlobalIntegrals.txt',ACCESS='append')
+ open(fileID, file='GlobalIntegrals.txt',POSITION='append')
endif
- write(fileID,'(100es24.16)') gh, gpv, gpe, ge, e_cor, abs(ke/e_cor), ke_sink_val+pe_sink_val,abs(ke/(ke_sink_val+pe_sink_val))
+ write(fileID,'(1i, 100es24.16)') timeIndex, timeIndex*dt, globalFluidThickness, globalPotentialVorticity, globalPotentialEnstrophy, &
+ globalEnergy, globalCoriolisEnergyTendency, globalKineticEnergyTendency+globalPotentialEnergyTendency, &
+ globalKineticEnergy, globalPotentialEnergy
close(fileID)
end if
Modified: trunk/mpas/src/core_sw/mpas_interface.F
===================================================================
--- trunk/mpas/src/core_sw/mpas_interface.F        2010-09-01 23:21:01 UTC (rev 492)
+++ trunk/mpas/src/core_sw/mpas_interface.F        2010-09-02 17:14:50 UTC (rev 493)
@@ -64,18 +64,20 @@
call timestep(domain, dt)
- if (mod(itimestep, config_stats_interval) == 0) then
- block_ptr => domain % blocklist
- if(associated(block_ptr % next)) then
- write(0,*) 'Error: computeGlobalDiagnostics assumes ',&
- 'that there is only one block per processor.'
+ if(config_stats_interval .gt. 0) then
+ if(mod(itimestep, config_stats_interval) == 0) then
+ block_ptr => domain % blocklist
+ if(associated(block_ptr % next)) then
+ write(0,*) 'Error: computeGlobalDiagnostics assumes ',&
+ 'that there is only one block per processor.'
+ end if
+
+ call timer_start("global_diagnostics")
+ call computeGlobalDiagnostics(domain % dminfo, &
+ block_ptr % time_levs(2) % state, block_ptr % mesh, &
+ itimestep, dt)
+ call timer_stop("global_diagnostics")
end if
-
- call timer_start("global_diagnostics")
- call computeGlobalDiagnostics(domain % dminfo, &
- block_ptr % time_levs(2) % state, block_ptr % mesh, &
- itimestep, dt)
- call timer_stop("global_diagnostics")
end if
end subroutine mpas_timestep
</font>
</pre>