<p><b>mpetersen@lanl.gov</b> 2012-01-17 15:23:38 -0700 (Tue, 17 Jan 2012)</p><p>Merge trunk to ale branch<br>
</p><hr noshade><pre><font color="gray">
Property changes on: branches/ocean_projects/ale_vert_coord
___________________________________________________________________
Modified: svn:mergeinfo
   - /branches/ocean_projects/imp_vert_mix_mrp:754-986
/branches/ocean_projects/performance:1120-1357
/branches/ocean_projects/split_explicit_mrp:1134-1138
/branches/ocean_projects/split_explicit_timestepping:1044-1097
/branches/ocean_projects/vert_adv_mrp:704-745
/branches/source_renaming:1082-1113
/branches/time_manager:924-962
   + /branches/cam_mpas_nh:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp:754-986
/branches/ocean_projects/performance:1120-1357
/branches/ocean_projects/split_explicit_mrp:1134-1138
/branches/ocean_projects/split_explicit_timestepping:1044-1097
/branches/ocean_projects/vert_adv_mrp:704-745
/branches/source_renaming:1082-1113
/branches/time_manager:924-962
/trunk/mpas:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/namelist.input.ocean
===================================================================
--- branches/ocean_projects/ale_vert_coord/namelist.input.ocean        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/namelist.input.ocean        2012-01-17 22:23:38 UTC (rev 1385)
@@ -45,8 +45,6 @@
    config_visc_vorticity_term = .true.
    config_h_tracer_eddy_diff2 = 1.0e4
    config_h_tracer_eddy_diff4 = 0.0
-   config_mom_decay      = .false.
-   config_mom_decay_time = 3600.0
 /
 &amp;vmix
    config_vert_visc_type  = 'rich'


Property changes on: branches/ocean_projects/ale_vert_coord/src
___________________________________________________________________
Modified: svn:mergeinfo
   - /branches/ocean_projects/imp_vert_mix_mrp/src:754-986
/branches/ocean_projects/performance/src:1120-1357
/branches/ocean_projects/rayleigh/src:1298-1311
/branches/ocean_projects/split_explicit_mrp/src:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src:1044-1097
/branches/ocean_projects/vert_adv_mrp/src:704-745
/branches/source_renaming/src:1082-1113
/branches/time_manager/src:924-962
/trunk/mpas/src:1224-1336
   + /branches/ocean_projects/imp_vert_mix_mrp/src:754-986
/branches/ocean_projects/performance/src:1120-1357
/branches/ocean_projects/rayleigh/src:1298-1311
/branches/ocean_projects/split_explicit_mrp/src:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src:1044-1097
/branches/ocean_projects/vert_adv_mrp/src:704-745
/branches/source_renaming/src:1082-1113
/branches/time_manager/src:924-962
/trunk/mpas/src:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_camrad_init.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_camrad_init.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_camrad_init.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -327,36 +327,6 @@
  DM_BCAST_MACRO(ln_ah2ow)
  DM_BCAST_MACRO(ln_eh2ow)
 
-!101 format(i4,7(1x,e19.12))
-!write(0,*) 'ah2onw'
-!do i_te = 1,21
-!   write(0,101) i_te,(ah2onw(1,1,1,i_te,i_rh),i_rh=1,n_rh)
-!enddo
-!write(0,*) 'eh2onw'
-!do i_te = 1,21
-!   write(0,101) i_te,(eh2onw(1,1,1,i_te,i_rh),i_rh=1,n_rh)
-!enddo
-!write(0,*) 'ah2ow'
-!do i_te = 1,21
-!   write(0,101) i_te,(ah2ow(1,1,1,i_te,i_rh),i_rh=1,n_rh)
-!enddo
-!write(0,*) 'cn_ah2ow'
-!do i_te = 1,21
-!   write(0,101) i_te,(cn_ah2ow(1,1,1,i_te,i_rh),i_rh=1,n_rh)
-!enddo
-!write(0,*) 'cn_eh2ow'
-!do i_te = 1,21
-!   write(0,101) i_te,(cn_eh2ow(1,1,1,i_te,i_rh),i_rh=1,n_rh)
-!enddo
-!write(0,*) 'ln_ah2ow'
-!do i_te = 1,21
-!   write(0,101) i_te,(ln_ah2ow(1,1,1,i_te,i_rh),i_rh=1,n_rh)
-!enddo
-!write(0,*) 'ln_eh2ow'
-!do i_te = 1,21
-!   write(0,101) i_te,(ln_eh2ow(1,1,1,i_te,i_rh),i_rh=1,n_rh)
-!enddo
-
  if(dminfo % my_proc_id == IO_NODE) close(cam_abs_unit)
 
 ! Set up table of H2O saturation vapor pressures for use in calculation effective path RH.  
@@ -612,14 +582,15 @@
  DM_BCAST_MACRO(kcb)
  DM_BCAST_MACRO(wcb)
  DM_BCAST_MACRO(gcb)
- DM_BCAST_MACRO(kvolc)
- DM_BCAST_MACRO(wvolc)
  DM_BCAST_MACRO(kdst)
  DM_BCAST_MACRO(wdst)
  DM_BCAST_MACRO(gdst)
  DM_BCAST_MACRO(kbg)
  DM_BCAST_MACRO(wbg)
  DM_BCAST_MACRO(gbg)
+ DM_BCAST_MACRO(kvolc)
+ DM_BCAST_MACRO(wvolc)
+ DM_BCAST_MACRO(gvolc)
 
  if(dminfo % my_proc_id == IO_NODE) close(cam_aer_unit)
 
@@ -676,7 +647,7 @@
     enddo
  enddo
 
-!write(0,*) '--- end subroutine aer_optics_initialize:'
+ write(0,*) '--- end subroutine aer_optics_initialize:'
 
  end subroutine aer_optics_initialize
 

Modified: branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_date_time.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_date_time.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_date_time.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -108,7 +108,7 @@
 
  write(0,*)
  write(0,*) '--- enter subroutine monthly_interp_to_date:'
- write(0,*) '--- current_date =',date_str
+ write(0,*) '--- current_date  = ',date_str
 
  write(day15,fmt='(I2.2)') 15
  do l = 1 , 12
@@ -141,8 +141,8 @@
              month2 = month1 + 1
           endif
           if(n == 1) then
-             write(0,*) '--- month1=',month1
-             write(0,*) '--- month2=',month2
+             write(0,*) '--- month1 =',month1
+             write(0,*) '--- month2 =',month2
           endif
           field_out(n) = ( field_in(month2,n) * (target_date - middle(l))    &amp;
                        +   field_in(month1,n) * (middle(l+1) - target_date)) &amp;

Modified: branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_driver.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_driver.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_driver.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -71,7 +71,7 @@
 
     !call to short wave radiation scheme:
     if(l_radtsw) then
-       call allocate_radiation_sw
+       call allocate_radiation_sw(xtime_s)
        call driver_radiation_sw(itimestep,block%mesh,block%state%time_levs(1)%state,  &amp;
                                block%diag_physics,block%sfc_input,block%tend_physics, &amp;
                                xtime_s)
@@ -79,7 +79,7 @@
 
     !call to long wave radiation scheme:
     if(l_radtlw) then
-       call allocate_radiation_lw
+       call allocate_radiation_lw(xtime_s)
        call driver_radiation_lw(itimestep,block%mesh,block%state%time_levs(1)%state,  &amp;
                                block%diag_physics,block%sfc_input,block%tend_physics, &amp;
                                xtime_s)
@@ -89,7 +89,7 @@
     if(config_radt_sw_scheme.ne.'off' .or. config_radt_lw_scheme.ne.'off') &amp;
        call deallocate_cloudiness
     if(config_radt_sw_scheme.ne.'off') call deallocate_radiation_sw
-    if(config_radt_lw_scheme.ne.'off') call deallocate_radiation_lw
+    if(config_radt_lw_scheme.ne.'off') call deallocate_radiation_lw(xtime_s)
 
     !call to surface-layer scheme:
     if(config_sfclayer_scheme .ne. 'off') then

Modified: branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_driver_radiation_lw.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_driver_radiation_lw.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_driver_radiation_lw.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -1,5 +1,6 @@
 !=============================================================================================
  module mpas_atmphys_driver_radiation_lw
+ use mpas_configure, only: config_do_restart
  use mpas_grid_types
  use mpas_timer
 
@@ -26,12 +27,17 @@
  contains
 
 !=============================================================================================
- subroutine allocate_radiation_lw
+ subroutine allocate_radiation_lw(xtime_s)
 !=============================================================================================
 
- if(.not.allocated(f_ice)        ) allocate(f_ice(ims:ime,kms:kme,jms:jme)                )
- if(.not.allocated(f_rain)       ) allocate(f_rain(ims:ime,kms:kme,jms:jme)               )
+!input arguments:
+ real(kind=RKIND),intent(in):: xtime_s
 
+!---------------------------------------------------------------------------------------------
+
+ if(.not.allocated(f_ice)        ) allocate(f_ice(ims:ime,kms:kme,jms:jme)        )
+ if(.not.allocated(f_rain)       ) allocate(f_rain(ims:ime,kms:kme,jms:jme)       )
+
  if(.not.allocated(sfc_emiss_p)  ) allocate(sfc_emiss_p(ims:ime,jms:jme)          )
  if(.not.allocated(snow_p)       ) allocate(snow_p(ims:ime,jms:jme)               )
  if(.not.allocated(tsk_p)        ) allocate(tsk_p(ims:ime,jms:jme)                )
@@ -93,12 +99,21 @@
        if(.not.allocated(aerosolcp_p)  ) &amp;
           allocate(aerosolcp_p(ims:ime,1:num_aerlevels,jms:jme,num_aerosols) )
 
-       if(.not.allocated(emstot_p) ) allocate(emstot_p(ims:ime,kms:kme,jms:jme) )
-       if(.not.allocated(abstot_p) ) &amp;
-          allocate(abstot_p(ims:ime,kms:kme,cam_abs_dim2,jms:jme) )
-       if(.not.allocated(absnxt_p) ) &amp;
-          allocate(absnxt_p(ims:ime,kms:kme,cam_abs_dim1,jms:jme) )
+       !allocate these arrays on the first time step, only:
+       if(xtime_s .lt. 1.e-12) then
 
+          write(0,*) 
+          write(0,*) '--- end subroutine allocate_radiation_lw:'
+          write(0,*) '--- allocate emstot_p,abstot_p,absnxt_p'
+
+          if(.not.allocated(emstot_p) ) allocate(emstot_p(ims:ime,kms:kme,jms:jme) )
+          if(.not.allocated(abstot_p) ) &amp;
+             allocate(abstot_p(ims:ime,kms:kme,cam_abs_dim2,jms:jme) )
+          if(.not.allocated(absnxt_p) ) &amp;
+             allocate(absnxt_p(ims:ime,kms:kme,cam_abs_dim1,jms:jme) )
+
+       endif
+
     case default
 
  end select radiation_lw_select
@@ -106,9 +121,14 @@
  end subroutine allocate_radiation_lw
 
 !=============================================================================================
- subroutine deallocate_radiation_lw
+ subroutine deallocate_radiation_lw(xtime_s)
 !=============================================================================================
 
+!input arguments:
+ real(kind=RKIND),intent(in):: xtime_s
+
+!---------------------------------------------------------------------------------------------
+
  if(allocated(f_ice)        ) deallocate(f_ice        )
  if(allocated(f_rain)       ) deallocate(f_rain       )
  if(allocated(sfc_emiss_p)  ) deallocate(sfc_emiss_p  )
@@ -169,10 +189,6 @@
        if(allocated(aerosolcn_p)  ) deallocate(aerosolcn_p  )
        if(allocated(aerosolcp_p)  ) deallocate(aerosolcp_p  )
 
-       if(allocated(abstot_p)     ) deallocate(abstot_p     )
-       if(allocated(absnxt_p)     ) deallocate(absnxt_p     )
-       if(allocated(emstot_p)     ) deallocate(emstot_p     )
-
     case default
 
  end select radiation_lw_select
@@ -180,7 +196,7 @@
  end subroutine deallocate_radiation_lw
 
 !=============================================================================================
- subroutine radiation_lw_from_MPAS(mesh,state,diag_physics,sfc_input)
+ subroutine radiation_lw_from_MPAS(mesh,state,diag_physics,sfc_input,xtime_s)
 !=============================================================================================
 
 !input arguments:
@@ -189,6 +205,8 @@
  type(sfc_input_type)   ,intent(in):: sfc_input
  type(diag_physics_type),intent(in):: diag_physics
 
+ real(kind=RKIND),intent(in):: xtime_s
+
 !---------------------------------------------------------------------------------------------
 
  do j = jts,jte
@@ -283,42 +301,34 @@
        enddo
        enddo
 
-       call mpas_timer_start(&quot;CAM lw: read arrays for infrared absorption&quot;)
-       !infrared absorption:
-       do j = jts,jte
-       do n = 1,cam_abs_dim1
-       do k = kts,kte
-       do i = its,ite
-          absnxt_p(i,k,n,j) = diag_physics % absnxt % array(k,n,i)
-       enddo
-       enddo
-       enddo
-       do n = 1,cam_abs_dim2
-       do k = kts,kte+1
-       do i = its,ite
-          abstot_p(i,k,n,j) = diag_physics % abstot % array(k,n,i)
-       enddo
-       enddo
-       enddo
-       do k = kts,kte+1
-       do i = its,ite
-          emstot_p(i,k,j) = diag_physics % emstot % array(k,i)
-       enddo
-       enddo
-       enddo
-       call mpas_timer_stop(&quot;CAM lw: read arrays for infrared absorption&quot;)
-!      write(0,*) '--- end radiation_lw_from_MPAS: doabsems=',doabsems
-!      do k = kts,kte+1
-!         write(0,102) k,(maxval(abstot_p(its:ite,k,n,jts:jte)),n=1,10)
-!      enddo
-!      write(0,*)
-!      do k = kts,kte
-!         write(0,102) k,(maxval(absnxt_p(its:ite,k,n,jts:jte)),n=1,cam_abs_dim1)
-!      enddo
-!      write(0,*)
-!      do k = kts,kte+1
-!         write(0,102) k,maxval(emstot_p(its:ite,k,jts:jte))
-!      enddo
+       call mpas_timer_start(&quot;CAM lw: fill arrays for infrared absorption&quot;)
+       if(xtime_s .lt. 1.e-12) then
+          write(0,*)
+          write(0,*) '--- radiation_lw_from_MPAS:'
+          write(0,*) '--- initialize emstot_p,abstot_p,absnxt_p'
+          !infrared absorption:
+          do j = jts,jte
+          do n = 1,cam_abs_dim1
+          do k = kts,kte
+          do i = its,ite
+             absnxt_p(i,k,n,j) = 0.
+          enddo
+          enddo
+          enddo
+          do n = 1,cam_abs_dim2
+          do k = kts,kte+1
+          do i = its,ite
+             abstot_p(i,k,n,j) = 0.
+          enddo
+          enddo
+          enddo
+          do k = kts,kte+1
+          do i = its,ite
+             emstot_p(i,k,j) = 0.
+          enddo
+          enddo
+          enddo
+       endif
 
        call mpas_timer_start(&quot;CAM lw: ozone and aerosols&quot;)
        !ozone mixing ratio:
@@ -367,13 +377,15 @@
  end subroutine radiation_lw_from_MPAS
 
 !=============================================================================================
- subroutine radiation_lw_to_MPAS(diag_physics,tend_physics)
+ subroutine radiation_lw_to_MPAS(diag_physics,tend_physics,xtime_s)
 !=============================================================================================
 
 !input arguments:
  type(diag_physics_type),intent(inout):: diag_physics
  type(tend_physics_type),intent(inout):: tend_physics
 
+ real(kind=RKIND),intent(in):: xtime_s
+
 !---------------------------------------------------------------------------------------------
 
  do j = jts,jte
@@ -407,48 +419,6 @@
  enddo
  enddo
 
- radiation_lw_select: select case (trim(radt_lw_scheme))
-
-    case(&quot;cam_lw&quot;)
-       !infrared absorption:
-       do j = jts,jte
-       do n = 1,cam_abs_dim1
-       do k = kts,kte
-       do i = its,ite
-          diag_physics % absnxt % array(k,n,i) = absnxt_p(i,k,n,j)
-       enddo
-       enddo
-       enddo
-       do n = 1,cam_abs_dim2
-       do k = kts,kte+1
-       do i = its,ite
-          diag_physics % abstot % array(k,n,i) = abstot_p(i,k,n,j)
-       enddo
-       enddo
-       enddo
-       do k = kts,kte+1
-       do i = its,ite
-          diag_physics % emstot % array(k,i) = emstot_p(i,k,j)
-       enddo
-       enddo
-       enddo
-!      write(0,*) '--- end radiation_lw_to_MPAS: doabsems=',doabsems
-!      do k = kts,kte+1
-!         write(0,102) k,(maxval(abstot_p(its:ite,k,n,jts:jte)),n=1,10)
-!      enddo
-!      write(0,*)
-!      do k = kts,kte
-!        write(0,102) k,(maxval(absnxt_p(its:ite,k,n,jts:jte)),n=1,cam_abs_dim1)
-!      enddo
-!      write(0,*)
-!      do k = kts,kte+1
-!         write(0,102) k,maxval(emstot_p(its:ite,k,jts:jte))
-!      enddo
-
-    case default
-
- end select radiation_lw_select
-
 !format:
  101 format(i3,2i6,12(1x,e15.8))
  102 format(i6,12(1x,e15.8))
@@ -517,7 +487,7 @@
  101 format(i8,12(1x,e15.8))
 
 !copy all MPAS arrays to rectangular grid:
- call radiation_lw_from_MPAS(mesh,state,diag_physics,sfc_input)
+ call radiation_lw_from_MPAS(mesh,state,diag_physics,sfc_input,xtime_s)
 
 !call to longwave radiation scheme:
  radiation_lw_select: select case (trim(radt_lw_scheme))
@@ -560,6 +530,8 @@
        radt = dt_radtlw/60.
 
        call mpas_timer_start(&quot;camrad&quot;)
+       write(0,*) '--- enter subroutine camrad_lw: doabsems=',doabsems
+       call mpas_timer_start(&quot;camrad&quot;)
        call camrad( dolw = .true. , dosw = .false. ,                                         &amp;
                 rthratenlw    = rthratenlw_p  , rthratensw    = rthratensw_p  ,              &amp;
                 swupt         = swupt_p       , swuptc        = swuptc_p      ,              &amp;
@@ -579,7 +551,7 @@
                 xlat          = xlat_p        , xlong         = xlon_p        ,              &amp;
                 t_phy         = t_p           , pi_phy        = pi_p          ,              &amp;
                 p_phy         = pres_p        , p8w           = pres2_p       ,              &amp;
-                z             = z_p           , dz8w          = dz_p          ,              &amp;            
+                z             = zmid_p        , dz8w          = dz_p          ,              &amp;            
                 rho_phy       = rho_p         , qv3d          = qv_p          ,              &amp; 
                 qc3d          = qc_p          , qr3d          = qr_p          ,              &amp;
                 qi3d          = qi_p          , qs3d          = qs_p          ,              &amp;
@@ -609,28 +581,27 @@
                 its = its , ite = ite , jts = jts , jte = jte , kts = kts , kte = kte        &amp;
                   )
        call mpas_timer_stop(&quot;camrad&quot;)
-       write(0,*) '--- end subroutince camrad lw: doabsems=',doabsems
        write(0,*) 'max lwupb      =',maxval(lwupb_p(its:ite,jms:jme))
-       write(0,*) 'min lwupb      =',minval(lwupb_p(its:ite,jms:jme))
        write(0,*) 'max lwupbc     =',maxval(lwupbc_p(its:ite,jms:jme))
-       write(0,*) 'min lwupbc     =',minval(lwupbc_p(its:ite,jms:jme))
        write(0,*) 'max lwupt      =',maxval(lwupt_p(its:ite,jms:jme))
-       write(0,*) 'min lwupt      =',minval(lwupt_p(its:ite,jms:jme))
        write(0,*) 'max lwuptc     =',maxval(lwuptc_p(its:ite,jms:jme))
-       write(0,*) 'min lwuptc     =',minval(lwuptc_p(its:ite,jms:jme))
        write(0,*) 'max rthratenlw =',maxval(rthratenlw_p(its:ite,kts:kte,jms:jme))
        write(0,*) 'min rthratenlw =',minval(rthratenlw_p(its:ite,kts:kte,jms:jme))
+       write(0,*) '--- end subroutine camrad lw: doabsems=',doabsems
 
     case default
 
  end select radiation_lw_select
 
 !copy all arrays back to MPAS geodesic grid:
- call radiation_lw_to_MPAS(diag_physics,tend_physics)
+ call radiation_lw_to_MPAS(diag_physics,tend_physics,xtime_s)
 
- write(0,*) '--- end subroutine driver_radiation_lw:'
+ write(0,*) '--- end subroutine driver_radiation_lw'
  call mpas_timer_stop(&quot;radiation_lw&quot;)
 
+!formats:
+ 200 format(i3,i3,8(1x,e15.8))
+
  end subroutine driver_radiation_lw
 
 !=============================================================================================

Modified: branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_driver_radiation_sw.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_driver_radiation_sw.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_driver_radiation_sw.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -27,12 +27,17 @@
  contains
 
 !=============================================================================================
- subroutine allocate_radiation_sw
+ subroutine allocate_radiation_sw(xtime_s)
 !=============================================================================================
 
- if(.not.allocated(f_ice)        ) allocate(f_ice(ims:ime,kms:kme,jms:jme)                )
- if(.not.allocated(f_rain)       ) allocate(f_rain(ims:ime,kms:kme,jms:jme)               )
+!input arguments:
+ real(kind=RKIND),intent(in):: xtime_s
 
+!---------------------------------------------------------------------------------------------
+
+ if(.not.allocated(f_ice)        ) allocate(f_ice(ims:ime,kms:kme,jms:jme)        )
+ if(.not.allocated(f_rain)       ) allocate(f_rain(ims:ime,kms:kme,jms:jme)       )
+
  if(.not.allocated(xlat_p)       ) allocate(xlat_p(ims:ime,jms:jme)               )
  if(.not.allocated(xlon_p)       ) allocate(xlon_p(ims:ime,jms:jme)               )
 
@@ -65,9 +70,6 @@
        if(.not.allocated(swupflxc_p)   ) allocate(swupflxc_p(ims:ime,kms:kme+1,jms:jme) )
 
     case(&quot;cam_sw&quot;)
-       if(.not.allocated(xlat_p)       ) allocate(xlat_p(ims:ime,jms:jme)               )
-       if(.not.allocated(xlon_p)       ) allocate(xlon_p(ims:ime,jms:jme)               )
-
        if(.not.allocated(glw_p)        ) allocate(glw_p(ims:ime,jms:jme)                )
        if(.not.allocated(lwcf_p)       ) allocate(lwcf_p(ims:ime,jms:jme)               )
        if(.not.allocated(lwdnb_p)      ) allocate(lwdnb_p(ims:ime,jms:jme)              )
@@ -98,12 +100,17 @@
        if(.not.allocated(aerosolcp_p)  ) &amp;
           allocate(aerosolcp_p(ims:ime,1:num_aerlevels,jms:jme,num_aerosols) )
 
-       if(.not.allocated(emstot_p) ) allocate(emstot_p(ims:ime,kms:kme,jms:jme) )
-       if(.not.allocated(abstot_p) ) &amp;
-          allocate(abstot_p(ims:ime,kms:kme,cam_abs_dim2,jms:jme) )
-       if(.not.allocated(absnxt_p) ) &amp;
-          allocate(absnxt_p(ims:ime,kms:kme,cam_abs_dim1,jms:jme) )
+       !allocate these arrays on the first time step, only:
+       if(xtime_s .lt. 1.e-12) then
 
+          if(.not.allocated(emstot_p) ) allocate(emstot_p(ims:ime,kms:kme,jms:jme) )
+          if(.not.allocated(abstot_p) ) &amp;
+             allocate(abstot_p(ims:ime,kms:kme,cam_abs_dim2,jms:jme) )
+          if(.not.allocated(absnxt_p) ) &amp;
+             allocate(absnxt_p(ims:ime,kms:kme,cam_abs_dim1,jms:jme) )
+
+       endif
+
     case default
 
  end select radiation_sw_select
@@ -176,10 +183,6 @@
        if(allocated(aerosolcn_p)  ) deallocate(aerosolcn_p  )
        if(allocated(aerosolcp_p)  ) deallocate(aerosolcp_p  )
 
-       if(allocated(abstot_p)     ) deallocate(abstot_p     )
-       if(allocated(absnxt_p)     ) deallocate(absnxt_p     )
-       if(allocated(emstot_p)     ) deallocate(emstot_p     )
-
     case default
 
  end select radiation_sw_select
@@ -187,7 +190,7 @@
  end subroutine deallocate_radiation_sw
 
 !=============================================================================================
- subroutine radiation_sw_from_MPAS(mesh,state,diag_physics,sfc_input)
+ subroutine radiation_sw_from_MPAS(mesh,state,diag_physics,sfc_input,xtime_s)
 !=============================================================================================
 
 !input arguments:
@@ -196,6 +199,8 @@
  type(diag_physics_type),intent(inout):: diag_physics
  type(sfc_input_type)   ,intent(inout):: sfc_input
 
+ real(kind=RKIND),intent(in):: xtime_s
+
 !---------------------------------------------------------------------------------------------
 
  do j = jts,jte
@@ -265,7 +270,6 @@
        enddo
 
     case(&quot;cam_sw&quot;)
-
        do j = jts,jte
        do i = its,ite
           sfc_emiss_p(i,j) = diag_physics % sfc_emiss % array(i)
@@ -292,39 +296,29 @@
        enddo
        enddo
        !infrared absorption:
-       do j = jts,jte
-       do n = 1,cam_abs_dim1
-       do k = kts,kte
-       do i = its,ite
-          absnxt_p(i,k,n,j) = diag_physics % absnxt % array(k,n,i)
-       enddo
-       enddo
-       enddo
-       do n = 1,cam_abs_dim2
-       do k = kts,kte+1
-       do i = its,ite
-          abstot_p(i,k,n,j) = diag_physics % abstot % array(k,n,i)
-       enddo
-       enddo
-       enddo
-       do k = kts,kte+1
-       do i = its,ite
-          emstot_p(i,k,j) = diag_physics % emstot % array(k,i)
-       enddo
-       enddo
-       enddo
-!      write(0,*) '--- end radiation_sw_from_MPAS: doabsems=',doabsems
-!      do k = kts,kte+1
-!         write(0,102) k,(maxval(abstot_p(its:ite,k,n,jts:jte)),n=1,10)
-!      enddo
-!      write(0,*)
-!      do k = kts,kte
-!         write(0,102) k,(maxval(absnxt_p(its:ite,k,n,jts:jte)),n=1,cam_abs_dim1)
-!      enddo
-!      write(0,*)
-!      do k = kts,kte+1
-!         write(0,102) k,maxval(emstot_p(its:ite,k,jts:jte))
-!      enddo
+       if(xtime_s .lt. 1.e-12) then
+          do j = jts,jte
+          do n = 1,cam_abs_dim1
+          do k = kts,kte
+          do i = its,ite
+             absnxt_p(i,k,n,j) = 0.
+          enddo
+          enddo
+          enddo
+          do n = 1,cam_abs_dim2
+          do k = kts,kte+1
+          do i = its,ite
+             abstot_p(i,k,n,j) = 0.
+          enddo
+          enddo
+          enddo
+          do k = kts,kte+1
+          do i = its,ite
+             emstot_p(i,k,j) = 0.
+          enddo
+          enddo
+          enddo
+       endif
        !ozone mixing ratio:
        do k = 1, num_oznlevels
           pin_p(k) = mesh % pin % array(k,1)
@@ -412,35 +406,6 @@
 
  enddo
 
- radiation_sw_select: select case (trim(radt_sw_scheme))
-    case(&quot;cam_sw&quot;)
-       !infrared absorption:
-       do j = jts,jte
-       do n = 1,cam_abs_dim1
-       do k = kts,kte
-       do i = its,ite
-          diag_physics % absnxt % array(k,n,i) = absnxt_p(i,k,n,j)
-       enddo
-       enddo
-       enddo
-       do n = 1,cam_abs_dim2
-       do k = kts,kte+1
-       do i = its,ite
-          diag_physics % abstot % array(k,n,i) = abstot_p(i,k,n,j)
-       enddo
-       enddo
-       enddo
-       do k = kts,kte+1
-       do i = its,ite
-          diag_physics % emstot % array(k,i) = emstot_p(i,k,j)
-       enddo
-       enddo
-       enddo
-
-    case default
-
- end select radiation_sw_select
-
  end subroutine radiation_sw_to_MPAS
 
 !=============================================================================================
@@ -456,7 +421,7 @@
 
 !---------------------------------------------------------------------------------------------
  write(0,*)
- write(0,*) '--- begin radiation_sw initialization:'
+ write(0,*) '--- enter radiation_sw initialization:'
 
 !call to shortwave radiation scheme:
  radiation_sw_select: select case (trim(radt_sw_scheme))
@@ -475,7 +440,7 @@
 
  end select radiation_sw_select
 
- write(0,*) '--- end radiation_sw initialization:'
+ write(0,*) '--- end radiation_sw initialization'
 
  end subroutine init_radiation_sw
 
@@ -517,7 +482,7 @@
   xtime_m = xtime_s/60.
 
 !copy all MPAS arrays to rectangular grid:
- call radiation_sw_from_MPAS(mesh,state,diag_physics,sfc_input)
+ call radiation_sw_from_MPAS(mesh,state,diag_physics,sfc_input,xtime_s)
 
 !... calculates solar declination:
 !call radconst(declin,solcon,julday,degrad,dpd)
@@ -531,6 +496,7 @@
  radiation_sw_select: select case (trim(radt_sw_scheme))
 
     case (&quot;rrtmg_sw&quot;)
