improve solve_mod performance

[dsm]

Currently solve_mod(x^2==41, 100) breaks this down into solving modulus 2^2 and 5^2, but it doesn't do further reductions, so it actually loops over 4 and 25 possibilities directly. This works well for small exponents:

sage: time solve_mod(x^2 == 41, 10^2) 
CPU times: user 0.01 s, sys: 0.00 s, total: 0.01 s Wall time: 0.01 s 
[(29,), (21,), (79,), (71,)] 

but performs poorly on larger ones, even if the primes involved are small:

sage: time solve_mod(x^2 == 41, 10^6)
CPU times: user 0.86 s, sys: 0.01 s, total: 0.87 s Wall time: 0.99 s 
[(703821,), (452429,), (47571,), (796179,), (203821,), (952429,), (547571,), (296179,)]

---

Apply only trac_11036-solve-mod.patch​ to the Sage library.

[kcrisman]

Agreed that this is a place we need better work. It's been fairly low priority, as you can see. Ideas? Hensel's Lemma?

[dsm]

This patch doesn't attempt to do anything sophisticated, it simply reaches for the low-hanging fruit:

(1) It modifies solve_mod_enumerate so that if the modulus is 5^100, instead of iterating over 5^100 possibilities, it finds solutions mod 5, then uses those to find solutions mod 5^2, etc.

(2) It modifies solve_mod to short-circuit in some cases.

(3) It also makes an unrelated change: currently solve_mod(x==1, 1) returns [()], but should probably return [(0,)].

The change has little effect for small moduli but wins easily on anything nontrivial:

#OLD 
sage: time solve_mod(x^2==41, 10^7)
CPU times: user 4.21 s, sys: 0.05 s, total: 4.27 s
Wall time: 4.33 s
[(3703821,), (1452429,), (3547571,), (1296179,), (8703821,), (6452429,), (8547571,), (6296179,)]

# NEW 
sage: time solve_mod(x^2==41, 10^7)
CPU times: user 0.04 s, sys: 0.00 s, total: 0.04 s
Wall time: 0.04 s
[(3703821,), (1452429,), (3547571,), (1296179,), (8703821,), (6452429,), (8547571,), (6296179,)]

It still behaves very badly for large primes, but a lot better than the current approach on powers of small primes. This includes my use case of powers of 10-- there was an OEIS sequence I was extending which I had to do manually because solve_mod was too slow, and it irritated me that Mma could do it quickly..

Beyond the doctests I've tested it for small-coefficient univariate polynomials of order <= 5 with modulus <= 12, and for lots of random multivariate polynomials. Tests were against a four-line brute force solver. Probably I've missed something, though, so more tests would be appreciated!

[mariah]

trac_11036_improve_solve_mod_perf.patch​ does not apply to sage-4.7.1.alpha2 on sknet/taurus (x86_64-Linux-nehalem-fc):

taurus% ./sage
----------------------------------------------------------------------
| Sage Version 4.7.1.alpha2, Release Date: 2011-06-07                |
| Type notebook() for the GUI, and license() for information.        |
----------------------------------------------------------------------
**********************************************************************
*                                                                    *
* Warning: this is a prerelease version, and it may be unstable.     *
*                                                                    *
**********************************************************************
sage: hg_sage.apply("/home/mariah/trac_11036_improve_solve_mod_perf.patch")
cd "/home/mariah/sage/sage-4.7.1.alpha2-x86_64-Linux-nehalem-fc-review-11036/devel/sage" && hg status
cd "/home/mariah/sage/sage-4.7.1.alpha2-x86_64-Linux-nehalem-fc-review-11036/devel/sage" && hg status
cd "/home/mariah/sage/sage-4.7.1.alpha2-x86_64-Linux-nehalem-fc-review-11036/devel/sage" && hg import   "/home/mariah/trac_11036_improve_solve_mod_perf.patch"
applying /home/mariah/trac_11036_improve_solve_mod_perf.patch
patching file sage/symbolic/relation.py
Hunk #3 FAILED at 735
Hunk #4 succeeded at 810 with fuzz 2 (offset 44 lines).
Hunk #5 FAILED at 821
Hunk #6 FAILED at 830
3 out of 6 hunks FAILED -- saving rejects to file sage/symbolic/relation.py.rej
abort: patch failed to apply
sage:

Perhaps it needs to be rebased?

[cremona]

I did a rebase but also have managed to mess up my 4.7.1.rc[01] builds while testing other tickets, so no guarantees...

The code seems to assume that the expressions is a polynomial (in one or more variables), so there should definitely be a more efficient way to do the enumeration of solutions modulo prime powers -- this should be a linear (algebra mod p) problem once the solutions mod p are known, right?

[mderickx]

removed a lot of bugs and refactored the patch

[spice]

Looks fine; can't find any bugs, tests work, speedups as advertised. Revised patch works with 4.7.2alpha2.

[dsm]

If there were things I missed -- either preexisting bugs I didn't notice or bugs I introduced -- then we should probably add doctests to make sure they stay fixed.

[mderickx]

the bugs I fixed where all in the existing doctests, I don't know how you got the timings because the code wasn't working at all