form_qr Interface

public interface form_qr

Module Procedures

private subroutine form_qr_no_pivot(r, tau, q, work, olwork, err)

Forms the full M-by-M orthogonal matrix from the elementary reflectors returned by the base QR factorization algorithm.

Arguments

Type IntentOptional Attributes Name
real(kind=real64), intent(inout), dimension(:,:) :: r

On input, an M-by-N matrix where the elements below the diagonal contain the elementary reflectors generated from the QR factorization. On and above the diagonal, the matrix contains the matrix . On output, the elements below the diagonal are zeroed such that the remaining matrix is simply the M-by-N matrix .

real(kind=real64), intent(in), dimension(:) :: tau

A MIN(M, N)-element array containing the scalar factors of each elementary reflector defined in .

real(kind=real64), intent(out), dimension(:,:) :: q

An M-by-M matrix where the full orthogonal matrix will be written. In the event that M > N, may be supplied as M-by-N, and therefore only return the useful submatrix as the factorization can be written as .

real(kind=real64), intent(out), optional, target, dimension(:) :: work

An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork.

integer(kind=int32), intent(out), optional :: olwork

An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work, and returns without performing any actual calculations.

class(errors), intent(inout), optional, target :: err

The error object to be updated.

private subroutine form_qr_no_pivot_cmplx(r, tau, q, work, olwork, err)

Forms the full M-by-M orthogonal matrix from the elementary reflectors returned by the base QR factorization algorithm.

Arguments

Type IntentOptional Attributes Name
complex(kind=real64), intent(inout), dimension(:,:) :: r

On input, an M-by-N matrix where the elements below the diagonal contain the elementary reflectors generated from the QR factorization. On and above the diagonal, the matrix contains the matrix . On output, the elements below the diagonal are zeroed such that the remaining matrix is simply the M-by-N matrix .

complex(kind=real64), intent(in), dimension(:) :: tau

A MIN(M, N)-element array containing the scalar factors of each elementary reflector defined in .

complex(kind=real64), intent(out), dimension(:,:) :: q

An M-by-M matrix where the full orthogonal matrix will be written. In the event that M > N, may be supplied as M-by-N, and therefore only return the useful submatrix as the factorization can be written as .

complex(kind=real64), intent(out), optional, target, dimension(:) :: work

An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork.

integer(kind=int32), intent(out), optional :: olwork

An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work, and returns without performing any actual calculations.

class(errors), intent(inout), optional, target :: err

The error object to be updated.

private subroutine form_qr_pivot(r, tau, pvt, q, p, work, olwork, err)

Forms the full M-by-M orthogonal matrix from the elementary reflectors returned by the base QR factorization algorithm.

Arguments

Type IntentOptional Attributes Name
real(kind=real64), intent(inout), dimension(:,:) :: r

On input, an M-by-N matrix where the elements below the diagonal contain the elementary reflectors generated from the QR factorization. On and above the diagonal, the matrix contains the matrix . On output, the elements below the diagonal are zeroed such that the remaining matrix is simply the M-by-N matrix .

real(kind=real64), intent(in), dimension(:) :: tau

A MIN(M, N)-element array containing the scalar factors of each elementary reflector defined in .

integer(kind=int32), intent(in), dimension(:) :: pvt

An N-element column pivot array as returned by the QR factorization.

real(kind=real64), intent(out), dimension(:,:) :: q

An M-by-M matrix where the full orthogonal matrix will be written. In the event that M > N, may be supplied as M-by-N, and therefore only return the useful submatrix as the factorization can be written as .

real(kind=real64), intent(out), dimension(:,:) :: p

An N-by-N matrix where the pivot matrix will be written.

real(kind=real64), intent(out), optional, target, dimension(:) :: work

An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork.

integer(kind=int32), intent(out), optional :: olwork

An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work, and returns without performing any actual calculations.

class(errors), intent(inout), optional, target :: err

The error object to be updated.

private module subroutine form_qr_pivot_cmplx(r, tau, pvt, q, p, work, olwork, err)

Forms the full M-by-M orthogonal matrix from the elementary reflectors returned by the base QR factorization algorithm.

Arguments

Type IntentOptional Attributes Name
complex(kind=real64), intent(inout), dimension(:,:) :: r

On input, an M-by-N matrix where the elements below the diagonal contain the elementary reflectors generated from the QR factorization. On and above the diagonal, the matrix contains the matrix . On output, the elements below the diagonal are zeroed such that the remaining matrix is simply the M-by-N matrix .

complex(kind=real64), intent(in), dimension(:) :: tau

A MIN(M, N)-element array containing the scalar factors of each elementary reflector defined in .

integer(kind=int32), intent(in), dimension(:) :: pvt

An N-element column pivot array as returned by the QR factorization.

complex(kind=real64), intent(out), dimension(:,:) :: q

An M-by-M matrix where the full orthogonal matrix will be written. In the event that M > N, may be supplied as M-by-N, and therefore only return the useful submatrix as the factorization can be written as .

complex(kind=real64), intent(out), dimension(:,:) :: p

An N-by-N matrix where the pivot matrix will be written.

complex(kind=real64), intent(out), optional, target, dimension(:) :: work

An optional input, that if provided, prevents any local memory allocation. If not provided, the memory required is allocated within. If provided, the length of the array must be at least olwork.

integer(kind=int32), intent(out), optional :: olwork

An optional output used to determine workspace size. If supplied, the routine determines the optimal size for work, and returns without performing any actual calculations.

class(errors), intent(inout), optional, target :: err

The error object to be updated.