+
        write(0,*) '--- enter subroutine rrtmg_swrad:'
        call rrtmg_swrad( &amp;
                 rthratensw = rthratensw_p , swupt     = swupt_p    , swuptc    = swuptc_p  , &amp;
@@ -562,6 +528,8 @@
        write(0,*) '--- exit subroutine rrtmg_swrad'
 
     case (&quot;cam_sw&quot;)
+
+       write(0,*) '--- enter subroutine camrad_sw:'
        call camrad( dolw = .false. , dosw = .true. ,                                         &amp;
                 rthratenlw    = rthratenlw_p  , rthratensw    = rthratensw_p  ,              &amp;
                 swupt         = swupt_p       , swuptc        = swuptc_p      ,              &amp;
@@ -581,7 +549,7 @@
                 xlat          = xlat_p        , xlong         = xlon_p        ,              &amp;
                 t_phy         = t_p           , pi_phy        = pi_p          ,              &amp;
                 p_phy         = pres_p        , p8w           = pres2_p       ,              &amp;
-                z             = z_p           , dz8w          = dz_p          ,              &amp;            
+                z             = zmid_p        , dz8w          = dz_p          ,              &amp;            
                 rho_phy       = rho_p         , qv3d          = qv_p          ,              &amp; 
                 qc3d          = qc_p          , qr3d          = qr_p          ,              &amp;
                 qi3d          = qi_p          , qs3d          = qs_p          ,              &amp;
@@ -610,19 +578,14 @@
                 ims = ims , ime = ime , jms = jms , jme = jme , kms = kms , kme = kme ,      &amp;
                 its = its , ite = ite , jts = jts , jte = jte , kts = kts , kte = kte        &amp;
                   )
-
-       write(0,*) '--- end subroutine camrad sw'
        write(0,*) 'doabsems =',doabsems
        write(0,*) 'max swupb      =',maxval(swupb_p(its:ite,jms:jme))
-       write(0,*) 'min swupb      =',minval(swupb_p(its:ite,jms:jme))
        write(0,*) 'max swupbc     =',maxval(swupbc_p(its:ite,jms:jme))
-       write(0,*) 'min swupbc     =',minval(swupbc_p(its:ite,jms:jme))
        write(0,*) 'max swupt      =',maxval(swupt_p(its:ite,jms:jme))
-       write(0,*) 'min swupt      =',minval(swupt_p(its:ite,jms:jme))
        write(0,*) 'max swuptc     =',maxval(swuptc_p(its:ite,jms:jme))
-       write(0,*) 'min swuptc     =',minval(swuptc_p(its:ite,jms:jme))
        write(0,*) 'max rthratensw =',maxval(rthratensw_p(its:ite,kts:kte,jms:jme))
        write(0,*) 'min rthratensw =',minval(rthratensw_p(its:ite,kts:kte,jms:jme))
+       write(0,*) '--- end subroutine camrad sw'
 
     case default
 
@@ -631,8 +594,11 @@
 !copy all arrays back to MPAS geodesic grid:
  call radiation_sw_to_MPAS(diag_physics,tend_physics)
 
- write(0,*) '--- end subroutine driver_radiation_sw:'
+ write(0,*) '--- end subroutine driver_radiation_sw'
 
+!formats:
+ 200 format(i3,i6,8(1x,e15.8))
+
  end subroutine driver_radiation_sw
 
 !=============================================================================================

Modified: branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_initialize_real.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_initialize_real.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_initialize_real.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -425,8 +425,8 @@
  nCells        = mesh % nCells
  nSoilLevels   = mesh % nSoilLevels
  nFGSoilLevels = mesh % nFGSoilLevels
- write(0,*) 'nSoilLevels  =',nSoilLevels
- write(0,*) 'nFGSoilLevels=',nFGSoilLevels
+ write(0,*) 'nSoilLevels   =',nSoilLevels
+ write(0,*) 'nFGSoilLevels =',nFGSoilLevels
 
  landmask =&gt; mesh % landmask % array
 
@@ -606,7 +606,7 @@
  endif
 
 !make sure that all the cells flagged as sea-ice cells are defined as ocean cells:
-!num_seaice_changes = 0
+ num_seaice_changes = 0
  do iCell = 1, nCellsSolve
     if((landmask(iCell).eq.1 .and. xice(iCell).gt.0.) .or. xice(iCell).gt.200.) then
        num_seaice_changes = num_seaice_changes + 1

Modified: branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_interface_nhyd.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_interface_nhyd.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_interface_nhyd.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -23,15 +23,18 @@
  subroutine allocate_forall_physics
 !=============================================================================================
 
- if(.not.allocated(psfc_p))  allocate(psfc_p(ims:ime,jms:jme)           )
- if(.not.allocated(ptop_p))  allocate(ptop_p(ims:ime,jms:jme)           )
+ if(.not.allocated(psfc_p) )   allocate(psfc_p(ims:ime,jms:jme)         )
+ if(.not.allocated(ptop_p) )   allocate(ptop_p(ims:ime,jms:jme)         )
 
  if(.not.allocated(u_p)    )  allocate(u_p(ims:ime,kms:kme,jms:jme)     )
  if(.not.allocated(v_p)    )  allocate(v_p(ims:ime,kms:kme,jms:jme)     )
+ if(.not.allocated(fzm_p)  )   allocate(fzm_p(ims:ime,kms:kme,jms:jme)  )
+ if(.not.allocated(fzp_p)  )   allocate(fzp_p(ims:ime,kms:kme,jms:jme)  )
  if(.not.allocated(zz_p)   )  allocate(zz_p(ims:ime,kms:kme,jms:jme)    )
  if(.not.allocated(pres_p) )  allocate(pres_p(ims:ime,kms:kme,jms:jme)  )
  if(.not.allocated(pi_p)   )  allocate(pi_p(ims:ime,kms:kme,jms:jme)    )
  if(.not.allocated(z_p)    )  allocate(z_p(ims:ime,kms:kme,jms:jme)     )
+ if(.not.allocated(zmid_p) )  allocate(zmid_p(ims:ime,kms:kme,jms:jme)  )
  if(.not.allocated(dz_p)   )  allocate(dz_p(ims:ime,kms:kme,jms:jme)    )
  if(.not.allocated(t_p)    )  allocate(t_p(ims:ime,kms:kme,jms:jme)     )
  if(.not.allocated(th_p)   )  allocate(th_p(ims:ime,kms:kme,jms:jme)    )
@@ -42,44 +45,7 @@
  if(.not.allocated(w_p)    )  allocate(w_p(ims:ime,kms:kme,jms:jme)     )
  if(.not.allocated(pres2_p))  allocate(pres2_p(ims:ime,kms:kme,jms:jme) )
  if(.not.allocated(t2_p)   )  allocate(t2_p(ims:ime,kms:kme,jms:jme)    )
-
- if(.not.allocated(pres_hyd_p) ) allocate(pres_hyd_p(ims:ime,kms:kme,jms:jme)  )
- if(.not.allocated(pres2_hyd_p)) allocate(pres2_hyd_p(ims:ime,kms:kme,jms:jme) )
  
- do j = jms,jme
- do i = ims,ime
-    psfc_p(i,j)   = 0.
-    ptop_p(i,j)   = 0.
- enddo
- enddo
-
- do j = jms,jme
- do k = kms,kme
- do i = ims,ime
-    u_p(i,k,j)     = 0.
-    v_p(i,k,j)     = 0.
-    w_p(i,k,j)     = 0.
-    pres_p(i,k,j)  = 0.
-    pi_p(i,k,j)    = 0.
-    z_p(i,k,j)     = 0.
-    dz_p(i,k,j)    = 0.
-    t_p(i,k,j)     = 0.
-    th_p(i,k,j)    = 0.
-    al_p(i,k,j)    = 0.
-    rho_p(i,k,j)   = 0.
-    rh_p(i,k,j)    = 0.
-
-    w_p(i,k,j)     = 0.
-    pres2_p(i,k,j) = 0.
-    t2_p(i,k,j)    = 0.
-
-    pres_hyd_p(i,k,j)  = 0.
-    pres2_hyd_p(i,k,j) = 0.
- enddo
- enddo
- enddo
-
-!allocate moist species (to be revisited!):
  if(.not.allocated(qv_p)   ) allocate(qv_p(ims:ime,kms:kme,jms:jme)     )
  if(.not.allocated(qc_p)   ) allocate(qc_p(ims:ime,kms:kme,jms:jme)     )
  if(.not.allocated(qr_p)   ) allocate(qr_p(ims:ime,kms:kme,jms:jme)     )
@@ -93,16 +59,18 @@
  subroutine deallocate_forall_physics
 !=============================================================================================
 
-!de-allocation of all physics arrays:
  if(allocated(psfc_p)  ) deallocate(psfc_p  )
  if(allocated(ptop_p)  ) deallocate(ptop_p  )
 
  if(allocated(u_p)     ) deallocate(u_p     )
  if(allocated(v_p)     ) deallocate(v_p     )
+ if(allocated(fzm_p)   ) deallocate(fzm_p   )
+ if(allocated(fzp_p)   ) deallocate(fzp_p   )
  if(allocated(zz_p)    ) deallocate(zz_p    )
  if(allocated(pres_p)  ) deallocate(pres_p  )
  if(allocated(pi_p)    ) deallocate(pi_p    )
  if(allocated(z_p)     ) deallocate(z_p     )
+ if(allocated(zmid_p)  ) deallocate(zmid_p  )
  if(allocated(dz_p)    ) deallocate(dz_p    )
  if(allocated(t_p)     ) deallocate(t_p     )
  if(allocated(th_p)    ) deallocate(th_p    )
@@ -114,9 +82,6 @@
  if(allocated(pres2_p) ) deallocate(pres2_p )
  if(allocated(t2_p)    ) deallocate(t2_p    )
 
- if(allocated(pres_hyd_p) ) deallocate(pres_hyd_p  )
- if(allocated(pres2_hyd_p)) deallocate(pres2_hyd_p )

  if(allocated(qv_p)    ) deallocate(qv_p    )
  if(allocated(qc_p)    ) deallocate(qc_p    )
  if(allocated(qr_p)    ) deallocate(qr_p    )
@@ -141,12 +106,13 @@
 
  real(kind=RKIND),dimension(:),pointer:: latCell,lonCell
  real(kind=RKIND),dimension(:),pointer  :: fzm,fzp,rdzw
+ real(kind=RKIND),dimension(:),pointer  :: sfc_pressure
  real(kind=RKIND),dimension(:,:),pointer:: zgrid
  real(kind=RKIND),dimension(:,:),pointer:: zz,exner,pressure_b,rtheta_p,rtheta_b
  real(kind=RKIND),dimension(:,:),pointer:: rho_zz,theta_m,qv,pressure_p,u,v,w
  real(kind=RKIND),dimension(:,:),pointer:: qvs,rh
 
- integer:: ip,iEdg
+ real(kind=RKIND):: rho1,rho2,tem1,tem2
 
 !---------------------------------------------------------------------------------------------
 
@@ -164,25 +130,50 @@
  latCell   =&gt; mesh % latCell % array
  lonCell   =&gt; mesh % lonCell % array
 
- fzm        =&gt; mesh % fzm % array
- fzp        =&gt; mesh % fzp % array
- rdzw       =&gt; mesh % rdzw % array
- zgrid      =&gt; mesh % zgrid % array
- zz         =&gt; mesh % zz % array
- exner      =&gt; diag % exner % array
- pressure_b =&gt; diag % pressure_base % array
- pressure_p =&gt; diag % pressure_p % array
- rtheta_p   =&gt; diag % rtheta_p % array
- rtheta_b   =&gt; diag % rtheta_base % array
+ fzm          =&gt; mesh % fzm % array
+ fzp          =&gt; mesh % fzp % array
+ rdzw         =&gt; mesh % rdzw % array
+ zgrid        =&gt; mesh % zgrid % array
+ zz           =&gt; mesh % zz % array
+ sfc_pressure =&gt; diag % surface_pressure % array
+ exner        =&gt; diag % exner % array
+ pressure_b   =&gt; diag % pressure_base % array
+ pressure_p   =&gt; diag % pressure_p % array
+ rtheta_p     =&gt; diag % rtheta_p % array
+ rtheta_b     =&gt; diag % rtheta_base % array
 
- rho_zz     =&gt; state % rho_zz % array
- theta_m    =&gt; state % theta_m % array
- qv         =&gt; state % scalars % array(state%index_qv,:,:)
+ rho_zz  =&gt; state % rho_zz % array
+ theta_m =&gt; state % theta_m % array
+ qv      =&gt; state % scalars % array(state%index_qv,:,:)
 
- w          =&gt; state % w % array
- u          =&gt; diag  % uReconstructZonal % array
- v          =&gt; diag  % uReconstructMeridional % array
+ w =&gt; state % w % array
+ u =&gt; diag  % uReconstructZonal % array
+ v =&gt; diag  % uReconstructMeridional % array
 
+!ldf (2012-01-06): updates the surface pressure as is done in subroutine microphysics_to_MPAS.
+!do j = jts,jte
+!do i = its,ite
+!   sfc_pressure(i) = 0.5*g*(zgrid(2,i)-zgrid(1,i)) &amp;
+!                   * (1.25 * rho_zz(1,i) * zz(1,i) * (1. + qv(1,i))  &amp;
+!                   -  0.25 * rho_zz(2,i) * zz(2,i) * (1. + qv(1,i)))
+!   sfc_pressure(i) = sfc_pressure(i) + pressure_p(1,i) + pressure_b(1,i)
+!enddo
+!enddo
+!ldf end.
+!ldf (2012-01-09): updates the surface pressure using zgrid.
+ do j = jts,jte
+ do i = its,ite
+    tem1 = zgrid(2,i)-zgrid(1,i)
+    tem2 = zgrid(3,i)-zgrid(2,i)
+    rho1 = rho_zz(1,i) * zz(1,i) * (1. + qv(1,i))
+    rho2 = rho_zz(2,i) * zz(2,i) * (1. + qv(2,i))
+    sfc_pressure(i) = 0.5*g*(zgrid(2,i)-zgrid(1,i)) &amp;
+                    * (rho1 + 0.5*(rho2-rho1)*tem1/(tem1+tem2))
+    sfc_pressure(i) = sfc_pressure(i) + pressure_p(1,i) + pressure_b(1,i)
+ enddo
+ enddo
+!ldf end.
+
 !copy sounding variables from the geodesic grid to the rectangular grid:
  do j = jts, jte
  do i = its, ite
@@ -215,11 +206,18 @@
     pi_p(i,k,j)   = exner(k,i)
     pres_p(i,k,j) = pressure_p(k,i) + pressure_b(k,i)
 
-    dz_p(i,k,j) = zgrid(k+1,i)-zgrid(k,i)
+    zmid_p(i,k,j) = 0.5*(zgrid(k+1,i)+zgrid(k,i))
+    dz_p(i,k,j)   = zgrid(k+1,i)-zgrid(k,i)
 
-    !arrays located at w points:
+ enddo
+ enddo
+ enddo
+
+ do j = jts, jte
+ do k = kts,kte+1
+ do i = its,ite
     w_p(i,k,j) = w(k,i)
-   
+    z_p(i,k,j) = zgrid(k,i)
  enddo
  enddo
  enddo
@@ -247,11 +245,23 @@
  enddo
 
 !interpolation of pressure and temperature from theta points to w points:
+!do j = jts,jte
+!do k = kts+1,kte
+!do i = its,ite
+!   t2_p(i,k,j)    = fzm(k)*t_p(i,k,j) + fzp(k)*t_p(i,k-1,j)
+!   pres2_p(i,k,j) = fzm(k)*pres_p(i,k,j) + fzp(k)*pres_p(i,k-1,j)
+!enddo
+!enddo
+!enddo
+!ldf(2011-01-10):
  do j = jts,jte
  do k = kts+1,kte
  do i = its,ite
-    t2_p(i,k,j)    = fzm(k)*t_p(i,k,j) + fzp(k)*t_p(i,k-1,j)
-    pres2_p(i,k,j) = fzm(k)*pres_p(i,k,j) + fzp(k)*pres_p(i,k-1,j)
+    tem1 = 1./(zgrid(k+1,i)-zgrid(k-1,i))
+    fzm_p(i,k,j) = (zgrid(k,i)-zgrid(k-1,i)) * tem1
+    fzp_p(i,k,j) = (zgrid(k+1,i)-zgrid(k,i)) * tem1
+    t2_p(i,k,j)    = fzm_p(i,k,j)*t_p(i,k,j) + fzp_p(i,k,j)*t_p(i,k-1,j)
+    pres2_p(i,k,j) = fzm_p(i,k,j)*pres_p(i,k,j) + fzp_p(i,k,j)*pres_p(i,k-1,j)
  enddo
  enddo
  enddo
@@ -311,17 +321,6 @@
  enddo
  enddo
 
-!calculation of the hydrostatic pressure at w points:
- do j = jts,jte
- do i = its,ite
-    pres2_hyd_p(i,1,j) = psfc_p(i,j)
-    do k = kts+1,kte+1
-       pres2_hyd_p(i,k,j) = pres2_hyd_p(i,k-1,j) &amp;
-                          - rho_p(i,k-1,j)*(1+qv_p(i,k-1,j))*g*dz_p(i,k-1,j)
-    enddo
- enddo
- enddo
-
 !formats: 
  201 format(3i8,10(1x,e15.8))
  202 format(2i6,10(1x,e15.8))
@@ -480,6 +479,13 @@
  real(kind=RKIND),dimension(:,:),pointer:: rho_zz,theta_m,pressure_p
  real(kind=RKIND),dimension(:,:),pointer:: rt_diabatic_tend
 
+!ldf(2011-11-12): surface pressure.
+ real(kind=RKIND):: rho1,rho2,tem1,tem2
+ real(kind=RKIND),dimension(:),pointer:: rdzw
+ real(kind=RKIND),dimension(:),pointer:: sfc_pressure
+ real(kind=RKIND),dimension(:,:),pointer:: zgrid
+!ldf end.
+
 !---------------------------------------------------------------------------------------------
 
  write(0,*)
@@ -487,6 +493,7 @@
 
 !initialization:
  zz         =&gt; mesh % zz % array
+ zgrid      =&gt; mesh % zgrid % array
  exner      =&gt; diag % exner % array
  exner_b    =&gt; diag % exner_base % array
  pressure_b =&gt; diag % pressure_base % array
@@ -499,6 +506,11 @@
 
  rt_diabatic_tend =&gt; tend % rt_diabatic_tend % array
 
+!ldf (2011-11-12): update surface pressure.
+ rdzw =&gt; mesh % rdzw % array
+ sfc_pressure =&gt; diag % surface_pressure % array
+!ldf end.
+
 !variables common to all cloud microphysics schemes:
 
  do j = jts, jte
@@ -529,7 +541,30 @@
  enddo
  enddo
  enddo

+
+!updates the surface pressure.
+!do j = jts,jte
+!do i = its,ite
+!   sfc_pressure(i) = 0.5*g*(zgrid(2,i)-zgrid(1,i)) &amp;
+!                   * (1.25 * rho_zz(1,i) * zz(1,i) * (1. + qv_p(i,1,j))  &amp;
+!                   -  0.25 * rho_zz(2,i) * zz(2,i) * (1. + qv_p(i,2,j)))
+!   sfc_pressure(i) = sfc_pressure(i) + pressure_p(1,i) + pressure_b(1,i)
+!enddo
+!enddo
+!ldf (2012-01-09):
+ do j = jts,jte
+ do i = its,ite
+    tem1 = zgrid(2,i)-zgrid(1,i)
+    tem2 = zgrid(3,i)-zgrid(2,i)
+    rho1 = rho_zz(1,i) * zz(1,i) * (1. + qv_p(i,1,j))
+    rho2 = rho_zz(2,i) * zz(2,i) * (1. + qv_P(i,2,j))
+    sfc_pressure(i) = 0.5*g*(zgrid(2,i)-zgrid(1,i)) &amp;
+                    * (rho1 + 0.5*(rho2-rho1)*tem1/(tem1+tem2))
+    sfc_pressure(i) = sfc_pressure(i) + pressure_p(1,i) + pressure_b(1,i)
+ enddo
+ enddo
+!ldf end.
+
 !variables specific to different cloud microphysics schemes:
 
  microp_select_init: select case(microp_scheme)

Modified: branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_lsm_noahinit.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_lsm_noahinit.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_lsm_noahinit.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -9,7 +9,7 @@
  use mpas_configure, only: restart =&gt; config_do_restart,  &amp;
                            mminlu  =&gt; input_landuse_data, &amp;
                            mminsl  =&gt; input_soil_data   , &amp;
-                           input_sfc_albedo =&gt; config_sfc_albedo
+                           input_sfc_albedo =&gt; config_sfc_snowalbedo
  use mpas_dmpar
  use mpas_grid_types
 

Modified: branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_todynamics.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_todynamics.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_todynamics.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -3,7 +3,7 @@
  use mpas_configure
  use mpas_grid_types
 
- use mpas_atmphys_constants, only: R_d,R_v
+ use mpas_atmphys_constants, only: R_d,R_v,degrad
 
  implicit none
  private
@@ -34,9 +34,9 @@
 
 !local variables:
 !----------------
- integer:: i,k,nCells,nCellsSolve,nEdges,nEdgesSolve,nVertLevels
+ integer:: i,iCell,k,n,nCells,nCellsSolve,nEdges,nEdgesSolve,nVertLevels
 
- real(kind=RKIND),dimension(:,:),pointer:: theta,theta_m,qv
+ real(kind=RKIND),dimension(:,:),pointer:: theta_m,qv
  real(kind=RKIND),dimension(:,:),pointer:: rthblten,rqvblten,rqcblten, &amp;
                                            rqiblten,rublten,rvblten
  real(kind=RKIND),dimension(:,:),pointer:: rthcuten,rqvcuten,rqccuten, &amp;
@@ -49,6 +49,10 @@
  real(kind=RKIND):: tem
  real(kind=RKIND),dimension(:,:),allocatable:: rublten_Edge
 
+!ldf (2011-12-16):
+ real(kind=RKIND),dimension(:,:),allocatable:: theta,tend_th
+!ldf end.
+
 !=============================================================================================
 !write(0,*)
 !write(0,*) '--- enter subroutine physics_add_tend:'
@@ -59,7 +63,7 @@
  nEdgesSolve = mesh % nEdgesSolve
  nVertLevels = mesh % nVertLevels
 
- theta   =&gt; diag % theta % array
+!theta   =&gt; diag % theta % array
  theta_m =&gt; state % theta_m % array
  qv      =&gt; state % scalars % array(state%index_qv,:,:)
 
@@ -84,6 +88,11 @@
  tend_theta   =&gt; tend % theta_m % array
  tend_scalars =&gt; tend % scalars % array
 
+!initialize the tendency for the potential temperature and all scalars due to PBL, convection,
+!and longwave and shortwave radiation:
+ allocate(theta(nVertLevels,nCellsSolve)  )
+ allocate(tend_th(nVertLevels,nCellsSolve))
+ tend_th      = 0.
  tend_scalars = 0.
 
 !add coupled tendencies due to PBL processes:
@@ -96,24 +105,25 @@
        tend_u(k,i)=tend_u(k,i)+rublten_Edge(k,i)*mass_edge(k,i)
     enddo
     enddo
-
     deallocate(rublten_Edge)
 
     do i = 1, nCellsSolve
     do k = 1, nVertLevels
-       tend_theta(k,i)=tend_theta(k,i)+rthblten(k,i)*mass(k,i)
+       tend_th(k,i) = tend_th(k,i)+rthblten(k,i)*mass(k,i)
        tend_scalars(tend%index_qv,k,i)=tend_scalars(tend%index_qv,k,i)+rqvblten(k,i)*mass(k,i)
        tend_scalars(tend%index_qc,k,i)=tend_scalars(tend%index_qc,k,i)+rqcblten(k,i)*mass(k,i)
        tend_scalars(tend%index_qi,k,i)=tend_scalars(tend%index_qi,k,i)+rqiblten(k,i)*mass(k,i)
     enddo
     enddo
  endif
+ write(0,*) 'max rthblten = ',maxval(rthblten(:,1:nCellsSolve))
+ write(0,*) 'min rthblten = ',minval(rthblten(:,1:nCellsSolve))
 
 !add coupled tendencies due to convection:
  if(config_conv_deep_scheme .ne. 'off') then
     do i = 1, nCellsSolve
     do k = 1, nVertLevels
-       tend_theta(k,i)=tend_theta(k,i)+rthcuten(k,i)*mass(k,i)
+       tend_th(k,i)=tend_th(k,i)+rthcuten(k,i)*mass(k,i)
        tend_scalars(tend%index_qv,k,i)=tend_scalars(tend%index_qv,k,i)+rqvcuten(k,i)*mass(k,i)
        tend_scalars(tend%index_qc,k,i)=tend_scalars(tend%index_qc,k,i)+rqccuten(k,i)*mass(k,i)
        tend_scalars(tend%index_qr,k,i)=tend_scalars(tend%index_qr,k,i)+rqrcuten(k,i)*mass(k,i)
@@ -122,64 +132,55 @@
     enddo
     enddo
  endif
+ write(0,*) 'max rthcuten = ',maxval(rthcuten(:,1:nCellsSolve))
+ write(0,*) 'min rthcuten = ',minval(rthcuten(:,1:nCellsSolve))
 
 !add coupled tendencies due to longwave radiation:
  if(config_radt_lw_scheme .ne. 'off') then
     do i = 1, nCellsSolve
     do k = 1, nVertLevels
-       tend_theta(k,i)=tend_theta(k,i)+rthratenlw(k,i)*mass(k,i)
+       tend_th(k,i)=tend_th(k,i)+rthratenlw(k,i)*mass(k,i)
     enddo
     enddo
  endif
+ write(0,*) 'max rthratenlw = ',maxval(rthratenlw(:,1:nCellsSolve))
+ write(0,*) 'min rthratenlw = ',minval(rthratenlw(:,1:nCellsSolve))
     
 !add coupled tendencies due to shortwave radiation:
  if(config_radt_sw_scheme .ne. 'off') then
     do i = 1, nCellsSolve
     do k = 1, nVertLevels
-       tend_theta(k,i)=tend_theta(k,i)+rthratensw(k,i)*mass(k,i)
+       tend_th(k,i)=tend_th(k,i)+rthratensw(k,i)*mass(k,i)
     enddo
     enddo
  endif
+ write(0,*) 'max rthratensw = ',maxval(rthratensw(:,1:nCellsSolve))
+ write(0,*) 'min rthratensw = ',minval(rthratensw(:,1:nCellsSolve))
 
 !if non-hydrostatic core, convert the tendency for the potential temperature to a
 !tendency for the modified potential temperature:
 #ifdef non_hydrostatic_core
     do i = 1, nCellsSolve
     do k = 1, nVertLevels
-
        theta(k,i) = theta_m(k,i) / (1. + R_v/R_d * qv(k,i))
-       tend_theta(k,i) = (1. + R_v/R_d * qv(k,i)) * tend_theta(k,i) &amp;
+       tend_th(k,i) = (1. + R_v/R_d * qv(k,i)) * tend_th(k,i) &amp;
                        + R_v/R_d * theta(k,i) * tend_scalars(tend%index_qv,k,i)
+       tend_theta(k,i) = tend_theta(k,i) + tend_th(k,i)
     enddo
-    enddo 
+    enddo
+#elif hydrostatic_core
+    do i = 1, nCellsSolve
+    do k = 1, nVertLevels
+       tend_theta(k,i) = tend_theta(k,i) + tend_th(k,i)
+    enddo
+    enddo
 #endif
+ deallocate(theta)
+ deallocate(tend_th)
 
-!write(0,*) 'max PBL tendencies:'
-!write(0,*) 'max rthblten=',maxval(rthblten(:,:))
-!write(0,*) 'max rqvblten=',maxval(rqvblten(:,:))
-!write(0,*) 'max rqcblten=',maxval(rqcblten(:,:))
-!write(0,*) 'max rqiblten=',maxval(rqiblten(:,:))
-!write(0,*) 'max rublten =',maxval(rublten(:,:))
-!write(0,*) 'max rvblten =',maxval(rvblten(:,:))
-!write(0,*)
-!write(0,*) 'max CU tendencies:'
-!write(0,*) 'max rthcuten=',maxval(rthcuten(:,:))
-!write(0,*) 'max rqvcuten=',maxval(rqvcuten(:,:))
-!write(0,*) 'max rqccuten=',maxval(rqccuten(:,:))
-!write(0,*) 'max rqrcuten=',maxval(rqrcuten(:,:))
-!write(0,*) 'max rqicuten=',maxval(rqicuten(:,:))
-!write(0,*) 'max rqscuten=',maxval(rqscuten(:,:))
-!write(0,*)
-!write(0,*) 'max tend_scalars:'
-!write(0,*) 'max tend qv=',maxval(tend_scalars(tend%index_qv,:,:))
-!write(0,*) 'max tend qc=',maxval(tend_scalars(tend%index_qc,:,:))
-!write(0,*) 'max tend qr=',maxval(tend_scalars(tend%index_qr,:,:))
-!write(0,*) 'max tend qi=',maxval(tend_scalars(tend%index_qi,:,:))
-!write(0,*) 'max tend qs=',maxval(tend_scalars(tend%index_qs,:,:))
-!write(0,*)
-
 !formats:
- 201 format(2i6,8(1x,e15.8))
+ 201 format(2i6,10(1x,e15.8))
+ 202 format(3i6,10(1x,e15.8))
 
  end subroutine  physics_addtend
 

Modified: branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_vars.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_vars.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_atmos_physics/mpas_atmphys_vars.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -32,7 +32,7 @@
  integer,public:: its,ite,jts,jte,kts,kte
  integer,public:: n_microp
 
- integer,public:: num_months          !number of months                                     [-]
+ integer,public:: num_months          !number of months                                    [-]
 
  real(kind=RKIND),public:: dt_dyn     !time-step for dynamics
  real(kind=RKIND),public:: dt_microp  !time-step for cloud microphysics parameterization.
@@ -42,15 +42,19 @@
  real(kind=RKIND),public:: xice_threshold
 
  real(kind=RKIND),dimension(:,:),allocatable:: &amp;
-    area_p             !grid cell area                                                     [m2]
+    area_p             !grid cell area                                                    [m2]
 
 !... arrays related to surface:
  real(kind=RKIND),dimension(:,:),allocatable:: &amp;
     psfc_p,           &amp;!surface pressure                                                  [Pa]
     ptop_p             !model-top pressure                                                [Pa]
 
+ real(kind=RKIND),dimension(:,:,:),allocatable:: &amp;
+    fzm_p,            &amp;!weight for interpolation to w points                               [-]
+    fzp_p              !weight for interpolation to w points                               [-]
+
+ real(kind=RKIND),dimension(:,:,:),allocatable:: &amp;
 !... arrays related to u- and v-velocities interpolated to theta points:
- real(kind=RKIND),dimension(:,:,:),allocatable:: &amp;
     u_p,              &amp;!u-velocity interpolated to theta points                          [m/s]
     v_p                !v-velocity interpolated to theta points                          [m/s]
     
@@ -60,6 +64,7 @@
     pres_p,           &amp;!pressure                                                          [Pa]
     pi_p,             &amp;!(p_phy/p0)**(r_d/cp)                                               [-]
     z_p,              &amp;!height of layer                                                    [m]
+    zmid_p,           &amp;!height of middle of layer                                          [m]
     dz_p,             &amp;!layer thickness                                                    [m]
     t_p,              &amp;!temperature                                                        [K]
     th_p,             &amp;!potential temperature                                              [K]
@@ -68,10 +73,6 @@
     rh_p               !relative humidity                                                  [-]
 
  real(kind=RKIND),dimension(:,:,:),allocatable:: &amp;
-    pres_hyd_p,       &amp;!hydrostatic pressure at theta points                              [Pa]
-    pres2_hyd_p        !hydrostatic pressure at w points                                  [Pa]
-
- real(kind=RKIND),dimension(:,:,:),allocatable:: &amp;
     qv_p,             &amp;!water vapor mixing ratio                                       [kg/kg]
     qc_p,             &amp;!cloud water mixing ratio                                       [kg/kg]
     qr_p,             &amp;!rain mixing ratio                                              [kg/kg]
@@ -323,6 +324,7 @@
 
  integer,public:: &amp;
     num_soils          !number of soil layers                                              [-]
+    
  integer,dimension(:,:),allocatable:: &amp;
     isltyp_p,         &amp;!dominant soil type category                                        [-]
     ivgtyp_p           !dominant vegetation category                                       [-]

Modified: branches/ocean_projects/ale_vert_coord/src/core_init_nhyd_atmos/Registry
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_init_nhyd_atmos/Registry        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_init_nhyd_atmos/Registry        2012-01-17 22:23:38 UTC (rev 1385)
@@ -11,6 +11,7 @@
 namelist integer   dimensions config_nsoillevels          4
 namelist integer   dimensions config_nfglevels            27
 namelist integer   dimensions config_nfgsoillevels        4
+namelist integer   dimensions config_months               12
 namelist character data_sources config_geog_data_path     /data3/mp/wrfhelp/WPS_GEOG/
 namelist character data_sources config_met_prefix         FILE
 namelist character data_sources config_sst_prefix         FILE
@@ -54,7 +55,7 @@
 dim nFGLevels namelist:config_nfglevels
 dim nFGSoilLevels namelist:config_nfgsoillevels
 dim nVertLevelsP1 nVertLevels+1
-dim nMonths 12
+dim nMonths namelist:config_months
 
 #
 # var  type  name_in_file  ( dims )  iro-  name_in_code super-array array_class
@@ -178,6 +179,7 @@
 var persistent real    sh2o ( nSoilLevels nCells Time ) 1 io sh2o fg - -
 var persistent real    smois ( nSoilLevels nCells Time ) 1 io smois fg - -
 var persistent real    tslb ( nSoilLevels nCells Time ) 1 io tslb fg - -
+var persistent real    smcrel ( nSoilLevels nCells Time ) 1 io smcrel fg - -
 var persistent real    tmn ( nCells Time ) 1 io tmn fg - -
 var persistent real    skintemp ( nCells Time ) 1 io skintemp fg - -
 var persistent real    sst ( nCells Time ) 1 iso sst fg - -
@@ -227,7 +229,7 @@
 
 var persistent real    cqw  ( nVertLevels nCells Time ) 1 - cqw diag - -
 
-var persistent real    surface_pressure    ( nCells Time    ) 1  o surface_pressure    diag - -
+var persistent real    surface_pressure    ( nCells Time    ) 1  io surface_pressure    diag - -
 
 #  coupled variables needed by the solver, but not output...
 var persistent real    ru ( nVertLevels nEdges Time ) 1 - ru diag - -

Modified: branches/ocean_projects/ale_vert_coord/src/core_init_nhyd_atmos/mpas_init_atm_test_cases.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_init_nhyd_atmos/mpas_init_atm_test_cases.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_init_nhyd_atmos/mpas_init_atm_test_cases.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -5,12 +5,16 @@
    use mpas_constants
    use mpas_dmpar
    use atm_advection
-   use mpas_sort
-   use mpas_timekeeping
-
    use mpas_atmphys_initialize_real
+   use mpas_RBF_interpolation
+   use mpas_vector_reconstruction
 
+   ! Added only clause to keep xlf90 from getting confused from the overloaded abs intrinsic in mpas_timekeeping
+   use mpas_timekeeping !, only: MPAS_Time_type, MPAS_TimeInterval_type, MPAS_Clock_type, &amp;
+                        !        mpas_set_time, mpas_set_timeInterval, mpas_get_time, operator(+), add_t_ti
+   
 
+
    contains
 
 
@@ -30,17 +34,9 @@
       integer :: i
       type (block_type), pointer :: block_ptr
 
-      if (config_test_case == 0) then
-         write(0,*) ' Using initial conditions from input file'
-         block_ptr =&gt; domain % blocklist
-         do while (associated(block_ptr))
-            do i=2,nTimeLevs
-               call mpas_copy_state(block_ptr % state % time_levs(i) % state, block_ptr % state % time_levs(1) % state)
-            end do
-            block_ptr =&gt; block_ptr % next
-         end do
 
-      else if ((config_test_case == 1) .or. (config_test_case == 2) .or. (config_test_case == 3)) then
+      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 '
          if (config_test_case == 1) write(0,*) ' no initial perturbation '
          if (config_test_case == 2) write(0,*) ' initial perturbation included '
@@ -49,11 +45,8 @@
          do while (associated(block_ptr))
             write(0,*) ' calling test case setup '
             call init_atm_test_case_jw(block_ptr % mesh, block_ptr % state % time_levs(1) % state, block_ptr % diag, config_test_case)
+            call decouple_variables(block_ptr % mesh, block_ptr % state % time_levs(1) % state, block_ptr % diag)
             write(0,*) ' returned from test case setup '
-            do i=2,nTimeLevs
-               call mpas_copy_state(block_ptr % state % time_levs(i) % state, block_ptr % state % time_levs(1) % state)
-            end do
-
             block_ptr =&gt; block_ptr % next
          end do
 
@@ -66,11 +59,8 @@
          do while (associated(block_ptr))
             write(0,*) ' calling test case setup '
             call init_atm_test_case_squall_line(domain % dminfo, block_ptr % mesh, block_ptr % state % time_levs(1) % state, block_ptr % diag, config_test_case)
+            call decouple_variables(block_ptr % mesh, block_ptr % state % time_levs(1) % state, block_ptr % diag)
             write(0,*) ' returned from test case setup '
-            do i=2,nTimeLevs
-               call mpas_copy_state(block_ptr % state % time_levs(i) % state, block_ptr % state % time_levs(1) % state)
-            end do
-
             block_ptr =&gt; block_ptr % next
          end do
 
@@ -81,11 +71,8 @@
          do while (associated(block_ptr))
             write(0,*) ' calling test case setup '
             call init_atm_test_case_mtn_wave(block_ptr % mesh, block_ptr % state % time_levs(1) % state, block_ptr % diag, config_test_case)
+            call decouple_variables(block_ptr % mesh, block_ptr % state % time_levs(1) % state, block_ptr % diag)
             write(0,*) ' returned from test case setup '
-            do i=2,nTimeLevs
-               call mpas_copy_state(block_ptr % state % time_levs(i) % state, block_ptr % state % time_levs(1) % state)
-            end do
-
             block_ptr =&gt; block_ptr % next
          end do
 
@@ -98,11 +85,6 @@
                                     block_ptr % diag, config_test_case, block_ptr % parinfo)
             if(config_physics_init) &amp;
                call physics_initialize_real(block_ptr % mesh, block_ptr % fg)
