Merge branch 'master' into 2642_implement_extends_procedure_and_class

changelog

@@ -310,6 +310,10 @@
 	117) PR #2803 for #2796. Fixes bug where PSyData nodes inserted within
 	LFRic loops did not have the associated variable declared.
 
+	118) PR #2792 for # 2766. Updates to all the examples and test
+	infrastructure now that the index ordering for operators in LFRic
+	has been changed.
+
 release 2.5.0 14th of February 2024
 
 	1) PR #2199 for #2189. Fix bugs with missing maps in enter data

doc/user_guide/dynamo0p3.rst

@@ -1989,11 +1989,10 @@ conventions, are:
       is a rank-3, ``real`` array. Its precision (kind) depends on how
       it is defined in the algorithm layer, see the
       :ref:`lfric-mixed-precision` section for more details. The
-      extents of the first two dimensions are the local degrees of
+      extent of the first dimension is ``<operator_name>"_ncell_3d"``,
+      and of the second and third dimension are the local degrees of
       freedom for the ``to`` and ``from`` function spaces,
-      respectively, and that of the third is
-      ``<operator_name>"_ncell_3d"``. The name of the operator is
-      ``"op_"<argument_position>``. Again the intent is determined
+      respectively. Again the intent is determined
       from the metadata (see :ref:`lfric-api-meta-args`).
 
 4) For each function space in the order they appear in the metadata arguments
@@ -2255,10 +2254,10 @@ as the number of DoFs for each of the dofmaps. The full set of rules is:
 5) For each argument in the ``meta_args`` metadata array:
 
    1) If it is a LMA operator, include a ``real``, 3-dimensional
-      array. The first two dimensions are the local degrees of freedom
-      for the ``to`` and ``from`` spaces, respectively. The third
-      dimension is ``ncell_3d``. The precision of the array depends on
-      how it is defined in the algorithm layer, see the
+      array. The first dimension is ``ncell_3d``. The second and third
+      dimension are the local degrees of freedom for the ``to`` and
+      ``from`` spaces, respectively.  The precision of the array depends
+      on how it is defined in the algorithm layer, see the
       :ref:`lfric-mixed-precision` section for more details;
 
    2) If it is a CMA operator, include a ``real``, 3-dimensional array

doc/user_guide/introduction_to_psykal.rst

@@ -313,7 +313,7 @@ the Built-ins since they are provided as part of the environment.
 
   module solver_mod
     ...
-    use matrix_vector_mm_mod, only: matrix_vector_kernel_mm_type
+    use matrix_vector_mm_kernel_mod, only: matrix_vector_kernel_mm_type
     ...
 
     subroutine jacobi_solver_algorithm(lhs, rhs, mm, mesh, n_iter)

examples/lfric/code/dg_matrix_vector_kernel_mod.F90

@@ -8,7 +8,7 @@
 ! -----------------------------------------------------------------------------
 ! BSD 3-Clause License
 !
-! Modifications copyright (c) 2017-2021, Science and Technology Facilities Council
+! Modifications copyright (c) 2017-2024, Science and Technology Facilities Council
 ! All rights reserved.
 !
 ! Redistribution and use in source and binary forms, with or without
@@ -70,38 +70,38 @@ module dg_matrix_vector_kernel_mod
 
     integer :: operates_on = CELL_COLUMN
   contains
-    procedure, nopass :: dg_matrix_vector_kernel_code
+    procedure, nopass :: dg_matrix_vector_code
   end type
 
   !----------------------------------------------------------------------------
   ! Contained functions/subroutines
   !----------------------------------------------------------------------------
-  public dg_matrix_vector_kernel_code
+  public dg_matrix_vector_code
 
 contains
 
 !> @brief Computes lhs = matrix*x for discontinuous function spaces
 !> @param[in] cell Horizontal cell index
-!! @param[in] nlayers Number of layers
-!! @param[in] ncell_3d Total number of cells
-!! @param[in] ndf1 Number of degrees of freedom per cell for the output field
-!! @param[in] undf1 Unique number of degrees of freedom  for the output field
-!! @param[in] map1 Dofmap for the cell at the base of the column for the
+!> @param[in] nlayers Number of layers
+!> @param[in,out] lhs Output lhs (A*x)
+!> @param[in] x Input data
+!> @param[in] ncell_3d Total number of cells
+!> @param[in] matrix Matrix values in LMA form
+!> @param[in] ndf1 Number of degrees of freedom per cell for the output field
+!> @param[in] undf1 Unique number of degrees of freedom  for the output field
+!> @param[in] map1 Dofmap for the cell at the base of the column for the
 !!            output field
-!! @param[in] map2 Dofmap for the cell at the base of the column for the
+!> @param[in] map2 Dofmap for the cell at the base of the column for the
 !!            input field
-!! @param[in] ndf2 Number of degrees of freedom per cell for the input field
-!! @param[in] undf2 Unique number of degrees of freedom for the input field
-!! @param[in] x input data
-!> @param[in,out] lhs Output lhs (A*x)
-!! @param[in] matrix Matrix values in LMA form
-subroutine dg_matrix_vector_kernel_code(cell,              &
-                                        nlayers,           &
-                                        lhs, x,            &
-                                        ncell_3d,          &
-                                        matrix,            &
-                                        ndf1, undf1, map1, &
-                                        ndf2, undf2, map2)
+!> @param[in] ndf2 Number of degrees of freedom per cell for the input field
+!> @param[in] undf2 Unique number of degrees of freedom for the input field
+subroutine dg_matrix_vector_code(cell,              &
+                                 nlayers,           &
+                                 lhs, x,            &
+                                 ncell_3d,          &
+                                 matrix,            &
+                                 ndf1, undf1, map1, &
+                                 ndf2, undf2, map2)
 
   implicit none
 
@@ -114,7 +114,7 @@ subroutine dg_matrix_vector_kernel_code(cell,              &
   integer(kind=i_def), dimension(ndf2), intent(in)    :: map2
   real(kind=r_def), dimension(undf2),              intent(in)    :: x
   real(kind=r_def), dimension(undf1),              intent(inout) :: lhs
-  real(kind=r_def), dimension(ndf1,ndf2,ncell_3d), intent(in)    :: matrix
+  real(kind=r_def), dimension(ncell_3d,ndf1,ndf2), intent(in)    :: matrix
 
   ! Internal variables
   integer(kind=i_def)               :: df, k, ik
@@ -126,11 +126,11 @@ subroutine dg_matrix_vector_kernel_code(cell,              &
       x_e(df) = x(map2(df)+k)
     end do
     ik = (cell-1)*nlayers + k + 1
-    lhs_e = matmul(matrix(:,:,ik),x_e)
+    lhs_e = matmul(matrix(ik,:,:),x_e)
     do df = 1,ndf1
        lhs(map1(df)+k) = lhs_e(df)
     end do
   end do
-end subroutine dg_matrix_vector_kernel_code
+end subroutine dg_matrix_vector_code
 
 end module dg_matrix_vector_kernel_mod

examples/lfric/code/matrix_vector_kernel_mod.F90

@@ -8,7 +8,7 @@
 ! -----------------------------------------------------------------------------
 ! BSD 3-Clause License
 !
-! Modifications copyright (c) 2019-2021, Science and Technology Facilities Council
+! Modifications copyright (c) 2019-2024, Science and Technology Facilities Council
 ! All rights reserved.
 !
 ! Redistribution and use in source and binary forms, with or without
@@ -39,6 +39,8 @@
 ! Modified by I. Kavcic, Met Office
 ! Modified by R. W. Ford, STFC Daresbury Lab
 !
+!> @brief Matrix-vector multiplication of LMA form of an operator
+!!        by a vector
 module matrix_vector_kernel_mod
 
   use argument_mod,            only : arg_type,                 &
@@ -66,13 +68,13 @@ module matrix_vector_kernel_mod
          /)
     integer :: operates_on = CELL_COLUMN
   contains
-    procedure, nopass :: matrix_vector_kernel_code
+    procedure, nopass :: matrix_vector_code
   end type
 
   !------------------------------------------------------------------------------
   ! Contained functions/subroutines
   !------------------------------------------------------------------------------
-  public matrix_vector_kernel_code
+  public matrix_vector_code
 
 contains
 
@@ -91,13 +93,13 @@ module matrix_vector_kernel_mod
 !! @param[in] undf2 Unique number of degrees of freedom for the input field
 !! @param[in] map2 Dofmap for the cell at the base of the column for the
 !!                 input field
-subroutine matrix_vector_kernel_code(cell,              &
-                                     nlayers,           &
-                                     lhs, x,            &
-                                     ncell_3d,          &
-                                     matrix,            &
-                                     ndf1, undf1, map1, &
-                                     ndf2, undf2, map2)
+subroutine matrix_vector_code(cell,              &
+                              nlayers,           &
+                              lhs, x,            &
+                              ncell_3d,          &
+                              matrix,            &
+                              ndf1, undf1, map1, &
+                              ndf2, undf2, map2)
 
   implicit none
 
@@ -109,24 +111,20 @@ subroutine matrix_vector_kernel_code(cell,              &
   integer(kind=i_def), dimension(ndf2), intent(in) :: map2
   real(kind=r_def), dimension(undf2),              intent(in)    :: x
   real(kind=r_def), dimension(undf1),              intent(inout) :: lhs
-  real(kind=r_def), dimension(ndf1,ndf2,ncell_3d), intent(in)    :: matrix
+  real(kind=r_def), dimension(ncell_3d,ndf1,ndf2), intent(in)    :: matrix
 
   ! Internal variables
-  integer(kind=i_def)               :: df, k, ik
-  real(kind=r_def), dimension(ndf2) :: x_e
-  real(kind=r_def), dimension(ndf1) :: lhs_e
+  integer(kind=i_def) :: df, df2, k, ik
 
-  do k = 0, nlayers-1
-    do df = 1, ndf2
-      x_e(df) = x(map2(df)+k)
-    end do
-    ik = (cell-1)*nlayers + k + 1
-    lhs_e = matmul(matrix(:,:,ik),x_e)
-    do df = 1,ndf1
-       lhs(map1(df)+k) = lhs(map1(df)+k) + lhs_e(df)
+  do df = 1, ndf1
+    do df2 = 1, ndf2
+      do k = 0, nlayers-1
+        ik = (cell-1)*nlayers + k + 1
+        lhs(map1(df)+k) = lhs(map1(df)+k) + matrix(ik,df,df2)*x(map2(df2)+k)
+      end do
     end do
   end do
 
-end subroutine matrix_vector_kernel_code
+end subroutine matrix_vector_code
 
 end module matrix_vector_kernel_mod

examples/lfric/eg11/scaled_matrix_vector_kernel_mod.F90

@@ -9,7 +9,7 @@
 !
 ! BSD 3-Clause License
 !
-! Modifications copyright (c) 2018-2021, Science and Technology Facilities Council
+! Modifications copyright (c) 2018-2024, Science and Technology Facilities Council
 ! All rights reserved.
 !
 ! Redistribution and use in source and binary forms, with or without
@@ -39,7 +39,7 @@
 ! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
 ! POSSIBILITY OF SUCH DAMAGE.
 ! -----------------------------------------------------------------------------
-! Modified I. Kavcic, Met Office
+! Modified by I. Kavcic, Met Office
 
 module scaled_matrix_vector_kernel_mod
 
@@ -133,23 +133,19 @@ subroutine scaled_matrix_vector_code(cell,              &
 
   real(kind=r_def), dimension(undf2),              intent(in)    :: x
   real(kind=r_def), dimension(undf1),              intent(inout) :: lhs
-  real(kind=r_def), dimension(ndf1,ndf2,ncell_3d), intent(in)    :: matrix
+  real(kind=r_def), dimension(ncell_3d,ndf1,ndf2), intent(in)    :: matrix
   real(kind=r_def), dimension(undf1),              intent(in)    :: y
   real(kind=r_def), dimension(undf1),              intent(in)    :: z
 
   ! Internal variables
-  integer(kind=i_def)               :: df, k, ik
-  real(kind=r_def), dimension(ndf2) :: x_e
-  real(kind=r_def), dimension(ndf1) :: lhs_e
-
-  do k = 0, nlayers-1
-    do df = 1, ndf2
-      x_e(df) = x(map2(df)+k)
-    end do
-    ik = (cell-1)*nlayers + k + 1
-    lhs_e = matmul(matrix(:,:,ik),x_e)
-    do df = 1,ndf1
-       lhs(map1(df)+k) = lhs(map1(df)+k) + lhs_e(df)*y(map1(df)+k)*z(map1(df)+k)
+  integer(kind=i_def) :: df, df2, k, ik
+
+  do df = 1, ndf1
+    do df2 = 1, ndf2
+      do k = 0, nlayers-1
+        ik = (cell-1)*nlayers + k + 1
+        lhs(map1(df)+k) = lhs(map1(df)+k) + matrix(ik,df,df2)*x(map2(df2)+k)*y(map1(df)+k)*z(map1(df)+k)
+      end do
     end do
   end do
 
@@ -193,19 +189,19 @@ subroutine opt_scaled_matrix_vector_code(cell,              &
   integer(kind=i_def), dimension(ndf2),  intent(in) :: map2
   real(kind=r_def), dimension(undf2),        intent(in)    :: x
   real(kind=r_def), dimension(undf1),        intent(inout) :: lhs
-  real(kind=r_def), dimension(6,1,ncell_3d), intent(in)    :: matrix
+  real(kind=r_def), dimension(ncell_3d,6,1), intent(in)    :: matrix
   real(kind=r_def), dimension(undf1),        intent(in)    :: y
   real(kind=r_def), dimension(undf1),        intent(in)    :: z
 
   ! Internal variables
-  integer(kind=i_def) :: k, ik, df
+  integer(kind=i_def) :: df, k, ik
 
   ! Hard wired optimisation for desired configuration
-  do df = 1,6
+  do df = 1, ndf1
     do k = 0, nlayers-1
       ik = (cell-1)*nlayers + k + 1
-      lhs(map1(df)+k) = lhs(map1(df)+k) + matrix(df,1,ik)*x(map2(1)+k)*y(map1(df)+k)*z(map1(df)+k)
-    end do
+      lhs(map1(df)+k) = lhs(map1(df)+k) + matrix(ik,df,1)*x(map2(1)+k)*y(map1(df)+k)*z(map1(df)+k)
+      end do
   end do
 
   ! Apply zero flux boundary conditions

examples/lfric/eg12/extract_kernel_with_transformations.py

@@ -65,7 +65,7 @@
 
 # Specify the Kernel name as it appears in the Kernel calls
 # (ending with "_code")
-KERNEL_NAME = "dg_matrix_vector_kernel_code"
+KERNEL_NAME = "dg_matrix_vector_code"
 # Specify the name of Invoke containing the Kernel call. If the name
 # does not correspond to Invoke names in the Algorithm file no Kernels
 # will be extracted.

examples/lfric/eg18/impose_min_flux_kernel_mod.f90

@@ -8,7 +8,7 @@
 ! ----------------------------------------------------------------------------
 ! BSD 3-Clause License
 !
-! Modifications copyright (c) 2021, Science and Technology Facilities Council
+! Modifications copyright (c) 2021-2024, Science and Technology Facilities Council
 ! All rights reserved.
 !
 ! Redistribution and use in source and binary forms, with or without
@@ -36,7 +36,8 @@
 ! OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 ! OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 ! -----------------------------------------------------------------------------
-! Modified: R. W. Ford, STFC Daresbury Lab
+! Modified by R. W. Ford, STFC Daresbury Lab
+! Modified by I. Kavcic, Met Office
 
 module impose_min_flux_kernel_mod
 
@@ -121,8 +122,8 @@ subroutine impose_min_flux_code(cell,              &
   integer(kind=i_def), dimension(ndf2), intent(in) :: map2
   real(kind=r_def), dimension(undf2),               intent(inout) :: flux
   real(kind=r_def), dimension(undf1),               intent(in)    :: field
-  real(kind=r_def), dimension(ndf1,ndf2,ncell_3d1), intent(in)    :: div
-  real(kind=r_def), dimension(ndf1,ndf1,ncell_3d2), intent(in)    :: div_multiplier
+  real(kind=r_def), dimension(ncell_3d1,ndf1,ndf2), intent(in)    :: div
+  real(kind=r_def), dimension(ncell_3d2,ndf1,ndf1), intent(in)    :: div_multiplier
   real(kind=r_def), intent(in)                                    :: field_min
   real(kind=r_def), intent(in)                                    :: dts
 
@@ -147,7 +148,7 @@ subroutine impose_min_flux_code(cell,              &
        flux_change_id = 0_i_def
 
        do df2 = 1, ndf2
-         inc = - dts*div(df1,df2,ik)*div_multiplier(df1,df1,ik)*cell_fluxes(df2)
+         inc = - dts*div(ik,df1,df2)*div_multiplier(ik,df1,df1)*cell_fluxes(df2)
          if ( inc < 0.0_r_def ) then
              inc_n = inc_n - inc
              flux_change_id(df2) = 1_i_def

examples/lfric/eg19/mixed_kernel_mod.f90

@@ -84,7 +84,7 @@ subroutine mixed_code_32(cell, nlayers, rscalar,    &
     real*4,              intent(inout), dimension(undf_w3) :: field_w3
     integer(kind=i_def), intent(in) :: cell
     integer(kind=i_def), intent(in) :: op_ncell_3d
-    real*4,              intent(in), dimension(ndf_w0,ndf_w0,op_ncell_3d) :: op
+    real*4,              intent(in), dimension(op_ncell_3d,ndf_w0,ndf_w0) :: op
 
     write(*,*) "32-bit example called"
 
@@ -104,7 +104,7 @@ subroutine mixed_code_64(cell, nlayers, rscalar,    &
     real*8,              intent(inout), dimension(undf_w3) :: field_w3
     integer(kind=i_def), intent(in) :: cell
     integer(kind=i_def), intent(in) :: op_ncell_3d
-    real*8,              intent(in), dimension(ndf_w0,ndf_w0,op_ncell_3d) :: op
+    real*8,              intent(in), dimension(op_ncell_3d,ndf_w0,ndf_w0) :: op
 
     write(*,*) "64-bit example called"