Increase Performance in Projective Morphism

[drose]

Increase Performance in Projective Morphism by replacing all evaluation of polynomials with fast_callable from sage.ext.fast_callable

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​5bafb1b

# Tue Oct 22 20:59:00 2013 -0400

[bhutz]

I haven't looked closely at functionality yet, but there are a few things that need to be fixed first.

1) There are several doctests that fail with this:

sage -t --long sage/schemes/elliptic_curves/lseries_ell.py  # 1 doctest failed
sage -t --long sage/schemes/plane_conics/con_field.py  # 3 doctests failed
sage -t --long sage/schemes/generic/morphism.py  # 3 doctests failed

It looks like most of these are because the __call__ in generic\morphism first coerces the input point into the domain and the __call__ in projective\morphism does not. Although the elliptic curve failure seems to be a precision issue.

2) The new functions must have appropriate documentation. See the developers guide for the appropriate format.

SCORE projective_morphism.py: 91.9% (34 of 37)

Missing documentation:

• line 182: def call(self, x, check=True)
• line 187: def _fast_eval(self, x, check=True)
• line 2625: def _fast_eval(self, x, check=True)

3) There are a number of Tab characters in projective_morphism.py

sage -t --long sage/schemes/projective/projective_morphism.py # Tab character found

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​52ec0b3

trac 15780. Change to projective morphism _call_ function.

[bhutz]

Testing functionality. Everywhere where the evaluation works it gives the right answer. However, there are several places where you can no longer do the evaluation which give a variety of errors. These all previously worked.

T.<z>=PowerSeriesRing(ZZ)
P.<x,y>=ProjectiveSpace(T,1)
H=End(P)
f=H([x^2+x*y,y^2])
Q=P(z,1)
f(Q)

T.<z>=LaurentSeriesRing(ZZ)
P.<x,y>=ProjectiveSpace(T,1)
H=End(P)
f=H([x^2+x*y,y^2])
Q=P(z,1)
f(Q)

T.<z>=PolynomialRing(Qp(7))
I=T.ideal(z^3)
P.<x,y>=ProjectiveSpace(T.quotient_ring(I),1)
H=End(P)
f=H([x^2+x*y,y^2])
Q=P(z^2,1)
f(Q)

T.<z>=PolynomialRing(CC)
I=T.ideal(z^3)
P.<x,y>=ProjectiveSpace(T.quotient_ring(I),1)
H=End(P)
f=H([x^2+x*y,y^2])
Q=P(z^2,1)
f(Q)

T.<z>=LaurentSeriesRing(CC)
R.<t>=PolynomialRing(T)
P.<x,y>=ProjectiveSpace(R,1)
H=End(P)
f=H([x^2+x*y,y^2])
F=f.dehomogenize(1)
Q=P(t^2,z)
f(Q)

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​5ca376b

trac 15780. Added documentation.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​888a2ab

trac 15780. Added documentation to projective morphism.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​87572ab

trac 15780. Added documentation.

[bhutz]

Ticket #16051 fixes the functionality issues. So what remains to be done is the list of issues back in comment 5.

[nbruin]

initialization of the fast_polys should probably be done lazily, unless you can show it is always quick. There are many reasons to create morphisms/rational maps that do not involve evaluating them at points. For instance, you might want to take direct or inverse images of subvarieties. For doing that one does not evaluate the defining polynomials, but one does computations with the polynomials. It would be a shame to waste time on computing fast_polys in those cases.

You'd also have to experiment to see if doing this really is an improvement: what you're basically doing is representing the polynomials by predetermined straight-line programs, and using those programs regardless of the ring over which the evaluation point is defined. Different rings may benefit from different evaluation strategies. So there's something to be said for deferring the choice of evaluation algorithm until the ring in which this needs to happen is known. Of course, in many common cases the ring for the point is the same as the base ring of the polynomials. Then fast_eval probably does pretty well (also because no intermediate coercions are required).

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​2993f75

Removed more camelCase variables.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​27663fe	Trac 16051: comparing fast_callable constants by value and type
​9eb0789	Trac #16051: changed type to parent for constant comparison
​1aaeab2	Merge branch 'u/bhutz/ticket/16051' of git://trac.sagemath.org/sage into ticket/15780

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​83f8b79

Fixed doctests.

[drose]