-
-            do i=2,nTimeLevs
-               call mpas_copy_state(block_ptr % state % time_levs(i) % state, block_ptr % state % time_levs(1) % state)
-            end do
-
             block_ptr =&gt; block_ptr % next
          end do
 
@@ -113,20 +95,35 @@
          do while (associated(block_ptr))
             call init_atm_test_case_sst(domain, domain % dminfo, block_ptr % mesh, block_ptr % fg, block_ptr % state % time_levs(1) % state, &amp;
                                     block_ptr % diag, config_test_case, block_ptr % parinfo)
-            do i=2,nTimeLevs
-               call mpas_copy_state(block_ptr % state % time_levs(i) % state, block_ptr % state % time_levs(1) % state)
-            end do
-
             block_ptr =&gt; block_ptr % next
          end do
 
       else
 
+         write(0,*) ' Only test cases 1, 2, 3, 4, 5, 6, 7, and 8 are currently supported for nonhydrostatic core '
+         stop
 
-         write(0,*) ' Only test case 1, 2, 3, 4, 5, 6, and 7 are currently supported for nonhydrostatic core '
-         stop
       end if
 
+      block_ptr =&gt; domain % blocklist
+      do while (associated(block_ptr))
+         do i=2,nTimeLevs
+            call mpas_copy_state(block_ptr % state % time_levs(i) % state, block_ptr % state % time_levs(1) % state)
+         end do
+         block_ptr =&gt; block_ptr % next
+      end do
+
+      !initialization of surface input variables technically not needed to run our current set of
+      !idealized test cases:
+      if (config_test_case &lt; 7)  then
+         block_ptr =&gt; domain % blocklist
+         do while (associated(block_ptr))
+            call physics_idealized_init(block_ptr % mesh, block_ptr % fg)
+            block_ptr =&gt; block_ptr % next
+         end do
+      endif
+
+
    end subroutine init_atm_setup_test_case
 
 !----------------------------------------------------------------------------------------------------------
@@ -141,6 +138,7 @@
       type (mesh_type), intent(inout) :: grid
       type (state_type), intent(inout) :: state
       type (diag_type), intent(inout) :: diag
+      !type (diag_physics_type), intent(inout) :: diag_physics
       integer, intent(in) :: test_case
 
       real (kind=RKIND), parameter :: u0 = 35.0
@@ -150,15 +148,22 @@
       real (kind=RKIND), parameter :: u_perturbation = 1., pert_radius = 0.1, latitude_pert = 40., longitude_pert = 20.
       real (kind=RKIND), parameter :: theta_c = pii/4.0
       real (kind=RKIND), parameter :: lambda_c = 3.0*pii/2.0
-      real (kind=RKIND), parameter :: rh_max = 0.4       ! Maximum relative humidity
       real (kind=RKIND), parameter :: k_x = 9.           ! Normal mode wave number
 
       real (kind=RKIND), dimension(:), pointer :: rdzw, dzu, rdzu, fzm, fzp
+      real (kind=RKIND), dimension(:), pointer :: surface_pressure
       real (kind=RKIND), dimension(:,:), pointer :: zgrid, zx, zz, hx
       real (kind=RKIND), dimension(:,:), pointer :: pressure, ppb, pb, rho_zz, rb, rr, tb, rtb, p, pp, dss, t, rt
       real (kind=RKIND), dimension(:,:,:), pointer :: zf, zf3, zb, zb3
-      real (kind=RKIND), dimension(:,:,:), pointer :: scalars
       real (kind=RKIND), dimension(:,:,:), pointer :: deriv_two
+      
+!.. initialization of moisture:
+      integer:: index_qv
+      real (kind=RKIND),parameter :: rh_max = 0.40 ! Maximum relative humidity
+!      real (kind=RKIND),parameter :: rh_max = 0.70 ! Maximum relative humidity
+      real (kind=RKIND),dimension(grid % nVertLevels, grid % nCells) :: qsat, relhum
+      real (kind=RKIND),dimension(:,:,:),pointer:: scalars
+!.. end initialization of moisture.
 
       integer :: iCell, iCell1, iCell2 , iEdge, vtx1, vtx2, ivtx, i, k, nz, nz1, itr, itrp, cell1, cell2, nCellsSolve
 
@@ -166,8 +171,8 @@
       integer :: eoe, j
       integer, dimension(:), pointer :: nEdgesOnEdge 
       integer, dimension(:,:), pointer :: edgesOnEdge, CellsOnEdge
-      real, dimension(:), pointer :: dvEdge, AreaCell 
-      real, dimension(:,:), pointer :: weightsOnEdge
+      real (kind=RKIND), dimension(:), pointer :: dvEdge, AreaCell 
+      real (kind=RKIND), dimension(:,:), pointer :: weightsOnEdge
 
       real (kind=RKIND) :: u, v, flux, fluxk, lat1, lat2, eta_v, r_pert, u_pert, lat_pert, lon_pert, r
 
@@ -176,8 +181,7 @@
 
       real (kind=RKIND) :: r_earth, etavs, ztemp, zd, zt, dz, gam, delt, str
 
-      real (kind=RKIND), dimension(grid % nVertLevels, grid % nCells) :: rel_hum, temperature, qv
-      real (kind=RKIND) :: ptmp, es, qvs, xnutr, znut, ptemp
+      real (kind=RKIND) :: es, qvs, xnutr, znut, ptemp 
       integer :: iter
 
       real (kind=RKIND), dimension(grid % nVertLevels + 1 ) :: hyai, hybi, znu, znw, znwc, znwv, hyam, hybm
@@ -185,20 +189,23 @@
 
       real (kind=RKIND), dimension(grid % nVertLevels + 1 ) :: sh, zw, ah
       real (kind=RKIND), dimension(grid % nVertLevels ) :: zu, dzw, rdzwp, rdzwm
-      real (kind=RKIND), dimension(grid % nVertLevels ) :: eta, etav, teta, ppi, tt
+      real (kind=RKIND), dimension(grid % nVertLevels ) :: eta, etav, teta, ppi, tt, temperature_1d
 
       real (kind=RKIND) :: d1, d2, d3, cf1, cf2, cf3, cof1, cof2, psurf
 
       !  storage for (lat,z) arrays for zonal velocity calculation
 
-      integer, parameter :: nlat=361
-      real (kind=RKIND), dimension(grid % nVertLevels + 1) :: zz_1d, zgrid_1d, hx_1d
+      logical, parameter :: rebalance = .true.
+      integer, parameter :: nlat=721
       real (kind=RKIND), dimension(grid % nVertLevels) :: flux_zonal
-      real (kind=RKIND), dimension(nlat, grid % nVertLevels) :: u_2d, etavs_2d
+      real (kind=RKIND), dimension(grid % nVertLevels + 1, nlat) :: zgrid_2d
+      real (kind=RKIND), dimension(grid % nVertLevels, nlat) :: u_2d, pp_2d, rho_2d, qv_2d, etavs_2d, zz_2d
+      real (kind=RKIND), dimension(grid % nVertLevels, nlat-1) :: zx_2d 
       real (kind=RKIND), dimension(nlat) :: lat_2d
-      real (kind=RKIND) :: dlat
+      real (kind=RKIND) :: dlat, hx_1d
       real (kind=RKIND) :: z_edge, z_edge3, d2fdx2_cell1, d2fdx2_cell2
 
+      logical, parameter :: moisture = .true.
       !
       ! Scale all distances and areas from a unit sphere to one with radius a
       !
@@ -259,13 +266,25 @@
       t =&gt; state % theta_m % array      
       rt =&gt; diag % rtheta_p % array
 
+      surface_pressure =&gt; diag % surface_pressure % array
+
+!.. initialization of moisture:
       scalars =&gt; state % scalars % array
+      !qsat    =&gt; diag_physics % qsat % array
+      !relhum  =&gt; diag_physics % relhum % array
+      scalars(:,:,:) = 0.0
+      qsat(:,:)      = 0.0
+      relhum(:,:)    = 0.0
+      qv_2d(:,:)     = 0.0
+!.. end initialization of moisture.
 
-      scalars(:,:,:) = 0.
+      surface_pressure(:) = 0.0
 
       call atm_initialize_advection_rk(grid) 
       call atm_initialize_deformation_weights(grid) 
 
+      index_qv = state % index_qv
+
       xnutr = 0.
       zd = 12000.
       znut = eta_t
@@ -293,7 +312,7 @@
 
       !     Metrics for hybrid coordinate and vertical stretching
 
-      str = 1.5
+      str = 1.8
       zt = 45000.
       dz = zt/float(nz1)
 
@@ -323,7 +342,7 @@
  
             ah(k) = 1.-cos(.5*pii*(k-1)*dz/zt)**6
 !            ah(k) = 0.
-            write(0,*) ' k, sh, zw, ah ',k,sh(k),zw(k),ah(k)                        
+            write(0,*) ' k, sh, zw, ah ',k,sh(k),zw(k),ah(k)
       end do
       do k=1,nz1
          dzw (k) = zw(k+1)-zw(k)
@@ -388,47 +407,43 @@
         end do
       enddo
 
-      do k=1,nz1
-        write(0,*) ' k, zgrid(k,1),hx(k,1) ',k,zgrid(k,1),hx(k,1)
-      enddo
+      !do k=1,nz1
+      !  write(0,*) ' k, zgrid(k,1),hx(k,1) ',k,zgrid(k,1),hx(k,1)
+      !enddo
 
-      do k=1,nz1
-        write(0,*) ' k, zx(k,1) ',k,zx(k,1)
-      enddo
+      !do k=1,nz1
+      !  write(0,*) ' k, zx(k,1) ',k,zx(k,1)
+      !enddo
 
       write(0,*) ' grid metrics setup complete '
 
-!**************  section for 2d (lat,z) calc for zonal velocity
+!**************  section for 2d (z,lat) calc for zonal velocity
 
       dlat = 0.5*pii/float(nlat-1)
       do i = 1,nlat
 
         lat_2d(i) = float(i-1)*dlat
-!        write(0,*) ' zonal setup, latitude = ',lat_2d(i)*180./pii
+        phi = lat_2d(i)
+        hx_1d    = u0/gravity*cos(etavs)**1.5                           &amp;
+                   *((-2.*sin(phi)**6                                   &amp;
+                         *(cos(phi)**2+1./3.)+10./63.)                  &amp;
+                         *(u0)*cos(etavs)**1.5                          &amp;
+                    +(1.6*cos(phi)**3                                   &amp;
+                         *(sin(phi)**2+2./3.)-pii/4.)*r_earth*omega_e)
 
-        do k=1,nz
-          phi = lat_2d(i)
-          hx_1d(k) = u0/gravity*cos(etavs)**1.5                            &amp;
-                      *((-2.*sin(phi)**6                                   &amp;
-                            *(cos(phi)**2+1./3.)+10./63.)                  &amp;
-                            *(u0)*cos(etavs)**1.5                          &amp;
-                       +(1.6*cos(phi)**3                                   &amp;
-                            *(sin(phi)**2+2./3.)-pii/4.)*r_earth*omega_e)
-        enddo
-
         do k=1,nz        
-          zgrid_1d(k) = (1.-ah(k))*(sh(k)*(zt-hx_1d(k))+hx_1d(k))  &amp;
+          zgrid_2d(k,i) = (1.-ah(k))*(sh(k)*(zt-hx_1d)+hx_1d)  &amp;
                          + ah(k) * sh(k)* zt        
         end do
         do k=1,nz1
-          zz_1d (k) = (zw(k+1)-zw(k))/(zgrid_1d(k+1)-zgrid_1d(k))
+          zz_2d (k,i) = (zw(k+1)-zw(k))/(zgrid_2d(k+1,i)-zgrid_2d(k,i))
         end do
 
         do k=1,nz1
-          ztemp    = .5*(zgrid_1d(k+1)+zgrid_1d(k))
+          ztemp    = .5*(zgrid_2d(k+1,i)+zgrid_2d(k,i))
           ppb(k,i) = p0*exp(-gravity*ztemp/(rgas*t0b)) 
           pb (k,i) = (ppb(k,i)/p0)**(rgas/cp)
-          rb (k,i) = ppb(k,i)/(rgas*t0b*zz_1d(k))
+          rb (k,i) = ppb(k,i)/(rgas*t0b*zz_2d(k,i))
           tb (k,i) = t0b/pb(k,i)
           rtb(k,i) = rb(k,i)*tb(k,i)
           p  (k,i) = pb(k,i)
@@ -448,40 +463,43 @@
               teta(k) = t0*eta(k)**(rgas*dtdz/gravity) + delta_t*(znut-eta(k))**5
             end if
           end do
-          ! phi = grid % latCell % array (i)
+
           phi = lat_2d (i)
           do k=1,nz1
-            tt(k) = 0.
-            tt(k) = teta(k)+.75*eta(k)*pii*u0/rgas*sin(etav(k))      &amp;
+            temperature_1d(k) = teta(k)+.75*eta(k)*pii*u0/rgas*sin(etav(k))      &amp;
                             *sqrt(cos(etav(k)))*                   &amp;
                               ((-2.*sin(phi)**6                    &amp;
                                    *(cos(phi)**2+1./3.)+10./63.)   &amp;
                                    *2.*u0*cos(etav(k))**1.5        &amp;
                               +(1.6*cos(phi)**3                    &amp;
-                                *(sin(phi)**2+2./3.)-pii/4.)*r_earth*omega_e)
+                                *(sin(phi)**2+2./3.)-pii/4.)*r_earth*omega_e)/(1.+0.61*qv_2d(k,i))
 
 
-            ztemp   = .5*(zgrid_1d(k)+zgrid_1d(k+1))
+            ztemp   = .5*(zgrid_2d(k,i)+zgrid_2d(k+1,i))
             ptemp   = ppb(k,i) + pp(k,i)
-            qv(k,i) = 0.
 
+            !get moisture 
+            if (moisture) then
+              qv_2d(k,i) = env_qv( ztemp, temperature_1d(k), ptemp, rh_max )
+            end if
+
+            tt(k) = temperature_1d(k)*(1.+1.61*qv_2d(k,i))
           end do
-                
+
           do itrp = 1,25
             do k=1,nz1                                
-              rr(k,i)  = (pp(k,i)/(rgas*zz_1d(k))  &amp;
-                          -rb(k,i)*(tt(k)-t0b))/tt(k)
+              rr(k,i)  = (pp(k,i)/(rgas*zz_2d(k,i)) - rb(k,i)*(tt(k)-t0b))/tt(k)
             end do
 
-            ppi(1) = p0-.5*dzw(1)*gravity                         &amp;
-                          *(1.25*(rr(1,i)+rb(1,i))*(1.+qv(1,i))   &amp;
-                            -.25*(rr(2,i)+rb(2,i))*(1.+qv(2,i)))
+            ppi(1) = p0-.5*dzw(1)*gravity                            &amp;
+                          *(1.25*(rr(1,i)+rb(1,i))*(1.+qv_2d(1,i))   &amp;
+                            -.25*(rr(2,i)+rb(2,i))*(1.+qv_2d(2,i)))
 
             ppi(1) = ppi(1)-ppb(1,i)
             do k=1,nz1-1
-              ppi(k+1) = ppi(k)-.5*dzu(k+1)*gravity*                     &amp;
-                            (rr(k  ,i)+(rr(k  ,i)+rb(k  ,i))*qv(k  ,i)   &amp;
-                            +rr(k+1,i)+(rr(k+1,i)+rb(k+1,i))*qv(k+1,i))
+              ppi(k+1) = ppi(k)-.5*dzu(k+1)*gravity*                        &amp;
+                            (rr(k  ,i)+(rr(k  ,i)+rb(k  ,i))*qv_2d(k  ,i)   &amp;
+                            +rr(k+1,i)+(rr(k+1,i)+rb(k+1,i))*qv_2d(k+1,i))
             end do
 
             do k=1,nz1
@@ -493,21 +511,28 @@
         end do  ! end outer iteration loop itr
 
         do k=1,nz1
-          etavs_2d(i,k) = (0.5*(ppb(k,i)+ppb(k,i)+pp(k,i)+pp(k,i))/p0 - 0.252)*pii/2.
-!          u_2d(i,k) = u0*(sin(2.*lat_2d(i))**2) *(cos(etavs_2d(i,k))**1.5)
-          u_2d(i,k) = u0*(sin(2.*lat_2d(i))**2) *(cos(etavs_2d(i,k))**1.5)*(rb(k,i)+rr(k,i))
+          rho_2d(k,i) = rr(k,i)+rb(k,i)
+          pp_2d(k,i) = pp(k,i)
+          etavs_2d(k,i) = ((ppb(k,i)+pp(k,i))/p0 - 0.252)*pii/2.
+          u_2d(k,i) = u0*(sin(2.*lat_2d(i))**2) *(cos(etavs_2d(k,i))**1.5)
         end do
 
       end do  ! end loop over latitudes for 2D zonal wind field calc
 
-!      do i=1,nlat
-!        do k=1,nz1
-!          u_2d(i,k) = u_2d(i,k) - u0*(sin(2.*lat_2d(i))**2) *(cos(etavs_2d(nlat/2,k))**1.5)
-!        end do
-!      end do
-!
-!      write(22,*) nz1,nlat,u_2d
+      !SHP-balance:: in case of rebalacing for geostrophic wind component
+      if (rebalance) then
 
+        do i=1,nlat-1
+          do k=1,nz1
+            zx_2d(k,i) = (zgrid_2d(k,i+1)-zgrid_2d(k,i))/(dlat*r_earth)
+          end do
+        end do
+
+        call init_atm_recompute_geostrophic_wind(u_2d,rho_2d,pp_2d,qv_2d,lat_2d,zz_2d,zx_2d,     &amp;
+                                        cf1,cf2,cf3,fzm,fzp,rdzw,nz1,nlat,dlat)
+
+      end if
+
 !******************************************************************      
 
 !
@@ -516,7 +541,6 @@
 !     reference sounding based on dry isothermal atmosphere
 !
       do i=1, grid % nCells
-        !write(0,*) ' thermodynamic setup, cell ',i
         do k=1,nz1
           ztemp    = .5*(zgrid(k+1,i)+zgrid(k,i))
           ppb(k,i) = p0*exp(-gravity*ztemp/(rgas*t0b)) 
@@ -529,12 +553,17 @@
           rr (k,i) = 0.
         end do
 
-        if(i == 1) then
-          do k=1,nz1
-            write(0,*) ' k, ppb, pb, rb, tb (k,1) ',k,ppb(k,1),pb(k,1),rb(k,1)*zz(k,1),tb(k,1)
-          enddo
-        end if
-!
+!       if(i == 1) then
+!         do k=1,nz1
+!           write(0,*) ' k, ppb, pb, rb, tb (k,1) ',k,ppb(k,1),pb(k,1),rb(k,1)*zz(k,1),tb(k,1)
+!         enddo
+!       end if
+
+      200 format(4i6,8(1x,e15.8))
+      201 format(3i6,8(1x,e15.8))
+      202 format(2i6,10(1x,e15.8))
+      203 format(i6,10(1x,e15.8))
+
 !     iterations to converge temperature as a function of pressure
 !
         do itr = 1,10
@@ -550,42 +579,60 @@
           end do
           phi = grid % latCell % array (i)
           do k=1,nz1
-            tt(k) = 0.
-            tt(k) = teta(k)+.75*eta(k)*pii*u0/rgas*sin(etav(k))      &amp;
+            temperature_1d(k) = teta(k)+.75*eta(k)*pii*u0/rgas*sin(etav(k))      &amp;
                             *sqrt(cos(etav(k)))*                   &amp;
                               ((-2.*sin(phi)**6                    &amp;
                                    *(cos(phi)**2+1./3.)+10./63.)   &amp;
                                    *2.*u0*cos(etav(k))**1.5        &amp;
                               +(1.6*cos(phi)**3                    &amp;
-                                *(sin(phi)**2+2./3.)-pii/4.)*r_earth*omega_e)
+                                *(sin(phi)**2+2./3.)-pii/4.)*r_earth*omega_e)/(1.+0.61*scalars(index_qv,k,i))
 
-
-            !write(0,*) ' k, tt(k) ',k,tt(k)
             ztemp   = .5*(zgrid(k,i)+zgrid(k+1,i))
             ptemp   = ppb(k,i) + pp(k,i)
-!            qv(k,i) = env_qv( ztemp, tt(k), ptemp, 0 )
-            qv(k,i) = 0.
 
+            !get moisture 
+            if (moisture) then

+                !scalars(index_qv,k,i) = env_qv( ztemp, temperature_1d(k), ptemp, rh_max )
+
+               if(ptemp &lt; 50000.) then
+                  relhum(k,i) = 0.0
+               elseif(ptemp &gt; p0) then
+                  relhum(k,i) = 1.0
+               else
+                  relhum(k,i) = (1.-((p0-ptemp)/50000.)**1.25)
+               endif
+               relhum(k,i) = min(rh_max,relhum(k,i))
+
+               !.. calculation of water vapor mixing ratio:
+               if (temperature_1d(k) &gt; 273.15) then
+                   es  = 1000.*0.6112*exp(17.67*(temperature_1d(k)-273.15)/(temperature_1d(k)-29.65))
+               else
+                   es  = 1000.*0.6112*exp(21.8745584*(temperature_1d(k)-273.15)/(temperature_1d(k)-7.66))
+               end if
+               qsat(k,i) = (287.04/461.6)*es/(ptemp-es)
+               if(relhum(k,i) .eq. 0.0) qsat(k,i) = 0.0
+               scalars(index_qv,k,i) = relhum(k,i)*qsat(k,i)
+            end if
+
+            tt(k) = temperature_1d(k)*(1.+1.61*scalars(index_qv,k,i))
+
           end do
-!          do k=2,nz1
-!            cqw(k,i) = 1./(1.+.5*(qv(k,i)+qv(k-1,i)))
-!          end do
                 
           do itrp = 1,25
             do k=1,nz1                                
-              rr(k,i)  = (pp(k,i)/(rgas*zz(k,i))  &amp;
-                          -rb(k,i)*(tt(k)-t0b))/tt(k)
+              rr(k,i)  = (pp(k,i)/(rgas*zz(k,i)) - rb(k,i)*(tt(k)-t0b))/tt(k)
             end do
 
             ppi(1) = p0-.5*dzw(1)*gravity                         &amp;
-                          *(1.25*(rr(1,i)+rb(1,i))*(1.+qv(1,i))   &amp;
-                            -.25*(rr(2,i)+rb(2,i))*(1.+qv(2,i)))
+                          *(1.25*(rr(1,i)+rb(1,i))*(1.+scalars(index_qv,1,i))   &amp;
+                            -.25*(rr(2,i)+rb(2,i))*(1.+scalars(index_qv,2,i)))
 
             ppi(1) = ppi(1)-ppb(1,i)
             do k=1,nz1-1
               ppi(k+1) = ppi(k)-.5*dzu(k+1)*gravity*                     &amp;
-                            (rr(k  ,i)+(rr(k  ,i)+rb(k  ,i))*qv(k  ,i)   &amp;
-                            +rr(k+1,i)+(rr(k+1,i)+rb(k+1,i))*qv(k+1,i))
+                            (rr(k  ,i)+(rr(k  ,i)+rb(k  ,i))*scalars(index_qv,k  ,i)   &amp;
+                            +rr(k+1,i)+(rr(k+1,i)+rb(k+1,i))*scalars(index_qv,k+1,i))
             end do
 
             do k=1,nz1
@@ -603,14 +650,25 @@
           rho_zz (k,i) = rb(k,i) + rr(k,i)
         end do
 
-        if(i == 1) then
-          do k=1,nz1
-            write(0,*) ' k, p, t, rt ',k,p(k,1),t(k,1),rt(k,1)
-          enddo
-        end if
+        !calculation of surface pressure:
+        surface_pressure(i) = 0.5*dzw(1)*gravity                                    &amp;
+                        * (1.25*(rr(1,i) + rb(1,i)) * (1. + scalars(index_qv,1,i))  &amp;
+                        -  0.25*(rr(2,i) + rb(2,i)) * (1. + scalars(index_qv,2,i)))
+        surface_pressure(i) = surface_pressure(i) + pp(1,i) + ppb(1,i)
 
       end do  ! end loop over cells
 
+      !write(0,*)
+      !write(0,*) '--- initialization of water vapor:'
+      !do iCell = 1, grid % nCells
+      !   if(iCell == 1 .or. iCell == grid % nCells) then
+      !      do k = nz1, 1, -1
+      !         write(0,202) iCell,k,t(k,iCell),relhum(k,iCell),qsat(k,iCell),scalars(index_qv,k,iCell)
+      !      enddo
+      !      write(0,*)
+      !   endif
+      !enddo
+
       lat_pert = latitude_pert*pii/180.
       lon_pert = longitude_pert*pii/180.
 
@@ -626,7 +684,7 @@
 
          if (config_test_case == 2) then
             r_pert = sphere_distance( grid % latEdge % array (iEdge), grid % lonEdge % array (iEdge), &amp;
-                                      lat_pert, lon_pert, 1.)/(pert_radius)
+                                      lat_pert, lon_pert, 1.0_RKIND)/(pert_radius)
             u_pert = u_perturbation*exp(-r_pert**2)*(lat2-lat1)*a/grid % dvEdge % array(iEdge)
 
          else if (config_test_case == 3) then
@@ -637,34 +695,30 @@
             u_pert = 0.0
          end if
 
-         call init_atm_calc_flux_zonal(u_2d,etavs_2d,lat_2d,flux_zonal,lat1,lat2,grid % dvEdge % array(iEdge),a,u0,nz1,nlat)
+         if (rebalance) then
 
-         do k=1,grid % nVertLevels
-!!           etavs = (0.5*(ppb(k,iCell1)+ppb(k,iCell2)+pp(k,iCell1)+pp(k,iCell2))/p0 - 0.252)*pii/2.
-!           etavs = (0.5*(ppb(k,1)+ppb(k,1)+pp(k,1)+pp(k,1))/p0 - 0.252)*pii/2.
-           etavs = (0.5*(ppb(k,440)+ppb(k,440)+pp(k,440)+pp(k,440))/p0 - 0.252)*pii/2.  ! 10262 mesh
-!           etavs = (0.5*(ppb(k,505)+ppb(k,505)+pp(k,505)+pp(k,505))/p0 - 0.252)*pii/2.  ! 40962 mesh
-  
-!           fluxk = u0*flux*(cos(etavs)**1.5)
+           call init_atm_calc_flux_zonal(u_2d,etavs_2d,lat_2d,flux_zonal,lat1,lat2,grid % dvEdge % array(iEdge),a,u0,nz1,nlat)
+           do k=1,grid % nVertLevels
+             fluxk = u0*flux_zonal(k)/(0.5*(rb(k,iCell1)+rb(k,iCell2)+rr(k,iCell1)+rr(k,iCell2)))
+             state % u % array(k,iEdge) = fluxk + u_pert
+           end do
 
-            fluxk = u0*flux_zonal(k)/(0.5*(rb(k,iCell1)+rb(k,iCell2)+rr(k,iCell1)+rr(k,iCell2)))
+         else 
 
-!           if(k.eq.18) then
-!              write(21,*) ' iEdge, u1, u2 ',iEdge,fluxk,u0*flux_zonal(k)
-!           end if
-!!           fluxk = u0*flux*(cos(znuv(k))**(1.5))
-!!           fluxk = u0 * cos(grid % angleEdge % array(iEdge)) * (sin(lat1+lat2)**2) *(cos(etavs)**1.5)
-           state % u % array(k,iEdge) = fluxk + u_pert
-         end do
+           do k=1,grid % nVertLevels
+             etavs = (0.5*(ppb(k,iCell1)+ppb(k,iCell2)+pp(k,iCell1)+pp(k,iCell2))/p0 - 0.252)*pii/2.
+             fluxk = u0*flux*(cos(etavs)**1.5)
+             state % u % array(k,iEdge) = fluxk + u_pert
+           end do
 
-         cell1 = grid % CellsOnEdge % array(1,i)
-         cell2 = grid % CellsOnEdge % array(2,i)
-         if(cell1 &lt;= nCellsSolve .or. cell2 &lt;= nCellsSolve) then
-           do k=1,nz1
-             diag % ru % array (k,iEdge)  = 0.5*(rho_zz(k,cell1)+rho_zz(k,cell2))*state % u % array (k,iEdge)
-           end do
          end if
 
+         cell1 = grid % CellsOnEdge % array(1,iEdge)
+         cell2 = grid % CellsOnEdge % array(2,iEdge)
+         do k=1,nz1
+            diag % ru % array (k,iEdge)  = 0.5*(rho_zz(k,cell1)+rho_zz(k,cell2))*state % u % array (k,iEdge)
+         end do
+
       !
       ! Generate rotated Coriolis field
       !
@@ -692,55 +746,49 @@
       do iEdge = 1,grid % nEdges
          cell1 = CellsOnEdge(1,iEdge)
          cell2 = CellsOnEdge(2,iEdge)
-         if (cell1 &lt;= nCellsSolve .or. cell2 &lt;= nCellsSolve ) then
 
-            do k = 1, grid%nVertLevels
+         do k = 1, grid%nVertLevels
 
-               if (config_theta_adv_order == 2) then
+            if (config_theta_adv_order == 2) then
 
-                  z_edge = (zgrid(k,cell1)+zgrid(k,cell2))/2.
+               z_edge = (zgrid(k,cell1)+zgrid(k,cell2))/2.
 
-               else if (config_theta_adv_order == 3 .or. config_theta_adv_order ==4) then !theta_adv_order == 3 or 4 
+            else if (config_theta_adv_order == 3 .or. config_theta_adv_order ==4) then !theta_adv_order == 3 or 4 
 
-                  d2fdx2_cell1 = deriv_two(1,1,iEdge) * zgrid(k,cell1)
-                  d2fdx2_cell2 = deriv_two(1,2,iEdge) * zgrid(k,cell2)
-                  do i=1, grid % nEdgesOnCell % array (cell1)
-                     if ( grid % CellsOnCell % array (i,cell1) &gt; 0)       &amp;
-                     d2fdx2_cell1 = d2fdx2_cell1 + deriv_two(i+1,1,iEdge) * zgrid(k,grid % CellsOnCell % array (i,cell1))
-                  end do
-                  do i=1, grid % nEdgesOnCell % array (cell2)
-                     if ( grid % CellsOnCell % array (i,cell2) &gt; 0)       &amp;
-                     d2fdx2_cell2 = d2fdx2_cell2 + deriv_two(i+1,2,iEdge) * zgrid(k,grid % CellsOnCell % array (i,cell2))
-                  end do
+               d2fdx2_cell1 = deriv_two(1,1,iEdge) * zgrid(k,cell1)
+               d2fdx2_cell2 = deriv_two(1,2,iEdge) * zgrid(k,cell2)
+               do i=1, grid % nEdgesOnCell % array (cell1)
+                  d2fdx2_cell1 = d2fdx2_cell1 + deriv_two(i+1,1,iEdge) * zgrid(k,grid % CellsOnCell % array (i,cell1))
+                  d2fdx2_cell2 = d2fdx2_cell2 + deriv_two(i+1,2,iEdge) * zgrid(k,grid % CellsOnCell % array (i,cell2))
+               end do
 
-                  z_edge =  0.5*(zgrid(k,cell1) + zgrid(k,cell2))         &amp;
-                                - (grid % dcEdge % array(iEdge) **2) * (d2fdx2_cell1 + d2fdx2_cell2) / 12.
+               z_edge =  0.5*(zgrid(k,cell1) + zgrid(k,cell2))         &amp;
+                             - (grid % dcEdge % array(iEdge) **2) * (d2fdx2_cell1 + d2fdx2_cell2) / 12.
 
-                  if (config_theta_adv_order == 3) then
-                     z_edge3 =  - (grid % dcEdge % array(iEdge) **2) * (d2fdx2_cell1 - d2fdx2_cell2) / 12.
-                  else
-                     z_edge3 = 0.
-                  end if
-
+               if (config_theta_adv_order == 3) then
+                  z_edge3 =  - (grid % dcEdge % array(iEdge) **2) * (d2fdx2_cell1 - d2fdx2_cell2) / 12.
+               else
+                  z_edge3 = 0.
                end if
 
-                  zb(k,1,iEdge) = (z_edge-zgrid(k,cell1))*dvEdge(iEdge)/AreaCell(cell1)
-                  zb(k,2,iEdge) = (z_edge-zgrid(k,cell2))*dvEdge(iEdge)/AreaCell(cell2)
-                  zb3(k,1,iEdge)=  z_edge3*dvEdge(iEdge)/AreaCell(cell1)
-                  zb3(k,2,iEdge)=  z_edge3*dvEdge(iEdge)/AreaCell(cell2)
+            end if
 
-                  if (k /= 1) then
-                     zf(k,1,iEdge) = ( fzm(k)*zz(k,cell1)+fzp(k)*zz(k-1,cell1))*zb(k,1,iEdge)
-                     zf(k,2,iEdge) = ( fzm(k)*zz(k,cell2)+fzp(k)*zz(k-1,cell2))*zb(k,2,iEdge)
-                     zf3(k,1,iEdge)= ( fzm(k)*zz(k,cell1)+fzp(k)*zz(k-1,cell1))*zb3(k,1,iEdge)
-                     zf3(k,2,iEdge)= ( fzm(k)*zz(k,cell2)+fzp(k)*zz(k-1,cell2))*zb3(k,2,iEdge)
-                  end if
+               zb(k,1,iEdge) = (z_edge-zgrid(k,cell1))*dvEdge(iEdge)/AreaCell(cell1)
+               zb(k,2,iEdge) = (z_edge-zgrid(k,cell2))*dvEdge(iEdge)/AreaCell(cell2)
+               zb3(k,1,iEdge)=  z_edge3*dvEdge(iEdge)/AreaCell(cell1)
+               zb3(k,2,iEdge)=  z_edge3*dvEdge(iEdge)/AreaCell(cell2)
 
-            end do
+               if (k /= 1) then
+                  zf(k,1,iEdge) = ( fzm(k)*zz(k,cell1)+fzp(k)*zz(k-1,cell1))*zb(k,1,iEdge)
+                  zf(k,2,iEdge) = ( fzm(k)*zz(k,cell2)+fzp(k)*zz(k-1,cell2))*zb(k,2,iEdge)
+                  zf3(k,1,iEdge)= ( fzm(k)*zz(k,cell1)+fzp(k)*zz(k-1,cell1))*zb3(k,1,iEdge)
+                  zf3(k,2,iEdge)= ( fzm(k)*zz(k,cell2)+fzp(k)*zz(k-1,cell2))*zb3(k,2,iEdge)
+               end if
 
-         end if
-       end do
+         end do
 
+      end do
+
       ! for including terrain
       diag % rw % array = 0.
       state % w % array = 0.
@@ -749,7 +797,6 @@
          cell1 = CellsOnEdge(1,iEdge)
          cell2 = CellsOnEdge(2,iEdge)
 
-         if (cell1 &lt;= nCellsSolve .or. cell2 &lt;= nCellsSolve ) then
          do k = 2, grid%nVertLevels
             flux =  (fzm(k)*diag % ru % array(k,iEdge)+fzp(k)*diag % ru % array(k-1,iEdge))
             diag % rw % array(k,cell2) = diag % rw % array(k,cell2) + zf(k,2,iEdge)*flux
@@ -757,13 +804,12 @@
 
             if (config_theta_adv_order ==3) then 
                diag % rw % array(k,cell2) = diag % rw % array(k,cell2)    &amp;
-                                            - sign(1.,diag % ru % array(k,iEdge))*config_coef_3rd_order*zf3(k,2,iEdge)*flux
+                                            - sign(1.0_RKIND,diag % ru % array(k,iEdge))*config_coef_3rd_order*zf3(k,2,iEdge)*flux
                diag % rw % array(k,cell1) = diag % rw % array(k,cell1)    &amp;
-                                            + sign(1.,diag % ru % array(k,iEdge))*config_coef_3rd_order*zf3(k,1,iEdge)*flux
+                                            + sign(1.0_RKIND,diag % ru % array(k,iEdge))*config_coef_3rd_order*zf3(k,1,iEdge)*flux
             end if
 
          end do
-         end if
 
       end do
 
@@ -783,11 +829,9 @@
       do iEdge = 1, grid%nEdges
          do i=1,nEdgesOnEdge(iEdge)
             eoe = edgesOnEdge(i,iEdge)
-            if (eoe &gt; 0) then
-               do k = 1, grid%nVertLevels
-                 diag % v % array(k,iEdge) = diag % v %array(k,iEdge) + weightsOnEdge(i,iEdge) * state % u % array(k, eoe)
-              end do
-            end if
+            do k = 1, grid%nVertLevels
+               diag % v % array(k,iEdge) = diag % v %array(k,iEdge) + weightsOnEdge(i,iEdge) * state % u % array(k, eoe)
+           end do
          end do
       end do
 
@@ -795,8 +839,8 @@
         psurf = (cf1*(ppb(1,i)+pp(1,i)) + cf2*(ppb(2,i)+pp(2,i)) + cf3*(ppb(3,i)+pp(3,i)))/100.
 
             psurf = (ppb(1,i)+pp(1,i)) + .5*dzw(1)*gravity        &amp;
-                          *(1.25*(rr(1,i)+rb(1,i))*(1.+qv(1,i))   &amp;
-                            -.25*(rr(2,i)+rb(2,i))*(1.+qv(2,i)))
+                          *(1.25*(rr(1,i)+rb(1,i))*(1.+scalars(index_qv,1,i))   &amp;
+                            -.25*(rr(2,i)+rb(2,i))*(1.+scalars(index_qv,2,i)))
 
         write(0,*) ' i, psurf, lat ',i,psurf,grid%latCell%array(i)*180./3.1415828
       enddo
@@ -805,7 +849,7 @@
       do iCell=1,grid%nCells
          do k=1,grid%nVertLevels
             diag % rho % array(k,iCell) = state % rho_zz % array(k,iCell) * zz(k,iCell)
-            diag % theta % array(k,iCell) = state % theta_m % array(k,iCell) / (1.0 + 1.61 * scalars(state % index_qv,k,iCell))
+            diag % theta % array(k,iCell) = state % theta_m % array(k,iCell) / (1.0 + 1.61 * scalars(index_qv,k,iCell))
          end do
       end do
 
@@ -815,7 +859,7 @@
 
    implicit none
    integer, intent(in) :: nz1,nlat
-   real (kind=RKIND), dimension(nlat,nz1), intent(in) :: u_2d,etavs_2d
+   real (kind=RKIND), dimension(nz1,nlat), intent(in) :: u_2d,etavs_2d
    real (kind=RKIND), dimension(nlat), intent(in) :: lat_2d
    real (kind=RKIND), dimension(nz1), intent(out) :: flux_zonal
    real (kind=RKIND), intent(in) :: lat1_in, lat2_in, dvEdge, a, u0
@@ -845,7 +889,7 @@
      w2 = 0.5*(db-da)**2
 
      do k=1,nz1
-       flux_zonal(k) = flux_zonal(k) + w1*u_2d(i,k) + w2*u_2d(i+1,k)
+       flux_zonal(k) = flux_zonal(k) + w1*u_2d(k,i) + w2*u_2d(k,i+1)
      end do
 
      end if
@@ -855,12 +899,105 @@
 !  renormalize for setting cell-face fluxes
 
    do k=1,nz1
-     flux_zonal(k) = sign(1.,lat2_in-lat1_in)*flux_zonal(k)*dlat*a/dvEdge/u0
+     flux_zonal(k) = sign(1.0_RKIND,lat2_in-lat1_in)*flux_zonal(k)*dlat*a/dvEdge/u0
    end do
      
    end subroutine init_atm_calc_flux_zonal
 
+   !SHP-balance
+   subroutine init_atm_recompute_geostrophic_wind(u_2d,rho_2d,pp_2d,qv_2d,lat_2d,zz_2d,zx_2d,     &amp;
+                                         cf1,cf2,cf3,fzm,fzp,rdzw,nz1,nlat,dlat)
 
+   implicit none
+   integer, intent(in) :: nz1,nlat
+   real (kind=RKIND), dimension(nz1,nlat), intent(inout) :: u_2d
+   real (kind=RKIND), dimension(nz1,nlat), intent(in) :: rho_2d, pp_2d, qv_2d, zz_2d
+   real (kind=RKIND), dimension(nz1,nlat-1), intent(in) :: zx_2d
+   real (kind=RKIND), dimension(nlat), intent(in) :: lat_2d
+   real (kind=RKIND), dimension(nz1), intent(in) :: fzm, fzp, rdzw
+   real (kind=RKIND), intent(in) :: cf1, cf2, cf3, dlat
+
+   !local variable
+   real (kind=RKIND), dimension(nz1,nlat-1) :: pgrad, ru, u
+   real (kind=RKIND), dimension(nlat-1) :: f
+   real (kind=RKIND), dimension(nz1+1)  :: dpzx
+
+   real (kind=RKIND), parameter :: omega_e = 7.29212e-05
+   real (kind=RKIND) :: rdx, qtot, r_earth, phi
+   integer :: k,i, itr
+
+   r_earth  = a
+   rdx = 1./(dlat*r_earth)
+
+   do i=1,nlat-1
+     do k=1,nz1
+       pgrad(k,i) = rdx*(pp_2d(k,i+1)/zz_2d(k,i+1)-pp_2d(k,i)/zz_2d(k,i))
+     end do
+
+     dpzx(:) = 0.
+
+     k=1
+     dpzx(k) = .5*zx_2d(k,i)*(cf1*(pp_2d(k  ,i+1)+pp_2d(k  ,i))        &amp;
+                             +cf2*(pp_2d(k+1,i+1)+pp_2d(k+1,i))        &amp;
+                             +cf3*(pp_2d(k+2,i+1)+pp_2d(k+2,i)))
+     do k=2,nz1
+        dpzx(k) = .5*zx_2d(k,i)*(fzm(k)*(pp_2d(k  ,i+1)+pp_2d(k  ,i))   &amp;
+                                +fzp(k)*(pp_2d(k-1,i+1)+pp_2d(k-1,i)))
+     end do
+
+     do k=1,nz1
+        pgrad(k,i) = pgrad(k,i) - rdzw(k)*(dpzx(k+1)-dpzx(k))
+     end do
+   end do
+
+
+   !initial value of v and rv -&gt; that is from analytic sln. 
+   do i=1,nlat-1
+      do k=1,nz1
+         u(k,i) = .5*(u_2d(k,i)+u_2d(k,i+1))
+         ru(k,i) = u(k,i)*(rho_2d(k,i)+rho_2d(k,i+1))*.5
+      end do
+   end do
+
+   write(0,*) &quot;MAX U wind before REBALANCING ----&gt;&quot;, maxval(abs(u))
+
+   !re-calculate geostrophic wind using iteration 
+   do itr=1,50
+   do i=1,nlat-1
+      phi = (lat_2d(i)+lat_2d(i+1))/2.
+      f(i) = 2.*omega_e*sin(phi)
+      do k=1,nz1
+         if (f(i).eq.0.) then
+           ru(k,i) = 0.
+         else
+           qtot = .5*(qv_2d(k,i)+qv_2d(k,i+1))
+           ru(k,i) = - ( 1./(1.+qtot)*pgrad(k,i) + tan(phi)/r_earth*u(k,i)*ru(k,i) )/f(i)
+         end if
+           u(k,i) = ru(k,i)*2./(rho_2d(k,i)+rho_2d(k,i+1))
+      end do
+   end do
+   end do
+
+   write(0,*) &quot;MAX U wind after REBALANCING ----&gt;&quot;, maxval(abs(u))
+
+   !update 2d ru
+   do i=2,nlat-1
+     do k=1,nz1
+       u_2d(k,i) = (ru(k,i-1)+ru(k,i))*.5
+     end do
+   end do
+
+   i=1
+   do k=1,nz1
+      u_2d(k,i) = (3.*u_2d(k,i+1)-u_2d(k,i+2))*.5
+   end do
+   i=nlat
+   do k=1,nz1
+      u_2d(k,i) = (3.*u_2d(k,i-1)-u_2d(k,i-2))*.5
+   end do
+
+
+   end subroutine init_atm_recompute_geostrophic_wind
 !----------------------------------------------------------------------------------------------------------
 
    subroutine init_atm_test_case_squall_line(dminfo, grid, state, diag, test_case)
@@ -885,7 +1022,7 @@
       integer :: eoe, j
       integer, dimension(:), pointer :: nEdgesOnEdge 
       integer, dimension(:,:), pointer :: edgesOnEdge
-      real, dimension(:,:), pointer :: weightsOnEdge
+      real (kind=RKIND), dimension(:,:), pointer :: weightsOnEdge
 
       integer :: iCell, iCell1, iCell2 , iEdge, vtx1, vtx2, ivtx, i, k, nz, nz1, itr, cell1, cell2, nCellsSolve
       integer :: index_qv
@@ -1174,7 +1311,7 @@
       call mpas_dmpar_bcast_real(dminfo, pibtop)
 
       ptopb = p0*pibtop**(1./rcp)
-      write(0,*) 'ptopb = ',.01*ptopb
+      write(6,*) 'ptopb = ',.01*ptopb
 
       do i=1, grid % nCells
          pb(nz1,i) = pibtop+.5*dzw(nz1)*gravity*(1.+qvb(nz1))/(cp*tb(nz1,i)*zz(nz1,i))
@@ -1189,6 +1326,7 @@
             rb (k,i) = pb(k,i)**(1./rcv)/((rgas/p0)*tb(k,i)*zz(k,i))
             rtb(k,i) = rb(k,i)*tb(k,i)
             rr (k,i) = p (k,i)**(1./rcv)/((rgas/p0)*t (k,i)*zz(k,i))-rb(k,i)
+            ppb(k,i) = p0*(zz(k,i)*rgas*rtb(k,i)/p0)**(cp/cv)
          end do
       end do
 
@@ -1200,7 +1338,7 @@
             temp     = p(k,i)*thi(k,i)
             pres     = p0*p(k,i)**(1./rcp)
             qvs      = 380.*exp(17.27*(temp-273.)/(temp-36.))/pres
-            scalars(index_qv,k,i) = amin1(0.014,rh(k,i)*qvs)
+            scalars(index_qv,k,i) = min(0.014_RKIND,rh(k,i)*qvs)
          end do
       end do
 
@@ -1283,7 +1421,7 @@
         ptop = p0*pitop**(1./rcp)
         write(0,*) 'ptop  = ',.01*ptop, .01*ptopb
 
-        call mpas_dmpar_bcast_real(dminfo, pitop)
+        call mpas_dmpar_bcast_real(dminfo, ptop)
 
         do i = 1, grid % nCells
 
@@ -1301,7 +1439,7 @@
           end do
           if (itr==1.and.i==1) then
           do k=1,nz1
-          print *, &quot;pp-check&quot;, pp(k,i) 
+          write(0,*) &quot;pp-check&quot;, pp(k,i) 
           end do
           end if
           do k=1,nz1
@@ -1389,7 +1527,7 @@
       do iCell=1,grid%nCells
          do k=1,grid%nVertLevels
             diag % rho % array(k,iCell) = state % rho_zz % array(k,iCell) * zz(k,iCell)
-            diag % theta % array(k,iCell) = state % theta_m % array(k,iCell) / (1.0 + 1.61 * scalars(state % index_qv,k,iCell))
+            diag % theta % array(k,iCell) = state % theta_m % array(k,iCell) / (1.0 + 1.61 * scalars(index_qv,k,iCell))
          end do
       end do
 
@@ -1422,8 +1560,8 @@
       integer :: eoe, j
       integer, dimension(:), pointer :: nEdgesOnEdge 
       integer, dimension(:,:), pointer :: edgesOnEdge, CellsOnEdge
-      real, dimension(:), pointer :: dvEdge, AreaCell, xCell, yCell 
-      real, dimension(:,:), pointer :: weightsOnEdge
+      real (kind=RKIND), dimension(:), pointer :: dvEdge, AreaCell, xCell, yCell 
+      real (kind=RKIND), dimension(:,:), pointer :: weightsOnEdge
 
       integer :: iCell, iCell1, iCell2 , iEdge, vtx1, vtx2, ivtx, i, k, nz, nz1, itr, itrp, cell1, cell2, nCellsSolve
       integer :: index_qv
@@ -1433,7 +1571,7 @@
       real (kind=RKIND) :: ztemp, zd, zt, dz, str
 
       real (kind=RKIND), dimension(grid % nVertLevels, grid % nCells) :: rh
-      real (kind=RKIND) :: ptmp, es, qvs, xnutr, ptemp
+      real (kind=RKIND) :: es, qvs, xnutr, ptemp
       integer :: iter
 
       real (kind=RKIND), dimension(grid % nVertLevels + 1 ) :: zc, zw, ah
@@ -1649,7 +1787,7 @@
 ! smoothing grid for the upper level &gt;&gt; but not propoer for parallel programing 
       dzmin=.7
       do k=2,nz1
-         sm = .25*min((zc(k)-zc(k-1))/dz,1.)
+         sm = .25*min((zc(k)-zc(k-1))/dz,1.0_RKIND)
          do i=1,grid % nCells
             hx(k,i) = hx(k-1,i)
          end do
@@ -1674,7 +1812,7 @@
             end if
 
          end do !end of iteration for smoothing
-99       print *,&quot;PASS-SHP&quot;
+99       write(0,*) &quot;PASS-SHP&quot;
       end do
 
       do iCell=1,grid % nCells
@@ -1826,7 +1964,7 @@
               temp   = p(k,i)*t(k,i)
               pres   = p0*p(k,i)**(1./rcp)
               qvs    = 380.*exp(17.27*(temp-273.)/(temp-36.))/pres
-              scalars(index_qv,k,i) = amin1(0.014,rh(k,i)*qvs)
+              scalars(index_qv,k,i) = min(0.014_RKIND,rh(k,i)*qvs)
            end do
                          
            do k=1,nz1
@@ -1845,7 +1983,7 @@
       write(0,*) ' *** sounding for the simulation ***'
       write(0,*) '    z       theta       pres         qv       rho_m        u        rr'
       do k=1,nz1
-         write(0,'(8(f14.9,2x))') .5*(zgrid(k,1)+zgrid(k+1,1))/1000.,   &amp;
+         write(6,'(8(f14.9,2x))') .5*(zgrid(k,1)+zgrid(k+1,1))/1000.,   &amp;
                        t(k,1)/(1.+1.61*scalars(index_qv,k,1)),        &amp;
                        .01*p0*p(k,1)**(1./rcp),                       &amp;
                        1000.*scalars(index_qv,k,1),                   &amp;
@@ -1949,9 +2087,9 @@
 
             if (config_theta_adv_order ==3) then
                diag % rw % array(k,cell2) = diag % rw % array(k,cell2)    &amp;
-                                            - sign(1.,ru(k,iEdge))*config_coef_3rd_order*zf3(k,2,iEdge)*flux
+                                            - sign(1.0_RKIND,ru(k,iEdge))*config_coef_3rd_order*zf3(k,2,iEdge)*flux
                diag % rw % array(k,cell1) = diag % rw % array(k,cell1)    &amp;
-                                            + sign(1.,ru(k,iEdge))*config_coef_3rd_order*zf3(k,1,iEdge)*flux
+                                            + sign(1.0_RKIND,ru(k,iEdge))*config_coef_3rd_order*zf3(k,1,iEdge)*flux
             end if
 
          end do
@@ -1999,7 +2137,7 @@
       do iCell=1,grid%nCells
          do k=1,grid%nVertLevels
             diag % rho % array(k,iCell) = state % rho_zz % array(k,iCell) * zz(k,iCell)
-            diag % theta % array(k,iCell) = state % theta_m % array(k,iCell) / (1.0 + 1.61 * scalars(state % index_qv,k,iCell))
+            diag % theta % array(k,iCell) = state % theta_m % array(k,iCell) / (1.0 + 1.61 * scalars(index_qv,k,iCell))
          end do
       end do
 
@@ -2069,7 +2207,7 @@
       real (kind=RKIND), dimension(:,:), pointer :: sorted_arr
 
       real(kind=RKIND), dimension(:), pointer :: hs, hs1
-      real(kind=RKIND) :: hm, zh, dzmin, dzmina, dzminf, sm
+      real(kind=RKIND) :: hm, zh, dzmin, dzmina, dzmina_global, dzminf, sm
       integer :: nsmterrain, kz, sfc_k
       logical :: hybrid, smooth
 
@@ -2520,12 +2658,12 @@
       grid % soiltemp % array(:) = 0.0
 
       call map_set(PROJ_LATLON, proj, &amp;
-                   latinc = 1.0, &amp;
-                   loninc = 1.0, &amp;
-                   knowni = 1.0, &amp;
-                   knownj = 1.0, &amp;
-                   lat1 = -89.5, &amp;
-                   lon1 = -179.5)
+                   latinc = 1.0_RKIND, &amp;
+                   loninc = 1.0_RKIND, &amp;
+                   knowni = 1.0_RKIND, &amp;
+                   knownj = 1.0_RKIND, &amp;
+                   lat1 = -89.5_RKIND, &amp;
+                   lon1 = -179.5_RKIND)
 
       write(fname,'(a,i5.5,a1,i5.5,a1,i5.5,a1,i5.5)') trim(config_geog_data_path)//'soiltemp_1deg/',1,'-',180,'.',1,'-',180
 write(0,*) trim(fname)
@@ -2568,8 +2706,8 @@
             end if
 if (y &lt; 1.0) y = 1.0
 if (y &gt; 179.0) y = 179.0
-!            grid % soiltemp % array(iCell) = interp_sequence(x, y, 1, soiltemp_1deg, 1, 360, 1, 180, 1, 1, -1.e30, interp_list, 1)
-            grid % soiltemp % array(iCell) = interp_sequence(x, y, 1, soiltemp_1deg, 1, 360, 1, 180, 1, 1, 0., interp_list, 1)
+!            grid % soiltemp % array(iCell) = interp_sequence(x, y, 1, soiltemp_1deg, 1, 360, 1, 180, 1, 1, -1.e30_RKIND, interp_list, 1)
+            grid % soiltemp % array(iCell) = interp_sequence(x, y, 1, soiltemp_1deg, 1, 360, 1, 180, 1, 1, 0.0_RKIND, interp_list, 1)
          else
             grid % soiltemp % array(iCell) = 0.0
          end if
@@ -2595,12 +2733,12 @@
       grid % snoalb % array(:) = 0.0
 
       call map_set(PROJ_LATLON, proj, &amp;
-                   latinc = 1.0, &amp;
-                   loninc = 1.0, &amp;
-                   knowni = 1.0, &amp;
-                   knownj = 1.0, &amp;
-                   lat1 = -89.5, &amp;
-                   lon1 = -179.5)
+                   latinc = 1.0_RKIND, &amp;
+                   loninc = 1.0_RKIND, &amp;
+                   knowni = 1.0_RKIND, &amp;
+                   knownj = 1.0_RKIND, &amp;
+                   lat1 = -89.5_RKIND, &amp;
+                   lon1 = -179.5_RKIND)
 
       write(fname,'(a,i5.5,a1,i5.5,a1,i5.5,a1,i5.5)') trim(config_geog_data_path)//'maxsnowalb/',1,'-',180,'.',1,'-',180
 write(0,*) trim(fname)
@@ -2643,8 +2781,8 @@
             end if
 if (y &lt; 1.0) y = 1.0
 if (y &gt; 179.0) y = 179.0
-!            grid % snoalb % array(iCell) = interp_sequence(x, y, 1, maxsnowalb, 1, 360, 1, 180, 1, 1, -1.e30, interp_list, 1)
-            grid % snoalb % array(iCell) = interp_sequence(x, y, 1, maxsnowalb, 1, 360, 1, 180, 1, 1, 0., interp_list, 1)
+!            grid % snoalb % array(iCell) = interp_sequence(x, y, 1, maxsnowalb, 1, 360, 1, 180, 1, 1, -1.e30_RKIND, interp_list, 1)
+            grid % snoalb % array(iCell) = interp_sequence(x, y, 1, maxsnowalb, 1, 360, 1, 180, 1, 1, 0.0_RKIND, interp_list, 1)
          else
             grid % snoalb % array(iCell) = 0.0
          end if
@@ -2673,12 +2811,12 @@
       grid % greenfrac % array(:,:) = 0.0
 
       call map_set(PROJ_LATLON, proj, &amp;
-                   latinc = 0.144, &amp;
-                   loninc = 0.144, &amp;
-                   knowni = 1.0, &amp;
-                   knownj = 1.0, &amp;
-                   lat1 = -89.928, &amp;
-                   lon1 = -179.928)
+                   latinc = 0.144_RKIND, &amp;
+                   loninc = 0.144_RKIND, &amp;
+                   knowni = 1.0_RKIND, &amp;
+                   knownj = 1.0_RKIND, &amp;
+                   lat1 = -89.928_RKIND, &amp;
+                   lon1 = -179.928_RKIND)
 
       write(fname,'(a,i5.5,a1,i5.5,a1,i5.5,a1,i5.5)') trim(config_geog_data_path)//'greenfrac/',1,'-',1250,'.',1,'-',1250
 write(0,*) trim(fname)
@@ -2715,7 +2853,7 @@
 if (y &lt; 1.0) y = 1.0
 if (y &gt; 1249.0) y = 1249.0
             do k=1,12
-               grid % greenfrac % array(k,iCell) = interp_sequence(x, y, k, vegfra, 1, 2500, 1, 1250, 1, 12, -1.e30, interp_list, 1)
+               grid % greenfrac % array(k,iCell) = interp_sequence(x, y, k, vegfra, 1, 2500, 1, 1250, 1, 12, -1.e30_RKIND, interp_list, 1)
             end do
          else
             grid % greenfrac % array(:,iCell) = 0.0
@@ -2744,12 +2882,12 @@
       grid % albedo12m % array(:,:) = 0.0
 
       call map_set(PROJ_LATLON, proj, &amp;
-                   latinc = 0.144, &amp;
-                   loninc = 0.144, &amp;
-                   knowni = 1.0, &amp;
-                   knownj = 1.0, &amp;
-                   lat1 = -89.928, &amp;
-                   lon1 = -179.928)
+                   latinc = 0.144_RKIND, &amp;
+                   loninc = 0.144_RKIND, &amp;
+                   knowni = 1.0_RKIND, &amp;
+                   knownj = 1.0_RKIND, &amp;
+                   lat1 = -89.928_RKIND, &amp;
+                   lon1 = -179.928_RKIND)
 
       write(fname,'(a,i5.5,a1,i5.5,a1,i5.5,a1,i5.5)') trim(config_geog_data_path)//'albedo_ncep/',1,'-',1250,'.',1,'-',1250
 write(0,*) trim(fname)
@@ -2786,7 +2924,7 @@
 if (y &lt; 1.0) y = 1.0
 if (y &gt; 1249.0) y = 1249.0
             do k=1,12
-               grid % albedo12m % array(k,iCell) = interp_sequence(x, y, k, vegfra, 1, 2500, 1, 1250, 1, 12, 0.0, interp_list, 1)
+               grid % albedo12m % array(k,iCell) = interp_sequence(x, y, k, vegfra, 1, 2500, 1, 1250, 1, 12, 0.0_RKIND, interp_list, 1)
             end do
          else
             grid % albedo12m % array(:,iCell) = 8.0
@@ -2822,12 +2960,12 @@
           
             if (field % iproj == PROJ_LATLON) then
                call map_set(PROJ_LATLON, proj, &amp;
-                            latinc = real(field % deltalat), &amp;
-                            loninc = real(field % deltalon), &amp;
-                            knowni = 1.0, &amp;
-                            knownj = 1.0, &amp;
-                            lat1 = real(field % startlat), &amp;
-                            lon1 = real(field % startlon))
+                            latinc = real(field % deltalat,RKIND), &amp;
+                            loninc = real(field % deltalon,RKIND), &amp;
+                            knowni = 1.0_RKIND, &amp;
+                            knownj = 1.0_RKIND, &amp;
+                            lat1 = real(field % startlat,RKIND), &amp;
+                            lon1 = real(field % startlon,RKIND))
             end if
 
 
@@ -2848,9 +2986,9 @@
                   call latlon_to_ij(proj, lat, lon, x, y)
                end if
                if (ndims == 1) then
-                  destField1d(i) = interp_sequence(x, y, 1, field % slab, 1, field % nx, 1, field % ny, 1, 1, -1.e30, interp_list, 1)
+                  destField1d(i) = interp_sequence(x, y, 1, field % slab, 1, field % nx, 1, field % ny, 1, 1, -1.e30_RKIND, interp_list, 1)
                else if (ndims == 2) then
-                  destField2d(k,i) = interp_sequence(x, y, 1, field % slab, 1, field % nx, 1, field % ny, 1, 1, -1.e30, interp_list, 1)
+                  destField2d(k,i) = interp_sequence(x, y, 1, field % slab, 1, field % nx, 1, field % ny, 1, 1, -1.e30_RKIND, interp_list, 1)
                end if
             end do
          end if
@@ -3014,9 +3152,9 @@
             hx(k,:) = hx(k-1,:)
             dzminf = zw(k)-zw(k-1)
 
-!            dzmin = max(.5,1.-.5*zw(k)/hm)
+!            dzmin = max(0.5_RKIND,1.-.5*zw(k)/hm)
 
-            sm = .05*min(.5*zw(k)/hm,1.)
+            sm = .05*min(0.5_RKIND*zw(k)/hm,1.0_RKIND)
           
             do i=1,50
                do iCell=1,grid %nCells
@@ -3045,10 +3183,11 @@
 
              !  dzmina = minval(hs(:)-hx(k-1,:))
                dzmina = minval(zw(k)+ah(k)*hs(:)-zw(k-1)-ah(k-1)*hx(k-1,:))
+               call mpas_dmpar_min_real(dminfo, dzmina, dzmina_global)
              !  write(0,*) ' k,i, dzmina, dzmin, zw(k)-zw(k-1) ', k,i, dzmina, dzmin, zw(k)-zw(k-1)
-               if (dzmina &gt;= dzmin*(zw(k)-zw(k-1))) then
+               if (dzmina_global &gt;= dzmin*(zw(k)-zw(k-1))) then
                   hx(k,:)=hs(:)
-                  dzminf = dzmina
+                  dzminf = dzmina_global
                else
                   exit
                end if
@@ -3169,6 +3308,13 @@
 
       if (config_met_interp) then
 
+      !ldf (2011-11-19): added initialization of the sea-surface temperature, seaice fraction, and
+      !seaice flag:
+       fg % sst % array    = 0.
+       fg % xice % array   = 0.
+       fg % seaice % array = 0.
+      !ldf end.
+
       !
       ! First, try to locate the LANDSEA field for use as an interpolation mask
       !
@@ -3289,18 +3435,18 @@
           
             if (field % iproj == PROJ_LATLON) then
                call map_set(PROJ_LATLON, proj, &amp;
-                            latinc = real(field % deltalat), &amp;
-                            loninc = real(field % deltalon), &amp;
-                            knowni = 1.0, &amp;
-                            knownj = 1.0, &amp;
-                            lat1 = real(field % startlat), &amp;
-                            lon1 = real(field % startlon))
+                            latinc = real(field % deltalat,RKIND), &amp;
+                            loninc = real(field % deltalon,RKIND), &amp;
+                            knowni = 1.0_RKIND, &amp;
+                            knownj = 1.0_RKIND, &amp;
+                            lat1 = real(field % startlat,RKIND), &amp;
+                            lon1 = real(field % startlon,RKIND))
             else if (field % iproj == PROJ_GAUSS) then
                call map_set(PROJ_GAUSS, proj, &amp;
                             nlat = nint(field % deltalat), &amp;
-                            loninc = real(field % deltalon), &amp;
-                            lat1 = real(field % startlat), &amp;
-                            lon1 = real(field % startlon))
+                            loninc = real(field % deltalon,RKIND), &amp;
+                            lat1 = real(field % startlat,RKIND), &amp;
+                            lon1 = real(field % startlon,RKIND))
 !                            nxmax = nint(360.0 / field % deltalon), &amp;
             end if
 
@@ -3686,7 +3832,8 @@
       ! Set SST based on SKINTEMP field if it wasn't found in input data
       if (minval(fg % sst % array) == 0.0 .and. maxval(fg % sst % array) == 0.0) then
          write(0,*) 'Setting SST from SKINTEMP'
