-
Notifications
You must be signed in to change notification settings - Fork 441
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
Adding Dynamic Mode Decomposition (DMD) into LAPACK #736
Conversation
Thanks @dbielich. Note that this code will stress a lot the SVD solvers, so this is good. Hopefully we won't find any issue in our SVD solvers, but running the testing of this DMD code will stress the SVD codes a lot. This is a good thing. |
@dbielich this implementation looks incredible. But, I would like to suggest that the interface should also allow the user to provide an additional parameter |
Hi @erichson, I think this is what the input parameter NRNK attempts to achieve. See below. It would be great if you can follow up. (Yes or no?) Cheers, @langou.
|
@langou fantastic! This will do the trick. |
Codecov ReportPatch and project coverage have no change.
Additional details and impacted files@@ Coverage Diff @@
## master #736 +/- ##
=========================================
Coverage 0.00% 0.00%
=========================================
Files 1908 1916 +8
Lines 186954 188509 +1555
=========================================
- Misses 186954 188509 +1555
☔ View full report in Codecov by Sentry. |
… overall addition is complete
Failure occurs in AppVeyor build and seems unrelated to the changes in this PR: 40/102 Test #41: LAPACK-xeigtstd_dec_in ...........***Failed 0.49 sec 98% tests passed, 2 tests failed out of 102 |
You are correct! Those failures are unrelated to your changes. Thanks for the clarification. |
(Sorry for the accidental closing - inadvertently stacked laptop on keyboard) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This is an amazing PR. Thank you @dbielich for contributing this!
Hi @dbielich, this patch will solve a couple of typos in your branch: diff --git a/.github/workflows/cmake.yml b/.github/workflows/cmake.yml
index 4c4b87f26..20cf92bde 100644
--- a/.github/workflows/cmake.yml
+++ b/.github/workflows/cmake.yml
@@ -138,11 +138,19 @@ jobs:
- name: Checkout LAPACK
uses: actions/checkout@8e5e7e5ab8b370d6c329ec480221332ada57f0ab # v3.5.2
+ - name: Install ninja-build tool
+ uses: seanmiddleditch/gha-setup-ninja@16b940825621068d98711680b6c3ff92201f8fc0 # v3
+
+ - name: Install Basics
+ run: |
+ sudo apt update
+ sudo apt install -y cmake gfortran
+
- name: Configure CMake
# Configure CMake in a 'build' subdirectory. `CMAKE_BUILD_TYPE` is only required if you are using a single-configuration generator such as make.
# See https://cmake.org/cmake/help/latest/variable/CMAKE_BUILD_TYPE.html?highlight=cmake_build_type
run: >
- cmake -B build
+ cmake -B build -G Ninja
-D CMAKE_BUILD_TYPE=${{env.BUILD_TYPE}}
-D CMAKE_INSTALL_PREFIX=${{github.workspace}}/lapack_install
-D CBLAS:BOOL=ON
@@ -152,7 +160,7 @@ jobs:
-D BUILD_SHARED_LIBS:BOOL=ON
- name: Install
- run: cmake --build build --target install -j2
+ run: cmake --build build --target install
- name: Coverage
run: |
diff --git a/LAPACKE/include/lapacke.h b/LAPACKE/include/lapacke.h
index bc4c5fa3b..9a9ab4753 100644
--- a/LAPACKE/include/lapacke.h
+++ b/LAPACKE/include/lapacke.h
@@ -5720,7 +5720,7 @@ lapack_int LAPACKE_sgedmd_work( int matrix_layout, char jobs, char jobz,
lapack_int ldz, float* res, float* b,
lapack_int ldb, float* w, lapack_int ldw,
float* s, lapack_int lds, float* work,
- lapack_int lwork, float* iwork,
+ lapack_int lwork, lapack_int* iwork,
lapack_int liwork );
lapack_int LAPACKE_dgedmd_work( int matrix_layout, char jobs, char jobz,
@@ -5731,7 +5731,7 @@ lapack_int LAPACKE_dgedmd_work( int matrix_layout, char jobs, char jobz,
lapack_int ldz, double* res, double* b,
lapack_int ldb, double* w, lapack_int ldw,
double* s, lapack_int lds, double* work,
- lapack_int lwork, double* iwork,
+ lapack_int lwork, lapack_int* iwork,
lapack_int liwork );
lapack_int LAPACKE_cgedmd_work( int matrix_layout, char jobs, char jobz,
@@ -5747,7 +5747,7 @@ lapack_int LAPACKE_cgedmd_work( int matrix_layout, char jobs, char jobz,
lapack_complex_float* w, lapack_int ldw,
lapack_complex_float* s, lapack_int lds,
lapack_complex_float* work, lapack_int lwork,
- lapack_complex_float* iwork, lapack_int liwork );
+ lapack_int* iwork, lapack_int liwork );
lapack_int LAPACKE_zgedmd_work( int matrix_layout, char jobs, char jobz,
char jobf, lapack_int whtsvd, lapack_int m,
@@ -5762,7 +5762,7 @@ lapack_int LAPACKE_zgedmd_work( int matrix_layout, char jobs, char jobz,
lapack_complex_double* w, lapack_int ldw,
lapack_complex_double* s, lapack_int lds,
lapack_complex_double* work, lapack_int lwork,
- lapack_complex_double* iwork, lapack_int liwork );
+ lapack_int* iwork, lapack_int liwork );
lapack_int LAPACKE_sgedmdq_work( int matrix_layout, char jobs, char jobz,
char jobr, char jobq, char jobt, char jobf,
@@ -5774,7 +5774,7 @@ lapack_int LAPACKE_sgedmdq_work( int matrix_layout, char jobs, char jobz,
lapack_int ldz, float* res, float* b,
lapack_int ldb, float* v, lapack_int ldv,
float* s, lapack_int lds, float* work,
- lapack_int lwork, float* iwork,
+ lapack_int lwork, lapack_int* iwork,
lapack_int liwork );
lapack_int LAPACKE_dgedmdq_work( int matrix_layout, char jobs, char jobz,
@@ -5787,7 +5787,7 @@ lapack_int LAPACKE_dgedmdq_work( int matrix_layout, char jobs, char jobz,
lapack_int ldz, double* res, double* b,
lapack_int ldb, double* v, lapack_int ldv,
double* s, lapack_int lds, double* work,
- lapack_int lwork, double* iwork,
+ lapack_int lwork, lapack_int* iwork,
lapack_int liwork );
lapack_int LAPACKE_cgedmdq_work( int matrix_layout, char jobs, char jobz,
@@ -5805,7 +5805,7 @@ lapack_int LAPACKE_cgedmdq_work( int matrix_layout, char jobs, char jobz,
lapack_complex_float* v, lapack_int ldv,
lapack_complex_float* s, lapack_int lds,
lapack_complex_float* work, lapack_int lwork,
- lapack_complex_float* iwork,
+ lapack_int* iwork,
lapack_int liwork );
lapack_int LAPACKE_zgedmdq_work( int matrix_layout, char jobs, char jobz,
@@ -5823,7 +5823,7 @@ lapack_int LAPACKE_zgedmdq_work( int matrix_layout, char jobs, char jobz,
lapack_complex_double* v, lapack_int ldv,
lapack_complex_double* s, lapack_int lds,
lapack_complex_double* work, lapack_int lwork,
- lapack_complex_double* iwork,
+ lapack_int* iwork,
lapack_int liwork );
lapack_int LAPACKE_sgesv_work( int matrix_layout, lapack_int n, lapack_int nrhs,
diff --git a/LAPACKE/src/lapacke_cgedmdq_work.c b/LAPACKE/src/lapacke_cgedmdq_work.c
index d0f161bdb..5bdbd3f56 100644
--- a/LAPACKE/src/lapacke_cgedmdq_work.c
+++ b/LAPACKE/src/lapacke_cgedmdq_work.c
@@ -54,7 +54,7 @@ lapack_int LAPACKE_cgedmdq_work( int matrix_layout, char jobs, char jobz,
if( matrix_layout == LAPACK_COL_MAJOR ) {
/* Call LAPACK function and adjust info */
LAPACK_cgedmdq( &jobs, &jobz, &jobr, &jobq, &jobt, &jobf, &whtsvd, &m,
- &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, tol, &k, reig,
+ &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, &tol, &k, reig,
imeig, z, &ldz, res, b, &ldb, v, &ldv, s, &lds,
work, &lwork, iwork, &liwork, &info );
if( info < 0 ) {
@@ -114,7 +114,7 @@ lapack_int LAPACKE_cgedmdq_work( int matrix_layout, char jobs, char jobz,
/* Query optimal working array(s) size if requested */
if( lwork == -1 || liwork == -1 ) {
LAPACK_cgedmdq( &jobs, &jobz, &jobr, &jobq, &jobt, &jobf, &whtsvd, &m,
- &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, tol, &k, reig,
+ &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, &tol, &k, reig,
imeig, z, &ldz, res, b, &ldb, v, &ldv, s, &lds,
work, &lwork, iwork, &liwork, &info );
return (info < 0) ? (info - 1) : info;
@@ -165,7 +165,7 @@ lapack_int LAPACKE_cgedmdq_work( int matrix_layout, char jobs, char jobz,
LAPACKE_cge_trans( matrix_layout, m, n, s, lds, s_t, lds_t );
/* Call LAPACK function and adjust info */
LAPACK_cgedmdq( &jobs, &jobz, &jobr, &jobq, &jobt, &jobf, &whtsvd, &m,
- &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, tol, &k, reig,
+ &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, &tol, &k, reig,
imeig, z, &ldz, res, b, &ldb, v, &ldv, s, &lds,
work, &lwork, iwork, &liwork, &info );
if( info < 0 ) {
diff --git a/LAPACKE/src/lapacke_dgedmdq_work.c b/LAPACKE/src/lapacke_dgedmdq_work.c
index 04f2669bf..51b2a66d8 100644
--- a/LAPACKE/src/lapacke_dgedmdq_work.c
+++ b/LAPACKE/src/lapacke_dgedmdq_work.c
@@ -49,7 +49,7 @@ lapack_int LAPACKE_dgedmdq_work( int matrix_layout, char jobs, char jobz,
if( matrix_layout == LAPACK_COL_MAJOR ) {
/* Call LAPACK function and adjust info */
LAPACK_dgedmdq( &jobs, &jobz, &jobr, &jobq, &jobt, &jobf, &whtsvd, &m,
- &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, tol, &k, reig,
+ &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, &tol, &k, reig,
imeig, z, &ldz, res, b, &ldb, v, &ldv, s, &lds,
work, &lwork, iwork, &liwork, &info );
if( info < 0 ) {
@@ -109,7 +109,7 @@ lapack_int LAPACKE_dgedmdq_work( int matrix_layout, char jobs, char jobz,
/* Query optimal working array(s) size if requested */
if( lwork == -1 || liwork == -1 ) {
LAPACK_dgedmdq( &jobs, &jobz, &jobr, &jobq, &jobt, &jobf, &whtsvd, &m,
- &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, tol, &k, reig,
+ &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, &tol, &k, reig,
imeig, z, &ldz, res, b, &ldb, v, &ldv, s, &lds,
work, &lwork, iwork, &liwork, &info );
return (info < 0) ? (info - 1) : info;
@@ -160,7 +160,7 @@ lapack_int LAPACKE_dgedmdq_work( int matrix_layout, char jobs, char jobz,
LAPACKE_dge_trans( matrix_layout, m, n, s, lds, s_t, lds_t );
/* Call LAPACK function and adjust info */
LAPACK_dgedmdq( &jobs, &jobz, &jobr, &jobq, &jobt, &jobf, &whtsvd, &m,
- &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, tol, &k, reig,
+ &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, &tol, &k, reig,
imeig, z, &ldz, res, b, &ldb, v, &ldv, s, &lds,
work, &lwork, iwork, &liwork, &info );
if( info < 0 ) {
diff --git a/LAPACKE/src/lapacke_sgedmdq_work.c b/LAPACKE/src/lapacke_sgedmdq_work.c
index a9c228a45..9039898d2 100644
--- a/LAPACKE/src/lapacke_sgedmdq_work.c
+++ b/LAPACKE/src/lapacke_sgedmdq_work.c
@@ -49,7 +49,7 @@ lapack_int LAPACKE_sgedmdq_work( int matrix_layout, char jobs, char jobz,
if( matrix_layout == LAPACK_COL_MAJOR ) {
/* Call LAPACK function and adjust info */
LAPACK_sgedmdq( &jobs, &jobz, &jobr, &jobq, &jobt, &jobf, &whtsvd, &m,
- &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, tol, &k, reig,
+ &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, &tol, &k, reig,
imeig, z, &ldz, res, b, &ldb, v, &ldv, s, &lds,
work, &lwork, iwork, &liwork, &info );
if( info < 0 ) {
@@ -109,7 +109,7 @@ lapack_int LAPACKE_sgedmdq_work( int matrix_layout, char jobs, char jobz,
/* Query optimal working array(s) size if requested */
if( lwork == -1 || liwork == -1 ) {
LAPACK_sgedmdq( &jobs, &jobz, &jobr, &jobq, &jobt, &jobf, &whtsvd, &m,
- &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, tol, &k, reig,
+ &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, &tol, &k, reig,
imeig, z, &ldz, res, b, &ldb, v, &ldv, s, &lds,
work, &lwork, iwork, &liwork, &info );
return (info < 0) ? (info - 1) : info;
@@ -160,7 +160,7 @@ lapack_int LAPACKE_sgedmdq_work( int matrix_layout, char jobs, char jobz,
LAPACKE_sge_trans( matrix_layout, m, n, s, lds, s_t, lds_t );
/* Call LAPACK function and adjust info */
LAPACK_sgedmdq( &jobs, &jobz, &jobr, &jobq, &jobt, &jobf, &whtsvd, &m,
- &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, tol, &k, reig,
+ &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, &tol, &k, reig,
imeig, z, &ldz, res, b, &ldb, v, &ldv, s, &lds,
work, &lwork, iwork, &liwork, &info );
if( info < 0 ) {
diff --git a/LAPACKE/src/lapacke_zgedmdq.c b/LAPACKE/src/lapacke_zgedmdq.c
index 9b4b429b0..3648ffdf2 100644
--- a/LAPACKE/src/lapacke_zgedmdq.c
+++ b/LAPACKE/src/lapacke_zgedmdq.c
@@ -100,7 +100,7 @@ lapack_int LAPACKE_zgedmdq( int matrix_layout, char jobs, char jobz, char jobr,
info = LAPACK_WORK_MEMORY_ERROR;
goto exit_level_0;
}
- iwork = (lapack_complex_double*)LAPACKE_malloc( sizeof(lapack_complex_double) * liwork );
+ iwork = (lapack_int*)LAPACKE_malloc( sizeof(lapack_int) * liwork );
if( iwork == NULL ) {
info = LAPACK_WORK_MEMORY_ERROR;
goto exit_level_1;
diff --git a/LAPACKE/src/lapacke_zgedmdq_work.c b/LAPACKE/src/lapacke_zgedmdq_work.c
index 5aa3a4829..9afceba07 100644
--- a/LAPACKE/src/lapacke_zgedmdq_work.c
+++ b/LAPACKE/src/lapacke_zgedmdq_work.c
@@ -54,7 +54,7 @@ lapack_int LAPACKE_zgedmdq_work( int matrix_layout, char jobs, char jobz,
if( matrix_layout == LAPACK_COL_MAJOR ) {
/* Call LAPACK function and adjust info */
LAPACK_zgedmdq( &jobs, &jobz, &jobr, &jobq, &jobt, &jobf, &whtsvd, &m,
- &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, tol, &k, reig,
+ &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, &tol, &k, reig,
imeig, z, &ldz, res, b, &ldb, v, &ldv, s, &lds,
work, &lwork, iwork, &liwork, &info );
if( info < 0 ) {
@@ -114,7 +114,7 @@ lapack_int LAPACKE_zgedmdq_work( int matrix_layout, char jobs, char jobz,
/* Query optimal working array(s) size if requested */
if( lwork == -1 || liwork == -1 ) {
LAPACK_zgedmdq( &jobs, &jobz, &jobr, &jobq, &jobt, &jobf, &whtsvd, &m,
- &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, tol, &k, reig,
+ &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, &tol, &k, reig,
imeig, z, &ldz, res, b, &ldb, v, &ldv, s, &lds,
work, &lwork, iwork, &liwork, &info );
return (info < 0) ? (info - 1) : info;
@@ -165,7 +165,7 @@ lapack_int LAPACKE_zgedmdq_work( int matrix_layout, char jobs, char jobz,
LAPACKE_zge_trans( matrix_layout, m, n, s, lds, s_t, lds_t );
/* Call LAPACK function and adjust info */
LAPACK_zgedmdq( &jobs, &jobz, &jobr, &jobq, &jobt, &jobf, &whtsvd, &m,
- &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, tol, &k, reig,
+ &n, f, &ldf, x, &ldx, y, &ldy, &nrnk, &tol, &k, reig,
imeig, z, &ldz, res, b, &ldb, v, &ldv, s, &lds,
work, &lwork, iwork, &liwork, &info );
if( info < 0 ) { I couldn't make the suggestions directly to your branch and I am not sure why. |
@weslleyspereira thanks for the patch to make that quick. I have committed the changes to the routines. All except for Is this a file to add or ignore? I think ignore which is why I had left it out, but want confirmation. |
Ah, yes. Sorry for sending this one as well. This should go in another discussion I think. |
The work related to this PR is developed and coded by Zlatko Drmac, Faculty of Science, University of Zagreb. In cooperation with AIMdyn Inc., Santa Barbara, CA. Thanks for the great work!
For a detailed description of the algorithm and implementation, please refer here.
In short, XGEDMD computes the Dynamic Mode Decomposition (DMD) for a pair of data snapshot matrices. For the input matrices X and Y such that Y = A*X with an unaccessible matrix A, XGEDMD computes a certain number of Ritz pairs of A using the standard Rayleigh-Ritz extraction from a subspace of range(X) that is determined using the leading left singular vectors of X. Optionally, XGEDMD returns the residuals of the computed Ritz pairs, the information needed for a refinement of the Ritz vectors, or the eigenvectors of the Exact DMD.
Current work needed for the PR to be approved: