Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

linalg: Schur decomposition #892

Merged
merged 27 commits into from
Jan 3, 2025
Merged

Conversation

perazz
Copy link
Member

@perazz perazz commented Nov 20, 2024

Compute the Schur decomposition of a real or complex square matrix: $A = Z T Z^H$.

Proposed implementation

  • call schur(A, T [, Z] [, eigvals] [, overwrite_a] [, storage] [, err]) : subroutine interface
  • call schur_space(A, lwork [, err]) query internal storage size for pre-allocation.

Key facts

  • The Schur decomposition is based on LAPACK's General Schur decomposition (*GEES) and eigenvalue sorting mechanisms.
  • The output matrix $T$ is upper-triangular for complex matrices and quasi-upper-triangular for real matrices, with possible $2 \times 2$ blocks on the diagonal.
  • Optionally, the user may request eigenvalues of $A$, which are the diagonal elements of $T$.
  • The function supports an overwrite_a flag (default: .false.). If .true., and storage is user-provided, the input matrix a will be overwritten, avoiding memory allocation.

Progress

  • base implementation
  • tests
  • documentation
  • submodule
  • examples
  • workspace size evaluation

Prior art

  • Numpy: not available
  • Scipy: T, Z = schur(a, output='real', lwork=None, overwrite_a=False, sort=None, check_finite=True)

cc: @fortran-lang/stdlib @jvdp1 @jalvesz @loiseaujc

@perazz perazz marked this pull request as ready for review November 20, 2024 13:39
Copy link
Member

@jvdp1 jvdp1 left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you @perazz . It is a nice addition. I only have a few minor comments. So, it is close to be ready to be merged IMO

doc/specs/stdlib_linalg.md Show resolved Hide resolved
doc/specs/stdlib_linalg.md Outdated Show resolved Hide resolved
doc/specs/stdlib_linalg.md Show resolved Hide resolved
doc/specs/stdlib_linalg.md Outdated Show resolved Hide resolved
example/linalg/example_schur_complex.f90 Outdated Show resolved Hide resolved
example/linalg/example_schur_eigvals.f90 Show resolved Hide resolved
example/linalg/example_schur_eigvals.f90 Outdated Show resolved Hide resolved
example/linalg/example_schur_real.f90 Show resolved Hide resolved
src/stdlib_linalg.fypp Show resolved Hide resolved
src/stdlib_linalg_schur.fypp Outdated Show resolved Hide resolved
@perazz
Copy link
Member Author

perazz commented Dec 12, 2024

Thanks for the reviews @jvdp1, looks much better now. Regarding the eigenvalues, I'm open to changing the API

Sorry I just realized we had the approach in place for eigs. So, I implemented the same: option for real eigenvalues request; it will raise an error if any have non-zero imaginary component

Copy link
Member

@jvdp1 jvdp1 left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you @perazz . LGTM. It is ready to be merged IMO

@perazz
Copy link
Member Author

perazz commented Dec 26, 2024

Thank you @jvdp1, let's wait another couple of days, and then merge if there are no further comments.

src/stdlib_linalg_schur.fypp Outdated Show resolved Hide resolved
src/stdlib_linalg_schur.fypp Show resolved Hide resolved
src/stdlib_linalg_schur.fypp Outdated Show resolved Hide resolved
@perazz perazz merged commit 94f201d into fortran-lang:master Jan 3, 2025
15 checks passed
@perazz perazz deleted the linalg_schur branch January 3, 2025 07:21
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants