stdlib

`stdlib_linalg`

:local:
:depth: 2

The stdlib linear algebra library provides high-level APIs for dealing with common linear algebra operations.

BLAS and LAPACK

Status

Experimental

Description

BLAS and LAPACK backends provide efficient low level implementations of many linear algebra algorithms, and are employed for non-trivial operators. A Modern Fortran version of the Reference-LAPACK 3.10.1 implementation is provided as a backend. Modern Fortran modules with full explicit typing features are provided after an automated conversion of the legacy codes:

[stdlib_linalg_blas(module)], [stdlib_linalg_lapack(module)] provide kind-agnostic interfaces to all functions.
Both libraries are available for 32- (sp), 64- (dp) and 128-bit (qp) real and complex numbers (the latter if available in the current build)
Free format, lower-case style
implicit none(type, external) applied to all procedures and modules
intent added and all pure procedures where possible
stdlib provides all procedures in two different flavors: (a) original BLAS/LAPACK names with a prefix stdlib_?<name> (ex: stdlib_dgemv, stdlib_sgemv); (b) A generic, kind agnostic <name>, i.e. gemv.
F77-style parameters removed, and all numeric constants have been generalized with KIND-dependent Fortran intrinsics.
preprocessor-based OpenMP directives retained. The single-source module structure hopefully allows for cross-procedural inlining which is otherwise impossible without link-time optimization.

When available, highly optimized libraries that take advantage of specialized processor instructions should be preferred over the stdlib implementation. Examples of such libraries are: OpenBLAS, MKL (TM), Accelerate, and ATLAS. In order to enable their usage, simply ensure that the following pre-processor macros are defined:

STDLIB_EXTERNAL_BLAS wraps all BLAS procedures (except for the 128-bit ones) to an external library
STDLIB_EXTERNAL_LAPACK wraps all LAPACK procedures (except for the 128-bit ones) to an external library

These can be enabled during the build process. For example, with CMake, one can enable these preprocessor directives using add_compile_definitions(STDLIB_EXTERNAL_BLAS STDLIB_EXTERNAL_LAPACK). The same is possible from the fpm branch, where the cpp preprocessor is enabled by default. For example, the macros can be added to the project’s manifest:

# Link against appropriate external BLAS and LAPACK libraries, if necessary
[build]
link = ["blas", "lapack"]  

[dependencies]
stdlib="*"

# Macros are only needed if using an external library
[preprocess]
[preprocess.cpp]
macros = ["STDLIB_EXTERNAL_BLAS", "STDLIB_EXTERNAL_LAPACK"]

or directly via compiler flags:

fpm build --flag "-DSTDLIB_EXTERNAL_BLAS -DSTDLIB_EXTERNAL_LAPACK -lblas -llapack".

Syntax

All procedures in the BLAS and LAPACK backends follow the standard interfaces from the Reference LAPACK. So, the online Users Guide should be consulted for the full API and descriptions of procedure arguments and their usage.

The stdlib implementation makes both kind-agnostic and specific procedure interfaces available via modules [stdlib_linalg_blas(module)] and [stdlib_linalg_lapack(module)]. Because all procedures start with a letter that indicates the base datatype, the stdlib generic interface drops the heading letter and contains all kind-dependent implementations. For example, the generic interface to the axpy function looks like:

!> AXPY: constant times a vector plus a vector.
interface axpy
    module procedure stdlib_saxpy
    module procedure stdlib_daxpy
    module procedure stdlib_qaxpy
    module procedure stdlib_caxpy
    module procedure stdlib_zaxpy
    module procedure stdlib_waxpy
end interface axpy

The generic interface is the endpoint for using an external library. Whenever the latter is used, references to the internal module procedures are replaced with interfaces to the external library, for example:

!> AXPY: constant times a vector plus a vector.
interface axpy
    pure subroutine caxpy(n,ca,cx,incx,cy,incy)
        import sp,dp,qp,ilp,lk 
        implicit none 
        complex(sp), intent(in) :: ca,cx(*)
        integer(ilp), intent(in) :: incx,incy,n
        complex(sp), intent(inout) :: cy(*)
    end subroutine caxpy
    ! [....]
    module procedure stdlib_qaxpy
end interface axpy

Note that the 128-bit functions are only provided by stdlib and always point to the internal implementation. Because 128-bit precision is identified as [stdlib_kinds(module):qp], initials for 128-bit procedures were labelled as q (quadruple-precision reals) and w (“wide” or quadruple-precision complex numbers). Extended precision ([stdlib_kinds(module):xdp]) calculations are labelled as x (extended-precision reals). and y (extended-precision complex numbers).

Example

```{literalinclude} ../../example/linalg/example_blas_gemv.f90 :language: fortran

```{literalinclude} ../../example/linalg/example_lapack_getrf.f90
:language: fortran

Licensing

The Fortran Standard Library is distributed under the MIT License. LAPACK and its contained BLAS are a freely-available software package. They are available from netlib via anonymous ftp and the World Wide Web. Thus, they can be included in commercial software packages (and have been). The license used for the BLAS and LAPACK backends is the modified BSD license.

The header of the LICENSE.txt file has as its licensing requirements:

Copyright (c) 1992-2013 The University of Tennessee and The University
                        of Tennessee Research Foundation.  All rights
                        reserved.
Copyright (c) 2000-2013 The University of California Berkeley. All
                        rights reserved.
Copyright (c) 2006-2013 The University of Colorado Denver.  All rights
                        reserved.

$COPYRIGHT$

Additional copyrights may follow

$HEADER$

Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are
met:

- Redistributions of source code must retain the above copyright
  notice, this list of conditions and the following disclaimer.

- Redistributions in binary form must reproduce the above copyright
  notice, this list of conditions and the following disclaimer listed
  in this license in the documentation and/or other materials
  provided with the distribution.

- Neither the name of the copyright holders nor the names of its
  contributors may be used to endorse or promote products derived from
  this software without specific prior written permission.

The copyright holders provide no reassurances that the source code
provided does not infringe any patent, copyright, or any other
intellectual property rights of third parties.  The copyright holders
disclaim any liability to any recipient for claims brought against
recipient by any third party for infringement of that parties
intellectual property rights.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, 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.

So the license for the LICENSE.txt code is compatible with the use of modified versions of the code in the Fortran Standard Library under the MIT license.
Credit for the BLAS, LAPACK libraries should be given to the LAPACK authors. According to the original license, we also changed the name of the routines and commented the changes made to the original.

`diag` - Create a diagonal array or extract the diagonal elements of an array

Status

Stable

Class

Pure function.

Description

Create a diagonal array or extract the diagonal elements of an array

Syntax

d = diag (a [, k])

Arguments

a: Shall be a rank-1 or or rank-2 array. If a is a rank-1 array (i.e. a vector) then diag returns a rank-2 array with the elements of a on the diagonal. If a is a rank-2 array (i.e. a matrix) then diag returns a rank-1 array of the diagonal elements.

k (optional): Shall be a scalar of type integer and specifies the diagonal. The default k = 0 represents the main diagonal, k > 0 are diagonals above the main diagonal, k < 0 are diagonals below the main diagonal.

Return value

Returns a diagonal array or a vector with the extracted diagonal elements.

Example

```{literalinclude} ../../example/linalg/example_diag1.f90 :language: fortran

```{literalinclude} ../../example/linalg/example_diag2.f90
:language: fortran

```{literalinclude} ../../example/linalg/example_diag3.f90 :language: fortran

```{literalinclude} ../../example/linalg/example_diag4.f90
:language: fortran

```{literalinclude} ../../example/linalg/example_diag5.f90 :language: fortran

## `eye` - Construct the identity matrix

### Status

Stable

### Class

Pure function.

### Description

Constructs the identity matrix.

### Syntax

`I = ` `eye` `(dim1 [, dim2] [, mold])`

### Arguments

- `dim1`: A scalar of type `integer`. This is an `intent(in)` argument and specifies the number of rows.
- `dim2`: A scalar of type `integer`. This is an optional `intent(in)` argument specifying the number of columns. If not provided, the matrix is square (`dim1 = dim2`).
- `mold`: A scalar of any supported `integer`, `real`, or `complex` type. This is an optional `intent(in)` argument. If provided, the returned identity matrix will have the same type and kind as `mold`. If not provided, the matrix will be of type `real(real64)` by default.

### Return value

Returns the identity matrix, with ones on the main diagonal and zeros elsewhere. 

- By default, the return value is of type `real(real64)`, which is recommended for arithmetic safety.
- If the `mold` argument is provided, the return value will match the type and kind of `mold`, allowing for arbitrary `integer`, `real`, or `complex` return types.

### Example

```fortran
!> Return default type (real64)
A = eye(2,2)/2             !! A == diag([0.5_dp, 0.5_dp])
!> Return 32-bit complex
A = eye(2,2, mold=(0.0,0.0))/2   !! A == diag([(0.5,0.5), (0.5,0.5)])

```{literalinclude} ../../example/linalg/example_eye1.f90 :language: fortran

```{literalinclude} ../../example/linalg/example_eye2.f90
:language: fortran

`trace` - Trace of a matrix

Status

Stable

Description

Trace of a matrix (rank-2 array)

Syntax

result = trace (A)

Arguments

A: Shall be a rank-2 array. If A is not square, then trace(A) will return the sum of diagonal values from the square sub-section of A.

Return value

Returns the trace of the matrix, i.e. the sum of diagonal elements.

Example

```{literalinclude} ../../example/linalg/example_trace.f90 :language: fortran

## `outer_product` - Computes the outer product of two vectors

### Status

Experimental

### Description

Computes the outer product of two vectors

### Syntax

`d = ` `outer_product` `(u, v)`

### Arguments

`u`: Shall be a rank-1 array

`v`: Shall be a rank-1 array

### Return value

Returns a rank-2 array equal to `u v^T` (where `u, v` are considered column vectors). The shape of the returned array is `[size(u), size(v)]`.

### Example

```{literalinclude} ../../example/linalg/example_outer_product.f90
:language: fortran

`kronecker_product` - Computes the Kronecker product of two rank-2 arrays

Status

Experimental

Description

Computes the Kronecker product of two rank-2 arrays

Syntax

C = kronecker_product (A, B)

Arguments

A: Shall be a rank-2 array with dimensions M1, N1

B: Shall be a rank-2 array with dimensions M2, N2

Return value

Returns a rank-2 array equal to A \otimes B. The shape of the returned array is [M1*M2, N1*N2].

Example

```{literalinclude} ../../example/linalg/example_kronecker_product.f90 :language: fortran

## `cross_product` - Computes the cross product of two vectors

### Status

Experimental

### Description

Computes the cross product of two vectors

### Syntax

`c = ` `cross_product` `(a, b)`

### Arguments

`a`: Shall be a rank-1 and size-3 array

`b`: Shall be a rank-1 and size-3 array

### Return value

Returns a rank-1 and size-3 array which is perpendicular to both `a` and `b`.

### Example

```{literalinclude} ../../example/linalg/example_cross_product.f90
:language: fortran

`is_square` - Checks if a matrix is square

Status

Experimental

Description

Checks if a matrix is square

Syntax

d = is_square (A)

Arguments

A: Shall be a rank-2 array

Return value

Returns a logical scalar that is .true. if the input matrix is square, and .false. otherwise.

Example

```{literalinclude} ../../example/linalg/example_is_square.f90 :language: fortran

## `is_diagonal` - Checks if a matrix is diagonal

### Status

Experimental

### Description

Checks if a matrix is diagonal

### Syntax

`d = ` `is_diagonal` `(A)`

### Arguments

`A`: Shall be a rank-2 array

### Return value

Returns a `logical` scalar that is `.true.` if the input matrix is diagonal, and `.false.` otherwise.
Note that nonsquare matrices may be diagonal, so long as `a_ij = 0` when `i /= j`.

### Example

```{literalinclude} ../../example/linalg/example_is_diagonal.f90
:language: fortran

`is_symmetric` - Checks if a matrix is symmetric

Status

Experimental

Description

Checks if a matrix is symmetric

Syntax

d = is_symmetric (A)

Arguments

A: Shall be a rank-2 array

Return value

Returns a logical scalar that is .true. if the input matrix is symmetric, and .false. otherwise.

Example

```{literalinclude} ../../example/linalg/example_is_symmetric.f90 :language: fortran

## `is_skew_symmetric` - Checks if a matrix is skew-symmetric

### Status

Experimental

### Description

Checks if a matrix is skew-symmetric

### Syntax

`d = ` `is_skew_symmetric` `(A)`

### Arguments

`A`: Shall be a rank-2 array

### Return value

Returns a `logical` scalar that is `.true.` if the input matrix is skew-symmetric, and `.false.` otherwise.

### Example

```{literalinclude} ../../example/linalg/example_is_skew_symmetric.f90
:language: fortran

`hermitian` - Compute the Hermitian version of a rank-2 matrix

Status

Experimental

Description

Compute the Hermitian version of a rank-2 matrix. For complex matrices, the function returns the conjugate transpose (conjg(transpose(a))). For real or integer matrices, the function returns the transpose (transpose(a)).

Syntax

h = hermitian (a)

Arguments

a: Shall be a rank-2 array of type integer, real, or complex. The input matrix a is not modified.

Return value

Returns a rank-2 array of the same shape and type as a. If a is of type complex, the Hermitian matrix is computed as conjg(transpose(a)). For real or integer types, it is equivalent to the intrinsic transpose(a).

Example

```{literalinclude} ../../example/linalg/example_hermitian.f90 :language: fortran

## `is_hermitian` - Checks if a matrix is Hermitian

### Status

Experimental

### Description

Checks if a matrix is Hermitian

### Syntax

`d = ` `is_hermitian` `(A)`

### Arguments

`A`: Shall be a rank-2 array

### Return value

Returns a `logical` scalar that is `.true.` if the input matrix is Hermitian, and `.false.` otherwise.

### Example

```{literalinclude} ../../example/linalg/example_is_hermitian.f90
:language: fortran

`is_triangular` - Checks if a matrix is triangular

Status

Experimental

Description

Checks if a matrix is triangular

Syntax

d = is_triangular (A,uplo)

Arguments

A: Shall be a rank-2 array

uplo: Shall be a single character from {'u','U','l','L'}

Return value

Returns a logical scalar that is .true. if the input matrix is the type of triangular specified by uplo (upper or lower), and .false. otherwise. Note that the definition of triangular used in this implementation allows nonsquare matrices to be triangular. Specifically, upper triangular matrices satisfy a_ij = 0 when j < i, and lower triangular matrices satisfy a_ij = 0 when j > i.

Example

```{literalinclude} ../../example/linalg/example_is_triangular.f90 :language: fortran

## `is_hessenberg` - Checks if a matrix is hessenberg

### Status

Experimental

### Description

Checks if a matrix is Hessenberg

### Syntax

`d = ` `is_hessenberg` `(A,uplo)`

### Arguments

`A`: Shall be a rank-2 array

`uplo`: Shall be a single character from `{'u','U','l','L'}`

### Return value

Returns a `logical` scalar that is `.true.` if the input matrix is the type of Hessenberg specified by `uplo` (upper or lower), and `.false.` otherwise.
Note that the definition of Hessenberg used in this implementation allows nonsquare matrices to be Hessenberg.
Specifically, upper Hessenberg matrices satisfy `a_ij = 0` when `j < i-1`, and lower Hessenberg matrices satisfy `a_ij = 0` when `j > i+1`.

### Example

```{literalinclude} ../../example/linalg/example_is_hessenberg.f90
:language: fortran

`solve` - Solves a linear matrix equation or a linear system of equations.

Status

Stable

Description

This function computes the solution to a linear matrix equation ( A \cdot x = b ), where ( A ) is a square, full-rank, real or complex matrix.

Result vector or array x returns the exact solution to within numerical precision, provided that the matrix is not ill-conditioned. An error is returned if the matrix is rank-deficient or singular to working precision. The solver is based on LAPACK’s *GESV backends.

Syntax

Pure interface:

x = solve (a, b)

Expert interface:

x = solve (a, b [, overwrite_a], err)

Arguments

a: Shall be a rank-2 real or complex square array containing the coefficient matrix. It is normally an intent(in) argument. If overwrite_a=.true., it is an intent(inout) argument and is destroyed by the call.

b: Shall be a rank-1 or rank-2 array of the same kind as a, containing the right-hand-side vector(s). It is an intent(in) argument.

overwrite_a (optional): Shall be an input logical flag. if .true., input matrix a will be used as temporary storage and overwritten. This avoids internal data allocation. This is an intent(in) argument.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument. The function is not pure if this argument is provided.

Return value

For a full-rank matrix, returns an array value that represents the solution to the linear system of equations.

Raises LINALG_ERROR if the matrix is singular to working precision. Raises LINALG_VALUE_ERROR if the matrix and rhs vectors have invalid/incompatible sizes. If err is not present, exceptions trigger an error stop.

Example

```{literalinclude} ../../example/linalg/example_solve1.f90 :language: fortran

```{literalinclude} ../../example/linalg/example_solve2.f90
:language: fortran

`solve_lu` - Solves a linear matrix equation or a linear system of equations (subroutine interface).

Status

Stable

Description

This subroutine computes the solution to a linear matrix equation ( A \cdot x = b ), where ( A ) is a square, full-rank, real or complex matrix.

Result vector or array x returns the exact solution to within numerical precision, provided that the matrix is not ill-conditioned. An error is returned if the matrix is rank-deficient or singular to working precision. If all optional arrays are provided by the user, no internal allocations take place. The solver is based on LAPACK’s *GESV backends.

Syntax

Simple (Pure) interface:

call solve_lu (a, b, x)

Expert (Pure) interface:

call solve_lu (a, b, x [, pivot, overwrite_a, err])

Arguments

b: Shall be a rank-1 or rank-2 array of the same kind as a, containing the right-hand-side vector(s). It is an intent(in) argument.

x: Shall be a rank-1 or rank-2 array of the same kind and size as b, that returns the solution(s) to the system. It is an intent(inout) argument, and must have the contiguous property.

pivot (optional): Shall be a rank-1 array of the same kind and matrix dimension as a, providing storage for the diagonal pivot indices. It is an intent(inout) arguments, and returns the diagonal pivot indices.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument.

Return value

For a full-rank matrix, returns an array value that represents the solution to the linear system of equations.

Example

```{literalinclude} ../../example/linalg/example_solve3.f90 :language: fortran

## `solve_chol` - Solves a linear system using Cholesky factorization (one-shot interface). 

### Status

Experimental

### Description

This subroutine computes the solution to a linear matrix equation \( A \cdot x = b \), where \( A \) is a symmetric (or Hermitian) positive definite matrix. It combines Cholesky factorization and the solve step in a single call.

Result vector or array `x` returns the exact solution to within numerical precision, provided that the matrix is positive definite.
An error is returned if the matrix is not positive definite or has incompatible sizes with the right-hand-side.
Use this routine for one-time solves. For repeated solves with the same matrix but different right-hand sides, use `cholesky` followed by `solve_lower_chol`/`solve_upper_chol` for better performance.
The solver is based on LAPACK's `*POSV` backends.

### Syntax

Simple (`Pure`) interface:

`call ` `solve_chol` `(a, b, x)`

Expert (`Pure`) interface:

`call ` `solve_chol` `(a, b, x [, lower, overwrite_a, err])`

### Arguments

`a`: Shall be a rank-2 `real` or `complex` square array containing the coefficient matrix. It is an `intent(inout)` argument. By default (`overwrite_a=.false.`) its contents are preserved; if `overwrite_a=.true.`, it is used as temporary storage and its contents are destroyed by the call.

`b`: Shall be a rank-1 or rank-2 array of the same kind as `a`, containing the right-hand-side vector(s). It is an `intent(in)` argument.

`x`: Shall be a rank-1 or rank-2 array of the same kind and size as `b`, that returns the solution(s) to the system. It is an `intent(inout)` argument, and must have the `contiguous` property. 

`lower` (optional): Shall be an input `logical` flag. If `.true.` (default), the lower triangular Cholesky factorization is computed. If `.false.`, the upper triangular factorization is computed. It is an `intent(in)` argument.

`overwrite_a` (optional): Shall be an input `logical` flag. If `.true.`, input matrix `a` will be used as temporary storage and overwritten. This avoids internal data allocation. This is an `intent(in)` argument.

`err` (optional): Shall be a `type(linalg_state_type)` value. This is an `intent(out)` argument.

### Return value

For a positive definite matrix, returns an array value that represents the solution to the linear system of equations.

Raises `LINALG_ERROR` if the matrix is not positive definite.
Raises `LINALG_VALUE_ERROR` if the matrix and rhs vectors have invalid/incompatible sizes.
If `err` is not present, exceptions trigger an `error stop`.

### Example

```{literalinclude} ../../example/linalg/example_solve_chol.f90
:language: fortran

`solve_lower_chol` - Solves a linear system using pre-computed lower Cholesky factor.

Status

Experimental

Description

This subroutine computes the solution to a linear matrix equation ( A \cdot x = b ), where ( A ) is a symmetric (or Hermitian) positive definite matrix that has been previously factorized using the Cholesky decomposition (via cholesky with lower=.true.).

Result vector or array x returns the exact solution to within numerical precision, provided that the factorization is correct. An error is returned if the matrix and right-hand-side have incompatible sizes. The solver is based on LAPACK’s *POTRS backends.

Syntax

call solve_lower_chol (l, b, x [, err])

Arguments

l: Shall be a rank-2 real or complex square array containing the lower Cholesky factor L (output of cholesky(..., lower=.true.)). It is an intent(in) argument.

b: Shall be a rank-1 or rank-2 array of the same kind as l, containing the right-hand-side vector(s). It is an intent(in) argument.

x: Shall be a rank-1 or rank-2 array of the same kind and size as b, that returns the solution(s) to the system. It is an intent(inout) argument, and must have the contiguous property.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument.

Return value

For a correctly factorized matrix, returns an array value that represents the solution to the linear system of equations.

Raises LINALG_VALUE_ERROR if the matrix and rhs vectors have invalid/incompatible sizes. If err is not present, exceptions trigger an error stop.

Example

```{literalinclude} ../../example/linalg/example_solve_lower_chol.f90 :language: fortran

## `solve_upper_chol` - Solves a linear system using pre-computed upper Cholesky factor. 

### Status

Experimental

### Description

This subroutine computes the solution to a linear matrix equation \( A \cdot x = b \), where \( A \) is a symmetric (or Hermitian) positive definite matrix that has been **previously factorized** using the Cholesky decomposition (via `cholesky` with `lower=.false.`).

Result vector or array `x` returns the exact solution to within numerical precision, provided that the factorization is correct.
An error is returned if the matrix and right-hand-side have incompatible sizes.
The solver is based on LAPACK's `*POTRS` backends.

### Syntax

`call ` `solve_upper_chol` `(u, b, x [, err])`

### Arguments

`u`: Shall be a rank-2 `real` or `complex` square array containing the **upper** Cholesky factor `U` (output of `cholesky(..., lower=.false.)`). It is an `intent(in)` argument.

`b`: Shall be a rank-1 or rank-2 array of the same kind as `u`, containing the right-hand-side vector(s). It is an `intent(in)` argument.

`x`: Shall be a rank-1 or rank-2 array of the same kind and size as `b`, that returns the solution(s) to the system. It is an `intent(inout)` argument, and must have the `contiguous` property. 

`err` (optional): Shall be a `type(linalg_state_type)` value. This is an `intent(out)` argument.

### Return value

For a correctly factorized matrix, returns an array value that represents the solution to the linear system of equations.

Raises `LINALG_VALUE_ERROR` if the matrix and rhs vectors have invalid/incompatible sizes.
If `err` is not present, exceptions trigger an `error stop`.

### Example

```{literalinclude} ../../example/linalg/example_solve_upper_chol.f90
:language: fortran

`lstsq` - Computes the least squares solution to a linear matrix equation.

Status

Stable

Description

This function computes the least-squares solution to a linear matrix equation ( A \cdot x = b ).

Result vector x returns the approximate solution that minimizes the 2-norm ( A \cdot x - b _2 ), i.e., it contains the least-squares solution to the problem. Matrix A may be full-rank, over-determined, or under-determined. The solver is based on LAPACK’s *GELSD backends.

Syntax

x = lstsq (a, b, [, cond, overwrite_a, rank, err])

Arguments

a: Shall be a rank-2 real or complex array containing the coefficient matrix. It is an intent(inout) argument.

b: Shall be a rank-1 or rank-2 array of the same kind as a, containing one or more right-hand-side vector(s), each in its leading dimension. It is an intent(in) argument.

cond (optional): Shall be a scalar real value cut-off threshold for rank evaluation: s_i >= cond*maxval(s), i=1:rank. Shall be a scalar, intent(in) argument.

overwrite_a (optional): Shall be an input logical flag. If .true., input matrix A will be used as temporary storage and overwritten. This avoids internal data allocation. This is an intent(in) argument.

rank (optional): Shall be an integer scalar value, that contains the rank of input matrix A. This is an intent(out) argument.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument.

Return value

Returns an array value of the same kind and rank as b, containing the solution(s) to the least squares system.

Raises LINALG_ERROR if the underlying Singular Value Decomposition process did not converge. Raises LINALG_VALUE_ERROR if the matrix and right-hand-side vector have invalid/incompatible sizes. Exceptions trigger an error stop.

Example

```{literalinclude} ../../example/linalg/example_lstsq1.f90 :language: fortran

## `solve_lstsq` - Compute the least squares solution to a linear matrix equation (subroutine interface). 

### Status

Stable

### Description

This subroutine computes the least-squares solution to a linear matrix equation \( A \cdot x = b \).

Result vector `x` returns the approximate solution that minimizes the 2-norm \( || A \cdot x - b ||_2 \), i.e., it contains the least-squares solution to the problem. Matrix `A` may be full-rank, over-determined, or under-determined. The solver is based on LAPACK's `*GELSD` backends.

### Syntax

`call ` `solve_lstsq` `(a, b, x, [, real_storage, int_storage, [cmpl_storage, ] cond, singvals, overwrite_a, rank, err])`

### Arguments

`a`: Shall be a rank-2 `real` or `complex` array containing the coefficient matrix. It is an `intent(inout)` argument.

`b`: Shall be a rank-1 or rank-2 array of the same kind as `a`, containing one or more right-hand-side vector(s), each in its leading dimension. It is an `intent(in)` argument. 

`x`: Shall be an array of same kind and rank as `b`, and leading dimension of at least `n`, containing the solution(s) to the least squares system. It is an `intent(inout)` argument.

`real_storage` (optional): Shall be a `real` rank-1 array of the same kind `a`, providing working storage for the solver. It minimum size can be determined with a call to `lstsq_space`. It is an `intent(inout)` argument.

`int_storage` (optional): Shall be an `integer` rank-1 array, providing working storage for the solver. It minimum size can be determined with a call to `lstsq_space`. It is an `intent(inout)` argument.

`cmpl_storage` (optional): For `complex` systems, it shall be a `complex` rank-1 array, providing working storage for the solver. It minimum size can be determined with a call to `lstsq_space`. It is an `intent(inout)` argument.

`cond` (optional): Shall be a scalar `real` value cut-off threshold for rank evaluation: `s_i >= cond*maxval(s), i=1:rank`. Shall be a scalar, `intent(in)` argument.

`singvals` (optional): Shall be a `real` rank-1 array of the same kind `a` and size at least `min(m,n)`, returning the list of singular values `s(i)>=cond*maxval(s)` from the internal SVD, in descending order of magnitude. It is an `intent(out)` argument.

`overwrite_a` (optional): Shall be an input `logical` flag. If `.true.`, input matrix `A` will be used as temporary storage and overwritten. This avoids internal data allocation. This is an `intent(in)` argument.

`rank` (optional): Shall be an `integer` scalar value, that contains the rank of input matrix `A`. This is an `intent(out)` argument.

`err` (optional): Shall be a `type(linalg_state_type)` value. This is an `intent(out)` argument.

### Return value

Returns an array value that represents the solution to the least squares system.

Raises `LINALG_ERROR` if the underlying Singular Value Decomposition process did not converge.
Raises `LINALG_VALUE_ERROR` if the matrix and right-hand-side vector have invalid/incompatible sizes.
Exceptions trigger an `error stop`.

### Example

```{literalinclude} ../../example/linalg/example_lstsq2.f90
:language: fortran

`lstsq_space` - Compute internal working space requirements for the least squares solver.

Status

Stable

Description

This subroutine computes the internal working space requirements for the least-squares solver, solve_lstsq .

Syntax

call lstsq_space (a, b, lrwork, liwork [, lcwork])

Arguments

a: Shall be a rank-2 real or complex array containing the linear system coefficient matrix. It is an intent(in) argument.

b: Shall be a rank-1 or rank-2 array of the same kind as a, containing the system’s right-hand-side vector(s). It is an intent(in) argument.

lrwork: Shall be an integer scalar, that returns the minimum array size required for the real working storage to this system.

liwork: Shall be an integer scalar, that returns the minimum array size required for the integer working storage to this system.

lcwork (complex a, b): For a complex system, shall be an integer scalar, that returns the minimum array size required for the complex working storage to this system.

`constrained_lstsq` - Compute the solution of the equality-constrained least-squares problem

Status

Experimental

Description

