Merge branch 'hotfix-v8.0.1'

This merge corrects several issues in the MPAS-Atmosphere model, specifically: * Fix an OpenMP error in the deallocation of an array (rthdynten) when neither the Grell-Freitas nor the Tiedtke/nTiedtke cumulus schemes are used. * Fix a compilation issue due to a missing include path for the physics_mmm directory when certain Fortran compilers are used. * Fix an issue in reading real-valued global attributes from input files with SMIOL when the input file is of a different precision than the compiled precision of MPAS. * Fix a memory leak for the recloud_p, reice_p, and resnow_p arrays in the deallocate_microphysics routine. * Correct the units and description for the rt_diabatic_tend, pv_vertex, pv_edge, and pv_cell variables in the atmosphere core's Registry.xml file.
MPAS-Dev · Jul 6, 2023 · cd652c3 · cd652c3
2 parents 1d6ec7c + 2bfdc6b
commit cd652c3
Show file tree

Hide file tree

Showing 13 changed files with 83 additions and 87 deletions.
diff --git a/README.md b/README.md
@@ -1,4 +1,4 @@
-MPAS-v8.0.0
+MPAS-v8.0.1
 ====
 
 The Model for Prediction Across Scales (MPAS) is a collaborative project for

diff --git a/src/core_atmosphere/Makefile b/src/core_atmosphere/Makefile
@@ -79,7 +79,7 @@ clean:
 	$(RM) $@ $*.mod
 ifeq "$(GEN_F90)" "true"
 	$(CPP) $(CPPFLAGS) $(PHYSICS) $(CPPINCLUDES) -I./inc $< > $*.f90
-	$(FC) $(FFLAGS) -c $*.f90 $(FCINCLUDES) -I../framework -I../operators -I./physics -I./dynamics -I./diagnostics -I./physics/physics_wrf -I../external/esmf_time_f90
+	$(FC) $(FFLAGS) -c $*.f90 $(FCINCLUDES) -I../framework -I../operators -I./physics -I./dynamics -I./diagnostics -I./physics/physics_wrf -I./physics/physics_mmm -I../external/esmf_time_f90
 else
-	$(FC) $(CPPFLAGS) $(PHYSICS) $(FFLAGS) -c $*.F $(CPPINCLUDES) $(FCINCLUDES) -I./inc -I../framework -I../operators -I./physics -I./dynamics -I./diagnostics -I./physics/physics_wrf -I../external/esmf_time_f90
+	$(FC) $(CPPFLAGS) $(PHYSICS) $(FFLAGS) -c $*.F $(CPPINCLUDES) $(FCINCLUDES) -I./inc -I../framework -I../operators -I./physics -I./dynamics -I./diagnostics -I./physics/physics_wrf -I./physics/physics_mmm -I../external/esmf_time_f90
 endif
diff --git a/src/core_atmosphere/Registry.xml b/src/core_atmosphere/Registry.xml
@@ -1,5 +1,5 @@
 <?xml version="1.0"?>
-<registry model="mpas" core="atmosphere" core_abbrev="atm" version="8.0.0">
+<registry model="mpas" core="atmosphere" core_abbrev="atm" version="8.0.1">
 
 <!-- **************************************************************************************** -->
 <!-- ************************************** Dimensions ************************************** -->
@@ -1634,20 +1634,20 @@
                 <var name="vorticity" type="real" dimensions="nVertLevels nVertices Time" units="s^{-1}"
                      description="Relative vorticity at vertices"/>
 
-                <var name="pv_edge" type="real" dimensions="nVertLevels nEdges Time" units="s^{-1} kg^{-1} m^3"
-                     description="absolute vorticity/rho_zz averaged to the cell edge from the vertices"/>
+                <var name="pv_edge" type="real" dimensions="nVertLevels nEdges Time" units="s^{-1}"
+                     description="absolute vertical vorticity averaged to the cell edge from the vertices"/>
 
                 <var name="rho_edge" type="real" dimensions="nVertLevels nEdges Time" units="kg m^{-3}"
                      description="rho_zz averaged from cell centers to the cell edge"/>
 
                 <var name="ke" type="real" dimensions="nVertLevels nCells Time" units="m^2 s^{-2}"
                      description="Kinetic energy at a cell center"/>
 
-                <var name="pv_vertex" type="real" dimensions="nVertLevels nVertices Time" units="s^{-1} kg^{-1} m^3"
-                     description="absolute vorticity/rho_zz at a vertex"/>
+                <var name="pv_vertex" type="real" dimensions="nVertLevels nVertices Time" units="s^{-1}"
+                     description="absolute vertical vorticity at a vertex"/>
 
