Merge branch 'fortran-lang:master' into pseudo_inverse

fortran-lang · Dec 18, 2024 · ce36d85 · ce36d85
2 parents 38bae4d + cc30d4c
commit ce36d85
Show file tree

Hide file tree

Showing 11 changed files with 652 additions and 21 deletions.
diff --git a/.github/workflows/CI.yml b/.github/workflows/CI.yml
@@ -19,7 +19,7 @@ jobs:
     strategy:
       fail-fast: false
       matrix:
-        os: [ubuntu-latest, macos-12]
+        os: [ubuntu-latest, macos-13]
         toolchain:
           - {compiler: gcc, version: 10}
           - {compiler: gcc, version: 11}
@@ -34,9 +34,11 @@ jobs:
             toolchain:
               - {compiler: gcc, version: 10}
         exclude:
-          - os: macos-12
+          - os: macos-13
+            toolchain: {compiler: intel-classic, version: '2021.9'}
+          - os: macos-13
             toolchain: {compiler: intel, version: '2024.1'}
-          - os: macos-12
+          - os: macos-13
             toolchain: {compiler: gcc, version: 13}
     env:
       BUILD_DIR: ${{ matrix.build == 'cmake' && 'build' || '.' }}

diff --git a/doc/specs/stdlib_ascii.md b/doc/specs/stdlib_ascii.md
@@ -38,7 +38,7 @@ Converts input character variable to all lowercase.
 
 #### Class
 
-Pure function.
+Elemental function.
 
 #### Argument
 
@@ -70,7 +70,7 @@ Converts input character variable to all uppercase.
 
 #### Class
 
-Pure function.
+Elemental function.
 
 #### Argument
 
@@ -107,7 +107,7 @@ or numeral present next to either of its 2 ends.
 
 #### Class
 
-Pure function.
+Elemental function.
 
 #### Argument
 
@@ -142,7 +142,7 @@ transformed to lowercase.
 
 #### Class
 
-Pure function.
+Elemental function.
 
 #### Argument
 
@@ -174,7 +174,7 @@ Reverses the order of all characters in the input character type.
 
 #### Class
 
-Pure function.
+Elemental function.
 
 #### Argument
 