This function computes the solution (x) of the equality-constrained linear least-squares problem \(\begin{aligned} \mathrm{minimize} & \quad \| Ax - b \|^2 \\ \mathrm{subject~to} & \quad Cx = d, \end{aligned}\) where (A) is an ( m \times n ) matrix (with (m \geq n)) and (C) a ( p \times n) matrix (with (p \leq n)). The solver is based on LAPACK’s *GLSE backends.

Syntax

x = constrained_lstsq (A, b, C, d[, overwrite_matrices, err])

Arguments

a: Shall be a rank-2 real or complex array used in the definition of the least-squares cost. It is an intent(inout) argument.

b: Shall be a rank-1 array of the same kind as a appearing in the definition of the least-squares cost. It is an intent(inout) argument.

c: Shall be a rank-2 real or complex array of the same kind as a defining the linear equality constraints. It is an intent(inout) argument.

d: Shall be a rank-1 array of the same kind as a appearing in the definition of the linear equality constraints.

overwrite_matrices (optional): Shall be an input logical flag. If .true., the input matrices and vectors will be overwritten during the computation of the solution. It is an intent(in) argument.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument.

Return value

Returns an array of the same kind as a containing the solution of the equality constrained least-squares problem.

Raises LINALG_ERROR if the underlying constrained least-squares solver did not converge. Raises LINALG_VALUE_ERROR if the matrices and vectors have invalid/incompatible dimensions. Exceptions trigger an error stop.

Example

```{literalinclude} ../../example/linalg/example_constrained_lstsq1.f90 :language: fortran

## `solve_constrained_lstsq` - Compute the solution of the equality-constrained least squares problem (subroutine interface) {#solve-constrained-lstsq}

### Status

Experimental

### Description

This subroutine computes the solution \(x\) of the equality-constrained linear least-squares problem
$$
\begin{aligned}
    \mathrm{minimize}   &   \quad \| Ax - b \|^2 \\
    \mathrm{subject~to} &   \quad   Cx = d,
\end{aligned}
$$
where \(A\) is an \( m \times n \) matrix (with \(m \geq n\)) and 
\(C\) a \( p \times n\) matrix (with \(p \leq n\)). The solver is based on LAPACK's `*GLSE` backends.

### Syntax


`call ` `solve_constrained_lstsq` `(a, b, c, d, x [, storage, overwrite_matrices, err])`

### Arguments

`a`: Shall be a rank-2 `real` or `complex` array used in the definition of the least-squares cost. It is an `intent(inout)` argument.

`b`: Shall be a rank-1 array of the same kind as `a` appearing in the definition of the least-squares cost. It is an `intent(inout)` argument.

`c`: Shall be a rank-2 `real` or `complex` array of the same kind as `a` defining the linear equality constraints. It is an `intent(inout)` argument.

`d`: Shall be a rank-1 array of the same kind as `a` appearing in the definition of the linear equality constraints.

`x`: Shall be a rank-1 array of the same kind as `a`. On exit, it contains the solution of the constrained least-squares problem. It is an `intent(out)` argument.

`storage` (optional): Shall be a rank-1 array of the same kind as `a` providing working storage for the solver. Its minimum size can be determined with a call to [stdlib_linalg(module):constrained_lstsq_space(interface)]. It is an `intent(out)` argument.

`overwrite_matrices` (optional): Shall be an input `logical` flag. If `.true.`, the input matrices and vectors will be overwritten during the computation of the solution. It is an `intent(in)` argument.

`err` (optional): Shall be a `type(linalg_state_type)` value. This is an `intent(out)` argument.

### Example

```{literalinclude} ../../example/linalg/example_constrained_lstsq2.f90
:language: fortran

`constrained_lstsq_space` - Compute internal workspace requirements for the constrained least-square solver

Status

Experimental

Description

This subroutine computes the internal workspace requirements for the constrained least-squares solver, [stdlib_linalg(module):solve_constrained_lstsq(interface)].

Syntax

call [stdlib_linalg(module):constrained_lstsq_space(interface)](a, c, lwork [, err])

Arguments

a: Shall be a rank-2 real or complex array used in the definition of the least-squares cost. It is an intent(in) argument.

c: Shall be a rank-2 real or complex array of the same kind as a defining the linear equality constraints. It is an intent(in) argument. lwork: Shall be an integer scalar returning the optimal size required for the workspace array to solve the constrained least-squares problem.

`weighted_lstsq` - Computes the weighted least squares solution to a linear matrix equation

Status

Experimental

Description

This function computes the weighted least-squares solution to a linear matrix equation ( A \cdot x = b ) where each observation has a different weight.

The solver minimizes the weighted 2-norm ( | D(Ax - b) |_2^2 ) where ( D = \mathrm{diag}(\sqrt{w}) ) is a diagonal matrix formed from the square roots of the weight vector. This is equivalent to transforming the problem to ordinary least-squares through row scaling. The solver is based on LAPACK’s *GELSD backends.

Syntax

x = weighted_lstsq (w, a, b [, cond, overwrite_a, rank, err])

Arguments

w: Shall be a rank-1 real array containing the weight vector. For complex a and b, w shall use the same real kind as the components of a. All weights must be positive. It is an intent(in) argument.

a: Shall be a rank-2 real or complex array containing the coefficient matrix. It is an intent(inout) argument.

b: Shall be a rank-1 array of the same kind as a, containing the right-hand-side vector. It is an intent(in) argument.

cond (optional): Shall be a scalar real value cut-off threshold for rank evaluation: singular values with s_i <= cond*maxval(s) are treated as zero, and the rank counts values with s_i > cond*maxval(s), i=1:rank. Shall be a scalar, intent(in) argument.

overwrite_a (optional): Shall be an input logical flag. If .true., input matrix a will be used as temporary storage and overwritten. This avoids internal data allocation. This is an intent(in) argument.

rank (optional): Shall be an integer scalar value, that contains the rank of input matrix a. This is an intent(out) argument.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument.

Return value

Returns an array value of the same kind as a, containing the weighted least-squares solution.

Raises LINALG_VALUE_ERROR if any weight is non-positive or if the matrix and weight/right-hand-side have incompatible sizes. Raises LINALG_ERROR if the underlying Singular Value Decomposition process did not converge. Exceptions trigger an error stop.

Example

```{literalinclude} ../../example/linalg/example_weighted_lstsq.f90 :language: fortran

## `solve_weighted_lstsq` - Compute the weighted least squares solution to a linear matrix equation (subroutine interface). {#solve-weighted-lstsq}

### Status

Experimental

### Description

This subroutine computes the weighted least-squares solution to a linear matrix equation \( A \cdot x = b \) where each observation has a different weight.

The solver minimizes the weighted 2-norm \( \| D(Ax - b) \|_2^2 \) where \( D = \mathrm{diag}(\sqrt{w}) \) is a diagonal matrix formed from the square roots of the weight vector. This is equivalent to transforming the problem to ordinary least-squares through row scaling. The solver is based on LAPACK's `*GELSD` backends.

### Syntax

`call ` `solve_weighted_lstsq` `(w, a, b, x [, cond, overwrite_a, rank, err])`

### Arguments

`w`: Shall be a rank-1 `real` array containing the weight vector. For complex `a` and `b`, `w` shall use the same real kind as the components of `a`. All weights must be positive. It is an `intent(in)` argument.

`a`: Shall be a rank-2 `real` or `complex` array containing the coefficient matrix. It is an `intent(inout)` argument.

`b`: Shall be a rank-1 array of the same kind as `a`, containing the right-hand-side vector. It is an `intent(in)` argument.

`x`: Shall be a rank-1 array of the same kind as `a`, and size of at least `n`, containing the solution to the weighted least squares system. It is an `intent(inout)` argument.

`cond` (optional): Shall be a scalar `real` value cut-off threshold for rank evaluation: singular values with `s_i <= cond*maxval(s)` are treated as zero, and the rank counts values with `s_i > cond*maxval(s), i=1:rank`. Shall be a scalar, `intent(in)` argument.

`overwrite_a` (optional): Shall be an input `logical` flag. If `.true.`, input matrix `a` will be used as temporary storage and overwritten. This avoids internal data allocation. This is an `intent(in)` argument.

`rank` (optional): Shall be an `integer` scalar value, that contains the rank of input matrix `a`. This is an `intent(out)` argument.

`err` (optional): Shall be a `type(linalg_state_type)` value. This is an `intent(out)` argument.

### Return value

Returns an array value of the same kind as `a`, containing the weighted least-squares solution.

Raises `LINALG_VALUE_ERROR` if any weight is non-positive or if the matrix and weight/right-hand-side have incompatible sizes.
Raises `LINALG_ERROR` if the underlying Singular Value Decomposition process did not converge.
Exceptions trigger an `error stop`.

### Example

```{literalinclude} ../../example/linalg/example_weighted_lstsq.f90
:language: fortran

`det` - Computes the determinant of a square matrix

Status

Stable

Description

This function computes the determinant of a real or complex square matrix.

This interface comes with a pure version det(a), and a non-pure version det(a,overwrite_a,err) that allows for more expert control.

Syntax

c = det (a [, overwrite_a, err])

Arguments

a: Shall be a rank-2 square array

overwrite_a (optional): Shall be an input logical flag. if .true., input matrix a will be used as temporary storage and overwritten. This avoids internal data allocation. This is an intent(in) argument.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument.

Return value

Returns a real scalar value of the same kind of a that represents the determinant of the matrix.

Raises LINALG_ERROR if the matrix is singular. Raises LINALG_VALUE_ERROR if the matrix is non-square. Exceptions are returned to the err argument if provided; an error stop is triggered otherwise.

Example

```{literalinclude} ../../example/linalg/example_determinant.f90 :language: fortran

## `.det.` - Determinant operator of a square matrix

### Status

Stable

### Description

This operator returns the determinant of a real square matrix.

This interface is equivalent to the `pure` version of determinant `det`. 

### Syntax

`c = ` `operator(.det.)` `a`

### Arguments

`a`: Shall be a rank-2 square array of any `real` or `complex` kinds. It is an `intent(in)` argument.

### Return value

Returns a real scalar value that represents the determinant of the matrix.

Raises `LINALG_ERROR` if the matrix is singular.
Raises `LINALG_VALUE_ERROR` if the matrix is non-square.
Exceptions trigger an `error stop`.

### Example

```{literalinclude} ../../example/linalg/example_determinant2.f90
:language: fortran

