Merge branch 'gmtb/develop' of https://github.com/NCAR/ccpp-physics i…

…nto FA-HWRF-V4_0a

physics/GFS_PBL_generic.F90

@@ -100,18 +100,19 @@ subroutine GFS_PBL_generic_pre_run (im, levs, nvdiff, ntrac,
           enddo
         elseif (imp_physics == imp_physics_thompson) then
   ! Thompson
-          ! DH* Thompson ntrw and ntsw?
           if(ltaerosol) then
             do k=1,levs
               do i=1,im
-                vdftra(i,k,1) = qgrs(i,k,ntqv)
-                vdftra(i,k,2) = qgrs(i,k,ntcw)
-                vdftra(i,k,3) = qgrs(i,k,ntiw)
-                vdftra(i,k,4) = qgrs(i,k,ntlnc)
-                vdftra(i,k,5) = qgrs(i,k,ntinc)
-                vdftra(i,k,6) = qgrs(i,k,ntoz)
-                vdftra(i,k,7) = qgrs(i,k,ntwa)
-                vdftra(i,k,8) = qgrs(i,k,ntia)
+                vdftra(i,k,1)  = qgrs(i,k,ntqv)
+                vdftra(i,k,2)  = qgrs(i,k,ntcw)
+                vdftra(i,k,3)  = qgrs(i,k,ntiw)
+                vdftra(i,k,4)  = qgrs(i,k,ntrw)
+                vdftra(i,k,5)  = qgrs(i,k,ntsw)
+                vdftra(i,k,6)  = qgrs(i,k,ntlnc)
+                vdftra(i,k,7)  = qgrs(i,k,ntinc)
+                vdftra(i,k,8)  = qgrs(i,k,ntoz)
+                vdftra(i,k,9)  = qgrs(i,k,ntwa)
+                vdftra(i,k,10) = qgrs(i,k,ntia)
               enddo
             enddo
           else
@@ -120,8 +121,10 @@ subroutine GFS_PBL_generic_pre_run (im, levs, nvdiff, ntrac,
                 vdftra(i,k,1) = qgrs(i,k,ntqv)
                 vdftra(i,k,2) = qgrs(i,k,ntcw)
                 vdftra(i,k,3) = qgrs(i,k,ntiw)
-                vdftra(i,k,4) = qgrs(i,k,ntinc)
-                vdftra(i,k,5) = qgrs(i,k,ntoz)
+                vdftra(i,k,4) = qgrs(i,k,ntrw)
+                vdftra(i,k,5) = qgrs(i,k,ntsw)
+                vdftra(i,k,6) = qgrs(i,k,ntinc)
+                vdftra(i,k,7) = qgrs(i,k,ntoz)
               enddo
             enddo
           endif
@@ -362,18 +365,19 @@ subroutine GFS_PBL_generic_post_run (im, levs, nvdiff, ntrac,
           enddo
         elseif (imp_physics == imp_physics_thompson) then
   ! Thompson
-          ! DH* - Thompson ntrw, ntsw?
           if(ltaerosol) then
             do k=1,levs
               do i=1,im
                 dqdt(i,k,ntqv)  = dvdftra(i,k,1)
                 dqdt(i,k,ntcw)  = dvdftra(i,k,2)
                 dqdt(i,k,ntiw)  = dvdftra(i,k,3)
-                dqdt(i,k,ntlnc) = dvdftra(i,k,4)
-                dqdt(i,k,ntinc) = dvdftra(i,k,5)
-                dqdt(i,k,ntoz)  = dvdftra(i,k,6)
-                dqdt(i,k,ntwa)  = dvdftra(i,k,7)
-                dqdt(i,k,ntia)  = dvdftra(i,k,8)
+                dqdt(i,k,ntrw)  = dvdftra(i,k,4)
+                dqdt(i,k,ntsw)  = dvdftra(i,k,5)
+                dqdt(i,k,ntlnc) = dvdftra(i,k,6)
+                dqdt(i,k,ntinc) = dvdftra(i,k,7)
+                dqdt(i,k,ntoz)  = dvdftra(i,k,8)
+                dqdt(i,k,ntwa)  = dvdftra(i,k,9)
+                dqdt(i,k,ntia)  = dvdftra(i,k,10)
               enddo
             enddo
           else
@@ -382,8 +386,10 @@ subroutine GFS_PBL_generic_post_run (im, levs, nvdiff, ntrac,
                 dqdt(i,k,ntqv)  = dvdftra(i,k,1)
                 dqdt(i,k,ntcw)  = dvdftra(i,k,2)
                 dqdt(i,k,ntiw)  = dvdftra(i,k,3)
-                dqdt(i,k,ntinc) = dvdftra(i,k,4)
-                dqdt(i,k,ntoz)  = dvdftra(i,k,5)
+                dqdt(i,k,ntrw)  = dvdftra(i,k,4)
+                dqdt(i,k,ntsw)  = dvdftra(i,k,5)
+                dqdt(i,k,ntinc) = dvdftra(i,k,6)
+                dqdt(i,k,ntoz)  = dvdftra(i,k,7)
               enddo
             enddo
           endif

physics/GFS_phys_time_vary.scm.F90

@@ -367,25 +367,19 @@ subroutine GFS_phys_time_vary_run (Grid, Statein, Model, Tbd, Sfcprop, Cldprop,
         sec_zero = nint(Model%fhzero*con_hr)
         if (sec_zero >= nint(max(Model%fhswr,Model%fhlwr))) then
           if (mod(Model%kdt,Model%nszero) == 1) then
-            do nb = 1,nblks
               call Diag%rad_zero  (Model)
               call Diag%phys_zero (Model)
         !!!!  THIS IS THE POINT AT WHICH DIAG%ZHOUR NEEDS TO BE UPDATED
-            enddo
           endif
         else
           if (mod(Model%kdt,Model%nszero) == 1) then
-            do nb = 1,nblks
               call Diag%phys_zero (Model)
         !!!!  THIS IS THE POINT AT WHICH DIAG%ZHOUR NEEDS TO BE UPDATED
-            enddo
           endif
           kdt_rad = nint(min(Model%fhswr,Model%fhlwr)/Model%dtp)
           if (mod(Model%kdt, kdt_rad) == 1) then
-            do nb = 1,nblks
               call Diag%rad_zero  (Model)
         !!!!  THIS IS THE POINT AT WHICH DIAG%ZHOUR NEEDS TO BE UPDATED
-            enddo
           endif
         endif
 

physics/cu_gf_deep.F90

@@ -3227,9 +3227,9 @@ subroutine cup_up_aa0(aa0,z,zu,dby,gamma_cup,t_cup,       &
         enddo
         do k=kts+1,ktf
           do i=its,itf
-           if(ierr(i).ne.0) exit
-           if(k.lt.kbcon(i)) exit
-           if(k.gt.ktop(i)) exit
+           if(ierr(i).ne.0) cycle
+           if(k.lt.kbcon(i)) cycle
+           if(k.gt.ktop(i)) cycle
            dz=z(i,k)-z(i,k-1)
            da=zu(i,k)*dz*(9.81/(1004.*( &
                   (t_cup(i,k)))))*dby(i,k-1)/ &
@@ -4379,8 +4379,8 @@ subroutine cup_up_aa1bl(aa0,t,tn,q,qo,dtime,  &
         enddo
         do i=its,itf
           do k=kts,kbcon(i)
-            if(ierr(i).ne.0 ) exit
-!           if(k.gt.kbcon(i)) exit
+            if(ierr(i).ne.0 ) cycle
+!           if(k.gt.kbcon(i)) cycle
 
             dz = (z_cup (i,k+1)-z_cup (i,k))*g
             da = dz*(tn(i,k)*(1.+0.608*qo(i,k))-t(i,k)*(1.+0.608*q(i,k)))/dtime

physics/cu_gf_driver.F90

@@ -71,8 +71,6 @@ end subroutine cu_gf_driver_finalize
 !! \section arg_table_cu_gf_driver_run Argument Table
 !! | local_name     | standard_name                                             | long_name                                           | units         | rank | type      |    kind   | intent | optional |
 !! |----------------|-----------------------------------------------------------|-----------------------------------------------------|---------------|------|-----------|-----------|--------|----------|
-!! | tottracer      | number_of_total_tracers                                   | number of total tracers                             | count         |    0 | integer   |           | in     | F        |
-!! | ntrac          | number_of_vertical_diffusion_tracers                      | number of tracers to diffuse vertically             | count         |    0 | integer   |           | in     | F        |
 !! | garea          | cell_area                                                 | grid cell area                                      | m2            |    1 | real      | kind_phys | in     | F        |
 !! | im             | horizontal_loop_extent                                    | horizontal loop extent                              | count         |    0 | integer   |           | in     | F        |
 !! | ix             | horizontal_dimension                                      | horizontal dimension                                | count         |    0 | integer   |           | in     | F        |
@@ -99,7 +97,8 @@ end subroutine cu_gf_driver_finalize
 !! | xland          | sea_land_ice_mask                                         | landmask: sea/land/ice=0/1/2                        | flag          |    1 | integer   |           | in     | F        |
 !! | hfx2           | kinematic_surface_upward_sensible_heat_flux               | kinematic surface upward sensible heat flux         | K m s-1       |    1 | real      | kind_phys | in     | F        |
 !! | qfx2           | kinematic_surface_upward_latent_heat_flux                 | kinematic surface upward latent heat flux           | kg kg-1 m s-1 |    1 | real      | kind_phys | in     | F        |
-!! | clw            | convective_transportable_tracers                          | cloud water and other convective trans. tracers     | kg kg-1       |    3 | real      | kind_phys | inout  | F        |
+!! | cliw           | ice_water_mixing_ratio_convective_transport_tracer             | moist (dry+vapor, no condensates) mixing ratio of ice water in the convectively transported tracer array   | kg kg-1   |    2 | real        | kind_phys | inout  | F        |
+!! | clcw           | cloud_condensed_water_mixing_ratio_convective_transport_tracer | moist (dry+vapor, no condensates) mixing ratio of cloud water in the convectively transported tracer array | kg kg-1   |    2 | real        | kind_phys | inout  | F        |
 !! | pbl            | atmosphere_boundary_layer_thickness                       | PBL thickness                                       | m             |    1 | real      | kind_phys | in     | F        |
 !! | ud_mf          | instantaneous_atmosphere_updraft_convective_mass_flux     | (updraft mass flux) * delt                          | kg m-2        |    2 | real      | kind_phys | out    | F        |
 !! | dd_mf          | instantaneous_atmosphere_downdraft_convective_mass_flux   | (downdraft mass flux) * delt                        | kg m-2        |    2 | real      | kind_phys | out    | F        |
@@ -112,11 +111,12 @@ end subroutine cu_gf_driver_finalize
 !!
 !>\section gen_gf_driver GSD GF Cumulus Scheme General Algorithm
 !> @{
-      subroutine cu_gf_driver_run(tottracer,ntrac,garea,im,ix,km,dt,cactiv, &
-               forcet,forceqv_spechum,phil,raincv,qv_spechum,t,cld1d,       &
-               us,vs,t2di,w,qv2di_spechum,p2di,psuri,                       &
-               hbot,htop,kcnv,xland,hfx2,qfx2,clw,                          &
-               pbl,ud_mf,dd_mf,dt_mf,cnvw_moist,cnvc,imfshalcnv,errmsg,errflg)
+      subroutine cu_gf_driver_run(garea,im,ix,km,dt,cactiv,           &
+               forcet,forceqv_spechum,phil,raincv,qv_spechum,t,cld1d, &
+               us,vs,t2di,w,qv2di_spechum,p2di,psuri,                 &
+               hbot,htop,kcnv,xland,hfx2,qfx2,cliw,clcw,              &
+               pbl,ud_mf,dd_mf,dt_mf,cnvw_moist,cnvc,imfshalcnv,      &
+               errmsg,errflg)
 !-------------------------------------------------------------
       implicit none
       integer, parameter :: maxiens=1
@@ -137,15 +137,15 @@ subroutine cu_gf_driver_run(tottracer,ntrac,garea,im,ix,km,dt,cactiv, &
       integer            :: ishallow_g3 ! depend on imfshalcnv
 !-------------------------------------------------------------
    integer      :: its,ite, jts,jte, kts,kte 
-   integer, intent(in   ) :: im,ix,km,ntrac,tottracer
+   integer, intent(in   ) :: im,ix,km
 
    real(kind=kind_phys),  dimension( ix , km ),     intent(in ) :: forcet,forceqv_spechum,w,phil
    real(kind=kind_phys),  dimension( ix , km ),     intent(inout ) :: t,us,vs
    real(kind=kind_phys),  dimension( ix )   :: rand_mom,rand_vmas
    real(kind=kind_phys),  dimension( ix,4 ) :: rand_clos
    real(kind=kind_phys),  dimension( ix , km, 11 ) :: gdc,gdc2
    real(kind=kind_phys),  dimension( ix , km ),     intent(out ) :: cnvw_moist,cnvc
-   real(kind=kind_phys),  dimension( ix , km,tottracer+2 ), intent(inout ) :: clw
+   real(kind=kind_phys),  dimension( ix , km ), intent(inout ) :: cliw, clcw
 
 !hj change from ix to im
    integer, dimension (im), intent(inout) :: hbot,htop,kcnv
@@ -260,10 +260,13 @@ subroutine cu_gf_driver_run(tottracer,ntrac,garea,im,ix,km,dt,cactiv, &
      rand_mom(:)    = 0.
      rand_vmas(:)   = 0.
      rand_clos(:,:) = 0.
+!
      its=1
      ite=im
+     itf=ite
      jts=1
      jte=1
+     jtf=jte
      kts=1
      kte=km
      ktf=kte-1
@@ -307,9 +310,6 @@ subroutine cu_gf_driver_run(tottracer,ntrac,garea,im,ix,km,dt,cactiv, &
    iend=ite
    tcrit=258.
 
-   itf=ite
-   ktf=kte-1
-   jtf=jte
    ztm=0.
    ztq=0.
    hfm=0.
@@ -818,9 +818,9 @@ subroutine cu_gf_driver_run(tottracer,ntrac,garea,im,ix,km,dt,cactiv, &
                dsubclwm=0.
                dsubclws=0.
                dp=100.*(p2d(i,k)-p2d(i,k+1))
-               if (clw(i,k,2) .gt. -999.0 .and. clw(i,k+1,2) .gt. -999.0 )then
-                  clwtot = clw(i,k,1) + clw(i,k,2)
-                  clwtot1= clw(i,k+1,1) + clw(i,k+1,2)
+               if (clcw(i,k) .gt. -999.0 .and. clcw(i,k+1) .gt. -999.0 )then
+                  clwtot = cliw(i,k) + clcw(i,k)
+                  clwtot1= cliw(i,k+1) + clcw(i,k+1)
                   dsubclw=((-edt(i)*zd(i,k+1)+zu(i,k+1))*clwtot1   &
                        -(-edt(i)*zd(i,k)  +zu(i,k))  *clwtot  )*g/dp
                   dsubclwm=((-edtm(i)*zdm(i,k+1)+zum(i,k+1))*clwtot1   &
@@ -835,11 +835,11 @@ subroutine cu_gf_driver_run(tottracer,ntrac,garea,im,ix,km,dt,cactiv, &
 !                       +dsubclw*xmb(i)+dsubclws*xmbs(i)+dsubclwm*xmbm(i) &
                       )
                tem1 = max(0.0, min(1.0, (tcr-t(i,k))*tcrf))
-               if (clw(i,k,2) .gt. -999.0) then
-                clw(i,k,1) = max(0.,clw(i,k,1) + tem * tem1)            ! ice
-                clw(i,k,2) = max(0.,clw(i,k,2) + tem *(1.0-tem1))       ! water
+               if (clcw(i,k) .gt. -999.0) then
+                cliw(i,k) = max(0.,cliw(i,k) + tem * tem1)            ! ice
+                clcw(i,k) = max(0.,clcw(i,k) + tem *(1.0-tem1))       ! water
               else
-                clw(i,k,1) = max(0.,clw(i,k,1) + tem)
+                cliw(i,k) = max(0.,cliw(i,k) + tem)
               endif
 
             enddo

physics/cu_ntiedtke.F90

@@ -202,6 +202,10 @@ subroutine cu_ntiedtke_run(pu,pv,pt,pqv,pqvf,ptf,clw,poz,pzz,prsl,prsi,pomg,   &
       real(kind=kind_phys), dimension( ix , km ),     intent(inout) :: pu, pv, pt, pqv
       real(kind=kind_phys), dimension( ix , km ),     intent(in )   :: poz, prsl, pomg, pqvf, ptf
       real(kind=kind_phys), dimension( ix , km+1 ),   intent(in )   :: pzz, prsi
+      ! DH* TODO - check dimensions of clw, ktrac+2 seems to be smaller
+      ! than the actual dimensions (ok as long as only indices 1 and 2
+      ! are accessed here, and as long as these contain what is expected);
+      ! better to expand into the cloud-ice and cloud-water components *DH
       real(kind=kind_phys), dimension( ix , km, ktrac+2 ),    intent(inout ) ::  clw
 
       integer, dimension( lq ),   intent(out)  :: kbot, ktop, kcnv

physics/module_mp_thompson.F90

@@ -414,7 +414,8 @@ SUBROUTINE thompson_init(nwfa2d, nifa2d, nwfa, nifa,  &
                           ids, ide, jds, jde, kds, kde,     &
                           ims, ime, jms, jme, kms, kme,     &
                           its, ite, jts, jte, kts, kte,     &
-                          mpicomm, mpirank, mpiroot, threads)
+                          mpicomm, mpirank, mpiroot,        &
+                          threads, errmsg, errflg)
 
       IMPLICIT NONE
 
@@ -428,6 +429,8 @@ SUBROUTINE thompson_init(nwfa2d, nifa2d, nwfa, nifa,  &
       REAL, DIMENSION(ims:ime,jms:jme),         OPTIONAL, INTENT(IN) :: nwfa2d, nifa2d
       INTEGER, INTENT(IN) :: mpicomm, mpirank, mpiroot
       INTEGER, INTENT(IN) :: threads
+      CHARACTER(len=*), INTENT(INOUT) :: errmsg
+      INTEGER,          INTENT(INOUT) :: errflg
 
 
       INTEGER:: i, j, k, l, m, n
@@ -881,7 +884,8 @@ SUBROUTINE thompson_init(nwfa2d, nifa2d, nwfa, nifa,  &
       ! doing it always ensures that the correct data is in the SIONlib
       ! file containing the precomputed tables *DH
       WRITE (*,*) '  calling table_ccnAct routine'
-      call table_ccnAct
+      call table_ccnAct(errmsg,errflg)
+      if (.not. errflg==0) return
 
 !>  - Call table_efrw() and table_efsw() to creat collision efficiency table 
 !! between rain/snow and cloud water
@@ -4613,15 +4617,17 @@ end subroutine table_dropEvap
 !! vertical velocity, temperature, lognormal mean aerosol radius, and
 !! hygroscopicity, kappa.  The data are read from external file and
 !! contain activated fraction of CCN for given conditions.
-      subroutine table_ccnAct
+      subroutine table_ccnAct(errmess,errflag)
 
       implicit none
 
+!..Error handling variables
+      CHARACTER(len=*), INTENT(INOUT) :: errmess
+      INTEGER,          INTENT(INOUT) :: errflag
 
 !..Local variables
       INTEGER:: iunit_mp_th1, i
       LOGICAL:: opened
-      CHARACTER*64 errmess
 
       iunit_mp_th1 = -1
         DO i = 20,99
@@ -4649,19 +4655,11 @@ subroutine table_ccnAct
       RETURN
  9009 CONTINUE
       WRITE( errmess , '(A,I2)' ) 'module_mp_thompson: error opening CCN_ACTIVATE.BIN on unit ',iunit_mp_th1
-      write(0,*) errmess
-      ! DH* TEMPORARY FIX 20181203
-      call sleep(5)
-      stop
-      ! *DH
+      errflag = 1
       RETURN
  9010 CONTINUE
       WRITE( errmess , '(A,I2)' ) 'module_mp_thompson: error reading CCN_ACTIVATE.BIN on unit ',iunit_mp_th1
-      write(0,*) errmess
-      ! DH* TEMPORARY FIX 20181203
-      call sleep(5)
-      stop
-      ! *DH
+      errflag = 1
       RETURN
 
       end subroutine table_ccnAct

physics/mp_thompson.F90

@@ -124,7 +124,7 @@ subroutine mp_thompson_init(ncol, nlev, is_aerosol_aware, &
                                ims=ims, ime=ime, jms=jms, jme=jme, kms=kms, kme=kme, &
                                its=its, ite=ite, jts=jts, jte=jte, kts=kts, kte=kte, &
                                mpicomm=mpicomm, mpirank=mpirank, mpiroot=mpiroot,    &
-                               threads=threads)
+                               threads=threads, errmsg=errmsg, errflg=errflg)
             if (errflg /= 0) return
          else if (is_aerosol_aware) then
             write(errmsg,fmt='(*(a))') 'Logic error in mp_thompson_init:',                    &
@@ -137,7 +137,7 @@ subroutine mp_thompson_init(ncol, nlev, is_aerosol_aware, &
                                ims=ims, ime=ime, jms=jms, jme=jme, kms=kms, kme=kme, &
                                its=its, ite=ite, jts=jts, jte=jte, kts=kts, kte=kte, &
                                mpicomm=mpicomm, mpirank=mpirank, mpiroot=mpiroot,    &
-                               threads=threads)
+                               threads=threads, errmsg=errmsg, errflg=errflg)
             if (errflg /= 0) return
          end if
 

physics/samfdeepcnv.f

@@ -159,6 +159,7 @@ subroutine samfdeepcnv_run (im,ix,km,cliq,cp,cvap,                &
       logical, intent(in)  :: do_ca
 
       integer, intent(inout)  :: kcnv(im)
+      ! DH* TODO - check dimensions of qtr, ntr+2 correct?  *DH
       real(kind=kind_phys), intent(inout) ::   qtr(ix,km,ntr+2),        &
      &   q1(ix,km), t1(ix,km),   u1(ix,km), v1(ix,km),                  &
      &   cnvw(ix,km),  cnvc(ix,km)

physics/samfshalcnv.f

@@ -116,6 +116,7 @@ subroutine samfshalcnv_run(im,ix,km,cliq,cp,cvap,                 &
      &   prslp(ix,km), garea(im), hpbl(im), dot(ix,km), phil(ix,km)
 !
       integer, intent(inout)  :: kcnv(im)
+      ! DH* TODO - check dimensions of qtr, ntr+2 correct?  *DH
       real(kind=kind_phys), intent(inout) ::   qtr(ix,km,ntr+2),        &
      &   q1(ix,km), t1(ix,km), u1(ix,km), v1(ix,km)
 !