<p><b>duda</b> 2010-10-19 15:17:20 -0600 (Tue, 19 Oct 2010)</p><p>BRANCH COMMIT<br>
<br>
1) Edit IRO field in non-hydrostatic Registry file to: reduce the number<br>
of fields written to restart files while preserving bit-for-bit restarts;<br>
and add fields to input and output files such that an output file created<br>
by running the model for 0 time steps can be used as an input file yeilding<br>
identical simulation results compared with those generated through a normal <br>
test case initialization.<br>
<br>
2) Add a "Test Case 0" option to simply use initial fields given in the input<br>
file.<br>
<br>
3) Remove unused variables qtot and h_s from the Registry file.<br>
<br>
4) Move halo exchanges for cqu and cqw to a more reasonable location (immediately<br>
after they are computed) in module_time_integration.F.<br>
<br>
<br>
M src/core_nhyd_atmos/module_test_cases.F<br>
M src/core_nhyd_atmos/Registry<br>
M src/core_nhyd_atmos/module_time_integration.F<br>
</p><hr noshade><pre><font color="gray">Modified: branches/atmos_nonhydrostatic/src/core_nhyd_atmos/Registry
===================================================================
--- branches/atmos_nonhydrostatic/src/core_nhyd_atmos/Registry        2010-10-19 21:01:25 UTC (rev 570)
+++ branches/atmos_nonhydrostatic/src/core_nhyd_atmos/Registry        2010-10-19 21:17:20 UTC (rev 571)
@@ -52,14 +52,13 @@
dim FIFTEEN 15
dim TWENTYONE 21
dim R3 3
-#dim nVertLevels nVertLevels
dim nVertLevels namelist:config_nvertlevels
dim nVertLevelsP1 nVertLevels+1
#
# var type name_in_file ( dims ) iro- name_in_code super-array array_class
#
-var persistent real xtime ( Time ) 2 ro xtime state - -
+var persistent real xtime ( Time ) 2 iro xtime state - -
# horizontal grid structure
@@ -97,9 +96,9 @@
var persistent real areaCell ( nCells ) 0 iro areaCell mesh - -
var persistent real areaTriangle ( nVertices ) 0 iro areaTriangle mesh - -
-var persistent real edgeNormalVectors ( R3 nEdges ) 0 ro edgeNormalVectors mesh - -
-var persistent real localVerticalUnitVectors ( R3 nCells ) 0 ro localVerticalUnitVectors mesh - -
-var persistent real cellTangentPlane ( R3 TWO nEdges ) 0 ro cellTangentPlane mesh - -
+var persistent real edgeNormalVectors ( R3 nEdges ) 0 iro edgeNormalVectors mesh - -
+var persistent real localVerticalUnitVectors ( R3 nCells ) 0 iro localVerticalUnitVectors mesh - -
+var persistent real cellTangentPlane ( R3 TWO nEdges ) 0 iro cellTangentPlane mesh - -
var persistent integer cellsOnCell ( maxEdges nCells ) 0 iro cellsOnCell mesh - -
var persistent integer verticesOnCell ( maxEdges nCells ) 0 iro verticesOnCell mesh - -
@@ -109,7 +108,6 @@
var persistent real kiteAreasOnVertex ( vertexDegree nVertices ) 0 iro kiteAreasOnVertex mesh - -
var persistent real fEdge ( nEdges ) 0 iro fEdge mesh - -
var persistent real fVertex ( nVertices ) 0 iro fVertex mesh - -
-var persistent real h_s ( nCells ) 0 iro h_s mesh - -
# some solver scalar coefficients
@@ -137,18 +135,18 @@
# coefficients for the vertical tridiagonal solve
# Note: these could be local but...
-var persistent real cofrz ( nVertLevels Time ) 1 r cofrz diag - -
-var persistent real cofwr ( nVertLevels nCells Time ) 1 r cofwr diag - -
-var persistent real cofwz ( nVertLevels nCells Time ) 1 r cofwz diag - -
-var persistent real coftz ( nVertLevelsP1 nCells Time ) 1 r coftz diag - -
-var persistent real cofwt ( nVertLevels nCells Time ) 1 r cofwt diag - -
-var persistent real a_tri ( nVertLevels nCells Time ) 1 r a_tri diag - -
-var persistent real alpha_tri ( nVertLevels nCells Time ) 1 r alpha_tri diag - -
-var persistent real gamma_tri ( nVertLevels nCells Time ) 1 r gamma_tri diag - -
+var persistent real cofrz ( nVertLevels Time ) 1 - cofrz diag - -
+var persistent real cofwr ( nVertLevels nCells Time ) 1 - cofwr diag - -
+var persistent real cofwz ( nVertLevels nCells Time ) 1 - cofwz diag - -
+var persistent real coftz ( nVertLevelsP1 nCells Time ) 1 - coftz diag - -
+var persistent real cofwt ( nVertLevels nCells Time ) 1 - cofwt diag - -
+var persistent real a_tri ( nVertLevels nCells Time ) 1 - a_tri diag - -
+var persistent real alpha_tri ( nVertLevels nCells Time ) 1 - alpha_tri diag - -
+var persistent real gamma_tri ( nVertLevels nCells Time ) 1 - gamma_tri diag - -
# W-Rayleigh-damping coefficient
-var persistent real dss ( nVertLevels nCells ) 0 ir dss mesh - -
+var persistent real dss ( nVertLevels nCells ) 0 iro dss mesh - -
# Prognostic variables: read from input, saved in restart, and written to output
var persistent real u ( nVertLevels nEdges Time ) 2 iro u state - -
@@ -160,17 +158,17 @@
var persistent real qr ( nVertLevels nCells Time ) 2 iro qr state scalars moist
# Tendency variables
-var persistent real tend_u ( nVertLevels nEdges Time ) 1 r u tend - -
-var persistent real tend_w ( nVertLevelsP1 nCells Time ) 1 r w tend - -
-var persistent real tend_rho ( nVertLevels nCells Time ) 1 r rho tend - -
-var persistent real tend_theta ( nVertLevels nCells Time ) 1 r theta tend - -
-var persistent real tend_qv ( nVertLevels nCells Time ) 1 r qv tend scalars moist
-var persistent real tend_qc ( nVertLevels nCells Time ) 1 r qc tend scalars moist
-var persistent real tend_qr ( nVertLevels nCells Time ) 1 r qr tend scalars moist
-var persistent real rt_diabatic_tend ( nVertLevels nCells Time ) 1 r rt_diabatic_tend tend - -
+var persistent real tend_u ( nVertLevels nEdges Time ) 1 - u tend - -
+var persistent real tend_w ( nVertLevelsP1 nCells Time ) 1 - w tend - -
+var persistent real tend_rho ( nVertLevels nCells Time ) 1 - rho tend - -
+var persistent real tend_theta ( nVertLevels nCells Time ) 1 - theta tend - -
+var persistent real tend_qv ( nVertLevels nCells Time ) 1 - qv tend scalars moist
+var persistent real tend_qc ( nVertLevels nCells Time ) 1 - qc tend scalars moist
+var persistent real tend_qr ( nVertLevels nCells Time ) 1 - qr tend scalars moist
+var persistent real rt_diabatic_tend ( nVertLevels nCells Time ) 1 - rt_diabatic_tend tend - -
# state variables diagnosed from prognostic state
-var persistent real pressure ( nVertLevels nCells Time ) 2 ro pressure state - -
+var persistent real pressure ( nVertLevels nCells Time ) 2 iro pressure state - -
var persistent real u_init ( nVertLevels ) 0 iro u_init mesh - -
var persistent real t_init ( nVertLevels nCells ) 0 iro t_init mesh - -
@@ -185,31 +183,30 @@
var persistent real ke ( nVertLevels nCells Time ) 1 ro ke diag - -
var persistent real pv_vertex ( nVertLevels nVertices Time ) 1 ro pv_vertex diag - -
var persistent real pv_cell ( nVertLevels nCells Time ) 1 ro pv_cell diag - -
-var persistent real uReconstructX ( nVertLevels nCells Time ) 1 ro uReconstructX diag - -
-var persistent real uReconstructY ( nVertLevels nCells Time ) 1 ro uReconstructY diag - -
-var persistent real uReconstructZ ( nVertLevels nCells Time ) 1 ro uReconstructZ diag - -
-var persistent real uReconstructZonal ( nVertLevels nCells Time ) 1 ro uReconstructZonal diag - -
-var persistent real uReconstructMeridional ( nVertLevels nCells Time ) 1 ro uReconstructMeridional diag - -
+var persistent real uReconstructX ( nVertLevels nCells Time ) 1 o uReconstructX diag - -
+var persistent real uReconstructY ( nVertLevels nCells Time ) 1 o uReconstructY diag - -
+var persistent real uReconstructZ ( nVertLevels nCells Time ) 1 o uReconstructZ diag - -
+var persistent real uReconstructZonal ( nVertLevels nCells Time ) 1 o uReconstructZonal diag - -
+var persistent real uReconstructMeridional ( nVertLevels nCells Time ) 1 o uReconstructMeridional diag - -
# Other diagnostic variables: neither read nor written to any files
var persistent real rv ( nVertLevels nEdges Time ) 1 r rv diag - -
var persistent real circulation ( nVertLevels nVertices Time ) 1 r circulation diag - -
-var persistent real gradPVt ( nVertLevels nEdges Time ) 1 r gradPVt diag - -
-var persistent real gradPVn ( nVertLevels nEdges Time ) 1 r gradPVn diag - -
+var persistent real gradPVt ( nVertLevels nEdges Time ) 1 - gradPVt diag - -
+var persistent real gradPVn ( nVertLevels nEdges Time ) 1 - gradPVn diag - -
var persistent real h_divergence ( nVertLevels nCells Time ) 1 ro h_divergence diag - -
-var persistent real exner ( nVertLevels nCells Time ) 1 r exner diag - -
-var persistent real exner_base ( nVertLevels nCells Time ) 1 or exner_base diag - -
-var persistent real rtheta_base ( nVertLevels nCells Time ) 1 or rtheta_base diag - -
-var persistent real pressure_base ( nVertLevels nCells Time ) 1 or pressure_base diag - -
-var persistent real rho_base ( nVertLevels nCells Time ) 1 or rho_base diag - -
-var persistent real theta_base ( nVertLevels nCells Time ) 1 or theta_base diag - -
+var persistent real exner ( nVertLevels nCells Time ) 1 iro exner diag - -
+var persistent real exner_base ( nVertLevels nCells Time ) 1 iro exner_base diag - -
+var persistent real rtheta_base ( nVertLevels nCells Time ) 1 iro rtheta_base diag - -
+var persistent real pressure_base ( nVertLevels nCells Time ) 1 iro pressure_base diag - -
+var persistent real rho_base ( nVertLevels nCells Time ) 1 iro rho_base diag - -
+var persistent real theta_base ( nVertLevels nCells Time ) 1 iro theta_base diag - -
-var persistent real ruAvg ( nVertLevels nEdges Time ) 1 r ruAvg diag - -
-var persistent real wwAvg ( nVertLevelsP1 nCells Time ) 1 r wwAvg diag - -
-var persistent real qtot ( nVertLevels nCells Time ) 1 r qtot diag - -
-var persistent real cqu ( nVertLevels nEdges Time ) 1 r cqu diag - -
+var persistent real ruAvg ( nVertLevels nEdges Time ) 1 - ruAvg diag - -
+var persistent real wwAvg ( nVertLevelsP1 nCells Time ) 1 - wwAvg diag - -
+var persistent real cqu ( nVertLevels nEdges Time ) 1 - cqu diag - -
var persistent real cqw ( nVertLevels nCells Time ) 1 r cqw diag - -
# coupled variables needed by the solver, but not output...
@@ -219,11 +216,11 @@
var persistent real ru_save ( nVertLevels nEdges Time ) 1 r ru_save diag - -
-var persistent real rw ( nVertLevelsP1 nCells Time ) 1 r rw diag - -
+var persistent real rw ( nVertLevelsP1 nCells Time ) 1 iro rw diag - -
var persistent real rw_p ( nVertLevelsP1 nCells Time ) 1 r rw_p diag - -
var persistent real rw_save ( nVertLevelsP1 nCells Time ) 1 r rw_save diag - -
-var persistent real rtheta_p ( nVertLevels nCells Time ) 1 r rtheta_p diag - -
+var persistent real rtheta_p ( nVertLevels nCells Time ) 1 iro rtheta_p diag - -
var persistent real rtheta_pp ( nVertLevels nCells Time ) 1 r rtheta_pp diag - -
var persistent real rtheta_p_save ( nVertLevels nCells Time ) 1 r rtheta_p_save diag - -
var persistent real rtheta_pp_old ( nVertLevels nCells Time ) 1 r rtheta_pp_old diag - -
@@ -237,14 +234,14 @@
var persistent real qr_old ( nVertLevels nCells Time ) 1 r rqr diag scalars_old moist_old
# Space needed for advection
-var persistent real deriv_two ( FIFTEEN TWO nEdges ) 0 ro deriv_two mesh - -
-var persistent integer advCells ( TWENTYONE nCells ) 0 r advCells mesh - -
+var persistent real deriv_two ( FIFTEEN TWO nEdges ) 0 iro deriv_two mesh - -
+var persistent integer advCells ( TWENTYONE nCells ) 0 iro advCells mesh - -
# Space needed for deformation calculation weights
-var persistent real defc_a ( maxEdges nCells ) 0 r defc_a mesh - -
-var persistent real defc_b ( maxEdges nCells ) 0 r defc_b mesh - -
-var persistent real kdiff ( nVertLevels nCells Time ) 2 r kdiff diag - -
+var persistent real defc_a ( maxEdges nCells ) 0 iro defc_a mesh - -
+var persistent real defc_b ( maxEdges nCells ) 0 iro defc_b mesh - -
+var persistent real kdiff ( nVertLevels nCells Time ) 2 - kdiff diag - -
# Arrays required for reconstruction of velocity field
-var persistent real coeffs_reconstruct ( R3 maxEdges nCells ) 0 r coeffs_reconstruct mesh - -
+var persistent real coeffs_reconstruct ( R3 maxEdges nCells ) 0 iro coeffs_reconstruct mesh - -
Modified: branches/atmos_nonhydrostatic/src/core_nhyd_atmos/module_test_cases.F
===================================================================
--- branches/atmos_nonhydrostatic/src/core_nhyd_atmos/module_test_cases.F        2010-10-19 21:01:25 UTC (rev 570)
+++ branches/atmos_nonhydrostatic/src/core_nhyd_atmos/module_test_cases.F        2010-10-19 21:17:20 UTC (rev 571)
@@ -27,8 +27,14 @@
type (block_type), pointer :: block_ptr
if (config_test_case == 0) then
- write(0,*) ' need nonhydrostatic test case configuration, error stop '
- stop
+ write(0,*) ' Using initial conditions from input file'
+ block_ptr => domain % blocklist
+ do while (associated(block_ptr))
+ do i=2,nTimeLevs
+ call copy_state(block_ptr % state % time_levs(i) % state, block_ptr % state % time_levs(1) % state)
+ end do
+ block_ptr => block_ptr % next
+ end do
else if ((config_test_case == 1) .or. (config_test_case == 2) .or. (config_test_case == 3)) then
write(0,*) ' Jablonowski and Williamson baroclinic wave test case '
Modified: branches/atmos_nonhydrostatic/src/core_nhyd_atmos/module_time_integration.F
===================================================================
--- branches/atmos_nonhydrostatic/src/core_nhyd_atmos/module_time_integration.F        2010-10-19 21:01:25 UTC (rev 570)
+++ branches/atmos_nonhydrostatic/src/core_nhyd_atmos/module_time_integration.F        2010-10-19 21:17:20 UTC (rev 571)
@@ -132,16 +132,8 @@
call dmpar_exch_halo_field2dReal(domain % dminfo, block % diag % ru % array(:,:), &
block % mesh % nVertLevels, block % mesh % nEdges, &
block % parinfo % edgesToSend, block % parinfo % edgesToRecv)
-
+! rw
call dmpar_exch_halo_field2dReal(domain % dminfo, &
- block % diag % cqu % array, &
- block % mesh % nVertLevels, block % mesh % nEdges, &
- block % parinfo % EdgesToSend, block % parinfo % EdgesToRecv)
- call dmpar_exch_halo_field2dReal(domain % dminfo, &
- block % diag % cqw % array, &
- block % mesh % nVertLevels, block % mesh % nCells, &
- block % parinfo % CellsToSend, block % parinfo % CellsToRecv)
- call dmpar_exch_halo_field2dReal(domain % dminfo, &
block % diag % rw % array, &
block % mesh % nVertLevels+1, block % mesh % nCells, &
block % parinfo % CellsToSend, block % parinfo % CellsToRecv)
@@ -174,7 +166,20 @@
block => block % next
end do
+ block => domain % blocklist
+ do while (associated(block))
+ call dmpar_exch_halo_field2dReal(domain % dminfo, &
+ block % diag % cqu % array, &
+ block % mesh % nVertLevels, block % mesh % nEdges, &
+ block % parinfo % EdgesToSend, block % parinfo % EdgesToRecv)
+ call dmpar_exch_halo_field2dReal(domain % dminfo, &
+ block % diag % cqw % array, &
+ block % mesh % nVertLevels, block % mesh % nCells, &
+ block % parinfo % CellsToSend, block % parinfo % CellsToRecv)
+ block => block % next
+ end do
+
if (debug) write(0,*) ' compute_dyn_tend '
block => domain % blocklist
do while (associated(block))
@@ -2972,7 +2977,7 @@
real (kind=RKIND) :: flux, vorticity_abs, h_vertex, workpv, r
integer :: nCells, nEdges, nVertices, nVertLevels
- real (kind=RKIND), dimension(:), pointer :: h_s, fVertex, fEdge, dvEdge, dcEdge, areaCell, areaTriangle
+ real (kind=RKIND), dimension(:), pointer :: fVertex, fEdge, dvEdge, dcEdge, areaCell, areaTriangle
real (kind=RKIND), dimension(:,:), pointer :: vh, weightsOnEdge, kiteAreasOnVertex, h_edge, h, u, v, tend_h, tend_u, &
circulation, vorticity, ke, pv_edge, pv_vertex, pv_cell, gradPVn, gradPVt, &
divergence
@@ -3009,7 +3014,6 @@
dvEdge => grid % dvEdge % array
areaCell => grid % areaCell % array
areaTriangle => grid % areaTriangle % array
- h_s => grid % h_s % array
fVertex => grid % fVertex % array
fEdge => grid % fEdge % array
</font>
</pre>