`qr` - Compute the QR factorization of a matrix

Status

Experimental

Description

This subroutine computes the QR factorization of a real or complex matrix: ( A = Q R ) where ( Q ) is orthonormal and ( R ) is upper-triangular. Matrix ( A ) has size [m,n], with ( m \ge n ).

The results are returned in output matrices ( Q ) and (R ), that have the same type and kind as ( A ). Given k = min(m,n), one can write ( A = ( Q_1 Q_2 ) \cdot ( \frac{R_1}{0}) ). Because the lower rows of ( R ) are zeros, a reduced problem ( A = Q_1 R_1 ) may be solved. The size of the input arguments determines what problem is solved: on full matrices (shape(Q)==[m,m], shape(R)==[m,n]), the full problem is solved. On reduced matrices (shape(Q)==[m,k], shape(R)==[k,n]), the reduced problem is solved.

Syntax

call qr (a, q, r [, pivots] [, overwrite_a] [, storage] [, err])

Arguments

a: Shall be a rank-2 real or complex array containing the coefficient matrix of size [m,n]. It is an intent(in) argument, if overwrite_a=.false.. Otherwise, it is an intent(inout) argument, and is destroyed upon return.

q: Shall be a rank-2 array of the same kind as a, containing the orthonormal matrix q. It is an intent(out) argument. It should have a shape equal to either [m,m] or [m,k], whether the full or the reduced problem is sought for.

r: Shall be a rank-2 array of the same kind as a, containing the upper triangular matrix r. It is an intent(out) argument. It should have a shape equal to either [m,n] or [k,n], whether the full or the reduced problem is sought for.

pivots (optional): Shall be an integer array of size n. If provided, QR factorization with column-pivoting is being computed. It is an intent(out) argument.

overwrite_a (optional): Shall be an input logical flag (default: .false.). If .true., input matrix a will be used as temporary storage and overwritten. This avoids internal data allocation. It is an intent(in) argument.

storage (optional): Shall be a rank-1 array of the same type and kind as a, providing working storage for the solver. Its minimum size can be determined with a call to qr_space. It is an intent(out) argument.

err (optional): Shall be a type(linalg_state_type) value. It is an intent(out) argument.

Return value

Returns the QR factorization matrices into the ( Q ) and ( R ) arguments and the optional pivots in pivots.

Raises LINALG_VALUE_ERROR if any of the matrices has invalid or unsuitable size for the full/reduced problem. Raises LINALG_ERROR on insufficient user storage space. If the state argument err is not present, exceptions trigger an error stop.

Example

```{literalinclude} ../../example/linalg/example_qr.f90 :language: fortran

```{literalinclude} ../../example/linalg/example_pivoting_qr.f90
:language: fortran

`qr_space` - Compute internal working space requirements for the QR factorization.

Status

Experimental

Description

This subroutine computes the internal working space requirements for the QR factorization, qr .

Syntax

call qr_space (a, lwork, [, pivoting] [, err])

Arguments

a: Shall be a rank-2 real or complex array containing the coefficient matrix. It is an intent(in) argument.

lwork: Shall be an integer scalar, that returns the minimum array size required for the working storage in qr to factorize a.

pivoting (optional): Shall a logical flag (default: .false.). If .true., on exit lwork is the optimal workspace size for the QR factorization with column pivoting. If .false., lwork is the optimal workspace size for the standard QR factorization. It is an intent(in) argument.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument.

Example

```{literalinclude} ../../example/linalg/example_qr_space.f90 :language: fortran

```{literalinclude} ../../example/linalg/example_pivoting_qr_space.f90
:language: fortran

`schur` - Compute the Schur decomposition of a matrix

Status

Experimental

Description

This subroutine computes the Schur decomposition of a real or complex matrix: ( A = Z T Z^H ), where ( Z ) is unitary (orthonormal) and ( T ) is upper-triangular (for complex matrices) or quasi-upper-triangular (for real matrices, with possible ( 2 \times 2 ) blocks on the diagonal). Matrix ( A ) has size [n,n].

The results are returned in output matrices ( T ) and ( Z ). Matrix ( T ) is the Schur form, and matrix ( Z ) is the unitary transformation matrix such that ( A = Z T Z^H ). If requested, the eigenvalues of ( T ) can also be returned as a complex array of size [n].

Syntax

call schur (a, t [, z,] [, eigvals] [, overwrite_a] [, storage] [, err])

Arguments

a: Shall be a rank-2 real or complex array containing the matrix to be decomposed. It is an intent(inout) argument if overwrite_a = .true.; otherwise, it is an intent(in) argument.
t: Shall be a rank-2 array of the same kind as a, containing the Schur form ( T ) of the matrix. It is an intent(out) argument and should have a shape of [n,n].
z: Shall be a rank-2 array of the same kind as a, containing the unitary matrix ( Z ). It is an intent(out) argument and is optional. If provided, it should have the shape [n,n].
eigvals (optional): Shall be a rank-1 complex or real array of the same kind as a, containing the eigenvalues of ( A ) (the diagonal elements of ( T )), or their real component only. The array must be of size [n]. If not provided, the eigenvalues are not returned. It is an intent(out) argument.
overwrite_a (optional): Shall be a logical flag (default: .false.). If .true., the input matrix a will be overwritten and destroyed upon return, avoiding internal data allocation. It is an intent(in) argument.
storage (optional): Shall be a rank-1 array of the same type and kind as a, providing working storage for the solver. Its minimum size can be determined with a call to schur_space. It is an intent(inout) argument.
err (optional): Shall be a type(linalg_state_type) value. It is an intent(out) argument. If not provided, exceptions trigger an error stop.

Return value

Returns the Schur decomposition matrices into the ( T ) and ( Z ) arguments. If eigvals is provided, it will also return the eigenvalues of the matrix ( A ).

Raises LINALG_VALUE_ERROR if any of the matrices have invalid or unsuitable size for the decomposition.
Raises LINALG_VALUE_ERROR if the real component is only requested, but the eigenvalues have non-trivial imaginary parts. Raises LINALG_ERROR on insufficient user storage space.
If the state argument err is not present, exceptions trigger an error stop.

Example

```{literalinclude} ../../example/linalg/example_schur_eigvals.f90 :language: fortran

---

## `schur_space` - Compute internal working space requirements for the Schur decomposition

### Status

Experimental

### Description

This subroutine computes the internal working space requirements for the Schur decomposition, `schur`.

### Syntax

`call ` `schur_space` `(a, lwork, [, err])`

### Arguments

- `a`: Shall be a rank-2 `real` or `complex` array containing the matrix to be decomposed. It is an `intent(in)` argument.

- `lwork`: Shall be an `integer` scalar that returns the minimum array size required for the working storage in `schur` to decompose `a`. It is an `intent(out)` argument.

- `err` (optional): Shall be a `type(linalg_state_type)` value. It is an `intent(out)` argument.

### Return value

Returns the required working storage size `lwork` for the Schur decomposition. This value can be used to pre-allocate a workspace array in case multiple Schur decompositions of the same matrix size are needed. If pre-allocated working arrays are provided, no internal memory allocations will take place during the decomposition.

## `eig` - Eigenvalues and Eigenvectors of a Square Matrix

### Status

Stable

### Description

This subroutine computes the solution to the eigenproblem \( A \cdot \bar{v} - \lambda \cdot \bar{v} \), 
where \( A \) is a square, full-rank, `real` or `complex` matrix, or to the generalized eigenproblem \( A \cdot \bar{v} - \lambda \cdot B \cdot \bar{v} \), 
where \( B \) is a square matrix with the same type, kind and size as \( A \).

Result array `lambda` returns the eigenvalues of \( A \). The user can request eigenvectors to be returned: if provided, on output `left` will contain the left eigenvectors, `right` the right eigenvectors of \( A \).
Both `left` and `right` are rank-2 arrays, where eigenvectors are stored as columns.
The solver is based on LAPACK's `*GEEV` (standard eigenproblem) and `*GGEV` (generalized eigenproblem) backends.

### Syntax

For the standard eigenproblem: 

`call ` `eig` `(a, lambda [, right] [,left] [,overwrite_a] [,err])`

For the generalized eigenproblem: 

`call ` `eig` `(a, b, lambda [, right] [, left] [, overwrite_a] [, overwrite_b] [, err])

### Arguments

`a` : `real` or `complex` square array containing the coefficient matrix. If `overwrite_a=.false.`, it is an `intent(in)` argument. Otherwise, it is an `intent(inout)` argument and is destroyed by the call. 

`b`: `real` or `complex` square array containing the second coefficient matrix. If `overwrite_b=.false.`, it is an `intent(in)` argument.  Otherwise, it is an `intent(inout)` argument and is destroyed by the call. 

`lambda`: Shall be a `complex` or `real` rank-1 array of the same kind as `a`, containing the eigenvalues, or their `real` component only. It is an `intent(out)` argument.

`right` (optional): Shall be a `complex` rank-2 array of the same size and kind as `a`, containing the right eigenvectors of `a`. It is an `intent(out)` argument.

`left` (optional): Shall be a `complex` rank-2 array of the same size and kind as `a`, containing the left eigenvectors of `a`. It is an `intent(out)` argument.

`overwrite_a` (optional): Shall be an input logical flag. if `.true.`, input matrix `a` will be used as temporary storage and overwritten. This avoids internal data allocation. This is an `intent(in)` argument.

`overwrite_b` (optional): Shall be an input logical flag. If `.true.`, input matrix `b` will be used as temporary storage and overwritten. This avoids internal data allocation. This is an `intent(in)` argument.

`err` (optional): Shall be a `type(linalg_state_type)` value. This is an `intent(out)` argument.