diff --git a/doc/specs/stdlib_linalg.md b/doc/specs/stdlib_linalg.md
@@ -1701,4 +1701,53 @@ If `err` is not present, exceptions trigger an `error stop`.
 {!example/linalg/example_norm.f90!}
 ```
 
+## `mnorm` - Computes the matrix norm of a generic-rank array.
+
+### Status
+
+Experimental
+
+### Description
+
+This function computes one of several matrix norms of `real` or `complex` array \( A \), depending on 
+the value of the `order` input argument. \( A \) must be an array of rank 2 or higher. For arrays of rank > 2,
+matrix norms are computed over specified dimensions.
+
+### Syntax
+
+`x = ` [[stdlib_linalg(module):mnorm(interface)]] `(a [, order, dim, err])`
+
+### Arguments
+
+`a`: Shall be a rank-n `real` or `complex` array containing the data, where n >= 2. It is an `intent(in)` argument.
+
+`order` (optional): Shall be an `integer` value or a `character` flag that specifies the norm type, as follows. It is an `intent(in)` argument. 
+
+| Integer input    | Character Input                 | Norm type                                                                   |
+|------------------|---------------------------------|-----------------------------------------------------------------------------|
+| `1`              | `'1'`                           | 1-norm (maximum column sum) \( \max_j \sum_i{ \left|a_{i,j}\right| } \)     |
+| `2`              | `'2'`                           | 2-norm (largest singular value)                                             |
+| (not prov.)      | `'Euclidean','Frobenius','Fro'` | Frobenius norm \( \sqrt{\sum_{i,j}{ \left|a_{i,j}\right|^2 }} \)            |
+| `huge(0)`        | `'inf', 'Inf', 'INF'`           | Infinity norm (maximum row sum) \( \max_i \sum_j{ \left|a_{i,j}\right| } \) |
+
+`dim` (optional): For arrays of rank > 2, shall be an integer array of size 2 specifying the dimensions over which to compute the matrix norm. Default value is `[1,2]`. It is an `intent(in)` argument.
+
+`err` (optional): Shall be a `type(linalg_state_type)` value. This is an `intent(out)` argument. 
+
+### Return value
+
+For rank-2 input arrays, the return value `x` is a scalar containing the matrix norm.
+For arrays of rank > 2, if the optional `dim` argument is present, `x` is a rank `n-2` array with the same shape as \( A \) except 
+for dimensions `dim(1)` and `dim(2)`, which are dropped. Each element of `x` contains the matrix norm of the corresponding submatrix of \( A \),
+evaluated over the specified dimensions only, with the given order.
+
+If an invalid norm type is provided, defaults to 1-norm and raises `LINALG_ERROR`.
+Raises `LINALG_VALUE_ERROR` if any of the arguments has an invalid size.
+If `err` is not present, exceptions trigger an `error stop`.
+
+### Example
+
+```fortran
+{!example/linalg/example_mnorm.f90!}
+```
 
diff --git a/example/linalg/CMakeLists.txt b/example/linalg/CMakeLists.txt
@@ -30,6 +30,7 @@ ADD_EXAMPLE(lapack_getrf)
 ADD_EXAMPLE(lstsq1)
 ADD_EXAMPLE(lstsq2)
 ADD_EXAMPLE(norm)
+ADD_EXAMPLE(mnorm)
 ADD_EXAMPLE(get_norm)
 ADD_EXAMPLE(solve1)
 ADD_EXAMPLE(solve2)

diff --git a/example/linalg/example_mnorm.f90 b/example/linalg/example_mnorm.f90
@@ -0,0 +1,26 @@
+program example_mnorm
+    use stdlib_linalg, only: mnorm
+    use stdlib_kinds, only: sp
+    implicit none
+    real(sp) :: a(3,3), na
+    real(sp) :: b(3,3,4), nb(4)  ! Array of 4 3x3 matrices
+
+    ! Initialize example matrix
+    a = reshape([1, 2, 3, 4, 5, 6, 7, 8, 9], [3, 3])
+
+    ! Compute Euclidean norm of single matrix
+    na = mnorm(a, 'Euclidean')
+    print *, "Euclidean norm of matrix a:", na
+
+    ! Initialize array of matrices
+    b(:,:,1) = a
+    b(:,:,2) = 2*a
+    b(:,:,3) = 3*a
+    b(:,:,4) = 4*a
+
+    ! Compute infinity norm of each 3x3 matrix in b
+    nb = mnorm(b, 'inf', dim=[1,2])
+
+    ! 18.0000000       36.0000000       54.0000000       72.0000000
+    print *, "Infinity norms of matrices in b:", nb
+end program example_mnorm
diff --git a/include/common.fypp b/include/common.fypp
@@ -157,6 +157,22 @@ $:"s" if cmplx=="c" else "d" if cmplx=="z" else "x" if cmplx=="y" else "q" if cm
 #{if rank > 0}#(${"0" + ",0" * (rank - 1)}$)#{endif}#
 #:enddef
 
+#! Generates an array rank suffix with a fixed integer size for all dimensions.
+#!
+#! Args:
+#!     rank (int): Rank of the variable
+#!     size (int): Size along each dimension
+#!
+#! Returns:
+#!     Array rank suffix string
+#! E.g.,
+#!   fixedranksuffix(3,4)
+#!   -> (4,4,4) 
+#!
+#:def fixedranksuffix(rank,size)
+#{if rank > 0}#(${str(size) + (","+str(size)) * (rank - 1)}$)#{endif}#
+#:enddef
+
 #! Joins stripped lines with given character string
 #!
 #! Args:
@@ -227,7 +243,7 @@ ${prefix + joinstr.join([line.strip() for line in txt.split("\n")]) + suffix}$
 #!   Array rank suffix string enclosed in braces
 #!
 #! E.g.,
-#!   select_subarray(5 , [(4, 'i'), (5, 'j')])}$
+#!   select_subarray(5 , [(4, 'i'), (5, 'j')])
 #!   -> (:, :, :, i, j)
 #!
 #:def select_subarray(rank, selectors)
@@ -327,6 +343,34 @@ ${prefix + joinstr.join([line.strip() for line in txt.split("\n")]) + suffix}$
     #:endcall
 #:enddef
 
+#!
+#! Generates a list of loop variables from an array
+#!
+#! Args:
+#!   varname(str): Name of the array variable to be used as prefix
+#!   n      (int): Number of loop variables to be created
+#!   offset (int): Optional index offset
+#!
+#! Returns:
+#!   Variable definition string 
+#!
+#! E.g.,
+#!  loop_array_variables('j', 5)
+#!   -> "j(1), j(2), j(3), j(4), j(5)
+#!
+#!  loop_array_variables('j', 5, 2)
+#!   -> "j(3), j(4), j(5), j(6), j(7)
+#!
+#:def loop_array_variables(varname, n, offset=0)
+  #:assert n > 0
+    #:call join_lines(joinstr=", ")
+      #:for i in range(1, n + 1)
+        ${varname}$(${i+offset}$)
+      #:endfor
+    #:endcall
+#:enddef
+
+
 #! Generates an array shape specifier from an N-D array size
 #!
 #! Args:

diff --git a/src/stdlib_ascii.fypp b/src/stdlib_ascii.fypp
@@ -70,35 +70,35 @@ module stdlib_ascii
 
     !> Returns a new character sequence which is the lower case 
     !> version of the input character sequence
-    !> This method is pure and returns a character sequence
+    !> This method is elemental and returns a character sequence
     interface to_lower
         module procedure :: to_lower
     end interface to_lower
 
     !> Returns a new character sequence which is the upper case
     !> version of the input character sequence
-    !> This method is pure and returns a character sequence
+    !> This method is elemental and returns a character sequence
     interface to_upper
         module procedure :: to_upper
     end interface to_upper
 
     !> Returns a new character sequence which is the title case
     !> version of the input character sequence
-    !> This method is pure and returns a character sequence
+    !> This method is elemental and returns a character sequence
     interface to_title
         module procedure :: to_title
     end interface to_title
 
     !> Returns a new character sequence which is the sentence case
     !> version of the input character sequence
-    !> This method is pure and returns a character sequence
+    !> This method is elemental and returns a character sequence
     interface to_sentence
         module procedure :: to_sentence
     end interface to_sentence
 
     !> Returns a new character sequence which is reverse of
     !> the input charater sequence
-    !> This method is pure and returns a character sequence
+    !> This method is elemental and returns a character sequence
     interface reverse
         module procedure :: reverse
     end interface reverse
@@ -220,7 +220,7 @@ contains
 
     !> Returns the corresponding lowercase letter, if `c` is an uppercase
     !> ASCII character, otherwise `c` itself.
-    pure function char_to_lower(c) result(t)
+    elemental function char_to_lower(c) result(t)
         character(len=1), intent(in) :: c !! A character.
         character(len=1)             :: t
         integer, parameter :: wp= iachar('a')-iachar('A'), BA=iachar('A'), BZ=iachar('Z')
@@ -234,7 +234,7 @@ contains
 
     !> Returns the corresponding uppercase letter, if `c` is a lowercase
     !> ASCII character, otherwise `c` itself.
-    pure function char_to_upper(c) result(t)
+    elemental function char_to_upper(c) result(t)
         character(len=1), intent(in) :: c !! A character.
         character(len=1)             :: t
         integer, parameter :: wp= iachar('a')-iachar('A'), la=iachar('a'), lz=iachar('z')
@@ -250,7 +250,7 @@ contains
     !> ([Specification](../page/specs/stdlib_ascii.html#to_lower))
     !>
     !> Version: experimental
-    pure function to_lower(string) result(lower_string)
+    elemental function to_lower(string) result(lower_string)
         character(len=*), intent(in) :: string
         character(len=len(string)) :: lower_string
         integer :: i
@@ -265,7 +265,7 @@ contains
     !> ([Specification](../page/specs/stdlib_ascii.html#to_upper))
     !>
     !> Version: experimental
-    pure function to_upper(string) result(upper_string)
+    elemental function to_upper(string) result(upper_string)
         character(len=*), intent(in) :: string
         character(len=len(string)) :: upper_string
         integer :: i
@@ -280,7 +280,7 @@ contains
     !> ([Specification](../page/specs/stdlib_ascii.html#to_title))
     !>
     !> Version: experimental
-    pure function to_title(string) result(title_string)
+    elemental function to_title(string) result(title_string)
         character(len=*), intent(in) :: string
         character(len=len(string)) :: title_string
         integer :: i
@@ -307,7 +307,7 @@ contains
     !> ([Specification](../page/specs/stdlib_ascii.html#to_sentence))
     !>
     !> Version: experimental
-    pure function to_sentence(string) result(sentence_string)
+    elemental function to_sentence(string) result(sentence_string)
         character(len=*), intent(in) :: string
         character(len=len(string)) :: sentence_string
         integer :: i, n
@@ -333,7 +333,7 @@ contains
     !> ([Specification](../page/specs/stdlib_ascii.html#reverse))
     !>
     !> Version: experimental
-    pure function reverse(string) result(reverse_string)
+    elemental function reverse(string) result(reverse_string)
         character(len=*), intent(in) :: string
         character(len=len(string)) :: reverse_string
         integer :: i, n