-         where (grid % landmask % array == 0) fg % sst % array = fg % skintemp % array
+         !where (grid % landmask % array == 0) fg % sst % array = fg % skintemp % array
+         fg % sst % array = fg % skintemp % array
       end if
 
       ! Set SNOWC (snow-cover flag) based on SNOW
@@ -3736,13 +3883,13 @@
           
                if (field % iproj == PROJ_PS) then
                   call map_set(PROJ_PS, proj, &amp;
-                               dx = real(field % dx * 1000.0), &amp;
-                               truelat1 = real(field % truelat1), &amp;
-                               stdlon = real(field % xlonc), &amp;
-                               knowni = real(field % nx / 2.0), &amp;
-                               knownj = real(field % ny / 2.0), &amp;
-                               lat1 = real(field % startlat), &amp;
-                               lon1 = real(field % startlon))
+                               dx = real(field % dx * 1000.0,RKIND), &amp;
+                               truelat1 = real(field % truelat1,RKIND), &amp;
+                               stdlon = real(field % xlonc,RKIND), &amp;
+                               knowni = real(field % nx / 2.0,RKIND), &amp;
+                               knownj = real(field % ny / 2.0,RKIND), &amp;
+                               lat1 = real(field % startlat,RKIND), &amp;
+                               lon1 = real(field % startlon,RKIND))
                end if
 
                if (index(field % field, 'SEAICE') /= 0) then
@@ -3948,6 +4095,20 @@
 
 
       !
+      ! Reconstruct zonal and meridional winds for diagnostic puposes:
+      !
+      call mpas_rbf_interp_initialize(grid)
+      call mpas_init_reconstruct(grid)
+      call mpas_reconstruct(grid, state % u % array,                        &amp;
+                            diag % uReconstructX % array,                   &amp;
+                            diag % uReconstructY % array,                   &amp;
+                            diag % uReconstructZ % array,                   &amp;
+                            diag % uReconstructZonal % array,               &amp;
+                            diag % uReconstructMeridional % array           &amp;
+                           )
+   
+
+      !
       ! Adjust surface pressure for difference in topography
       !
       do sfc_k=1,config_nfglevels
@@ -4083,9 +4244,9 @@
 
             if (config_theta_adv_order ==3) then 
                diag % rw % array(k,cell2) = diag % rw % array(k,cell2)    &amp;
-                                            - sign(1.,diag % ru % array(k,iEdge))*config_coef_3rd_order*zf3(k,2,iEdge)*flux
+                                            - sign(1.0_RKIND,diag % ru % array(k,iEdge))*config_coef_3rd_order*zf3(k,2,iEdge)*flux
                diag % rw % array(k,cell1) = diag % rw % array(k,cell1)    &amp;
-                                            + sign(1.,diag % ru % array(k,iEdge))*config_coef_3rd_order*zf3(k,1,iEdge)*flux
+                                            + sign(1.0_RKIND,diag % ru % array(k,iEdge))*config_coef_3rd_order*zf3(k,1,iEdge)*flux
             end if
 
          end do
@@ -4106,7 +4267,8 @@
       end if    ! config_met_interp
 
      
-      ! Calculate surface pressure
+      ! Calculate surface pressure (This is an ad-hoc calculation. The actual surface pressure is actually re-calculated at
+      !the top of the subroutine MPAS_to_physics in ../core_atmos_physics/mpas_atmphys_interface_nhyd.F
       do iCell=1,grid%nCells
          diag % surface_pressure % array(iCell) = 0.5*gravity/rdzw(1)                                              &amp;
                                                 * (1.25* rho_zz(1,iCell) * (1. + scalars(state % index_qv, 1, iCell))  &amp;
@@ -4214,18 +4376,18 @@
              
                if (field % iproj == PROJ_LATLON) then
                   call map_set(PROJ_LATLON, proj, &amp;
-                               latinc = real(field % deltalat), &amp;
-                               loninc = real(field % deltalon), &amp;
-                               knowni = 1.0, &amp;
-                               knownj = 1.0, &amp;
-                               lat1 = real(field % startlat), &amp;
-                               lon1 = real(field % startlon))
+                               latinc = real(field % deltalat,RKIND), &amp;
+                               loninc = real(field % deltalon,RKIND), &amp;
+                               knowni = 1.0_RKIND, &amp;
+                               knownj = 1.0_RKIND, &amp;
+                               lat1 = real(field % startlat,RKIND), &amp;
+                               lon1 = real(field % startlon,RKIND))
                else if (field % iproj == PROJ_GAUSS) then
                   call map_set(PROJ_GAUSS, proj, &amp;
                                nlat = nint(field % deltalat), &amp;
-                               loninc = real(field % deltalon), &amp;
-                               lat1 = real(field % startlat), &amp;
-                               lon1 = real(field % startlon))
+                               loninc = real(field % deltalon,RKIND), &amp;
+                               lat1 = real(field % startlat,RKIND), &amp;
+                               lon1 = real(field % startlon,RKIND))
 !                               nxmax = nint(360.0 / field % deltalon), &amp;
                end if
    
@@ -4246,7 +4408,7 @@
                      lon = lon - 360.0
                      call latlon_to_ij(proj, lat, lon, x, y)
                   end if
-                  fg % sst % array(iCell) = interp_sequence(x, y, 1, slab_r8, 1, field % nx, 1, field % ny, 1, 1, -1.e30, interp_list, 1)
+                  fg % sst % array(iCell) = interp_sequence(x, y, 1, slab_r8, 1, field % nx, 1, field % ny, 1, 1, -1.e30_RKIND, interp_list, 1)
                end do
 
                deallocate(slab_r8)
@@ -4271,18 +4433,18 @@
              
                if (field % iproj == PROJ_LATLON) then
                   call map_set(PROJ_LATLON, proj, &amp;
-                               latinc = real(field % deltalat), &amp;
-                               loninc = real(field % deltalon), &amp;
-                               knowni = 1.0, &amp;
-                               knownj = 1.0, &amp;
-                               lat1 = real(field % startlat), &amp;
-                               lon1 = real(field % startlon))
+                               latinc = real(field % deltalat,RKIND), &amp;
+                               loninc = real(field % deltalon,RKIND), &amp;
+                               knowni = 1.0_RKIND, &amp;
+                               knownj = 1.0_RKIND, &amp;
+                               lat1 = real(field % startlat,RKIND), &amp;
+                               lon1 = real(field % startlon,RKIND))
                else if (field % iproj == PROJ_GAUSS) then
                   call map_set(PROJ_GAUSS, proj, &amp;
                                nlat = nint(field % deltalat), &amp;
-                               loninc = real(field % deltalon), &amp;
-                               lat1 = real(field % startlat), &amp;
-                               lon1 = real(field % startlon))
+                               loninc = real(field % deltalon,RKIND), &amp;
+                               lat1 = real(field % startlat,RKIND), &amp;
+                               lon1 = real(field % startlon,RKIND))
 !                               nxmax = nint(360.0 / field % deltalon), &amp;
                end if
    
@@ -4303,7 +4465,7 @@
                      lon = lon - 360.0
                      call latlon_to_ij(proj, lat, lon, x, y)
                   end if
-                  fg % xice % array(iCell) = interp_sequence(x, y, 1, slab_r8, 1, field % nx, 1, field % ny, 1, 1, -1.e30, interp_list, 1)
+                  fg % xice % array(iCell) = interp_sequence(x, y, 1, slab_r8, 1, field % nx, 1, field % ny, 1, 1, -1.e30_RKIND, interp_list, 1)
 
                end do
 
@@ -4388,26 +4550,8 @@
    end function four_pt
 #endif
 
+!----------------------------------------------------------------------------------------------------------
 
-   real function sphere_distance(lat1, lon1, lat2, lon2, radius)
-   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-   ! Compute the great-circle distance between (lat1, lon1) and (lat2, lon2) on a
-   !   sphere with given radius.
-   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-
-      implicit none
-
-      real (kind=RKIND), intent(in) :: lat1, lon1, lat2, lon2, radius
-
-      real (kind=RKIND) :: arg1
-
-      arg1 = sqrt( sin(0.5*(lat2-lat1))**2 +  &amp;
-                   cos(lat1)*cos(lat2)*sin(0.5*(lon2-lon1))**2 )
-      sphere_distance = 2.*radius*asin(arg1)
-
-   end function sphere_distance
-
-
    integer function nearest_cell(target_lat, target_lon, &amp;
                                  start_cell, &amp;
                                  nCells, maxEdges, nEdgesOnCell, cellsOnCell, latCell, lonCell)
@@ -4432,12 +4576,12 @@
 
       do while (nearest_cell /= current_cell)
          current_cell = nearest_cell
-         current_distance = sphere_distance(latCell(current_cell), lonCell(current_cell), target_lat, target_lon, 1.0)
+         current_distance = sphere_distance(latCell(current_cell), lonCell(current_cell), target_lat, target_lon, 1.0_RKIND)
          nearest_cell = current_cell
          nearest_distance = current_distance
          do i = 1, nEdgesOnCell(current_cell)
             iCell = cellsOnCell(i,current_cell)
-            d = sphere_distance(latCell(iCell), lonCell(iCell), target_lat, target_lon, 1.0)
+            d = sphere_distance(latCell(iCell), lonCell(iCell), target_lat, target_lon, 1.0_RKIND)
             if (d &lt; nearest_distance) then
                nearest_cell = iCell
                nearest_distance = d
@@ -4475,13 +4619,13 @@
 
       do while (nearest_edge /= current_edge)
          current_edge = nearest_edge
-         current_distance = sphere_distance(latEdge(current_edge), lonEdge(current_edge), target_lat, target_lon, 1.0)
+         current_distance = sphere_distance(latEdge(current_edge), lonEdge(current_edge), target_lat, target_lon, 1.0_RKIND)
          nearest_edge = current_edge
          nearest_distance = current_distance
          cell1 = cellsOnEdge(1,current_edge)
          cell2 = cellsOnEdge(2,current_edge)
-         cell1_dist = sphere_distance(latCell(cell1), lonCell(cell1), target_lat, target_lon, 1.0)
-         cell2_dist = sphere_distance(latCell(cell2), lonCell(cell2), target_lat, target_lon, 1.0)
+         cell1_dist = sphere_distance(latCell(cell1), lonCell(cell1), target_lat, target_lon, 1.0_RKIND)
+         cell2_dist = sphere_distance(latCell(cell2), lonCell(cell2), target_lat, target_lon, 1.0_RKIND)
          if (cell1_dist &lt; cell2_dist) then
             iCell = cell1
          else
@@ -4489,7 +4633,7 @@
          end if
          do i = 1, nEdgesOnCell(iCell)
             iEdge = edgesOnCell(i,iCell)
-            d = sphere_distance(latEdge(iEdge), lonEdge(iEdge), target_lat, target_lon, 1.0)
+            d = sphere_distance(latEdge(iEdge), lonEdge(iEdge), target_lat, target_lon, 1.0_RKIND)
             if (d &lt; nearest_distance) then
                nearest_edge = iEdge
                nearest_distance = d
@@ -4500,7 +4644,7 @@
    end function nearest_edge
 
 
-   real function vertical_interp(target_z, nz, zf, order, extrap, surface_val, sealev_val)
+   real (kind=RKIND) function vertical_interp(target_z, nz, zf, order, extrap, surface_val, sealev_val)
 
       implicit none
 
@@ -4602,4 +4746,179 @@
 
    end subroutine init_atm_check_read_error
 
+
+!----------------------------------------------------------------------------------------------------------
+
+   real (kind=RKIND) function sphere_distance(lat1, lon1, lat2, lon2, radius)
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+   ! Compute the great-circle distance between (lat1, lon1) and (lat2, lon2) on a
+   !   sphere with given radius.
+   !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+      implicit none
+
+      real (kind=RKIND), intent(in) :: lat1, lon1, lat2, lon2, radius
+
+      real (kind=RKIND) :: arg1
+
+      arg1 = sqrt( sin(0.5*(lat2-lat1))**2 +  &amp;
+                   cos(lat1)*cos(lat2)*sin(0.5*(lon2-lon1))**2 )
+      sphere_distance = 2.*radius*asin(arg1)
+
+   end function sphere_distance
+
+!--------------------------------------------------------------------
+
+   real (kind=RKIND) function env_qv( z, temperature, pressure, rh_max )
+
+      implicit none
+      real (kind=RKIND) :: z, temperature, pressure, ztr, es, qvs, p0, rh_max
+
+      p0 = 100000.
+
+!      ztr = 5000.
+!
+!      if(z .gt. ztr) then
+!         env_qv = 0.
+!      else
+!         if(z.lt.2000.) then
+!            env_qv = .5
+!         else
+!            env_qv = .5*(1.-(z-2000.)/(ztr-2000.))
+!         end if
+!      end if
+
+       if (pressure .lt. 50000. ) then
+           env_qv = 0.0
+       else
+           env_qv = (1.-((p0-pressure)/50000.)**1.25)
+       end if
+
+       env_qv = min(rh_max,env_qv)
+
+! env_qv is the relative humidity, turn it into mixing ratio
+       if (temperature .gt. 273.15) then
+           es  = 1000.*0.6112*exp(17.67*(temperature-273.15)/(temperature-29.65))
+       else
+           es  = 1000.*0.6112*exp(21.8745584*(temperature-273.16)/(temperature-7.66))
+       end if
+       qvs = (287.04/461.6)*es/(pressure-es)
+
+       ! qvs =  380.*exp(17.27*(temperature-273.)/(temperature-36.))/pressure
+
+        env_qv = env_qv*qvs
+
+   end function env_qv
+
+
+   subroutine physics_idealized_init(mesh, fg)
+   
+      implicit none
+      
+      !input and output arguments:
+      type(mesh_type),intent(inout):: mesh
+      type (fg_type), intent(inout) :: fg
+      
+      !local variables:
+      integer:: iCell,iMonth,iSoil
+      
+      !---------------------------------------------------------------------------------------------
+      
+      !initialization of surface input variables that are not needed if we run the current set of
+      !idealized test cases:
+      
+      
+      do iCell = 1, mesh % nCells
+      
+         !terrain,soil type, and vegetation:
+         mesh % ter      % array(iCell) = 0.
+         fg % xice       % array(iCell) = 0.
+         mesh % landmask % array(iCell) = 0
+         mesh % lu_index  % array(iCell) = 0
+         mesh % soilcat_top % array(iCell) = 0
+         mesh % shdmin   % array(iCell) = 0.
+         mesh % shdmax   % array(iCell) = 0.
+         fg % vegfra   % array(iCell) = 0.
+         fg % sfc_albbck % array(iCell) = 0.
+         fg % xland % array(iCell) = 0.
+         fg % seaice % array(iCell) = 0.
+      
+         !snow coverage:
+         fg % snow     % array(iCell) = 0.
+         fg % snowc    % array(iCell) = 0.
+         mesh % snoalb % array(iCell) = 0.08
+         fg % snowh % array(iCell) = 0.
+      
+         !surface and sea-surface temperatures:
+         fg % skintemp % array(iCell) = 288.0
+         fg % sst      % array(iCell) = 288.0
+      
+         !soil layers:
+         fg % tmn % array(iCell) = 288.0
+         do iSoil = 1, mesh % nSoilLevels
+            fg % tslb % array(iSoil,iCell)   = 288.0
+            fg % smcrel % array(iSoil,iCell) =   0.0
+            fg % sh2o   % array(iSoil,iCell) =   0.0
+            fg % smois  % array(iSoil,iCell) =   0.0
+            fg % dzs    % array(iSoil,iCell) =   0.0
+         enddo
+      
+         !monthly climatological surface albedo and greeness fraction:
+         do iMonth = 1, mesh % nMonths
+            mesh % albedo12m % array(iMonth,iCell) = 0.08
+            mesh % greenfrac % array(iMonth,iCell) = 0.
+         enddo
+      
+      enddo
+   
+   end subroutine physics_idealized_init
+   
+   
+   subroutine decouple_variables(grid, state, diag)
+
+      implicit none
+
+      type (mesh_type), intent(in) :: grid
+      type (state_type), intent(inout) :: state
+      type (diag_type), intent(inout) :: diag
+
+      integer :: iCell, iEdge, k
+
+      integer, dimension(:,:), pointer :: cellsOnEdge
+      real (kind=RKIND), dimension(:), pointer :: fzm, fzp, rdzw
+      real (kind=RKIND), dimension(:,:), pointer :: zz, pp, ppb
+      real (kind=RKIND), dimension(:,:,:), pointer :: scalars
+
+      cellsOnEdge =&gt; grid % CellsOnEdge % array
+      fzp =&gt; grid % fzm % array
+      fzp =&gt; grid % fzp % array
+      rdzw =&gt; grid % rdzw % array
+      zz =&gt; grid % zz % array
+
+      pp  =&gt; diag % pressure_p % array
+      ppb =&gt; diag % pressure_base % array
+
+      scalars =&gt; state % scalars % array

+     
+      ! Compute surface pressure
+      do iCell=1,grid%nCells
+         diag % surface_pressure % array(iCell) = 0.5*gravity/rdzw(1)                                        &amp;
+                                                * (1.25* state % rho_zz % array(1,iCell) * (1. + scalars(state % index_qv, 1, iCell))  &amp;
+                                                -  0.25* state % rho_zz % array(2,iCell) * (1. + scalars(state % index_qv, 2, iCell)))
+         diag % surface_pressure % array(iCell) = diag % surface_pressure % array(iCell) + pp(1,iCell) + ppb(1,iCell)
+      end do
+
+
+      ! Compute rho and theta from rho_zz and theta_m
+      do iCell=1,grid%nCells
+         do k=1,grid%nVertLevels
+            diag % rho % array(k,iCell) = state % rho_zz % array(k,iCell) * zz(k,iCell)
+            diag % theta % array(k,iCell) = state % theta_m % array(k,iCell) / (1.0 + 1.61 * scalars(state % index_qv,k,iCell))
+         end do
+      end do
+
+   end subroutine decouple_variables
+
+
 end module init_atm_test_cases

Modified: branches/ocean_projects/ale_vert_coord/src/core_nhyd_atmos/Registry
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_nhyd_atmos/Registry        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_nhyd_atmos/Registry        2012-01-17 22:23:38 UTC (rev 1385)
@@ -329,6 +329,7 @@
 #... NAMELIST VARIABLES ADDED FOR PHYSICS CONFIGURATION:
 namelist logical   physics  config_frac_seaice          false
 namelist logical   physics  config_sfc_albedo           false
+namelist logical   physics  config_sfc_snowalbedo       false
 namelist logical   physics  config_sst_update           false
 namelist logical   physics  config_sstdiurn_update      false
 namelist logical   physics  config_deepsoiltemp_update  false

Modified: branches/ocean_projects/ale_vert_coord/src/core_nhyd_atmos/mpas_atm_time_integration.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_nhyd_atmos/mpas_atm_time_integration.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_nhyd_atmos/mpas_atm_time_integration.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -145,10 +145,6 @@
          call mpas_dmpar_exch_halo_field2d_real(domain % dminfo, block % diag % rtheta_p % array(:,:), &amp;
                                           block % mesh % nVertLevels, block % mesh % nCells, &amp;
                                           block % parinfo % cellsToSend, block % parinfo % cellsToRecv)
-!surface_pressure
-         call mpas_dmpar_exch_halo_field1d_real(domain % dminfo, block % diag % surface_pressure % array(:), &amp;
-                                          block % mesh % nCells, &amp;
-                                          block % parinfo % cellsToSend, block % parinfo % cellsToRecv)
          block =&gt; block % next
       end do
 
@@ -192,7 +188,7 @@
         block =&gt; domain % blocklist
         do while (associated(block))
            call physics_addtend( domain % dminfo , block % parinfo % cellsToSend, block % parinfo % cellsToRecv, &amp;
-                        block % mesh , block % state % time_levs(2) % state, block % diag, block % tend, &amp;
+                        block % mesh , block % state % time_levs(1) % state, block % diag, block % tend, &amp;
                         block % tend_physics , block % state % time_levs(2) % state % rho_zz % array(:,:), &amp;
                         block % diag % rho_edge % array(:,:) )
            block =&gt; block % next
@@ -1010,7 +1006,7 @@
                                                     rho_pp, rho_zz, rho_base, ruAvg, ru_save, ru_p, u, ru, &amp;
                                                     exner, exner_base, rtheta_base, pressure_p,         &amp;
                                                     zz, theta_m, pressure_b, qvapor
-      real (kind=RKIND), dimension(:), pointer :: fzm, fzp, dvEdge, AreaCell, rdzw, surface_pressure
+      real (kind=RKIND), dimension(:), pointer :: fzm, fzp, dvEdge, AreaCell
       real (kind=RKIND), dimension(:,:,:), pointer :: zb, zb3 
       integer, dimension(:,:), pointer :: cellsOnEdge
 
@@ -1053,9 +1049,7 @@
 
        pressure_p =&gt; diag % pressure_p % array
        pressure_b =&gt; diag % pressure_base % array
-       surface_pressure =&gt; diag % surface_pressure % array
 
-       rdzw =&gt; grid % rdzw % array
        zz =&gt; grid % zz % array
        zb =&gt; grid % zb % array
        zb3 =&gt; grid % zb3 % array
@@ -1130,13 +1124,6 @@
                                                           * (exner(k,iCell)-exner_base(k,iCell)))
           end do
 
-         !calculation of the surface pressure:
-         surface_pressure(iCell) = 0.5*gravity/rdzw(1)                                              &amp;
-                       * (1.25* rho_zz(1,iCell) * (1. + qvapor(1,iCell))  &amp;
-                       -  0.25* rho_zz(2,iCell) * (1. + qvapor(2,iCell)))
-         surface_pressure(iCell) = surface_pressure(iCell) + pressure_p(1,iCell) + pressure_b(1,iCell)
-
-
       end do
 
       ! recover time-averaged ruAvg on all edges of owned cells (for upcoming scalar transport).  
@@ -1872,7 +1859,8 @@
       real (kind=RKIND) :: r_earth
       real (kind=RKIND), dimension(:,:), pointer :: ur_cell, vr_cell
 