### Return value

Raises `LINALG_ERROR` if the calculation did not converge.
Raises `LINALG_VALUE_ERROR` if any matrix or arrays have invalid/incompatible sizes.
Raises `LINALG_VALUE_ERROR` if the `real` component is only requested, but the eigenvalues have non-trivial imaginary parts.
If `err` is not present, exceptions trigger an `error stop`.

### Example

```{literalinclude} ../../example/linalg/example_eig.f90
:language: fortran

`eigh` - Eigenvalues and Eigenvectors of a Real symmetric or Complex Hermitian Square Matrix

Status

Stable

Description

This subroutine computes the solution to the eigendecomposition ( A \cdot \bar{v} - \lambda \cdot \bar{v} ), where ( A ) is a square, full-rank, real symmetric ( A = A^T ) or complex Hermitian ( A = A^H ) matrix.

Result array lambda returns the real eigenvalues of ( A ). The user can request the orthogonal eigenvectors to be returned: on output vectors may contain the matrix of eigenvectors, returned as a column.

Normally, only the lower triangular part of ( A ) is accessed. On input, logical flag upper_a allows the user to request what triangular part of the matrix should be used.

The solver is based on LAPACK’s *SYEV and *HEEV backends.

Syntax

call eigh (a, lambda [, vectors] [, upper_a] [, overwrite_a] [,err])

Arguments

a : real or complex square array containing the coefficient matrix. It is an intent(in) argument. If overwrite_a=.true., it is an intent(inout) argument and is destroyed by the call.

lambda: Shall be a complex rank-1 array of the same precision as a, containing the eigenvalues. It is an intent(out) argument.

vectors (optional): Shall be a rank-2 array of the same type, size and kind as a, containing the eigenvectors of a. It is an intent(out) argument.

upper_a (optional): Shall be an input logical flag. If .true., the upper triangular part of a will be accessed. Otherwise, the lower triangular part will be accessed. It is an intent(in) argument.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument.

Return value

Raises LINALG_ERROR if the calculation did not converge. Raises LINALG_VALUE_ERROR if any matrix or arrays have invalid/incompatible sizes. If err is not present, exceptions trigger an error stop.

Example

```{literalinclude} ../../example/linalg/example_eigh.f90 :language: fortran

## `eigvals` - Eigenvalues of a Square Matrix

### Status

Stable

### Description

This function computes the eigenvalues for either a standard or generalized eigenproblem:

- **Standard eigenproblem**: \( A \cdot \bar{v} - \lambda \cdot \bar{v} \), where \( A \) is a square, full-rank `real` or `complex` matrix.
- **Generalized eigenproblem**: \( A \cdot \bar{v} - \lambda \cdot B \cdot \bar{v} \), where \( B \) is a square matrix with the same type and kind as \( A \).

The eigenvalues are stored in the result array `lambda`, which is `complex` (even for real input matrices).  
The solver uses LAPACK's `*GEEV` and `*GGEV` backends for the standard and generalized problems, respectively.

### Syntax

For the standard eigenproblem:

`lambda = ` `eigvals` `(a [, err])`

For the generalized eigenproblem:

`lambda = ` `eigvals` `(a, b [, err])`

### Arguments

`a`:  
Shall be a `real` or `complex` square array containing the coefficient matrix. It is an `intent(in)` argument.

`b` (optional):  
Shall be a `real` or `complex` square array containing the second coefficient matrix for the generalized problem. It is an `intent(in)` argument.

`err` (optional):  
Shall be a `type(linalg_state_type)` value. This is an `intent(out)` argument.

### Return Value

Returns a `complex` rank-1 array containing the eigenvalues of the problem.  

Raises `LINALG_ERROR` if the calculation did not converge.  
Raises `LINALG_VALUE_ERROR` if any matrix or arrays have invalid/incompatible sizes.  
If `err` is not present, exceptions trigger an `error stop`.

### Example

```{literalinclude} ../../example/linalg/example_eigvals.f90
:language: fortran

`eigvalsh` - Eigenvalues of a Real Symmetric or Complex Hermitian Square Matrix

Status

Stable

Description

This function returns the eigenvalues to matrix ( A ): a where ( A ) is a square, full-rank, real symmetric ( A = A^T ) or complex Hermitian ( A = A^H ) matrix. The eigenvalues are solutions to the eigenproblem ( A \cdot \bar{v} - \lambda \cdot \bar{v} ).

Result array lambda is real, and returns the eigenvalues of ( A ). The solver is based on LAPACK’s *SYEV and *HEEV backends.

Syntax

lambda = eigvalsh (a, [, upper_a] [,err])

Arguments

a : real or complex square array containing the coefficient matrix. It is an intent(in) argument.

upper_a (optional): Shall be an input logical flag. If .true., the upper triangular part of a will be used accessed. Otherwise, the lower triangular part will be accessed (default). It is an intent(in) argument.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument.

Return value

Returns a real array containing the eigenvalues of a.

Example

```{literalinclude} ../../example/linalg/example_eigvalsh.f90 :language: fortran

## `svd` - Compute the singular value decomposition of a rank-2 array (matrix).

### Status

Stable

### Description

This subroutine computes the singular value decomposition of a `real` or `complex` rank-2 array (matrix) \( A = U \cdot S \cdot \V^T \).
The solver is based on LAPACK's `*GESDD` backends.

Result vector `s` returns the array of singular values on the diagonal of \( S \). 
If requested, `u` contains the left singular vectors, as columns of \( U \).
If requested, `vt` contains the right singular vectors, as rows of \( V^T \).
 
### Syntax

`call ` `svd` `(a, s, [, u, vt, overwrite_a, full_matrices, err])`

### Class

Subroutine

### Arguments

`a`: Shall be a rank-2 `real` or `complex` array containing the coefficient matrix of size `[m,n]`. It is an `intent(inout)` argument, but returns unchanged unless `overwrite_a=.true.`.

`s`: Shall be a rank-1 `real` array, returning the list of `k = min(m,n)` singular values. It is an `intent(out)` argument. 

`u` (optional): Shall be a rank-2 array of same kind as `a`, returning the left singular vectors of `a` as columns. Its size should be `[m,m]` unless `full_matrices=.false.`, in which case, it can be `[m,min(m,n)]`. It is an `intent(out)` argument.

`vt` (optional): Shall be a rank-2 array of same kind as `a`, returning the right singular vectors of `a` as rows. Its size should be `[n,n]` unless `full_matrices=.false.`, in which case, it can be `[min(m,n),n]`. It is an `intent(out)` argument.

`overwrite_a` (optional): Shall be an input `logical` flag. If `.true.`, input matrix `A` will be used as temporary storage and overwritten. This avoids internal data allocation. By default, `overwrite_a=.false.`. It is an `intent(in)` argument.

`full_matrices` (optional): Shall be an input `logical` flag. If `.true.` (default), matrices `u` and `vt` shall be full-sized. Otherwise, their secondary dimension can be resized to `min(m,n)`. See `u`, `v` for details.

`err` (optional): Shall be a `type(linalg_state_type)` value. This is an `intent(out)` argument.

### Return values

Returns an array `s` that contains the list of singular values of matrix `a`.
If requested, returns a rank-2 array `u` that contains the left singular vectors of `a` along its columns.
If requested, returns a rank-2 array `vt` that contains the right singular vectors of `a` along its rows.

Raises `LINALG_ERROR` if the underlying Singular Value Decomposition process did not converge.
Raises `LINALG_VALUE_ERROR` if the matrix or any of the output arrays invalid/incompatible sizes.
Exceptions trigger an `error stop`, unless argument `err` is present.

### Example

```{literalinclude} ../../example/linalg/example_svd.f90
:language: fortran
:encoding: iso-8859-1

`svdvals` - Compute the singular values of a rank-2 array (matrix).

Status

Stable

Description

This subroutine computes the singular values of a real or complex rank-2 array (matrix) from its singular value decomposition ( A = U \cdot S \cdot \V^T ). The solver is based on LAPACK’s *GESDD backends.

Result vector s returns the array of singular values on the diagonal of ( S ).

Syntax

s = svdvals (a [, err])

Arguments

a: Shall be a rank-2 real or complex array containing the coefficient matrix of size [m,n]. It is an intent(in) argument.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument.

Return values

Returns an array s that contains the list of singular values of matrix a.

Raises LINALG_ERROR if the underlying Singular Value Decomposition process did not converge. Raises LINALG_VALUE_ERROR if the matrix or any of the output arrays invalid/incompatible sizes. Exceptions trigger an error stop, unless argument err is present.

Example

```{literalinclude} ../../example/linalg/example_svdvals.f90 :language: fortran

## `cholesky` - Compute the Cholesky factorization of a rank-2 square array (matrix)

### Status

Experimental

### Description

This subroutine computes the Cholesky factorization of a `real` or `complex` rank-2 square array (matrix), 
\( A = L \cdot L^T \), or \( A = U^T \cdot U \). \( A \) is symmetric or complex Hermitian, and \( L \), 
\( U \) are lower- or upper-triangular, respectively. 
The solver is based on LAPACK's `*POTRF` backends.
 
### Syntax

`call ` `cholesky` `(a, c, lower [, other_zeroed] [, err])`

### Class
Subroutine

### Arguments

`a`: Shall be a rank-2 square `real` or `complex` array containing the coefficient matrix of size `[n,n]`. It is an `intent(inout)` argument, but returns unchanged if the argument `c` is present.

`c` (optional): Shall be a rank-2 square `real` or `complex` of the same size and kind as `a`. It is an `intent(out)` argument, that returns the triangular Cholesky matrix `L` or `U`.

`lower`: Shall be an input `logical` flag. If `.true.`, the lower triangular decomposition \( A = L \cdot L^T \) will be performed. If `.false.`, the upper decomposition \( A = U^T \cdot U \) will be performed.

`other_zeroed` (optional): Shall be an input `logical` flag. If `.true.`, the unused part of the output matrix will contain zeroes. Otherwise, it will not be accessed. This saves cpu time. By default, `other_zeroed=.true.`. It is an `intent(in)` argument.

`err` (optional): Shall be a `type(linalg_state_type)` value. It is an `intent(out)` argument.

### Return values

The factorized matrix is returned in-place overwriting `a` if no other arguments are provided.
Otherwise, it can be provided as a second argument `c`. In this case, `a` is not overwritten.
The `logical` flag `lower` determines whether the lower- or the upper-triangular factorization should be performed. 

Results are returned on the applicable triangular region of the output matrix, while the unused triangular region 
is filled by zeroes by default. Optional argument `other_zeroed`, if `.false.` allows the expert user to avoid zeroing the unused part; 
however, in this case, the unused region of the matrix is not accessed and will usually contain invalid values. 

Raises `LINALG_ERROR` if the underlying process did not converge.
Raises `LINALG_VALUE_ERROR` if the matrix or any of the output arrays invalid/incompatible sizes.
Exceptions trigger an `error stop`, unless argument `err` is present.

### Example

```{literalinclude} ../../example/linalg/example_cholesky.f90
:language: fortran

`chol` - Compute the Cholesky factorization of a rank-2 square array (matrix)

Status

Experimental

Description

This is a pure functional interface for the Cholesky factorization of a real or complex rank-2 square array (matrix) computed as ( A = L \cdot L^T ), or ( A = U^T \cdot U ). ( A ) is symmetric or complex Hermitian, and ( L ), ( U ) are lower- or upper-triangular, respectively. The solver is based on LAPACK’s *POTRF backends.

Result matrix c has the same size and kind as a, and returns the lower or upper triangular factor.

Syntax

c = chol (a, lower [, other_zeroed])

Arguments

a: Shall be a rank-2 square real or complex array containing the coefficient matrix of size [n,n]. It is an intent(inout) argument, but returns unchanged if argument c is present.

lower: Shall be an input logical flag. If .true., the lower triangular decomposition ( A = L \cdot L^T ) will be performed. If .false., the upper decomposition ( A = U^T \cdot U ) will be performed.

other_zeroed (optional): Shall be an input logical flag. If .true., the unused part of the output matrix will contain zeroes. Otherwise, it will not be accessed. This saves cpu time. By default, other_zeroed=.true.. It is an intent(in) argument.

Return values

Returns a rank-2 array c of the same size and kind as a, that contains the triangular Cholesky matrix L or U.

Raises LINALG_ERROR if the underlying process did not converge. Raises LINALG_VALUE_ERROR if the matrix or any of the output arrays invalid/incompatible sizes. Exceptions trigger an error stop, unless argument err is present.

Example

```{literalinclude} ../../example/linalg/example_chol.f90 :language: fortran

## `.inv.` - Inverse operator of a square matrix

### Status

Stable

### Description

This operator returns the inverse of a `real` or `complex` square matrix \( A \).
The inverse \( A^{-1} \) is defined such that \( A \cdot A^{-1} = A^{-1} \cdot A = I_n \).

This interface is equivalent to the function  `inv`. 

### Syntax

`b = ` `operator(.inv.)` `a`

### Arguments

`a`: Shall be a rank-2 square array of any `real` or `complex` kinds. It is an `intent(in)` argument.

### Return value

Returns a rank-2 square array with the same type, kind and rank as `a`, that contains the inverse of `a`.

If an exception occurred on input errors, or singular matrix, `NaN`s will be returned.
For fine-grained error control in case of singular matrices prefer the `subroutine` and the `function`
interfaces.

### Example

```{literalinclude} ../../example/linalg/example_inverse_operator.f90
:language: fortran

`invert` - Inversion of a square matrix

Status

Stable

Description

This subroutine inverts a square real or complex matrix in-place. The inverse ( A^{-1} ) is defined such that ( A \cdot A^{-1} = A^{-1} \cdot A = I_n ).

On return, the input matrix a is replaced by its inverse. The solver is based on LAPACK’s *GETRF and *GETRI backends.

Syntax

call invert (a, [,inva] [, pivot] [, err])

Arguments

a: Shall be a rank-2, square, real or complex array containing the coefficient matrix. If inva is provided, it is an intent(in) argument. If inva is not provided, it is an intent(inout) argument: on output, it is replaced by the inverse of a.

inva (optional): Shall be a rank-2, square, real or complex array with the same size, and kind as a. On output, it contains the inverse of a.

pivot (optional): Shall be a rank-1 array of the same kind and matrix dimension as a, that contains the diagonal pivot indices on return. It is an intent(inout) argument.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument.

Return value

Computes the inverse of the matrix ( A ), (A^{-1}\, and returns it either in ( A ) or in another matrix.

Raises LINALG_ERROR if the matrix is singular or has invalid size. Raises LINALG_VALUE_ERROR if inva and a do not have the same size. If err is not present, exceptions trigger an error stop.

Example

```{literalinclude} ../../example/linalg/example_inverse_inplace.f90 :language: fortran

```{literalinclude} ../../example/linalg/example_inverse_subroutine.f90
:language: fortran

`inv` - Inverse of a square matrix.

Status

Stable

Description

This function returns the inverse of a square real or complex matrix in-place. The inverse, ( A^{-1} ), is defined such that ( A \cdot A^{-1} = A^{-1} \cdot A = I_n ).

The solver is based on LAPACK’s *GETRF and *GETRI backends.

Syntax

b inv (a, [, err])

Arguments

a: Shall be a rank-2, square, real or complex array containing the coefficient matrix. It is an intent(inout) argument.

err (optional): Shall be a type(linalg_state_type) value. It is an intent(out) argument.

Return value

Returns an array value of the same type, kind and rank as a, that contains the inverse matrix (A^{-1}).

Raises LINALG_ERROR if the matrix is singular or has invalid size. If err is not present, exceptions trigger an error stop.

Example

```{literalinclude} ../../example/linalg/example_inverse_function.f90 :language: fortran

## `pinv` - Moore-Penrose pseudo-inverse of a matrix

### Status

Experimental

### Description

This function computes the Moore-Penrose pseudo-inverse of a `real` or `complex` matrix.  
The pseudo-inverse, \( A^{+} \), generalizes the matrix inverse and satisfies the conditions:
- \( A \cdot A^{+} \cdot A = A \)
- \( A^{+} \cdot A \cdot A^{+} = A^{+} \)
- \( (A \cdot A^{+})^T = A \cdot A^{+} \)
- \( (A^{+} \cdot A)^T = A^{+} \cdot A \)

The computation is based on singular value decomposition (SVD). Singular values below a relative 
tolerance threshold \( \text{rtol} \cdot \sigma_{\max} \), where \( \sigma_{\max} \) is the largest 
singular value, are treated as zero.

### Syntax

`b =` `pinv` `(a, [, rtol, err])`

### Arguments

`a`: Shall be a rank-2, `real` or `complex` array of shape `[m, n]` containing the coefficient matrix. 
It is an `intent(in)` argument.

`rtol` (optional): Shall be a scalar `real` value specifying the relative tolerance for singular value cutoff.  
If `rtol` is not provided, the default relative tolerance is \( \text{rtol} = \text{max}(m, n) \cdot \epsilon \),  
where \( \epsilon \) is the machine precision for the element type of `a`. It is an `intent(in)` argument.

`err` (optional): Shall be a `type(linalg_state_type)` value. It is an `intent(out)` argument.

### Return value

Returns an array value of the same type, kind, and rank as `a` with shape `[n, m]`, that contains the pseudo-inverse matrix \( A^{+} \).

Raises `LINALG_ERROR` if the underlying SVD did not converge.
Raises `LINALG_VALUE_ERROR` if `a` has invalid size.
If `err` is not present, exceptions trigger an `error stop`.

### Example

```{literalinclude} ../../example/linalg/example_pseudoinverse.f90
:language: fortran

`pseudoinvert` - Moore-Penrose pseudo-inverse of a matrix

Status

Experimental

Description

This subroutine computes the Moore-Penrose pseudo-inverse of a real or complex matrix. The pseudo-inverse ( A^{+} ) is a generalization of the matrix inverse and satisfies the following properties:

( A \cdot A^{+} \cdot A = A )
( A^{+} \cdot A \cdot A^{+} = A^{+} )
( (A \cdot A^{+})^T = A \cdot A^{+} )
( (A^{+} \cdot A)^T = A^{+} \cdot A )

The computation is based on singular value decomposition (SVD). Singular values below a relative tolerance threshold ( \text{rtol} \cdot \sigma_{\max} ), where ( \sigma_{\max} ) is the largest singular value, are treated as zero.

On return, matrix pinva [n, m] will store the pseudo-inverse of a [m, n].

Syntax

call pseudoinvert (a, pinva [, rtol] [, err])

Arguments

a: Shall be a rank-2, real or complex array containing the coefficient matrix.
It is an intent(in) argument.

pinva: Shall be a rank-2 array of the same kind as a, and size equal to that of transpose(a).
On output, it contains the Moore-Penrose pseudo-inverse of a.

rtol (optional): Shall be a scalar real value specifying the relative tolerance for singular value cutoff.
If not provided, the default threshold is ( \text{max}(m, n) \cdot \epsilon ), where ( \epsilon ) is the machine precision for the element type of a.

err (optional): Shall be a type(linalg_state_type) value. It is an intent(out) argument.

Return value

Computes the Moore-Penrose pseudo-inverse of the matrix ( A ), ( A^{+} ), and returns it in matrix pinva.

Raises LINALG_ERROR if the underlying SVD did not converge. Raises LINALG_VALUE_ERROR if pinva and a have degenerate or incompatible sizes. If err is not present, exceptions trigger an error stop.

Example

```{literalinclude} ../../example/linalg/example_pseudoinverse.f90 :language: fortran

## `.pinv.` - Moore-Penrose Pseudo-Inverse operator

### Status

Experimental

### Description