I used timeit on a random degree 3 map on P^{3 and it took 717 microseconds per loop. This is very fast and doesn't really warrant lazy initialization.}

[bhutz]

Had some trouble with the rebase, so had to force the push on the branch. I think everything is still correct.

Although there are a couple doctest failures that are still present

generic/morphism.py projective/projective_morphism.py (from the parallel code #15920)

---

Last 10 new commits:

​1ef4bf6	trac 15780. Change to projective morphism _call_ function.
​b707ba2	trac 15780. Added documentation.
​aaade14	trac 15780. Added documentation to projective morphism.
​309f4ed	trac 15780. Added documentation.
​e351fcc	Trac 16051: comparing fast_callable constants by value and type
​79cbdf3	Trac #16051: changed type to parent for constant comparison
​c9cf64e	Removed camelCase variables.
​4ffd64c	Removed more camelCase variables.
​f7c887e	Fixed doctests.
​6be7361	rebase to 6.2.beta7 and remove whitespace

[nbruin]

Replying to drose:

I used timeit on a random degree 3 map on P³ and it took 717 microseconds per loop. This is very fast and doesn't really warrant lazy initialization.

Whether 717 ms is a lot depends on what you compare it to. I'm finding (on 6.0):

sage: P5=ProjectiveSpace(Rationals(),5)
sage: R=P2.coordinate_ring()
sage: mon=[R({tuple(a):1}) for a in WeightedIntegerVectors(4,[1,1,1,1,1,1])]
sage: L=[sum(random_sublist(mon,0.1)) for j in [0..5]]
sage: %timeit P2.hom(L,P2)
10000 loops, best of 3: 63.3 us per loop
sage: %timeit [fast_callable(l) for l in L]
100 loops, best of 3: 1.94 ms per loop

as you see, creating the fast callables is *much* slower than creating the map. If you increase the polynomials in L (by increasing the degree of the monomials, for instance), the difference gets worse. I'd say that's enough time penalty to hold off on it unless you actually need the fast_callables.

[bhutz]

Yes, it is definitely slower than creating the map, but < 1ms didn't seem to be worth bothering over. Note that you wrote 717ms not the actual 717us.

We could put a try/except in call and init them in the except. Any idea how much overhead a try block adds? As far as I can tell in a couple tests it doesn't actually add an execution time. If that's the case, then there is no reason not to do it lazily.

[nbruin]

Replying to bhutz:

We could put a try/except in call and init them in the except. Any idea how much overhead a try block adds? As far as I can tell in a couple tests it doesn't actually add an execution time. If that's the case, then there is no reason not to do it lazily.

Python is designed to let "try" be very fast. If an exception occurs there is a bit larger penalty but in this case it wouldn't be measurable compared to the fast callable.

The recommended way to store the fast_callable would be by accessing them through a cached_method or a lazy_attribute. See their docstrings for examples. Just don't go sticking cached_method decorators on methods arbitrarily. That's how we get memory leaks :-).

[drose]

New commits:

​d4c86a3	Fixed error message in __call__.
​a82afa0	Added lazy_attribute to fastpolys.
​6c0149c	Merge branch 'develop' into ticket/15780
​8215514	Fixed parallelizaion using p_iter_fork.
​9d93ea5	Minor changes to fix initial bugs in fix.
​3fd6fbd	Merge branch 'u/drose/ticket/16168' of git://trac.sagemath.org/sage into ticket/15780
​bf3c19f	Changed code to comply with style.
​d0489ec	Merge branch 'u/drose/ticket/16168' of git://trac.sagemath.org/sage into ticket/15780
​2a4342e	Fixed lazy attribute and merged 16168.

[bhutz]

Functionality all tests out fine. A couple minor issues that need to be addressed

1) _fast_eval does not use the parameter 'check'. Simply remove (from both _fast_eval functions)

2) @lazy_attribute function needs the documentation string + an example

2b) in the @lazy_attribute. Put the 'prime' and 'degree' inside the 'if' block so they are only computed if they are needed.

3) There is an extra blank line in _fast_eval for finite fields that can be removed

4) I'd like to see a few more doctests in _fast_eval. Maybe

a) Powerseries ring

b) quotient ring

You can pull them from comment 8 above.

[nbruin]

• justify constant 2**27 as bound on "no overflow in float" (you should use doubles).
• for the sake of defensive programming, would you use

  max(abs(c.lift()) for c in coefficients)
  

instead? max(GF(17)).lift() happens to work at present, but it doesn't really make mathematical sense.

Also, if you really want to scrape the bottom on this, you should use balanced representatives, also on the coordinates when you fill them in. That would give you a larger range where you can use float.

What you'd really need is a cdef ulong domain. Then you could just evaluate with 64-bit integers, which would give you a much better bound than 2^27.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​7484355

Fixed issues reported on trac server comment 28 and 29.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​9be3bba

Changed doctest for fast_eval.

[bhutz]

The previous issues all now test out and my concerns were all addressed. In response to Nils, the 27 is now replaced with the system float precision which is 53 on my 64-bit machine. I'm not overly concerned about aggressively optimizing the float usage for speed increase since if the prime is large enough to overrun the height bound for 'reasonable' maps, then the algorithms would never finish anyway (which is what the speed increase is supposed to help). Maps that overrun the bound due to large degree/moderate prime are theoretically fine in the algorithm, but won't get this speed increase and will do the fast_eval over ZZ (or whatever the base ring is).

If you have any other comments on these changes let me know before I mark this as positive.

[git]

Branch pushed to git repo; I updated commit sha1 and set ticket back to needs_review. New commits:

​c858d54	Switched parallelization to p_iter_fork.
​0a1cf6d	Added documentation on new keyword ncpus.
​be46079	moved import to header
​857d08c	Merge branch 'u/bhutz/ticket/16168' of git://trac.sagemath.org/sage into ticket/15780

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​81fba1b

Fixed import statements.

[bhutz]

The 16168 changes had not been merged in. Now they are.

[bhutz]

This passes all tests and the changes look fine.

[bhutz]

need to clean-up the history

[bhutz]

never mind. I think this is ok. I was concerned about commit 6c0149c, but it's where it is updated to 6.2.beta8