Sparse3 #1845
Sparse3 #1845
Conversation
|
This looks great so far. Obviously, the blocker for merging is that there are tests that crash and that the documentation doesn't build cleanly. A few nits to pick:
BTW, a great way to test the sparse vectors and matrices is to wrap them for generics like |
I may need some help with this, since the tests do not crash for me (I would not push otherwise)...are there settings I need to try when I configure/build flint to see these issues? I can try to get some other platforms to test on, I'm on an M3 Mac.
Happy with the latter, but the former was intended to convey that I'm multiplying by a
Sorry, that was copied over from a previous version, yes that will change.
Sounds good, will do.
That's fine.
Sorry, those should have all been replaced with
Ok, is it alright if it just returns
Thanks, I'll dig out all of those.
Will do, once things are in a good shape for the rings they definitely should work for. :-) |
| int gr_lil_mat_solve_block_wiedemann(gr_ptr x, const gr_lil_mat_t M, gr_srcptr b, slong block_size, flint_rand_t state, gr_ctx_t ctx) | ||
|
|
||
| Solve `Mx = b` for a sparse matrix `M` and *M->c* long column vector `b`, using either Lanczos' or Wiedemann's algorithm. | ||
| Both are randomized algorithms which use a given random state machine, and thus may fail without provably no solution |
There was a problem hiding this comment.
"fail without ... no solution" is very hard to read, sort of triple negation
| Both are randomized algorithms which use a given random state machine, and thus may fail without provably no solution | ||
| (returning `GR_UNABLE`). Both have block variants which are better for large matrices, and take an extra parameter of | ||
| `block_size`. The (block) Wiedemann algorithm requires the given matrix to be square, but not symmetric: the Lanczos | ||
| algorithm requires both, but assumes the given matrix may be neither, so instead solves `M^TMx = M^Tb`. Thus, it may |
There was a problem hiding this comment.
both -> both square and symmetric
|
"without having probably no solution"? "Without probably having no
solution"? Happy for other suggestions.
…On Thu, Apr 4, 2024, 7:42 AM Dima Pasechnik ***@***.***> wrote:
***@***.**** commented on this pull request.
------------------------------
In doc/source/gr_sparse_mat.rst
<#1845 (comment)>:
> +.. function:: int gr_csr_mat_mul_mat(gr_mat_t C, const gr_csr_mat_t A, const gr_mat_t B, gr_ctx_t ctx)
+ int gr_lil_mat_mul_mat(gr_mat_t C, const gr_lil_mat_t A, const gr_mat_t B, gr_ctx_t ctx)
+
+ Set *C* equal to `A \cdot B`, i.e., perform right multiplication by *B*.
+
+
+Solving, nullvector, and nullspace computation
+--------------------------------------------------------------------------------
+
+.. function:: int gr_lil_mat_solve_lanczos(gr_ptr x, const gr_lil_mat_t M, gr_srcptr b, flint_rand_t state, gr_ctx_t ctx)
+ int gr_lil_mat_solve_block_lanczos(gr_ptr x, const gr_lil_mat_t M, gr_srcptr b, slong block_size, flint_rand_t state, gr_ctx_t ctx)
+ int gr_lil_mat_solve_wiedemann(gr_ptr x, const gr_lil_mat_t M, gr_srcptr b, gr_ctx_t ctx)
+ int gr_lil_mat_solve_block_wiedemann(gr_ptr x, const gr_lil_mat_t M, gr_srcptr b, slong block_size, flint_rand_t state, gr_ctx_t ctx)
+
+ Solve `Mx = b` for a sparse matrix `M` and *M->c* long column vector `b`, using either Lanczos' or Wiedemann's algorithm.
+ Both are randomized algorithms which use a given random state machine, and thus may fail without provably no solution
"fail without ... no solution" is very hard to read, sort of triple
negation
—
Reply to this email directly, view it on GitHub
<#1845 (review)>,
or unsubscribe
<https://github.com/notifications/unsubscribe-auth/AAKDUIIBMP53TNWV6T4E3ZLY3VRDXAVCNFSM6AAAAABE3JHKI6VHI2DSMVQWIX3LMV43YUDVNRWFEZLROVSXG5CSMV3GSZLXHMYTSOBQGEZTQNBYHA>
.
You are receiving this because you authored the thread.Message ID:
***@***.***>
|
|
Let's be more precise here. Randomized algorithms are of two types - Las Vegas and Monte Carlo. https://en.wikipedia.org/wiki/Randomized_algorithm What type do we have here? |
|
Neither, at least in general:
For a floating point represented field (which I don't think we have?), I
think it's essentially a Monte Carlo algorithm, in that it will produce an
epsilon approximation of the correct answer.
For exact (eg finite) fields, it will give the correct answer or just fail,
mostly independent of the random input it uses and dependent instead on
whether the matrix is well conditioned...the randomness is mostly so that,
for things like kernel finding, it will eventually get the whole null
space. But there are no guarantees whatsoever.
…On Thu, Apr 4, 2024, 11:11 AM Dima Pasechnik ***@***.***> wrote:
Let's be more precise here. Randomized algorithms are of two types - Las
Vegas and Monte Carlo. https://en.wikipedia.org/wiki/Randomized_algorithm
What type do we have here?
—
Reply to this email directly, view it on GitHub
<#1845 (comment)>, or
unsubscribe
<https://github.com/notifications/unsubscribe-auth/AAKDUILFSDE4CL5F6D3Y2JLY3WJURAVCNFSM6AAAAABE3JHKI6VHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDAMZXHA3TIMZWGQ>
.
You are receiving this because you authored the thread.Message ID:
***@***.***>
|
|
an |
|
In the macro ...
a_nz_idx = a_nnz-1; \
b_nz_idx = b_nnz-1; \
dest_nz_idx = new_nnz-1; \
while (a_nz_idx >= -1 && b_nz_idx >= -1 && status == GR_SUCCESS) \
{ \
if (a_nz_idx == -1 && b_nz_idx == -1) break; \
a_ind = (A_VEC)->inds[a_nz_idx]; \
b_ind = (B_VEC)->inds[b_nz_idx]; \
if (b_nz_idx == -1 || (a_nz_idx >= 0 && a_ind > b_ind)) \
...This could appear when a row of a sparse matrix is empty, which leads to the crash of some tests of sparse matrix. |
| #include "test_helpers.h" | ||
| #include "gr_sparse_mat.h" | ||
|
|
||
| #define CHECK_TEST(x, name) { if (GR_SUCCESS != (x)) { flint_printf("FAIL %s\n", (name)); flint_abort(); } } |
There was a problem hiding this comment.
Do not abort in the test function. Instead, let the main function for the test (i.e. main.c) take care of any aborts. So instead, just return some non-zero value to the main function and it will abort appropriately.
If you can, use TEST_FUNCTION_FAIL, but this has to be put into the top-level function in this file, i.e. under TEST_FUNCTION_START(my_function, state).
There was a problem hiding this comment.
Also, this definition is in each file, which leads to a lot of code duplication.
| @@ -0,0 +1,27 @@ | |||
| #include <string.h> | |||
There was a problem hiding this comment.
A lot of copyright claims are missing. It is useful to know who wrote the code, so one knows who to contact if something needs to be fixed.
| // flint_printf("\n\ncoo_mat = "); status |= gr_coo_mat_print_nz(coo_mat, ctx); flint_printf("\nnnz = %d\n", coo_mat->nnz); | ||
| status |= gr_lil_mat_set_coo_mat(lil_mat, coo_mat, ctx); | ||
| // flint_printf("\n\nlil_mat = "); status |= gr_lil_mat_print_nz(lil_mat, ctx); flint_printf("\nnnz = %d\n", lil_mat->nnz); | ||
| status |= gr_lil_mat_set(lil_mat2, lil_mat, ctx); | ||
| // flint_printf("\n\nlil_mat = "); status |= gr_lil_mat_print_nz(lil_mat2, ctx); flint_printf("\nnnz = %d\n", lil_mat2->nnz); | ||
| status |= gr_csr_mat_set_lil_mat(csr_mat, lil_mat2, ctx); | ||
| // flint_printf("\n\ncsr_mat = "); status |= gr_csr_mat_print_nz(csr_mat, ctx); flint_printf("\nnnz = %d\n", csr_mat->nnz); | ||
| status |= gr_mat_set_csr_mat(dmat, csr_mat, ctx); | ||
| // flint_printf("\n\nmat = "); status |= gr_mat_print(dmat, ctx); | ||
| status |= gr_lil_mat_set_mat(lil_mat2, dmat, ctx); | ||
| // flint_printf("\n\nlil_mat = "); status |= gr_lil_mat_print_nz(lil_mat2, ctx); flint_printf("\nnnz = %d\n", lil_mat2->nnz); | ||
| status |= gr_coo_mat_set_lil_mat(coo_mat2, lil_mat2, ctx); | ||
| // flint_printf("\n\ncoo_mat = "); status |= gr_coo_mat_print_nz(coo_mat2, ctx); flint_printf("\nnnz = %d\n", coo_mat2->nnz); |
There was a problem hiding this comment.
Is it beneficial to have theses print statement left in the code?
|
Is there any I/O-functionality? |
Good catch. I pushed a fix and the Those "Solved 0" look suspect? |
|
Yes, I'm going to finish and check in the rref and lu stuff this weekend,
and will figure out what is going on with block lanzos. Thanks for catching
the zero length bug!
…On Sat, Apr 27, 2024, 5:33 AM Fredrik Johansson ***@***.***> wrote:
In the macro GR_SPV_RFL_TEMPLATE of gr_sparse_vec/arith.c, if a_nnz or
b_nnz is zero, the inital a_nz_idx or a_nz_idx is -1 and it's impossible
to visit (A_VEC)->inds[a_nz_idx] (or B) in this case:
...
a_nz_idx = a_nnz-1; \
b_nz_idx = b_nnz-1; \
dest_nz_idx = new_nnz-1; \
while (a_nz_idx >= -1 && b_nz_idx >= -1 && status == GR_SUCCESS) \
{ \
if (a_nz_idx == -1 && b_nz_idx == -1) break; \
a_ind = (A_VEC)->inds[a_nz_idx]; \
b_ind = (B_VEC)->inds[b_nz_idx]; \
if (b_nz_idx == -1 || (a_nz_idx >= 0 && a_ind > b_ind)) \
...
This could appear when a row of a sparse matrix is empty, which leads to
the crash of some tests of sparse matrix.
Good catch. I pushed a fix and the sparse_mat tests now pass as follows:
gr_sparse_mat_init...
gr_sparse_mat_init 0.00 (PASS)
gr_sparse_mat_conversion...
gr_sparse_mat_conversion 0.22 (PASS)
gr_sparse_mat_randtest...
gr_sparse_mat_randtest 0.11 (PASS)
gr_sparse_mat_arith...
gr_sparse_mat_arith 1.79 (PASS)
gr_sparse_mat_mul...
gr_sparse_mat_mul 0.06 (PASS)
gr_sparse_mat_solve...
solving Ax = b........................PASS
Solved 0 with rref
Solved 0 with lu
Solved 100 with Lanczos
Solved 0 with block Lanczos
Found no solution for 100/100 examples
Solved 100 with Wiedemann
Solved 100 with block Wiedemann
gr_sparse_mat_solve 2.73 (PASS)
All tests passed for gr_sparse_mat.
Those "Solved 0" look suspect?
—
Reply to this email directly, view it on GitHub
<#1845 (comment)>, or
unsubscribe
<https://github.com/notifications/unsubscribe-auth/AAKDUIPI6ALPI6FYSNDLQZ3Y7OLKDAVCNFSM6AAAAABE3JHKI6VHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDAOBQGU2DSMBRHE>
.
You are receiving this because you authored the thread.Message ID:
***@***.***>
|
|
Only write and print functions, nothing binary so far: happy to take
proposals, I don't know what gr does.
…On Sat, Apr 27, 2024, 5:19 AM Albin Ahlbäck ***@***.***> wrote:
Is there any I/O-functionality?
—
Reply to this email directly, view it on GitHub
<#1845 (comment)>, or
unsubscribe
<https://github.com/notifications/unsubscribe-auth/AAKDUINVCKIOXOEH24DEZUDY7OJUVAVCNFSM6AAAAABE3JHKI6VHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDAOBQGUZDAMJTGU>
.
You are receiving this because you authored the thread.Message ID:
***@***.***>
|
|
Only while I'm still debugging/adding requested functionality: they will be
removed before the merge.
…On Sat, Apr 27, 2024, 5:17 AM Albin Ahlbäck ***@***.***> wrote:
***@***.**** commented on this pull request.
------------------------------
In src/gr_sparse_mat/test/t-conversion.c
<#1845 (comment)>:
> + // flint_printf("\n\ncoo_mat = "); status |= gr_coo_mat_print_nz(coo_mat, ctx); flint_printf("\nnnz = %d\n", coo_mat->nnz);
+ status |= gr_lil_mat_set_coo_mat(lil_mat, coo_mat, ctx);
+ // flint_printf("\n\nlil_mat = "); status |= gr_lil_mat_print_nz(lil_mat, ctx); flint_printf("\nnnz = %d\n", lil_mat->nnz);
+ status |= gr_lil_mat_set(lil_mat2, lil_mat, ctx);
+ // flint_printf("\n\nlil_mat = "); status |= gr_lil_mat_print_nz(lil_mat2, ctx); flint_printf("\nnnz = %d\n", lil_mat2->nnz);
+ status |= gr_csr_mat_set_lil_mat(csr_mat, lil_mat2, ctx);
+ // flint_printf("\n\ncsr_mat = "); status |= gr_csr_mat_print_nz(csr_mat, ctx); flint_printf("\nnnz = %d\n", csr_mat->nnz);
+ status |= gr_mat_set_csr_mat(dmat, csr_mat, ctx);
+ // flint_printf("\n\nmat = "); status |= gr_mat_print(dmat, ctx);
+ status |= gr_lil_mat_set_mat(lil_mat2, dmat, ctx);
+ // flint_printf("\n\nlil_mat = "); status |= gr_lil_mat_print_nz(lil_mat2, ctx); flint_printf("\nnnz = %d\n", lil_mat2->nnz);
+ status |= gr_coo_mat_set_lil_mat(coo_mat2, lil_mat2, ctx);
+ // flint_printf("\n\ncoo_mat = "); status |= gr_coo_mat_print_nz(coo_mat2, ctx); flint_printf("\nnnz = %d\n", coo_mat2->nnz);
Is it beneficial to have theses print statement left in the code?
—
Reply to this email directly, view it on GitHub
<#1845 (review)>,
or unsubscribe
<https://github.com/notifications/unsubscribe-auth/AAKDUIIUB3TPV4643BDVA5LY7OJOTAVCNFSM6AAAAABE3JHKI6VHI2DSMVQWIX3LMV43YUDVNRWFEZLROVSXG5CSMV3GSZLXHMZDAMRWGY2DANJYGA>
.
You are receiving this because you authored the thread.Message ID:
***@***.***>
|
Ah, I see now. What's the current print format? |
|
Primarily comma-separated list of column: value for each row, since sparse
row form is the core implementation, though I added a coo presentation too.
…On Sat, Apr 27, 2024, 7:32 AM Albin Ahlbäck ***@***.***> wrote:
Only write and print functions, nothing binary so far: happy to take
proposals, I don't know what gr does.
Ah, I see now. What's the current print format?
—
Reply to this email directly, view it on GitHub
<#1845 (comment)>, or
unsubscribe
<https://github.com/notifications/unsubscribe-auth/AAKDUILYJO66PEKQB3Z4OBTY7OZGTAVCNFSM6AAAAABE3JHKI6VHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDAOBQG44TEMJYGY>
.
You are receiving this because you authored the thread.Message ID:
***@***.***>
|
Just the start, but converting the sparse2 stuff to the FLINT 3 generic ring system.