This operator returns the Moore-Penrose pseudo-inverse of a `real` or `complex` matrix \( A \).
The pseudo-inverse \( A^{+} \) is computed using Singular Value Decomposition (SVD), and singular values 
below a given threshold are treated as zero.

This interface is equivalent to the function `pinv`.

### Syntax

`b = ` `operator(.pinv.)` `a`

### Arguments

`a`: Shall be a rank-2 array of any `real` or `complex` kinds, with arbitrary dimensions \( m \times n \). It is an `intent(in)` argument.

### Return value

Returns a rank-2 array with the same type, kind, and rank as `a`, that contains the Moore-Penrose pseudo-inverse of `a`.

If an exception occurs, or if the input matrix is degenerate (e.g., rank-deficient), the returned matrix will contain `NaN`s.
For more detailed error handling, it is recommended to use the `subroutine` or `function` interfaces.

### Example

```{literalinclude} ../../example/linalg/example_pseudoinverse.f90
:language: fortran

`get_norm` - Computes the vector norm of a generic-rank array.

Status

Experimental

Description

This pure subroutine interface computes one of several vector norms of real or complex array ( A ), depending on the value of the order input argument. ( A ) may be an array of any rank.

Result nrm returns a real, scalar norm value for the whole array; if dim is specified, nrm is a rank n-1 array with the same shape as (A ) and dimension dim dropped, containing all norms evaluated along dim.

Syntax

call get_norm (a, nrm, order, [, dim, err])

Arguments

a: Shall be a rank-n real or complex array containing the data. It is an intent(in) argument.

nrm: if dim is absent, shall be a scalar with the norm evaluated over all the elements of the array. Otherwise, an array of rank n-1, and a shape similar to that of a with dimension dim dropped.

order: 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
`-huge(0)`	`'-inf', '-Inf'`	Minimum absolute value ( \min_i{ \left	a_i\right	} )
`1`	`'1'`	1-norm ( \sum_i{ \left	a_i\right	} )
`2`	`'2'`	Euclidean norm ( \sqrt{\sum_i{ a_i^2 }} )
`>=3`	`'3','4',...`	p-norm ( \left( \sum_i{ \left	a_i\right	^p }\right) ^{1/p} )
`huge(0)`	`'inf', 'Inf'`	Maximum absolute value ( \max_i{ \left	a_i\right	} )

dim (optional): Shall be a scalar integer value with a value in the range from 1 to n, where n is the rank of the array. It is an intent(in) argument.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument. If err is not present, the function is pure.

Return value

By default, the return value nrm is a scalar, and contains the norm as evaluated over all elements of the generic-rank array ( A ). If the optional dim argument is present, nrm is a rank n-1 array with the same shape as ( A ) except for dimension dim, that is collapsed. Each element of nrm contains the 1D norm of the elements of ( A ), evaluated along dimension dim only.

Raises LINALG_ERROR if the requested norm type is invalid. Raises LINALG_VALUE_ERROR if any of the arguments has an invalid size. If err is not present, exceptions trigger an error stop.

Example

```{literalinclude} ../../example/linalg/example_get_norm.f90 :language: fortran

## `norm` - Computes the vector norm of a generic-rank array.

### Status

Experimental

### Description

This function computes one of several vector norms of `real` or `complex` array \( A \), depending on 
the value of the `order` input argument. \( A \) may be an array of any rank. 

### Syntax

`x = ` `norm` `(a, order, [, dim, err])`

### Arguments

`a`: Shall be a rank-n `real` or `complex` array containing the data. It is an `intent(in)` argument.

`order`: 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                                               |
|------------------|------------------|---------------------------------------------------------|
| `-huge(0)`       | `'-inf', '-Inf'` | Minimum absolute value \( \min_i{ \left|a_i\right| } \) |
| `1`              | `'1'`            | 1-norm \( \sum_i{ \left|a_i\right| } \)                 |
| `2`              | `'2'`            | Euclidean norm \( \sqrt{\sum_i{ a_i^2 }} \)             |
| `>=3`            | `'3','4',...`    | p-norm \( \left( \sum_i{ \left|a_i\right|^p }\right) ^{1/p} \) |
| `huge(0)`        | `'inf', 'Inf'`   | Maximum absolute value \( \max_i{ \left|a_i\right| } \) |

`dim` (optional): Shall be a scalar `integer` value with a value in the range from `1` to `n`, where `n` is the rank of the array. It is an `intent(in)` argument.

`err` (optional): Shall be a `type(linalg_state_type)` value. This is an `intent(out)` argument. If `err` is not present, the function is `pure`. 

### Return value

By default, the return value `x` is a scalar, and contains the norm as evaluated over all elements of the generic-rank array \( A \). 
If the optional `dim` argument is present, `x` is a rank `n-1` array with the same shape as \( A \) except 
for dimension `dim`, that is dropped. Each element of `x` contains the 1D norm of the elements of \( A \), 
evaluated along dimension `dim` only.

Raises `LINALG_ERROR` if the requested norm type is invalid.
Raises `LINALG_VALUE_ERROR` if any of the arguments has an invalid size.
If `err` is not present, exceptions trigger an `error stop`.

### Example

```{literalinclude} ../../example/linalg/example_norm.f90
:language: fortran

`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 = mnorm (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

```{literalinclude} ../../example/linalg/example_mnorm.f90 :language: fortran

## `expm` - Computes the matrix exponential {#expm}

### Status

Experimental

### Description

Given a matrix \(A\), this function computes its matrix exponential \(E = \exp(A)\) using a Pade approximation.

### Syntax

`E = ` `expm` `(a [, order])`

### Arguments

`a`: Shall be a rank-2 `real` or `complex` array containing the data. It is an `intent(in)` argument.

`order` (optional): Shall be a non-negative `integer` value specifying the order of the Pade approximation. By default `order=10`. It is an `intent(in)` argument. 

### Return value

The returned array `E` contains the Pade approximation of \(\exp(A)\).

If `A` is non-square or `order` is negative, it raises a `LINALG_VALUE_ERROR`.

### Example

```{literalinclude} ../../example/linalg/example_expm.f90
:language: fortran

`matrix_exp` - Computes the matrix exponential

Status

Experimental

Description

Given a matrix (A), this function computes its matrix exponential (E = \exp(A)) using a Pade approximation.

Syntax

call matrix_exp (a [, e, order, err])

Arguments

a: Shall be a rank-2 real or complex array containing the data. If e is not passed, it is an intent(inout) argument and is overwritten on exit by the matrix exponential. If e is passed, it is an intent(in) argument and is left unchanged.

e (optional): Shall be a rank-2 real or complex array with the same dimensions as a. It is an intent(out) argument. On exit, it contains the matrix exponential of a.

order (optional): Shall be a non-negative integer value specifying the order of the Pade approximation. By default order=10. It is an intent(in) argument.

err (optional): Shall be a type(linalg_state_type) value. This is an intent(out) argument.

Return value

The returned array A (in-place) or E (out-of-place) contains the Pade approximation of (\exp(A)).

If A is non-square or order is negative, it raises a LINALG_VALUE_ERROR. If err is not present, exceptions trigger an error stop.

Example

{literalinclude} ../../example/linalg/example_matrix_exp.f90 :language: fortran

This site is open source. Improve this page.

stdlib

stdlib_linalg

BLAS and LAPACK

Status

Description

Syntax

Example

Licensing

diag - Create a diagonal array or extract the diagonal elements of an array

Status

Class

Description

Syntax

Arguments

Return value

Example

trace - Trace of a matrix

Status

Description

Syntax

Arguments

Return value

Example

kronecker_product - Computes the Kronecker product of two rank-2 arrays

Status

Description

Syntax

Arguments

Return value

Example

is_square - Checks if a matrix is square

Status

Description

Syntax

Arguments

Return value

Example

is_symmetric - Checks if a matrix is symmetric

Status

Description

Syntax

Arguments

Return value

Example

hermitian - Compute the Hermitian version of a rank-2 matrix

Status

Description

Syntax

Arguments

Return value

Example

is_triangular - Checks if a matrix is triangular

Status

Description

Syntax

Arguments

Return value

Example

solve - Solves a linear matrix equation or a linear system of equations.

Status

Description

Syntax

Arguments

Return value

Example

solve_lu - Solves a linear matrix equation or a linear system of equations (subroutine interface).

Status

Description

Syntax

Arguments

Return value

Example

solve_lower_chol - Solves a linear system using pre-computed lower Cholesky factor.

Status

Description

Syntax

Arguments

Return value

Example

lstsq - Computes the least squares solution to a linear matrix equation.

`stdlib_linalg`

`diag` - Create a diagonal array or extract the diagonal elements of an array

`trace` - Trace of a matrix

`kronecker_product` - Computes the Kronecker product of two rank-2 arrays

`is_square` - Checks if a matrix is square

`is_symmetric` - Checks if a matrix is symmetric

`hermitian` - Compute the Hermitian version of a rank-2 matrix

`is_triangular` - Checks if a matrix is triangular

`solve` - Solves a linear matrix equation or a linear system of equations.

`solve_lu` - Solves a linear matrix equation or a linear system of equations (subroutine interface).

`solve_lower_chol` - Solves a linear system using pre-computed lower Cholesky factor.

`lstsq` - Computes the least squares solution to a linear matrix equation.

`lstsq_space` - Compute internal working space requirements for the least squares solver.

`constrained_lstsq` - Compute the solution of the equality-constrained least-squares problem

`constrained_lstsq_space` - Compute internal workspace requirements for the constrained least-square solver

`weighted_lstsq` - Computes the weighted least squares solution to a linear matrix equation

`det` - Computes the determinant of a square matrix

`qr` - Compute the QR factorization of a matrix

`qr_space` - Compute internal working space requirements for the QR factorization.

`schur` - Compute the Schur decomposition of a matrix

`eigh` - Eigenvalues and Eigenvectors of a Real symmetric or Complex Hermitian Square Matrix

`eigvalsh` - Eigenvalues of a Real Symmetric or Complex Hermitian Square Matrix