<p><b>mpetersen@lanl.gov</b> 2011-05-17 10:38:11 -0600 (Tue, 17 May 2011)</p><p>Check-in with debugging print statements.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/ocean_projects/timesplitting_mrp/src/core_ocean/Registry
===================================================================
--- branches/ocean_projects/timesplitting_mrp/src/core_ocean/Registry 2011-05-16 19:34:10 UTC (rev 836)
+++ branches/ocean_projects/timesplitting_mrp/src/core_ocean/Registry 2011-05-17 16:38:11 UTC (rev 837)
@@ -177,9 +177,9 @@
# state variables for Higdon timesplitting
var persistent real uBtr ( nEdges Time ) 2 o uBtr state - -
-var persistent real eta ( nCells Time ) 2 o eta state - -
+var persistent real ssh ( nCells Time ) 2 o ssh state - -
var persistent real uBtrSubcycle ( nEdges Time ) 2 o uBtrSubcycle state - -
-var persistent real etaSubcycle ( nCells Time ) 2 o etaSubcycle state - -
+var persistent real sshSubcycle ( nCells Time ) 2 o sshSubcycle state - -
var persistent real FBtr ( nEdges Time ) 1 o FBtr state - -
var persistent real GBtrForcing ( nEdges Time ) 1 o GBtrForcing state - -
var persistent real uBcl ( nVertLevels nEdges Time ) 2 o uBcl state - -
Modified: branches/ocean_projects/timesplitting_mrp/src/core_ocean/module_mpas_core.F
===================================================================
--- branches/ocean_projects/timesplitting_mrp/src/core_ocean/module_mpas_core.F 2011-05-16 19:34:10 UTC (rev 836)
+++ branches/ocean_projects/timesplitting_mrp/src/core_ocean/module_mpas_core.F 2011-05-17 16:38:11 UTC (rev 837)
@@ -109,8 +109,12 @@
:block % mesh % nVertLevels,iCell) = 0.0
! mrp changed to 0
! :block % mesh % nVertLevels,iCell) = -1e34
+
+! mrp 110516 temp, added just to test for conservation of tracers
+ block % state % time_levs(1) % state % tracers % array(3,:,iCell) = 1.0
+
end do
-
+
if (.not. config_do_restart) then
block % state % time_levs(1) % state % xtime % scalar = 0.0
else
Modified: branches/ocean_projects/timesplitting_mrp/src/core_ocean/module_time_integration.F
===================================================================
--- branches/ocean_projects/timesplitting_mrp/src/core_ocean/module_time_integration.F 2011-05-16 19:34:10 UTC (rev 836)
+++ branches/ocean_projects/timesplitting_mrp/src/core_ocean/module_time_integration.F 2011-05-17 16:38:11 UTC (rev 837)
@@ -381,7 +381,7 @@
real (kind=RKIND), dimension(:), allocatable:: A,C,uTemp, G
real (kind=RKIND), dimension(:,:), allocatable:: tracersTemp
-print *, '1'
+!print *, '1'
!
@@ -391,17 +391,17 @@
do while (associated(block))
! printing:
- nCells = block % mesh % nCells
- nEdges = block % mesh % nEdges
- nVertLevels = block % mesh % nVertLevels
-print *, 'uold',minval(block % state % time_levs(1) % state % u % array(:,1:nEdges)),&
- maxval(block % state % time_levs(1) % state % u % array(:,1:nEdges))
-print *, 'hold',minval(block % state % time_levs(1) % state % h % array(:,1:nCells)),&
- maxval(block % state % time_levs(1) % state % h % array(:,1:nCells))
-print *, 'Told',minval(block % state % time_levs(1) % state % tracers % array(1,:,1:nCells)),&
- maxval(block % state % time_levs(1) % state % tracers % array(1,:,1:nCells))
-print *, 'Sold',minval(block % state % time_levs(1) % state % tracers % array(2,:,1:nCells)),&
- maxval(block % state % time_levs(1) % state % tracers % array(2,:,1:nCells))
+! nCells = block % mesh % nCells
+! nEdges = block % mesh % nEdges
+! nVertLevels = block % mesh % nVertLevels
+!print *, 'uold',minval(block % state % time_levs(1) % state % u % array(:,1:nEdges)),&
+! maxval(block % state % time_levs(1) % state % u % array(:,1:nEdges))
+!print *, 'hold',minval(block % state % time_levs(1) % state % h % array(:,1:nCells)),&
+! maxval(block % state % time_levs(1) % state % h % array(:,1:nCells))
+!print *, 'Told',minval(block % state % time_levs(1) % state % tracers % array(1,:,1:nCells)),&
+! maxval(block % state % time_levs(1) % state % tracers % array(1,:,1:nCells))
+!print *, 'Sold',minval(block % state % time_levs(1) % state % tracers % array(2,:,1:nCells)),&
+! maxval(block % state % time_levs(1) % state % tracers % array(2,:,1:nCells))
! printing end
do iEdge=1,block % mesh % nEdges
@@ -442,11 +442,11 @@
end do
end do
-!mrp 110512 need to add w, p, and eta here?
+!mrp 110512 need to add w, p, and ssh here?
block => block % next
end do
-print *, '2'
+!print *, '2'
!
@@ -473,7 +473,7 @@
block => block % next
end do
-print *, '3'
+!print *, '3'
!
! Stage 1: Baroclinic velocity (3D) prediction, explicit with long timestep
@@ -487,7 +487,7 @@
call compute_vertical_mix_coefficients(block % state % time_levs(2) % state, block % diagnostics, block % mesh)
end if
call compute_tend_u(block % tend, block % state % time_levs(2) % state , block % diagnostics, block % mesh)
-! mrp soon: for split version: add g/rho grad eta to tend_u
+! mrp soon: for split version: add g/rho grad ssh to tend_u
call enforce_boundaryEdge(block % tend, block % mesh)
block => block % next
end do
@@ -501,7 +501,7 @@
! call compute_uPerp(uBclNew,uBclPerp)
! need subroutine and variable just for perp, say uBclPerp - may already be there.
-print *, '4'
+!print *, '4'
block => domain % blocklist
do while (associated(block))
@@ -510,7 +510,7 @@
do iEdge=1,block % mesh % nEdges
G = 0.0 ! could put this after with G(maxleveledgetop+1:nvertlevels)=0
do k=1,block % mesh % maxLevelEdgeTop % array(iEdge)
- !{\bf G}_k = -f{\bf u}_{k,n+1}^{'\;\perp} + {\bf T}({\bf u}_k^*,w_k^*,p_k^*) +\frac{g}{\rho_0}</font>
<font color="blue">abla \eta^*
+ !{\bf G}_k = -f{\bf u}_{k,n+1}^{'\;\perp} + {\bf T}({\bf u}_k^*,w_k^*,p_k^*) +\frac{g}{\rho_0}</font>
<font color="gray">abla \ssh^*
! This soon, with coriolis:
! G(k) = -fEdge(iEdge)*uBclPerp(k,iEdge) + tend_u(k,iEdge)
@@ -547,7 +547,7 @@
end do
enddo ! config_n_bcl_iter
-print *, '5'
+!print *, '5'
@@ -566,7 +566,7 @@
end do ! block
-print *, '6'
+!print *, '6'
!
! Stage 3: Tracer, density, pressure, vertical velocity prediction
@@ -590,10 +590,15 @@
call compute_scalar_tend(block % tend, block % state % time_levs(2) % state , block % diagnostics, block % mesh)
+! printing:
+ nCells = block % mesh % nCells
+print *, 'tr_tend',minval(block % tend % tracers % array(3,1,1:nCells)),&
+ maxval(block % tend % tracers % array(3,1,1:nCells))
+
block => block % next
end do
-print *, '7'
+!print *, '7'
! --- update halos for prognostic variables
@@ -607,7 +612,7 @@
block % parinfo % cellsToSend, block % parinfo % cellsToRecv)
block => block % next
end do
-print *, '8'
+!print *, '8'
block => domain % blocklist
@@ -662,16 +667,16 @@
nCells = block % mesh % nCells
nEdges = block % mesh % nEdges
nVertLevels = block % mesh % nVertLevels
-print *, 'uold ',minval(block % state % time_levs(1) % state % u % array(:,1:nEdges)),&
- maxval(block % state % time_levs(1) % state % u % array(:,1:nEdges))
-print *, 'hold ',minval(block % state % time_levs(1) % state % h % array(:,1:nCells)),&
- maxval(block % state % time_levs(1) % state % h % array(:,1:nCells))
-print *, 'hold1 ',minval(block % state % time_levs(1) % state % h % array(1,1:nCells)),&
- maxval(block % state % time_levs(1) % state % h % array(1,1:nCells))
-print *, 'Told ',minval(block % state % time_levs(1) % state % tracers % array(1,:,1:nCells)),&
- maxval(block % state % time_levs(1) % state % tracers % array(1,:,1:nCells))
-print *, 'Sold ',minval(block % state % time_levs(1) % state % tracers % array(2,:,1:nCells)),&
- maxval(block % state % time_levs(1) % state % tracers % array(2,:,1:nCells))
+!print *, 'uold ',minval(block % state % time_levs(1) % state % u % array(:,1:nEdges)),&
+! maxval(block % state % time_levs(1) % state % u % array(:,1:nEdges))
+!print *, 'hold ',minval(block % state % time_levs(1) % state % h % array(:,1:nCells)),&
+! maxval(block % state % time_levs(1) % state % h % array(:,1:nCells))
+!print *, 'hold1 ',minval(block % state % time_levs(1) % state % h % array(1,1:nCells)),&
+! maxval(block % state % time_levs(1) % state % h % array(1,1:nCells))
+!print *, 'Told ',minval(block % state % time_levs(1) % state % tracers % array(1,:,1:nCells)),&
+! maxval(block % state % time_levs(1) % state % tracers % array(1,:,1:nCells))
+!print *, 'Sold ',minval(block % state % time_levs(1) % state % tracers % array(2,:,1:nCells)),&
+! maxval(block % state % time_levs(1) % state % tracers % array(2,:,1:nCells))
print *, 'unew ',minval(block % state % time_levs(2) % state % u % array(:,1:nEdges)),&
maxval(block % state % time_levs(2) % state % u % array(:,1:nEdges))
print *, 'hnew ',minval(block % state % time_levs(2) % state % h % array(:,1:nCells)),&
@@ -682,7 +687,11 @@
maxval(block % state % time_levs(2) % state % tracers % array(1,:,1:nCells))
print *, 'Snew ',minval(block % state % time_levs(2) % state % tracers % array(2,:,1:nCells)),&
maxval(block % state % time_levs(2) % state % tracers % array(2,:,1:nCells))
- block % state % time_levs(2) % state % tracers % array(3,1,1:nCells) = domain % dminfo % my_proc_id
+print *, 'Tr1Old1',minval(block % state % time_levs(1) % state % tracers % array(3,1,1:nCells)),&
+ maxval(block % state % time_levs(1) % state % tracers % array(3,1,1:nCells))
+print *, 'Tr1New1',minval(block % state % time_levs(2) % state % tracers % array(3,1,1:nCells)),&
+ maxval(block % state % time_levs(2) % state % tracers % array(3,1,1:nCells))
+! block % state % time_levs(2) % state % tracers % array(3,1,1:nCells) = domain % dminfo % my_proc_id
! printing end
u => block % state % time_levs(2) % state % u % array
@@ -927,7 +936,7 @@
end do
end do
do iCell=1,nCells
- tend_h(1,iCell) = tend_h(1,iCell) / areaCell(iCell)
+ tend_h(1,iCell) = tend_h(1,iCell) / areaCell(iCell)
end do
endif ! config_vert_grid_type
@@ -1024,6 +1033,7 @@
areaCell => grid % areaCell % array
areaTriangle => grid % areaTriangle % array
h_s => grid % h_s % array
+! mrp 110516 cleanup fvertex fedge not used in this subroutine
fVertex => grid % fVertex % array
fEdge => grid % fEdge % array
zMidZLevel => grid % zMidZLevel % array
@@ -1045,9 +1055,32 @@
!
! velocity tendency: start accumulating tendency terms
!
+ ! mrp 110516 efficiency: could remove next line and have first tend_u operation not be additive
tend_u(:,:) = 0.0
!
+ ! velocity tendency: nonlinear Coriolis term and grad of kinetic energy
+ !
+ do iEdge=1,grid % nEdgesSolve
+ cell1 = cellsOnEdge(1,iEdge)
+ cell2 = cellsOnEdge(2,iEdge)
+
+ do k=1,maxLevelEdgeTop(iEdge)
+
+ q = 0.0
+ do j = 1,nEdgesOnEdge(iEdge)
+ eoe = edgesOnEdge(j,iEdge)
+ workpv = 0.5 * (pv_edge(k,iEdge) + pv_edge(k,eoe))
+ q = q + weightsOnEdge(j,iEdge) * u(k,eoe) * workpv * h_edge(k,eoe)
+ end do
+ tend_u(k,iEdge) = tend_u(k,iEdge) &
+ + q &
+ - ( ke(k,cell2) - ke(k,cell1) ) / dcEdge(iEdge)
+
+ end do
+ end do
+
+ !
! velocity tendency: vertical advection term -w du/dz
!
if (config_vert_grid_type.eq.'zlevel') then
@@ -1069,7 +1102,8 @@
! Average w*du/dz from vertical interface to vertical middle of cell
do k=1,maxLevelEdgeTop(iEdge)
- tend_u(k,iEdge) = - 0.5 * (w_dudzTopEdge(k) + w_dudzTopEdge(k+1))
+ tend_u(k,iEdge) = tend_u(k,iEdge) &
+ - 0.5 * (w_dudzTopEdge(k) + w_dudzTopEdge(k+1))
enddo
enddo
deallocate(w_dudzTopEdge)
@@ -1224,28 +1258,6 @@
end if
!
- ! velocity tendency: nonlinear Coriolis term and grad of kinetic energy
- !
- do iEdge=1,grid % nEdgesSolve
- cell1 = cellsOnEdge(1,iEdge)
- cell2 = cellsOnEdge(2,iEdge)
-
- do k=1,maxLevelEdgeTop(iEdge)
-
- q = 0.0
- do j = 1,nEdgesOnEdge(iEdge)
- eoe = edgesOnEdge(j,iEdge)
- workpv = 0.5 * (pv_edge(k,iEdge) + pv_edge(k,eoe))
- q = q + weightsOnEdge(j,iEdge) * u(k,eoe) * workpv * h_edge(k,eoe)
- end do
- tend_u(k,iEdge) = tend_u(k,iEdge) &
- + q &
- - ( ke(k,cell2) - ke(k,cell1) ) / dcEdge(iEdge)
-
- end do
- end do
-
- !
! velocity tendency: forcing and bottom drag
!
! mrp 101115 note: in order to include flux boundary conditions, we will need to
@@ -1816,6 +1828,13 @@
end if
+! mrp 110516 printing
+!print *, 'tend_tr 1',minval(tend_tr(3,1,1:nCells)),&
+! maxval(tend_tr(3,1,1:nCells))
+!print *, 'tracer 1',minval(tracers(3,1,1:nCells)),&
+! maxval(tracers(3,1,1:nCells))
+! mrp 110516 printing end
+
!
! tracer tendency: vertical diffusion h d/dz( \kappa_v d\phi/dz)
!
@@ -1842,11 +1861,18 @@
+ fluxVertTop(iTracer,k) - fluxVertTop(iTracer,k+1)
enddo
enddo
+!print '(a,50e12.2)', 'fluxVertTop',fluxVertTop(3,1:maxLevelCell(iCell)+1)
+!print '(a,50e12.2)', 'tend_tr ',tend_tr(3,1,1:maxLevelCell(iCell))
enddo ! iCell loop
deallocate(fluxVertTop)
endif
+! mrp 110516 printing
+!print *, 'tend_tr 2',minval(tend_tr(3,1,1:nCells)),&
+! maxval(tend_tr(3,1,1:nCells))
+! mrp 110516 printing end
+
!
! add restoring to T and S in top model layer
!
</font>
</pre>