-      real (kind=RKIND), parameter :: c_s = 0.25
+      real (kind=RKIND), parameter :: c_s = 0.125
+!     real (kind=RKIND), parameter :: c_s = 0.25
       real (kind=RKIND), dimension( grid % nVertLevels ) :: d_diag, d_off_diag, flux_arr
       real (kind=RKIND), dimension(:,:), pointer :: defc_a, defc_b
       logical :: delsq_horiz_mixing, newpx


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean
___________________________________________________________________
Modified: svn:mergeinfo
   - /branches/cam_mpas_nh/src/core_ocean:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean:754-986
/branches/ocean_projects/performance/src/core_ocean:1120-1357
/branches/ocean_projects/rayleigh/src/core_ocean:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean:704-745
/branches/source_renaming/src/core_ocean:1082-1113
/branches/time_manager/src/core_ocean:924-962
/trunk/mpas/src/core_ocean:1224-1336
   + /branches/cam_mpas_nh/src/core_ocean:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean:754-986
/branches/ocean_projects/performance/src/core_ocean:1120-1357
/branches/ocean_projects/rayleigh/src/core_ocean:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean:704-745
/branches/source_renaming/src/core_ocean:1082-1113
/branches/time_manager/src/core_ocean:924-962
/trunk/mpas/src/core_ocean:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/Registry
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/Registry:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/Registry:754-986
/branches/ocean_projects/performance/src/core_ocean/Registry:1120-1372
/branches/ocean_projects/rayleigh/src/core_ocean/Registry:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/Registry:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/Registry:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/Registry:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/Registry:704-745
/branches/source_renaming/src/core_ocean/Registry:1082-1113
/branches/time_manager/src/core_ocean/Registry:924-962
/trunk/mpas/src/core_ocean/Registry:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_advection.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_advection.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_advection.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_advection.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_advection.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_advection.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_advection.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_advection.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_advection.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_advection.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_advection.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_advection.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_equation_of_state.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_equation_of_state.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_equation_of_state.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_equation_of_state.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_equation_of_state.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_equation_of_state.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_equation_of_state.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_equation_of_state.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_equation_of_state.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_equation_of_state.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_equation_of_state.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_equation_of_state.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_equation_of_state_jm.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_equation_of_state_jm.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_equation_of_state_jm.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_equation_of_state_jm.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_equation_of_state_jm.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_equation_of_state_jm.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_equation_of_state_jm.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_equation_of_state_jm.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_equation_of_state_jm.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_equation_of_state_jm.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_equation_of_state_jm.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_equation_of_state_jm.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_equation_of_state_linear.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_equation_of_state_linear.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_equation_of_state_linear.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_equation_of_state_linear.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_equation_of_state_linear.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_equation_of_state_linear.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_equation_of_state_linear.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_equation_of_state_linear.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_equation_of_state_linear.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_equation_of_state_linear.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_equation_of_state_linear.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_equation_of_state_linear.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_global_diagnostics.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_global_diagnostics.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_global_diagnostics.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_global_diagnostics.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_global_diagnostics.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_global_diagnostics.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_global_diagnostics.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_global_diagnostics.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_global_diagnostics.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_global_diagnostics.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_global_diagnostics.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_global_diagnostics.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_mpas_core.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_mpas_core.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_mpas_core.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_mpas_core.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_mpas_core.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_mpas_core.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_mpas_core.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_mpas_core.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_mpas_core.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_mpas_core.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_mpas_core.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_mpas_core.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_restoring.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_restoring.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_restoring.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_restoring.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_restoring.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_restoring.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_restoring.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_restoring.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_restoring.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_restoring.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_restoring.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_restoring.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tendency.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tendency.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tendency.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -788,6 +788,7 @@
       ! For zlevel, calculate in-situ density
       if (config_vert_grid_type.ne.'isopycnal') then
          call mpas_timer_start(&quot;equation of state&quot;, .false., diagEOSTimer)
+         call mpas_timer_start(&quot;equation of state&quot;, .false., diagEOSTimer)
          call ocn_equation_of_state_rho(s, grid, 0, 'relative', err)
       ! mrp 110324 In order to visualize rhoDisplaced, include the following
          call ocn_equation_of_state_rho(s, grid, 1, 'relative', err)


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tendency.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tendency.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tendency.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tendency.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tendency.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tendency.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tendency.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tendency.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tendency.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tendency.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tendency.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tendency.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_thick_hadv.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_thick_hadv.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_thick_hadv.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_thick_hadv.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_thick_hadv.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_thick_hadv.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_thick_hadv.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_thick_hadv.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_thick_hadv.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_thick_hadv.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_thick_hadv.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_thick_hadv.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_thick_vadv.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_thick_vadv.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_thick_vadv.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_thick_vadv.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_thick_vadv.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_thick_vadv.F:1134-1138
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_thick_vadv.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_thick_vadv.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_thick_vadv.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_thick_vadv.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_thick_vadv.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_thick_vadv.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_time_integration.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_time_integration.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_time_integration.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_time_integration.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_time_integration.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_time_integration.F:1134-1138
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_time_integration.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_time_integration.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_time_integration.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_time_integration.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_time_integration.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_time_integration.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_time_integration_rk4.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_time_integration_rk4.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_time_integration_rk4.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_time_integration_rk4.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_time_integration_rk4.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_time_integration_rk4.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_time_integration_rk4.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_time_integration_rk4.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_time_integration_rk4.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_time_integration_rk4.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_time_integration_rk4.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_time_integration_rk4.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_time_integration_split.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_time_integration_split.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_time_integration_split.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_time_integration_split.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_time_integration_split.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_time_integration_split.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_time_integration_split.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_time_integration_split.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_time_integration_split.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_time_integration_split.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_time_integration_split.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_time_integration_split.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_hadv.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_hadv.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_hadv.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_hadv.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_hadv.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_hadv.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_hadv.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_hadv.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_hadv.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_hadv.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_hadv.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_hadv.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_hadv2.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_hadv2.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_hadv2.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_hadv2.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_hadv2.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_hadv2.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_hadv2.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_hadv2.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_hadv2.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_hadv2.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_hadv2.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_hadv2.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_hadv3.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_hadv3.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_hadv3.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_hadv3.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_hadv3.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_hadv3.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_hadv3.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_hadv3.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_hadv3.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_hadv3.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_hadv3.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_hadv3.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_hadv4.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_hadv4.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_hadv4.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_hadv4.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_hadv4.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_hadv4.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_hadv4.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_hadv4.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_hadv4.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_hadv4.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_hadv4.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_hadv4.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_hmix.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_hmix.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_hmix.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_hmix.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_hmix.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_hmix.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_hmix.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_hmix.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_hmix.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_hmix.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_hmix.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_hmix.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_hmix_del2.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_hmix_del2.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_hmix_del2.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_hmix_del2.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_hmix_del2.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_hmix_del2.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_hmix_del2.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_hmix_del2.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_hmix_del2.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_hmix_del2.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_hmix_del2.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_hmix_del2.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_hmix_del4.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_hmix_del4.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_hmix_del4.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -44,7 +44,7 @@
    !
    !--------------------------------------------------------------------
 
-   integer :: del4On
+   logical :: del4On
 
    real (kind=RKIND) :: eddyDiff4
 
@@ -110,12 +110,10 @@
       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
+      integer, dimension(:,:), pointer :: edgeMask, cellsOnEdge
 
-      real (kind=RKIND) :: invAreaCell1, invAreaCell2, r, tracer_turb_flux, flux
+      real (kind=RKIND) :: invAreaCell1, invAreaCell2, tracer_turb_flux, flux, invdcEdge, r_tmp1, r_tmp2
 
       real (kind=RKIND), dimension(:,:,:), allocatable :: delsq_tracer
 
@@ -132,7 +130,7 @@
 
       err = 0
 
-      if (del4On == 0) return
+      if (.not.del4On) return
 
       nEdges = grid % nEdges
       nCells = grid % nCells
@@ -141,7 +139,6 @@
 
       maxLevelEdgeTop =&gt; grid % maxLevelEdgeTop % array
       maxLevelCell =&gt; grid % maxLevelCell % array
-      boundaryEdge =&gt; grid % boundaryEdge % array
       cellsOnEdge =&gt; grid % cellsOnEdge % array
 
       dcEdge =&gt; grid % dcEdge % array
@@ -149,67 +146,61 @@
       areaCell =&gt; grid % areaCell % array
       meshScalingDel4 =&gt; grid % meshScalingDel4 % array
 
-      allocate(boundaryMask(nVertLevels, nEdges+1))
-      boundaryMask = 1.0
-      where(boundaryEdge.eq.1) boundaryMask=0.0
+      edgeMask =&gt; grid % edgeMask % array
 
       allocate(delsq_tracer(num_tracers,nVertLevels, nCells+1))
 
-      delsq_tracer(:,:,:) = 0.
+      delsq_tracer(:,:,:) = 0.0
 
       ! first del2: div(h </font>
<font color="red">abla \phi) at cell center
       do iEdge=1,nEdges
          cell1 = cellsOnEdge(1,iEdge)
          cell2 = cellsOnEdge(2,iEdge)
 
+         invdcEdge = 1.0 / dcEdge(iEdge)
+
+         invAreaCell1 = 1.0 / areaCell(cell1)
+         invAreaCell2 = 1.0 / areaCell(cell2)
+
          do k=1,maxLevelEdgeTop(iEdge)
-           do iTracer=1,num_tracers
-              delsq_tracer(iTracer,k,cell1) = delsq_tracer(iTracer,k,cell1) &amp;
-                 + dvEdge(iEdge)*h_edge(k,iEdge) &amp;
-                   *(tracers(iTracer,k,cell2) - tracers(iTracer,k,cell1)) &amp;
-                   /dcEdge(iEdge) * boundaryMask(k,iEdge)
-              delsq_tracer(iTracer,k,cell2) = delsq_tracer(iTracer,k,cell2) &amp;
-                 - dvEdge(iEdge)*h_edge(k,iEdge) &amp;
-                 *(tracers(iTracer,k,cell2) - tracers(iTracer,k,cell1)) &amp;
-                 /dcEdge(iEdge) * boundaryMask(k,iEdge)
+           do iTracer=1,num_tracers * edgeMask(k, iEdge)
+
+              r_tmp1 = dvEdge(iEdge) * h_edge(k,iEdge) * invdcEdge
+
+              r_tmp2 = r_tmp1 * tracers(iTracer,k,cell2)
+              r_tmp1 = r_tmp1 * tracers(iTracer,k,cell1)
+
+              delsq_tracer(iTracer,k,cell1) = delsq_tracer(iTracer,k,cell1) + (r_tmp2 - r_tmp1) * invAreaCell1
+              delsq_tracer(iTracer,k,cell2) = delsq_tracer(iTracer,k,cell2) - (r_tmp2 - r_tmp1) * invAreaCell2
            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)
 
+         invdcEdge = 1.0 / dcEdge(iEdge)
+
          do k=1,maxLevelEdgeTop(iEdge)
-            do iTracer=1,num_tracers
+            do iTracer=1,num_tracers * edgeMask(k,iEdge)
                tracer_turb_flux = meshScalingDel4(iEdge) * eddyDiff4 &amp;
-                  *(  delsq_tracer(iTracer,k,cell2)  &amp;
-                    - delsq_tracer(iTracer,k,cell1))/dcEdge(iEdge)
+                  * (delsq_tracer(iTracer,k,cell2) - delsq_tracer(iTracer,k,cell1)) &amp;
+                  * invdcEdge
+
                flux = dvEdge (iEdge) * tracer_turb_flux
 
-               tend(iTracer,k,cell1) = tend(iTracer,k,cell1) &amp; 
-                  - flux * invAreaCell1 * boundaryMask(k,iEdge)
-               tend(iTracer,k,cell2) = tend(iTracer,k,cell2) &amp;
-                  + flux * invAreaCell2 * boundaryMask(k,iEdge)
-
+               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(delsq_tracer)
-      deallocate(boundaryMask)
    !--------------------------------------------------------------------
 
    end subroutine ocn_tracer_hmix_del4_tend!}}}
@@ -241,10 +232,10 @@
       integer, intent(out) :: err !&lt; Output: error flag
 
       err = 0
-      del4on = 0
+      del4on = .false.
 
       if ( config_h_tracer_eddy_diff4 &gt; 0.0 ) then
-          del4On = 1
+          del4On = .true.
           eddyDiff4 = config_h_tracer_eddy_diff4
       endif
 


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_hmix_del4.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_hmix_del4.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_hmix_del4.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_hmix_del4.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_hmix_del4.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_hmix_del4.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_hmix_del4.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_hmix_del4.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_hmix_del4.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_hmix_del4.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_hmix_del4.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_hmix_del4.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_vadv.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_vadv.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_vadv.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_vadv.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_vadv.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_vadv.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_vadv.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_vadv.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_vadv.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_vadv.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_vadv.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_spline.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_spline.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_spline.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -133,6 +133,7 @@
       call mpas_timer_start(&quot;spline 3&quot;, .false., spline3_timer)
       call ocn_tracer_vadv_spline3_tend(grid, h, wTop, tracers, tend, err2)
       call mpas_timer_stop(&quot;spline 3&quot;, spline3_timer)
+      call mpas_timer_stop(&quot;spline 3&quot;, spline3_timer)
 
       err = ior(err1, err2)
 


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_spline.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_vadv_spline.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_vadv_spline.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_vadv_spline.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_vadv_spline.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_vadv_spline.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_vadv_spline.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_vadv_spline.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_vadv_spline.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_vadv_spline.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_vadv_spline.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_vadv_spline.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_spline2.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_vadv_spline2.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_vadv_spline2.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_vadv_spline2.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_vadv_spline2.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_vadv_spline2.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_vadv_spline2.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_vadv_spline2.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_vadv_spline2.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_vadv_spline2.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_vadv_spline2.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_vadv_spline2.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_spline3.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_vadv_spline3.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_vadv_spline3.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_vadv_spline3.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_vadv_spline3.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_vadv_spline3.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_vadv_spline3.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_vadv_spline3.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_vadv_spline3.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_vadv_spline3.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_vadv_spline3.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_vadv_spline3.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -49,9 +49,8 @@
 
    type (timer_node), pointer :: stencil2_timer, stencil3_timer, stencil4_timer
 
-   integer :: stencilOn
+   logical :: stencilOn
 
-
 !***********************************************************************
 
 contains
@@ -127,7 +126,7 @@
 
       err = 0
 
-      if(stencilOn == 0) return
+      if(.not.stencilOn) return
 
       call mpas_timer_start(&quot;stencil 2&quot;, .false., stencil2_timer)
       call ocn_tracer_vadv_stencil2_tend(grid, wTop, tracers, tend, err1)


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_vadv_stencil.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -121,7 +121,7 @@
 
       err = 0
 
-      if(stencil2On == 0) return
+      if(.not.stencil2On) return
 
       nCells = grid % nCells
       nCellsSolve = grid % nCellsSolve


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_vadv_stencil2.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -123,7 +123,7 @@
 
       err = 0
 
-      if(stencil3On == 0) return
+      if(.not.stencil3On) return
 
       nCells = grid % nCells
       nCellsSolve = grid % nCellsSolve


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_vadv_stencil3.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -123,7 +123,7 @@
 
       err = 0
 
-      if(Stencil4On == 0) return
+      if(.not.Stencil4On) return
 
       nCells = grid % nCells
       nCellsSolve = grid % nCellsSolve


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_tracer_vadv_stencil4.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_coriolis.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vel_coriolis.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vel_coriolis.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vel_coriolis.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vel_coriolis.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vel_coriolis.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vel_coriolis.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vel_coriolis.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vel_coriolis.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vel_coriolis.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vel_coriolis.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vel_coriolis.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_forcing.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vel_forcing.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vel_forcing.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vel_forcing.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vel_forcing.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vel_forcing.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vel_forcing.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vel_forcing.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vel_forcing.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vel_forcing.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vel_forcing.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vel_forcing.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_forcing_bottomdrag.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vel_forcing_bottomdrag.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vel_forcing_bottomdrag.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vel_forcing_bottomdrag.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vel_forcing_bottomdrag.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vel_forcing_bottomdrag.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vel_forcing_bottomdrag.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vel_forcing_bottomdrag.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vel_forcing_bottomdrag.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vel_forcing_bottomdrag.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vel_forcing_bottomdrag.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vel_forcing_bottomdrag.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_forcing_windstress.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_forcing_windstress.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_forcing_windstress.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -169,7 +169,7 @@
 
       integer, intent(out) :: err !&lt; Output: error flag
 
-      windStressOn = 1
+      windStressOn = .true.
       rho_ref = 1000.0
 
       err = 0


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_forcing_windstress.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vel_forcing_windstress.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vel_forcing_windstress.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vel_forcing_windstress.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vel_forcing_windstress.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vel_forcing_windstress.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vel_forcing_windstress.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vel_forcing_windstress.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vel_forcing_windstress.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vel_forcing_windstress.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vel_forcing_windstress.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vel_forcing_windstress.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_hmix.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vel_hmix.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vel_hmix.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vel_hmix.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vel_hmix.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vel_hmix.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vel_hmix.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vel_hmix.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vel_hmix.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vel_hmix.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vel_hmix.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vel_hmix.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_hmix_del2.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_hmix_del2.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_hmix_del2.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -42,7 +42,7 @@
    !
    !--------------------------------------------------------------------
 
-   integer ::  hmixDel2On  !&lt; integer flag to determine whether del2 chosen
+   logical ::  hmixDel2On  !&lt; integer flag to determine whether del2 chosen
 
    real (kind=RKIND) :: &amp;
       eddyVisc2,        &amp;!&lt; base eddy diffusivity for Laplacian


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_hmix_del2.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vel_hmix_del2.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vel_hmix_del2.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vel_hmix_del2.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vel_hmix_del2.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vel_hmix_del2.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vel_hmix_del2.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vel_hmix_del2.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vel_hmix_del2.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vel_hmix_del2.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vel_hmix_del2.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vel_hmix_del2.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_hmix_del4.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_hmix_del4.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_hmix_del4.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -42,8 +42,7 @@
    !
    !--------------------------------------------------------------------
 
-   integer :: &amp;
-      hmixDel4On       !&lt; local flag to determine whether del4 chosen
+   logical :: hmixDel4On       !&lt; local flag to determine whether del4 chosen
 
    real (kind=RKIND) :: &amp;
       eddyVisc4,        &amp;!&lt; base eddy diffusivity for biharmonic


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_hmix_del4.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vel_hmix_del4.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vel_hmix_del4.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vel_hmix_del4.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vel_hmix_del4.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vel_hmix_del4.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vel_hmix_del4.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vel_hmix_del4.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vel_hmix_del4.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vel_hmix_del4.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vel_hmix_del4.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vel_hmix_del4.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_pressure_grad.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_pressure_grad.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_pressure_grad.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -156,7 +156,6 @@
 
       endif
 
-
    !--------------------------------------------------------------------
 
    end subroutine ocn_vel_pressure_grad_tend!}}}


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_pressure_grad.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vel_pressure_grad.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vel_pressure_grad.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vel_pressure_grad.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vel_pressure_grad.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vel_pressure_grad.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vel_pressure_grad.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vel_pressure_grad.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vel_pressure_grad.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vel_pressure_grad.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vel_pressure_grad.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vel_pressure_grad.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vel_vadv.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vel_vadv.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vel_vadv.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vel_vadv.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vel_vadv.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vel_vadv.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vel_vadv.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vel_vadv.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vel_vadv.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vel_vadv.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vel_vadv.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vel_vadv.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -592,8 +592,8 @@
       implicitOn = 0
 
       if(config_implicit_vertical_mix) then
-          explicitOn = 0
-          implicitOn = 1
+          explicitOn = .false.
+          implicitOn = .true.
       end if
 
       call ocn_vmix_coefs_const_init(err1)


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vmix.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vmix.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vmix.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vmix.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vmix.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vmix.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vmix.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vmix.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vmix.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vmix.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vmix.F:1224-1384


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix_coefs_const.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vmix_coefs_const.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vmix_coefs_const.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vmix_coefs_const.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vmix_coefs_const.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vmix_coefs_const.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vmix_coefs_const.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vmix_coefs_const.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vmix_coefs_const.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vmix_coefs_const.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vmix_coefs_const.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vmix_coefs_const.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix_coefs_rich.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix_coefs_rich.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix_coefs_rich.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -444,7 +444,7 @@
 
       err = 0
 
-      if((richViscOn == 0) .and. (richDiffOn == 0)) return
+      if((.not.richViscOn) .and. (.not.richDiffOn)) return
 
       nVertLevels = grid % nVertLevels
       nCells = grid % nCells


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix_coefs_rich.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vmix_coefs_rich.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vmix_coefs_rich.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vmix_coefs_rich.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vmix_coefs_rich.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vmix_coefs_rich.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vmix_coefs_rich.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vmix_coefs_rich.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vmix_coefs_rich.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vmix_coefs_rich.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vmix_coefs_rich.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vmix_coefs_rich.F:1224-1384

Modified: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F
===================================================================
--- branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F        2012-01-17 21:44:00 UTC (rev 1384)
+++ branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F        2012-01-17 22:23:38 UTC (rev 1385)
@@ -308,7 +308,7 @@
           tanhDiffOn = 1
       endif
 
-      if((tanhViscOn == 1) .or. (tanhDiffOn == 1)) then
+      if((tanhViscOn) .or. (tanhDiffOn)) then
          if (config_vert_grid_type.ne.'zlevel') then
             write(0,*) 'Abort: config_vert_diff_type.eq.tanh may only', &amp;
                        ' use config_vert_grid_type of zlevel at this time'


Property changes on: branches/ocean_projects/ale_vert_coord/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F
___________________________________________________________________
Added: svn:mergeinfo
   + /branches/cam_mpas_nh/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F:1260-1270
/branches/ocean_projects/imp_vert_mix_mrp/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F:754-986
/branches/ocean_projects/performance/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F:1120-1371
/branches/ocean_projects/rayleigh/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F:1298-1311
/branches/ocean_projects/split_explicit_mrp/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F:1133-1175
/branches/ocean_projects/split_explicit_timestepping/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F:1044-1097
/branches/ocean_projects/time_averaging/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F:1271-1305
/branches/ocean_projects/vert_adv_mrp/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F:704-745
/branches/source_renaming/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F:1082-1113
/branches/time_manager/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F:924-962
/trunk/mpas/src/core_ocean/mpas_ocn_vmix_coefs_tanh.F:1224-1384

</font>
</pre>