Opened 7 years ago
Closed 7 years ago
#16797 closed enhancement (fixed)
is_difference_matrix
Reported by:  ncohen  Owned by:  

Priority:  major  Milestone:  sage6.4 
Component:  combinatorial designs  Keywords:  
Cc:  vdelecroix  Merged in:  
Authors:  Nathann Cohen  Reviewers:  Vincent Delecroix 
Report Upstream:  N/A  Work issues:  
Branch:  dee4275 (Commits, GitHub, GitLab)  Commit:  dee4275d3b2af94bd625bbdd64c8b165d8414e1a 
Dependencies:  Stopgaps: 
Description
Spent a lifetime on that because I had not noticed that I was first enumerating elements using a list L, then enumerating the same elements using {x:i for i,x in enumerate(L)}
. And of course that was not the same ordering T_T
Nathann
Change History (42)
comment:1 Changed 7 years ago by
 Branch set to u/ncohen/16797
 Status changed from new to needs_review
comment:2 Changed 7 years ago by
 Commit set to e76b5ab0beef31356d6b640c523f380a92f7efbb
comment:3 Changed 7 years ago by
 Commit changed from e76b5ab0beef31356d6b640c523f380a92f7efbb to 8f10d9412cbb9d713dd27ec72acb97c14f3d479e
Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:
8f10d94  trac #16797: is_difference_matrix

comment:4 followup: ↓ 5 Changed 7 years ago by
Hi Nathann,
Why do you use this strange convention about rows/cols? From the Handbook chapter 17, there are k rows and lambda G columns...
Vincent
comment:5 in reply to: ↑ 4 Changed 7 years ago by
Why do you use this strange convention about rows/cols? From the Handbook chapter 17, there are k rows and lambda G columns...
Because in the OA it is the columns that must be orthogonal ...
Nathann
comment:6 Changed 7 years ago by
Hum... from the documentation of OA_from_quasi_difference_matrix
Let `G` be a group of order `g`. A *difference matrix* `M` is a `k \times g` matrix with entries from `G` such that for any `1\leq i < j < k` the set `\{d_{il}d_{jl}:1\leq l \leq g\}` is equal to `G`.
Vincent
comment:7 Changed 7 years ago by
I think that it will be easier to work with k columns. Look at the OA code, we spend our lifetime creating rows for every set of a PBD, or developping rows according to some groups... For instance look at the doc you just wrote for the helper function : you have one row in your matrix, and you let this H hyperplane act on it to generate new rows.
If you want to add columns instead it becomes a mess, you have to extend each column instead of just adding a row to the matrix.
By the way, many cals to OA_from_quasi_difference_matrix
are already preceded by a M=zip(*M) # transpose M
. Actually, I am looking at database.py and really a LOT of code would be simplified by changing that: look how often you have code like that to create the matrices:
for <something>: M[0].append(e) M[1].append(f) ....
If all that was the transpose instead, we would just have to write:
for <something>: M.append([e,f,...])
What do you think ?
Nathann
comment:8 Changed 7 years ago by
I am not against having the difference matrix the other way around (actually I do not care). But then:
 rewrite the documentation appropriately and warn the reader that the convention is different from the Handbook,
 modify all functions that use difference matrix
Vincent
comment:9 Changed 7 years ago by
 Status changed from needs_review to needs_work
comment:10 Changed 7 years ago by
Okay, then I will do this after the current patches have been reviewed. Otherwise we will have to deal with countless conflicts again. I write this to follow the current convention in Sage and will change it later, along with the implemented code.
Nathann
comment:11 Changed 7 years ago by
 Commit changed from 8f10d9412cbb9d713dd27ec72acb97c14f3d479e to aa097f61b7982ddd88cc5347f20da76098912d7a
Branch pushed to git repo; I updated commit sha1. New commits:
aa097f6  trac #16797: Fit with the current definition until we can change everything at once

comment:12 Changed 7 years ago by
 Status changed from needs_work to needs_review
comment:13 Changed 7 years ago by
I added a small check before M = zip(*M)
because this operation takes the min of the length of the lists in M:
sage: zip([1],[2],[3,4,5,6]) [(1, 2, 3)]
Nathann
comment:14 followup: ↓ 15 Changed 7 years ago by
 Status changed from needs_review to needs_info
Hello,
1) Why do you transpose the matrix to check the definition? You would avoid many multiplications in x_minus_y[M_c[ii*K+i]][M_c[ii*K+j]]
.
2) I do not like the order of the argument. For orthogonal arrays it is OA, k, n, but here it is G, k, M... I would prefer G,M,k.
3) We can win a lot of time with a C array for x_minus_y
and seen_values
. It is at u/vdelecroix/16797
Vincent
comment:15 in reply to: ↑ 14 Changed 7 years ago by
Yo !
1) Why do you transpose the matrix to check the definition?
When we will change the definition to have k columns instead of k rows we will just have to remove the block of lines before the transposition.
You would avoid many multiplications in
x_minus_y[M_c[ii*K+i]][M_c[ii*K+j]]
.
I am not sure that we pay these multiplications. Because of your experience with is_sum_of_square
and my useless optimization in is_orthogonal_array
2) I do not like the order of the argument. For orthogonal arrays it is OA, k, n, but here it is G, k, M... I would prefer G,M,k.
If you prefer ...
3) We can win a lot of time with a C array for
x_minus_y
andseen_values
. It is atu/vdelecroix/16797
I was not sure that it would make a difference, thus I kept the Cython version. Do you have timings ?
Nathann
comment:16 Changed 7 years ago by
old
sage: from sage.combinat.designs.designs_pyx import is_difference_matrix sage: q = 5**3 sage: F = GF(q,'x') sage: M = [[x*y for y in F] for x in F] sage: timeit("is_difference_matrix(F,q,M)") 5 loops, best of 3: 92 ms per loop
new
sage: from sage.combinat.designs.designs_pyx import is_difference_matrix sage: q = 5**3 sage: F = GF(q,'x') sage: M = [[x*y for y in F] for x in F] sage: timeit("is_difference_matrix(F,q,M)") 25 loops, best of 3: 11.9 ms per loop
Vincent
comment:17 Changed 7 years ago by
Hmmmmmmmmm.. I wonder how much come from the C version of x_minus_y instead of Cython. I'll do a test.
Nathann
comment:18 Changed 7 years ago by
HMmmmmm O_o
Well, with the Cython version of x_minus_y (at public/16797) your example takes 15ms, but what I get when I try with your branch is a segfault
sage: sage: from sage.combinat.designs.designs_pyx import is_difference_matrix sage: sage: q = 5**3 sage: sage: F = GF(q,'x') sage: sage: M = [[x*y for y in F] for x in F] sage: sage: timeit("is_difference_matrix(F,q,M)")  /home/ncohen/.Sage/local/lib/libcsage.so(print_backtrace+0x31)[0x7fca2ebe7af1] /home/ncohen/.Sage/local/lib/libcsage.so(sigdie+0x1e)[0x7fca2ebe7c85] /home/ncohen/.Sage/local/lib/libcsage.so(sage_signal_handler+0x1cb)[0x7fca2ebe7545] /lib/x86_64linuxgnu/libpthread.so.0(+0xf8d0)[0x7fca33d788d0] ...
Nathann
comment:19 Changed 7 years ago by
Stupid error in the malloc. I guess that on your computer sizeof(int)
and sizeof(int *)
are different.
comment:20 Changed 7 years ago by
(updated on my branch)
comment:21 Changed 7 years ago by
Ahaah Okay, 9.12ms : you win.
God it's hard to know when to optimize :p
Nathann
comment:22 Changed 7 years ago by
Hmmm... The way you manage memory, you don't free any of the x_minus_y if there is a problem during the allocation of G_seen. By the way you could spare yourself a loop of alloc/free if you do not allow a vector of length G_card for each entry of x_minus_y but rather alloc a big table of size G_card**2
and then let int ** x_minus_y
point toward parts of it.
Nathann
comment:23 Changed 7 years ago by
Oh, and could you also give more info when an element does not appear exactly lambda times ? The way you wrote the code you can say which element it is and how many times it appeared.
Nathann
comment:24 Changed 7 years ago by
I am on it!!
Vincent
comment:25 Changed 7 years ago by
done at u/vdelecroix/16797
.
Vincent
comment:26 Changed 7 years ago by
You still don't deallocate x_minus_y when there is a problem with other mallocs.
comment:27 Changed 7 years ago by
Okay I'll do it.
Nathann
comment:28 followup: ↓ 31 Changed 7 years ago by
Hey I just reread your comment above: if we want to copy is_orthogonal_array's input we should make M be the first parameter, shouldn't we ?
Nathann
comment:29 Changed 7 years ago by
 Commit changed from aa097f61b7982ddd88cc5347f20da76098912d7a to b78c3475cc243b0f50cdd5c5c7bc6e9285bb97fd
Branch pushed to git repo; I updated commit sha1. New commits:
5d21b8c  trac #16797: int * and ** instead of Python list

8f29bd1  trac #16797: change int to int * in a malloc

811e3e9  trac #16797: better malloc + better error msg

1b8c2e3  trac #16797: More compact mallocs

b78c347  trac #16797: Reorder the arguments to copy is_orthogonal_array

comment:30 Changed 7 years ago by
 Status changed from needs_info to needs_review
comment:31 in reply to: ↑ 28 ; followup: ↓ 32 Changed 7 years ago by
Replying to ncohen:
Hey I just reread your comment above: if we want to copy is_orthogonal_array's input we should make M be the first parameter, shouldn't we ?
I don't know. For difference family it is D,G
and not G,D
.
I like the rest... do what you want with the input and set to positive review.
(I am going for a bike run.... back in 3 hours)
Vincent
comment:32 in reply to: ↑ 31 Changed 7 years ago by
 Reviewers set to Vincent Delecroix
 Status changed from needs_review to positive_review
Yoooooo !
I like the rest... do what you want with the input and set to positive review.
(I am going for a bike run.... back in 3 hours)
Okayyyy, have fun ! And thanks for this patch !
Nathann
comment:33 Changed 7 years ago by
 Status changed from positive_review to needs_work
I have a problem
sage: B [[0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 1, 2, 3, 4, 2, 3, 4, 0, 1], [0, 2, 4, 1, 3, 3, 0, 2, 4, 1], [0, 3, 1, 4, 2, 3, 1, 4, 2, 0], [0, 4, 3, 2, 1, 2, 1, 0, 4, 3], [0, 1, 4, 4, 1, 0, 2, 3, 3, 2], [1, 0, 1, 4, 4, 2, 0, 2, 3, 3], [4, 1, 0, 1, 4, 3, 2, 0, 2, 3]] sage: sage: is_difference_matrix(B,G,3,2,verbose=True) The matrix has 8 columns and lambda.G=2.5=10 False
There are 8 rows and 10 columns...
Vincent
comment:34 Changed 7 years ago by
All right, it is because M_nrows
is initialized to soon. I will add the example in the doc.
Vincent
comment:35 Changed 7 years ago by
See at u/vdelecroix/16797
.
Vincent
comment:36 Changed 7 years ago by
 Status changed from needs_work to needs_review
comment:37 Changed 7 years ago by
Thaaanks ! Then it can get in.
By the way, I love "the future". Back in the old phone modem days, I wouldn't have thought that I'd be reviewing a trac ticket from a laptop in a car on the highway using a phone as a modem :PPP
Nathann
comment:38 Changed 7 years ago by
 Status changed from needs_review to positive_review
comment:39 Changed 7 years ago by
 Status changed from positive_review to needs_work
Merge conflict, please fix
comment:40 Changed 7 years ago by
 Commit changed from b78c3475cc243b0f50cdd5c5c7bc6e9285bb97fd to dee4275d3b2af94bd625bbdd64c8b165d8414e1a
Branch pushed to git repo; I updated commit sha1. New commits:
dee4275  trac #16797: Merged with 6.4.beta0

comment:41 Changed 7 years ago by
 Status changed from needs_work to positive_review
comment:42 Changed 7 years ago by
 Branch changed from u/ncohen/16797 to dee4275d3b2af94bd625bbdd64c8b165d8414e1a
 Resolution set to fixed
 Status changed from positive_review to closed
Branch pushed to git repo; I updated commit sha1. New commits:
trac #16797: is_difference_matrix