-                <var name="pv_cell" type="real" dimensions="nVertLevels nCells Time" units="s^{-1} kg^{-1} m^3"
-                     description="absolute vorticity/rho_zz averaged to the cell center from the vertices"/>
+                <var name="pv_cell" type="real" dimensions="nVertLevels nCells Time" units="s^{-1}"
+                     description="absolute vertical vorticity averaged to the cell center from the vertices"/>
 
                 <var name="dtheta_dt_mp" type="real" dimensions="nVertLevels nCells Time" units="K s^{-1}"
                      description="Potential temperature heating rate from microphysics"/>
@@ -1791,8 +1791,8 @@
                 <var name="tend_theta" name_in_code="theta_m" type="real" dimensions="nVertLevels nCells Time" units="kg K m^{-3} s^{-1}"
                      description="tendency of coupled potential temperature rho*theta_m/zz from dynamics and physics, updated each RK step"/>
 
-                <var name="rt_diabatic_tend" type="real" dimensions="nVertLevels nCells Time" units="kg K s^{-1}"
-                     description="Tendency of coupled potential temperature from physics"/>
+                <var name="rt_diabatic_tend" type="real" dimensions="nVertLevels nCells Time" units="K s^{-1}"
+                     description="Tendency of modified potential temperature due to cloud microphysics"/>
 
                 <var name="euler_tend_u" name_in_code="u_euler" type="real" dimensions="nVertLevels nEdges Time" units="m s^{-2}"
                      description="Tendency of u from dynamics"/>
@@ -2962,6 +2962,9 @@
                      description="Total dry air tendency from physics"
                      persistence="scratch" />
 
+                <var name="rthdynten" type="real" dimensions="nVertLevels nCells Time" units="K s^{-1}"
+                     description="tendency of temperature due to horizontal and vertical advections"/>
+
                 <!-- Scratch variables used when propagating cell-centered winds to edges -->
                 <var name="tend_uzonal" type="real" dimensions="nVertLevels nCells Time" units="m s^{-1} s^{-1}"
                      description="Total cell-centered zonal wind tendency from physics"
@@ -3007,10 +3010,6 @@
                      description="tendency of water vapor due to horizontal and vertical advections"
                      packages="cu_grell_freitas_in;cu_tiedtke_in"/>
 
-                <var name="rthdynten" type="real" dimensions="nVertLevels nCells Time" units="K s^{-1}"
-                     description="tendency of temperature due to horizontal and vertical advections"
-                     packages="cu_grell_freitas_in;cu_tiedtke_in"/>
-
                 <var name="rucuten" type="real" dimensions="nVertLevels nCells Time" units="m s^{-1} s^{-1}"
                      description="tendency of zonal wind due to cumulus convection"
                      packages="cu_grell_freitas_in;cu_tiedtke_in"/>

diff --git a/src/core_atmosphere/dynamics/Makefile b/src/core_atmosphere/dynamics/Makefile
@@ -20,7 +20,7 @@ clean:
 	$(RM) $@ $*.mod
 ifeq "$(GEN_F90)" "true"
 	$(CPP) $(CPPFLAGS) $(PHYSICS) $(CPPINCLUDES) $< > $*.f90
-	$(FC) $(FFLAGS) -c $*.f90 $(FCINCLUDES) -I.. -I../../framework -I../../operators -I../physics -I../physics/physics_wrf -I../../external/esmf_time_f90
+	$(FC) $(FFLAGS) -c $*.f90 $(FCINCLUDES) -I.. -I../../framework -I../../operators -I../physics -I../physics/physics_wrf -I../physics/physics_mmm -I../../external/esmf_time_f90
 else
-	$(FC) $(CPPFLAGS) $(PHYSICS) $(FFLAGS) -c $*.F $(CPPINCLUDES) $(FCINCLUDES) -I.. -I../../framework -I../../operators -I../physics -I../physics/physics_wrf -I../../external/esmf_time_f90
+	$(FC) $(CPPFLAGS) $(PHYSICS) $(FFLAGS) -c $*.F $(CPPINCLUDES) $(FCINCLUDES) -I.. -I../../framework -I../../operators -I../physics -I../physics/physics_wrf -I../physics/physics_mmm -I../../external/esmf_time_f90
 endif
diff --git a/src/core_atmosphere/dynamics/mpas_atm_time_integration.F b/src/core_atmosphere/dynamics/mpas_atm_time_integration.F
@@ -1142,9 +1142,7 @@ subroutine atm_srk3(domain, dt, itimestep, exchange_halo_group)
          call mpas_pool_get_array(tend_physics, 'rqvdynten', rqvdynten)
          call mpas_pool_get_array(state, 'theta_m', theta_m, 2)
 
-         if (associated(tend_physics)) then
-            call mpas_pool_get_array(tend_physics, 'rthdynten', rthdynten)
-         end if
+         call mpas_pool_get_array(tend_physics, 'rthdynten', rthdynten)
 
          !NOTE: The calculation of the tendency due to horizontal and vertical advection for the water vapor mixing ratio
          !requires that the subroutine atm_advance_scalars_mono was called on the third Runge Kutta step, so that a halo
@@ -3888,7 +3886,7 @@ subroutine atm_compute_dyn_tend(tend, tend_physics, state, diag, mesh, configs,
 
       real (kind=RKIND), dimension(:,:), pointer :: rr_save
 
-      real (kind=RKIND), dimension(:,:), pointer :: rthdynten => null()
+      real (kind=RKIND), dimension(:,:), pointer :: rthdynten
 
       real (kind=RKIND), dimension(:,:,:), pointer :: scalars
 
@@ -3931,8 +3929,6 @@ subroutine atm_compute_dyn_tend(tend, tend_physics, state, diag, mesh, configs,
       real (kind=RKIND), pointer :: config_rayleigh_damp_u_timescale_days
       integer, pointer :: config_number_rayleigh_damp_u_levels, config_number_cam_damping_levels
 
-      logical :: inactive_rthdynten
-
 
       call mpas_pool_get_config(mesh, 'sphere_radius', r_earth)
       call mpas_pool_get_config(configs, 'config_coef_3rd_order', coef_3rd_order)
@@ -3982,9 +3978,7 @@ subroutine atm_compute_dyn_tend(tend, tend_physics, state, diag, mesh, configs,
       call mpas_pool_get_array(diag, 'h_divergence', h_divergence)
       call mpas_pool_get_array(diag, 'exner', exner)
 
-      if (associated(tend_physics)) then
-         call mpas_pool_get_array(tend_physics, 'rthdynten', rthdynten)
-      end if
+      call mpas_pool_get_array(tend_physics, 'rthdynten', rthdynten)
 
       call mpas_pool_get_array(mesh, 'weightsOnEdge', weightsOnEdge)
       call mpas_pool_get_array(mesh, 'cellsOnEdge', cellsOnEdge)
@@ -4061,18 +4055,6 @@ subroutine atm_compute_dyn_tend(tend, tend_physics, state, diag, mesh, configs,
       call mpas_pool_get_array(mesh, 'cf2', cf2)
       call mpas_pool_get_array(mesh, 'cf3', cf3)
 
-      !
-      ! rthdynten is currently associated with packages, and if those packages
-      ! are not active at run-time, we need to produce an rthdynten array for
-      ! use in the atm_compute_dyn_tend_work routine
-      !
-      inactive_rthdynten = .false.
-      if (.not. associated(rthdynten)) then
-         allocate(rthdynten(nVertLevels,nCells+1))
-         rthdynten(:,nCells+1) = 0.0_RKIND
-         inactive_rthdynten = .true.
-      end if
-
       call atm_compute_dyn_tend_work(nCells, nEdges, nVertices, nVertLevels, &
          nCellsSolve, nEdgesSolve, vertexDegree, maxEdges, maxEdges2, num_scalars, moist_start, moist_end, &
          fEdge, dvEdge, dcEdge, invDcEdge, invDvEdge, invAreaCell, invAreaTriangle, meshScalingDel2, meshScalingDel4, &
@@ -4094,10 +4076,6 @@ subroutine atm_compute_dyn_tend(tend, tend_physics, state, diag, mesh, configs,
          cellStart, cellEnd, vertexStart, vertexEnd, edgeStart, edgeEnd, &
          cellSolveStart, cellSolveEnd, vertexSolveStart, vertexSolveEnd, edgeSolveStart, edgeSolveEnd)
 
-      if (inactive_rthdynten) then
-         deallocate(rthdynten)
-      end if
-
    end subroutine atm_compute_dyn_tend
 
 

diff --git a/src/core_atmosphere/physics/mpas_atmphys_driver_microphysics.F b/src/core_atmosphere/physics/mpas_atmphys_driver_microphysics.F
@@ -216,9 +216,9 @@ subroutine deallocate_microphysics(configs)
        if(allocated(graupelncv_p) ) deallocate(graupelncv_p )
 
        !cloud water,cloud ice,and snow effective radii:
-       if(.not.allocated(recloud_p) ) allocate(recloud_p(ims:ime,kms:kme,jms:jme) )
-       if(.not.allocated(reice_p)   ) allocate(reice_p(ims:ime,kms:kme,jms:jme)   )
-       if(.not.allocated(resnow_p)  ) allocate(resnow_p(ims:ime,kms:kme,jms:jme)  )
+       if(allocated(recloud_p) ) deallocate(recloud_p )
+       if(allocated(reice_p)   ) deallocate(reice_p   )
+       if(allocated(resnow_p)  ) deallocate(resnow_p  )
 
     microp2_select: select case(microp_scheme)
 

diff --git a/src/core_init_atmosphere/Registry.xml b/src/core_init_atmosphere/Registry.xml
@@ -1,5 +1,5 @@
 <?xml version="1.0"?>
-<registry model="mpas" core="init_atmosphere" core_abbrev="init_atm" version="8.0.0">
+<registry model="mpas" core="init_atmosphere" core_abbrev="init_atm" version="8.0.1">
 
 <!-- **************************************************************************************** -->
 <!-- ************************************** Dimensions ************************************** -->

diff --git a/src/core_landice/Registry.xml b/src/core_landice/Registry.xml
@@ -1,5 +1,5 @@
 <?xml version="1.0"?>
-<registry model="mpas" core="landice" core_abbrev="li" version="8.0.0">
+<registry model="mpas" core="landice" core_abbrev="li" version="8.0.1">
 
 
 <!-- ======================================================================= -->

diff --git a/src/core_ocean/Registry.xml b/src/core_ocean/Registry.xml
@@ -1,5 +1,5 @@
 <?xml version="1.0"?>
-<registry model="mpas" core="ocean" core_abbrev="ocn" version="8.0.0">
+<registry model="mpas" core="ocean" core_abbrev="ocn" version="8.0.1">
 
 	<dims>
 		<dim name="nCells" units="unitless"

diff --git a/src/core_seaice/Registry.xml b/src/core_seaice/Registry.xml
@@ -1,5 +1,5 @@
 <?xml version="1.0"?>
-<registry model="mpas" core="seaice" core_abbrev="seaice" version="8.0.0">
+<registry model="mpas" core="seaice" core_abbrev="seaice" version="8.0.1">
 
 	<dims>
 		<dim name="nCells"

diff --git a/src/core_sw/Registry.xml b/src/core_sw/Registry.xml
@@ -1,5 +1,5 @@
 <?xml version="1.0"?>
-<registry model="mpas" core="sw" core_abbrev="sw" version="8.0.0">
+<registry model="mpas" core="sw" core_abbrev="sw" version="8.0.1">
 	<dims>
 		<dim name="nCells"/>
 		<dim name="nEdges"/>

diff --git a/src/core_test/Registry.xml b/src/core_test/Registry.xml
@@ -1,5 +1,5 @@
 <?xml version="1.0"?>
-<registry model="mpas" core="test" core_abbrev="test" version="8.0.0">
+<registry model="mpas" core="test" core_abbrev="test" version="8.0.1">
 	<dims>
 		<dim name="nCells"/>
 		<dim name="nEdges"/>

diff --git a/src/framework/mpas_io.F b/src/framework/mpas_io.F
@@ -4603,92 +4603,111 @@ subroutine MPAS_io_get_att_real0d(handle, attName, attValue, fieldname, precisio
          local_precision = MPAS_IO_NATIVE_PRECISION
       end if
 
-      ! Query attribute value
 #ifdef MPAS_PIO_SUPPORT
+      ! Query attribute value
       pio_ierr = PIO_inq_att(handle % pio_file, varid, attName, xtype, len)
       if (pio_ierr /= PIO_noerr) then
          if (present(ierr)) ierr = MPAS_IO_ERR_PIO
          return
       end if
-#endif
 
       if ((local_precision == MPAS_IO_SINGLE_PRECISION) .and. &
           (MPAS_IO_NATIVE_PRECISION /= MPAS_IO_SINGLE_PRECISION)) then
 
-#ifdef MPAS_PIO_SUPPORT
          if (xtype /= PIO_REAL) then
             if (present(ierr)) ierr=MPAS_IO_ERR_WRONG_ATT_TYPE
             return
          end if
          pio_ierr = PIO_get_att(handle % pio_file, varid, attName, singleVal)
-#endif
-
-#ifdef MPAS_SMIOL_SUPPORT
-         if (present(fieldname)) then
-            local_ierr = SMIOLf_inquire_att(handle % smiol_file, trim(fieldname), trim(attName), singleVal)
-         else
-            local_ierr = SMIOLf_inquire_att(handle % smiol_file, '', trim(attName), singleVal)
-         end if
-#endif
 
          attValue = real(singleVal,RKIND)
 
       else if ((local_precision == MPAS_IO_DOUBLE_PRECISION) .and. &
           (MPAS_IO_NATIVE_PRECISION /= MPAS_IO_DOUBLE_PRECISION)) then
 
-#ifdef MPAS_PIO_SUPPORT
          if (xtype /= PIO_DOUBLE) then
             if (present(ierr)) ierr=MPAS_IO_ERR_WRONG_ATT_TYPE
             return
          end if
          pio_ierr = PIO_get_att(handle % pio_file, varid, attName, doubleVal)
-#endif
-
-#ifdef MPAS_SMIOL_SUPPORT
-         if (present(fieldname)) then
-            local_ierr = SMIOLf_inquire_att(handle % smiol_file, trim(fieldname), trim(attName), doubleVal)
-         else
-            local_ierr = SMIOLf_inquire_att(handle % smiol_file, '', trim(attName), doubleVal)
-         end if
-#endif
 
          attValue = real(doubleVal,RKIND)
 
       else
 
-#ifdef MPAS_PIO_SUPPORT
          if (xtype /= PIO_DOUBLE .and. xtype /= PIO_REAL) then
             if (present(ierr)) ierr=MPAS_IO_ERR_WRONG_ATT_TYPE
             return
          end if
          pio_ierr = PIO_get_att(handle % pio_file, varid, attName, attValue)
-#endif
-
-#ifdef MPAS_SMIOL_SUPPORT
-         if (present(fieldname)) then
-            local_ierr = SMIOLf_inquire_att(handle % smiol_file, trim(fieldname), trim(attName), attValue)
-         else
-            local_ierr = SMIOLf_inquire_att(handle % smiol_file, '', trim(attName), attValue)
-         end if
-#endif
 
       end if
-#ifdef MPAS_PIO_SUPPORT
       if (pio_ierr /= PIO_noerr) then
          if (present(ierr)) ierr = MPAS_IO_ERR_PIO
          return
       end if
 #endif
 
 #ifdef MPAS_SMIOL_SUPPORT
+      !
+      ! Try to read the attribute in the MPAS native precision
+      !
+      if (present(fieldname)) then
+         local_ierr = SMIOLf_inquire_att(handle % smiol_file, trim(fieldname), &
+                                         trim(attName), attValue)
+      else
+         local_ierr = SMIOLf_inquire_att(handle % smiol_file, '', &
+                                         trim(attName), attValue)
+      end if
+
+      !
+      ! If that fails, perhaps the attribute is in a different precision from
+      ! the native MPAS precision
+      !
+      if (local_ierr == SMIOL_WRONG_ARG_TYPE) then
+         if (MPAS_IO_NATIVE_PRECISION == MPAS_IO_DOUBLE_PRECISION) then
+
+            !
+            ! Try again, but read a single-precision value
+            !
+            if (present(fieldname)) then
+               local_ierr = SMIOLf_inquire_att(handle % smiol_file, trim(fieldname), &
+                                               trim(attName), singleVal)
+            else
+               local_ierr = SMIOLf_inquire_att(handle % smiol_file, '', &
+                                               trim(attName), singleVal)
+            end if
+            attValue = real(singleVal,RKIND)
+
+         else
+
+            !
+            ! Try again, but read a double-precision value
+            !
+            if (present(fieldname)) then
+               local_ierr = SMIOLf_inquire_att(handle % smiol_file, trim(fieldname), &
+                                               trim(attName), doubleVal)
+            else
+               local_ierr = SMIOLf_inquire_att(handle % smiol_file, '', &
+                                               trim(attName), doubleVal)
+            end if
+            attValue = real(doubleVal,RKIND)
+
+         end if
+      end if
+
+      !
+      ! If all of the above were unsuccessful, set attValue to a fill value
+      ! and return an error
+      !
       if (local_ierr /= SMIOL_SUCCESS) then
+         attValue = MPAS_REAL_FILLVAL
          if (local_ierr == SMIOL_WRONG_ARG_TYPE) then
             if (present(ierr)) ierr = MPAS_IO_ERR_WRONG_ATT_TYPE
-            return
          else
             if (present(ierr)) ierr = MPAS_IO_ERR_PIO
-            return
          end if
+         return
       end if
 #endif