Merge pull request #930 from SamuelTrahanNOAA/sing_prec_from_main

32-bit physics in FV3_RAP
NCAR · Jul 19, 2022 · 12c115e · 12c115e
2 parents dd49119 + 948471f
commit 12c115e
Show file tree

Hide file tree

Showing 21 changed files with 146 additions and 106 deletions.
diff --git a/physics/GFS_suite_interstitial_3.F90 b/physics/GFS_suite_interstitial_3.F90
@@ -148,11 +148,19 @@ subroutine GFS_suite_interstitial_3_run (otsptflag,                 &
           do k=1,levs
             do i=1,im
               kk = max(10,kpbl(i))
+#ifdef SINGLE_PREC
+              if (k < kk) then
+                tem    = rhcbot - (rhcbot-rhcpbl) * (one-prslk(i,k)) / max(one-prslk(i,kk),1e-7)
+              else
+                tem    = rhcpbl - (rhcpbl-rhctop) * (prslk(i,kk)-prslk(i,k)) / max(prslk(i,kk),1e-7)
+              endif
+#else
               if (k < kk) then
                 tem    = rhcbot - (rhcbot-rhcpbl) * (one-prslk(i,k)) / (one-prslk(i,kk))
               else
                 tem    = rhcpbl - (rhcpbl-rhctop) * (prslk(i,kk)-prslk(i,k)) / prslk(i,kk)
               endif
+#endif
               tem      = rhcmax * work1(i) + tem * work2(i)
               rhc(i,k) = max(zero, min(one,tem))
             enddo

diff --git a/physics/GFS_time_vary_pre.fv3.F90 b/physics/GFS_time_vary_pre.fv3.F90
@@ -70,7 +70,7 @@ subroutine GFS_time_vary_pre_timestep_init (jdat, idat, dtp, nsswr,
                   nslwr, nhfrad, idate, debug, me, master, nscyc, sec, phour, zhour, fhour,      &
                   kdt, julian, yearlen, ipt, lprnt, lssav, lsswr, lslwr, solhr, errmsg, errflg)
 
-        use machine,               only: kind_phys
+        use machine,               only: kind_phys, kind_dbl_prec, kind_sngl_prec
 
         implicit none
 
@@ -92,8 +92,10 @@ subroutine GFS_time_vary_pre_timestep_init (jdat, idat, dtp, nsswr,
 
         real(kind=kind_phys), parameter :: con_24  =   24.0_kind_phys
         real(kind=kind_phys), parameter :: con_hr  = 3600.0_kind_phys
-        real(kind=kind_phys) :: rinc(5)
+        real(kind=kind_sngl_prec) :: rinc4(5)
+        real(kind=kind_dbl_prec)  :: rinc8(5)
 
+        integer :: w3kindreal,w3kindint
         integer ::  iw3jdn
         integer :: jd0, jd1
         real    :: fjd
@@ -111,9 +113,19 @@ subroutine GFS_time_vary_pre_timestep_init (jdat, idat, dtp, nsswr,
 
         !--- jdat is being updated directly inside of FV3GFS_cap.F90
         !--- update calendars and triggers
-        rinc(1:5)   = 0
-        call w3difdat(jdat,idat,4,rinc)
-        sec = rinc(4)
+        call w3kind(w3kindreal,w3kindint)
+        if (w3kindreal == 8) then
+           rinc8(1:5) = 0
+           call w3difdat(jdat,idat,4,rinc8)
+           sec = rinc8(4)
+        else if (w3kindreal == 4) then
+           rinc4(1:5) = 0
+           call w3difdat(jdat,idat,4,rinc4)
+           sec = rinc4(4)
+        else
+           write(0,*)' FATAL ERROR: Invalid w3kindreal'
+           call abort
+        endif
         phour = sec/con_hr
         !--- set current bucket hour
         zhour = phour

diff --git a/physics/GFS_time_vary_pre.scm.F90 b/physics/GFS_time_vary_pre.scm.F90
@@ -69,7 +69,7 @@ subroutine GFS_time_vary_pre_timestep_init (jdat, idat, dtp, nsswr, &
         nslwr, idate, debug, me, master, nscyc, sec, phour, zhour, fhour, kdt,   &
         julian, yearlen, ipt, lprnt, lssav, lsswr, lslwr, solhr, errmsg, errflg)
 
-        use machine,               only: kind_phys
+        use machine,               only: kind_phys, kind_dbl_prec, kind_sngl_prec
 
         implicit none
 
@@ -91,8 +91,10 @@ subroutine GFS_time_vary_pre_timestep_init (jdat, idat, dtp, nsswr, &
 
         real(kind=kind_phys), parameter :: con_24  =   24.0_kind_phys
         real(kind=kind_phys), parameter :: con_hr  = 3600.0_kind_phys
-        real(kind=kind_phys) :: rinc(5)
-
+        real(kind=kind_sngl_prec) :: rinc4(5)
+        real(kind=kind_dbl_prec)  :: rinc8(5)
+
+        integer :: w3kindreal,w3kindint
         integer ::  iw3jdn      
         integer :: jd0, jd1
         real    :: fjd
@@ -112,9 +114,19 @@ subroutine GFS_time_vary_pre_timestep_init (jdat, idat, dtp, nsswr, &
         !--- jdat is being updated directly inside of the time integration
         !--- loop of scm.F90
         !--- update calendars and triggers
-        rinc(1:5)   = 0
-        call w3difdat(jdat,idat,4,rinc)
-        sec = rinc(4)
+        call w3kind(w3kindreal,w3kindint)
+        if (w3kindreal == 8) then
+           rinc8(1:5) = 0
+           call w3difdat(jdat,idat,4,rinc8)
+           sec = rinc8(4)
+        else if (w3kindreal == 4) then
+           rinc4(1:5) = 0
+           call w3difdat(jdat,idat,4,rinc4)
+           sec = rina4c(4)
+        else
+           write(0,*)' FATAL ERROR: Invalid w3kindreal'
+           call abort
+        endif
         phour = sec/con_hr
         !--- set current bucket hour
         zhour = phour

diff --git a/physics/aer_cloud.F b/physics/aer_cloud.F
@@ -3477,7 +3477,8 @@ SUBROUTINE EMPIRICAL_PARAM_PHILLIPS(SI, SIW, SW, D_grid_dust,
      & D_grid_bio, n_grid_bio, ijstop_bio, A_solo, n_iw, DSH,
      & Nhet_dep,Nhet_dhf,fdust_dep,P_ice, T_ice)
       implicit none
-      real, intent(IN):: SI, SIW, SW, A_solo,P_ice, T_ice
+      real, intent(IN):: SI, SIW, SW, A_solo
+      real*8, intent(IN):: P_ice, T_ice
       real, dimension(:), intent(IN):: D_grid_dust, n_grid_dust,
      & D_grid_soot, n_grid_soot, D_grid_bio, n_grid_bio
       integer, intent(IN):: ijstop_dust, ijstop_soot, ijstop_bio
@@ -3488,7 +3489,7 @@ SUBROUTINE EMPIRICAL_PARAM_PHILLIPS(SI, SIW, SW, D_grid_dust,
      & num_ic_solo_imm
 
       real, intent (inout) :: DSH, n_iw
-      real, intent (out) :: Nhet_dep,Nhet_dhf,fdust_dep
+      real*8, intent (out) :: Nhet_dep,Nhet_dhf,fdust_dep
 
       real :: dn_in_dust, dn_in_soot, dn_in_bio, dn_in_solo, dNall,
      & dNaux, naux, SS_w, dH_frac_dust, dH_frac_soot, dH_frac_solo, aux,

diff --git a/physics/aerinterp.F90 b/physics/aerinterp.F90
@@ -110,7 +110,8 @@ SUBROUTINE read_aerdataf ( me, master, iflip, idate, FHOUR, errmsg, errflg)
       integer      :: i, j, k, n, ii, imon, klev, n1, n2
       logical      :: file_exist
       integer  IDAT(8),JDAT(8)
-      real(kind=kind_phys) RINC(5), rjday
+      real(kind=kind_phys) rjday
+      real(8) RINC(5)
       integer jdow, jdoy, jday
       real(4) rinc4(5)
       integer w3kindreal,w3kindint      
@@ -244,8 +245,9 @@ SUBROUTINE aerinterpol( me,master,nthrds,npts,IDATE,FHOUR,iflip, jindx1,jindx2,
       real(kind=kind_phys) aerout(npts,lev,ntrcaer)
       real(kind=kind_phys) aerpm(npts,levsaer,ntrcaer)
       real(kind=kind_phys) prsl(npts,lev), aerpres(npts,levsaer)
-      real(kind=kind_phys) RINC(5), rjday
+      real(kind=kind_phys) rjday
       integer jdow, jdoy, jday
+      real(8) RINC(5)
       real(4) rinc4(5)
       integer w3kindreal,w3kindint
 

diff --git a/physics/cires_tauamf_data.F90 b/physics/cires_tauamf_data.F90
@@ -176,8 +176,9 @@ subroutine gfs_idate_calendar(idate, fhour, ddd, fddd)
 !
 !locals
 !
-      real(kind=kind_phys) :: rinc(5), rjday
+      real(kind=kind_phys) :: rjday
       integer              :: jdow, jdoy, jday
+      real(8)              :: rinc(5)
       real(4)              :: rinc4(5)
       integer              :: w3kindreal, w3kindint
 

diff --git a/physics/date_def.f b/physics/date_def.f
@@ -1,13 +1,13 @@
       module date_def
-      use machine,   ONLY: kind_evod
+      use machine,   ONLY: kind_phys
       implicit none
 
 !jw      integer idate(4)
-!jw      real(kind=kind_evod) fhour,shour,thour,z00
-       real(kind=kind_evod) shour,thour,z00
-       real(kind=kind_evod),target :: fhour, zhour
+!jw      real(kind=kind_phys) fhour,shour,thour,z00
+       real(kind=kind_phys) shour,thour,z00
+       real(kind=kind_phys),target :: fhour, zhour
        integer,target :: idate(4),idate7(7)
 !
-      REAL(KIND=KIND_EVOD) ,ALLOCATABLE :: spdmax(:)
+      REAL(KIND=KIND_PHYS) ,ALLOCATABLE :: spdmax(:)
 
       end module date_def
diff --git a/physics/flake.F90 b/physics/flake.F90
@@ -35,23 +35,9 @@ MODULE data_parameters
 !
 ! Description:
 !  Global parameters for the program are defined.
+!  Actually, scratch that. We'll import them from machine.F instead.
 !
-
-IMPLICIT NONE
-
-!=======================================================================
-! Global (i.e. public) Declarations:
-! Parameters for the Program:
-
-  INTEGER, PARAMETER       ::                                         &
-       ireals    = SELECTED_REAL_KIND (12,200),                       &
-                     ! number of desired significant digits for
-                     ! real variables
-                     ! corresponds to 8 byte real variables
-
-       iintegers = KIND  (1)
-                     ! kind-type parameter of the integer values
-                     ! corresponds to the default integers
+  use machine, only: ireals=>kind_phys, iintegers=>kind_INTEGER
 
 !=======================================================================
 

diff --git a/physics/h2ointerp.f90 b/physics/h2ointerp.f90
@@ -144,8 +144,9 @@ subroutine h2ointerpol(me,npts,idate,fhour,jindx1,jindx2,h2oplout,ddy)
 !
       real(kind=kind_phys) ddy(npts)
       real(kind=kind_phys) h2oplout(npts,levh2o,h2o_coeff)
-      real(kind=kind_phys) rinc(5), rjday
+      real(kind=kind_phys) rjday
       integer              jdow, jdoy, jday
+      real(8)              rinc(5)
       real(4)              rinc4(5)
       integer              w3kindreal, w3kindint
 !

diff --git a/physics/iccninterp.F90 b/physics/iccninterp.F90
@@ -143,8 +143,9 @@ SUBROUTINE ciinterpol(me,npts,IDATE,FHOUR,jindx1,jindx2,ddy, &
       real(kind=kind_phys) ciplout(npts,lev),cipm(npts,kcipl)
       real(kind=kind_phys) ccnout(npts,lev),ccnpm(npts,kcipl)
       real(kind=kind_phys) cipres(npts,kcipl), prsl(npts,lev)
-      real(kind=kind_phys) RINC(5), rjday
+      real(kind=kind_phys) rjday
       integer jdow, jdoy, jday
+      real(8) RINC(5)
       real(4) rinc4(5)
       integer w3kindreal,w3kindint
 !

diff --git a/physics/m_micro.meta b/physics/m_micro.meta
@@ -1,7 +1,7 @@
 [ccpp-table-properties]
   name = m_micro
   type = scheme
-  dependencies = aer_cloud.F,aerclm_def.F,cldmacro.F,cldwat2m_micro.F,machine.F,micro_mg_utils.F90,micro_mg2_0.F90,micro_mg3_0.F90,physcons.F90,wv_saturation.F
+  dependencies = aer_cloud.F,aerclm_def.F,cldmacro.F,cldwat2m_micro.F,machine.F,micro_mg_utils.F90,micro_mg2_0.F90,micro_mg3_0.F90,physcons.F90,wv_saturation.F,machine.F
 
 ########################################################################
 [ccpp-arg-table]

diff --git a/physics/machine.F b/physics/machine.F
@@ -5,32 +5,38 @@ module machine
 !!
 
       implicit none
-
-      integer, parameter :: kind_io4  = 4, kind_io8  = 8 , kind_ior = 8 &
-     &,                     kind_evod = 8, kind_dbl_prec = 8            &
-     &,                     kind_sngl_prec = 4, kind_INTEGER = 4        &
-     &,                     kind_LOGICAL = 4
 
-#ifdef SINGLE_PREC
-      integer, parameter :: kind_rad  = kind_sngl_prec                  &
-     &,                     kind_phys = kind_sngl_prec                  &
-     &,                     kind_grid = kind_dbl_prec                   &! atmos_cubed_sphere requres kind_grid=8
-     &,                     kind_REAL = kind_sngl_prec                   ! used in cmp_comm
+      integer, parameter :: kind_sngl_prec = 4                          &
+     &,                     kind_dbl_prec = 8                           &
+# ifdef __PGI
+     &,                     kind_qdt_prec = 8                           &
+# else
+     &,                     kind_qdt_prec = 16                          &
+# endif
+     &,                     kind_integer = 4                            &
+     &,                     kind_logical = 4
+     &,                     kind_io4  = kind_sngl_prec                  &
+     &,                     kind_ior = kind_dbl_prec                    &
+     &,                     kind_grid = kind_dbl_prec
+
+! Physics single precision flag
+#ifndef SINGLE_PREC
+      integer, parameter :: kind_phys = kind_dbl_prec
 #else
-      integer, parameter :: kind_rad  = kind_dbl_prec                   &
-     &,                     kind_phys = kind_dbl_prec                   &
-     &,                     kind_grid = kind_dbl_prec                   &! atmos_cubed_sphere requres kind_grid=8
-     &,                     kind_REAL = kind_dbl_prec                    ! used in cmp_comm
+      integer, parameter :: kind_phys = kind_sngl_prec
 #endif
 
+      integer, parameter :: kind_io8 = kind_phys                        &! Note kind_io8 is not always 8 bytes
+
+! Dynamics single precision flag
 #ifdef OVERLOAD_R4
-      integer, parameter :: kind_dyn  = 4
+      integer, parameter :: kind_dyn  = kind_sngl_prec
 #else
-      integer, parameter :: kind_dyn  = 8
+      integer, parameter :: kind_dyn  = kind_dbl_prec
 #endif
 
 !
-      real(kind=kind_evod), parameter :: mprec = 1.e-12           ! machine precision to restrict dep
-      real(kind=kind_evod), parameter :: grib_undef = 9.99e20     ! grib undefine value
+      real(kind=kind_phys), parameter :: mprec = 1.e-12           ! machine precision to restrict dep
+      real(kind=kind_phys), parameter :: grib_undef = 9.99e20     ! grib undefine value
 !
       end module machine
diff --git a/physics/module_MP_FER_HIRES.F90 b/physics/module_MP_FER_HIRES.F90
@@ -2249,9 +2249,10 @@ REAL FUNCTION DEPOSIT (PP,Tdum,WVdum,RHgrd,I,J,L)   !-- Debug 20120111
 !--- Also uses the Asai (1965) algorithm, but uses a different target
 !      vapor pressure for the adjustment
 !
+        use machine, only: HIGH_PRES => kind_dbl_prec
       IMPLICIT NONE      
 !
-      INTEGER, PARAMETER :: HIGH_PRES=Selected_Real_Kind(15)
+      !INTEGER, PARAMETER :: HIGH_PRES=Selected_Real_Kind(15)
       REAL (KIND=HIGH_PRES), PARAMETER :: RHLIMIT=.001,                 &
      & RHLIMIT1=-RHLIMIT
       REAL (KIND=HIGH_PRES) :: DEP, SSAT

diff --git a/physics/module_nst_parameters.f90 b/physics/module_nst_parameters.f90
@@ -9,8 +9,7 @@
 !! history:
 !!  20210305: X.Li, reduce z_w_max from 30 m to 20 m
 module module_nst_parameters
-  use machine, only :  kind_phys &
-       ,kind_rad ! for astronomy (date) calculations
+  use machine, only :  kind_phys
   !
   ! air constants and coefficients from the atmospehric model
   use physcons, only: &

diff --git a/physics/module_sf_noahmplsm.f90 b/physics/module_sf_noahmplsm.f90
@@ -10823,6 +10823,7 @@ real*8 function zolri(ri,za,z0,zt,zol1,psi_opt)
       real (kind=kind_phys) :: x1,x2,fx1,fx2
       integer :: n
       integer, parameter :: nmax = 20
+      real(kind=kind_phys) zolri_iteration
       !real, dimension(nmax):: zlhux
 !     real  :: zolri2
 
@@ -10855,7 +10856,9 @@ real*8 function zolri(ri,za,z0,zt,zol1,psi_opt)
 
       if (n==nmax .and. abs(x1 - x2) >= 0.01) then
          !if convergence fails, use approximate values:
-         call li_etal_2010(zolri, ri, za/z0, z0/zt)
+         zolri_iteration= zolri
+         call li_etal_2010(zolri_iteration, ri, za/z0, z0/zt)
+         zolri = zolri_iteration
          !zlhux(n)=zolri
          !print*,"iter fail, n=",n," ri=",ri," z0=",z0
       else
@@ -10921,7 +10924,7 @@ real*8 function zolrib(ri,za,z0,zt,logz0,logzt,zol1,psi_opt)
       integer :: n
       integer, parameter :: nmax = 20
       real (kind=kind_phys), dimension(nmax):: zlhux
-      real (kind=kind_phys) :: psit2,psix2
+      real (kind=kind_phys) :: psit2,psix2,zolrib_iteration
 
 !     real    :: psim_unstable, psim_stable
 !     real    :: psih_unstable, psih_stable
@@ -10972,14 +10975,18 @@ real*8 function zolrib(ri,za,z0,zt,logz0,logzt,zol1,psi_opt)
       if (n==nmax .and. abs(zolold - zolrib) > 0.01 ) then
          !print*,"iter fail, n=",n," ri=",ri," z/l=",zolri
          !if convergence fails, use approximate values:
-         call li_etal_2010(zolrib, ri, za/z0, z0/zt)
+         zolrib_iteration = zolrib
+         call li_etal_2010(zolrib_iteration, ri, za/z0, z0/zt)
+         zolrib = zolrib_iteration
          zlhux(n)=zolrib
          !print*,"failed, n=",n," ri=",ri," z0=",z0
          !print*,"z/l=",zlhux(1:nmax)
       else
          !if(zolrib*ri .lt. 0.) then
          !   !print*,"end: wrong quadrants: z/l=",zolrib," ri=",ri
+         !   !phys_temp = zolrib
          !   !call li_etal_2010(zolrib, ri, za/z0, z0/zt)
+         !   !zolrib = phys_temp
          !endif
          !print*,"success,n=",n," ri=",ri," z0=",z0
       endif

diff --git a/physics/mp_nssl.F90 b/physics/mp_nssl.F90
@@ -6,7 +6,7 @@
 !! This module contains the front end to NSSL microphysics scheme.
 module mp_nssl
 
-    use machine, only : kind_phys, kind_real
+    use machine, only : kind_phys
     use module_mp_nssl_2mom, only : nssl_2mom_init, nssl_2mom_driver
 
     implicit none

diff --git a/physics/ozinterp.f90 b/physics/ozinterp.f90
@@ -147,8 +147,9 @@ SUBROUTINE ozinterpol(me,npts,IDATE,FHOUR,jindx1,jindx2,ozplout,ddy)
 !
       real(kind=kind_phys) DDY(npts)
       real(kind=kind_phys) ozplout(npts,levozp,oz_coeff)
-      real(kind=kind_phys) RINC(5), rjday
+      real(kind=kind_phys) rjday
       integer jdow, jdoy, jday
+      real(8) rinc(5)
       real(4) rinc4(5)
       integer w3kindreal,w3kindint
 !