Prototype class EigenSystem for solving eigenvalue problems#768
Conversation
Use Apple Accelerate/veclib DGEEV (LAPACK interface) to prototype an eigenvalue solver for a general system. Provide Python binding using pybind11 and create basic tests in Python. In the Python wrapper, set `EigenSystem` to `None` when it is not built.
| namespace python | ||
| { | ||
|
|
||
| #ifdef __APPLE__ |
There was a problem hiding this comment.
This is a prototype. Fence for code that only builds with macOS/Apple.
| #ifdef __APPLE__ | ||
| WrapEigenSystem::commit(mod, "EigenSystem", "Eigen problem solver"); | ||
| #else // __APPLE__ | ||
| mod.attr("EigenSystem") = pybind11::none(); |
There was a problem hiding this comment.
Other than macOS/Apple, assign None to the symbol EigenSystem.
| * This header is only available on Apple platforms. | ||
| */ | ||
|
|
||
| #ifndef __APPLE__ |
There was a problem hiding this comment.
Only builds on macOS/Apple.
|
|
||
| public: | ||
|
|
||
| using value_type = double; |
There was a problem hiding this comment.
In this prototype, EigenSystem hard-codes to use float64 (double). It will be extended to be a template and support float32, complex64, and complex128 in the future.
|
|
||
| static std::string format_shape(array_type const & arr); | ||
|
|
||
| SimpleArray<value_type> const & m_matrix; |
There was a problem hiding this comment.
The lifecycle may not be correct with the const &. I am still considering to fix it with this or the next PR.
| np.testing.assert_allclose(A_np @ vr[:, j], wr[j] * vr[:, j], | ||
| rtol=1e-10, atol=1e-12) | ||
|
|
||
| def test_2x2_rotation_complex_conjugate_pair(self): |
There was a problem hiding this comment.
Test the rotation matrix that only has imaginary eigenvalues.
| np.testing.assert_allclose(A_np @ v_conj, lam_conj * v_conj, | ||
| rtol=1e-12, atol=1e-14) | ||
|
|
||
| def test_left_eigenvectors_satisfy_left_equation(self): |
| wr[j] * vl[:, j], | ||
| rtol=1e-10, atol=1e-12) | ||
|
|
||
| def test_rejects_non_square_and_non_2d_inputs(self): |
| ValueError, r"must be a square 2D SimpleArray"): | ||
| mm.EigenSystem(A_3d) | ||
|
|
||
| def test_accessors_before_run_are_inert(self): |
| wr = solver.wr.ndarray | ||
| self.assertTrue(np.all(np.isfinite(wr))) | ||
|
|
||
| def test_matrix_property_survives_input_gc(self): |
There was a problem hiding this comment.
matrix property lifecycle.
|
@KHLee529 @tigercosmos Please help take a look at the prototype. |
1 participant
Use Apple Accelerate/veclib DGEEV (LAPACK interface) to prototype an eigenvalue solver for a general linear system. The prototype uses double-precision floating point (float64) only.
Provide Python binding using pybind11 and create basic tests in Python. In the Python wrapper, set
EigenSystemtoNonewhen it is not built.All platforms has LAPACK, but properly set up the build system for all platforms will use some efforts. I start with Apple because it is the most popular platform among the contributors, and Apple Silicon is the most accessible platform for GPU and other interesting accelerating instructions.
Future work: