1950.sagews

time ecm.factor((2^197 + 1)/3)

[197002597249, 1348959352853811313, 251951573867253012259144010843] Time: CPU 0.16 s, Wall: 12.18 s

time factor((2^197 + 1)/3)

197002597249 * 1348959352853811313 * 251951573867253012259144010843 Time: CPU 2.54 s, Wall: 2.86 s

sage has a built-in "elliptic curve factoring" program ecm.  We can do some ecm factoring by hand.  We first create an RSF modulus $N$ that we want to factor and then produce a "reasonably random'' point on a random elliptic curve mod $N$.

N=next_prime(1667)*next_prime(2234); N; R=Integers(N); R.is_field = lambda : True

3733553

Unfortunately, we have to temporarily trick sage into thinking that $R$ is a field (which is isn't).

a=R.random_element(); b=R.random_element(); c=R.random_element();
d=b^2-a^3-c*a; [a,b,c,d]; E=EllipticCurve([c,d]); E

[994039, 2552712, 898226, 758912] Elliptic Curve defined by y^2 = x^3 + 898226*x + 758912 over Ring of integers modulo 3733553

Note that we have chosen $d$ in the coefficient of the elliptic curve so that $(a,b)$ is a point on the curve.  It is hard to construct random points on given curves because you don't know whether the $y^2$ for an $x$ that you choose is actually a square.

P=E([a,b]); P

(994039 : 2552712 : 1)

Note that the third coordinate of $P$ is $1$; $P$ is `"completely affine.''  However, some multiple of $P$ is likely to be the point at $\infty$ modulo one of the prime factors of $N$ but not the other factor of $N$.  If this occurs, we can factor $N$.

factorial(7)*P

Traceback (most recent call last): File "<stdin>", line 1, in <module> File "_sage_input_41.py", line 9, in <module> exec compile(ur'open("___code___.py","w").write("# -*- coding: utf-8 -*-\n" + _support_.preparse_worksheet_cell(base64.b64decode("ZmFjdG9yaWFsKDcpKlA="),globals())+"\n"); execfile(os.path.abspath("___code___.py"))' + '\n', '', 'single') File "", line 1, in <module> File "/private/var/folders/ck/ckRR1gdE2RWjAk+1YtnJGU+++TI/-Tmp-/tmprSSXiw/___code___.py", line 3, in <module> exec compile(ur'factorial(_sage_const_7 )*P' + '\n', '', 'single') File "", line 1, in <module> File "element.pyx", line 1398, in sage.structure.element.RingElement.__mul__ (sage/structure/element.c:11337) File "coerce.pyx", line 709, in sage.structure.coerce.CoercionModel_cache_maps.bin_op (sage/structure/coerce.c:6103) File "coerce_actions.pyx", line 453, in sage.structure.coerce_actions.IntegerMulAction._call_ (sage/structure/coerce_actions.c:6280) File "coerce_actions.pyx", line 518, in sage.structure.coerce_actions.fast_mul_long (sage/structure/coerce_actions.c:6967) File "element.pyx", line 927, in sage.structure.element.ModuleElement.__iadd__ (sage/structure/element.c:7710) File "element.pyx", line 919, in sage.structure.element.ModuleElement._add_ (sage/structure/element.c:7574) File "/Users/kribet/Desktop/sage/local/lib/python2.6/site-packages/sage/schemes/elliptic_curves/ell_point.py", line 625, in _add_ raise ZeroDivisionError, "Inverse of %s does not exist"%(x1-x2) ZeroDivisionError: Inverse of 265371 does not exist

When you see a ZeroDivisionError you know that something is up.

gcd(265371,N)

1669

N/_

2237

R = Integers(1669); E= EllipticCurve(R,[898226, 758912]); P = E([994039, 2552712])

P.order()

70

R = Integers(2237); E= EllipticCurve(R,[898226, 758912]); P = E([994039, 2552712])

P.order()

2213

Stein's code for a simple ecm routine is in his "Elementary number theory" book (available for download), p. 134.  Here it comes.

def ecm(N,B=10^3,trials=10):
    m = prod([p^int(math.log(B)/math.log(p)) for p in prime_range(B+1)])
    R = Integers(N)
    R.is_field = lambda: True
    for _ in range(trials):
        while True:
            a = R.random_element()
            if gcd(4*a.lift()^3 + 27, N) == 1: break
        try:
            m * EllipticCurve([a, 1])([0,1])
        except ZeroDivisionError, msg:
            # msg: "Inverse of <int> does not exist"
            return gcd(Integer(str(msg).split()[2]), N)
    return 1

ecm(12345644525857247889249058724507)

2011

Sato-Tate:

#%auto
def sato_tate(E, N):
    return [acos(E.ap(p)/(2*sqrt(p))) for p in prime_range(N+1) if N%p != 0]

def dist(v, b, left=float(0), right=float(pi)):
    """
    We divide the interval between left (default: 0) and 
    right (default: pi) up into b bins.
   
    For each number in v (which must left and right), 
    we find which bin it lies in and add this to a counter.
    This function then returns the bins and the number of
    elements of v that lie in each one. 

    ALGORITHM: To find the index of the bin that a given 
    number x lies in, we multiply x by b/length and take the 
    floor. 
    """
    length = right - left
    normalize = float(b/length)
    vals = {}
    d = dict([(i,0) for i in range(b)])
    for x in v:
        n = int(normalize*(float(x)-left))
        d[n] += 1
    return d, len(v)

#%auto
def graph(d, b, num=5000, left=float(0), right=float(pi)):
    s = Graphics()
    left = float(left); right = float(right)
    length = right - left
    w = length/b
    k = 0
    for i, n in d.iteritems():
        k += n
        # ith bin has n objects in it.
        s += polygon([(w*i+left,0), (w*(i+1)+left,0), \
                     (w*(i+1)+left, n/(num*w)), (w*i+left, n/(num*w))],\
                     rgbcolor=(0,0,0.5))
    return s

def sin2():
    PI = float(pi)
    return plot(lambda x: (2/PI)*math.sin(x)^2, 0,pi, plot_points=200, \
        rgbcolor=(0.3,0.1,0.1), thickness=2)

#%auto
def graph_ellcurve(E, b=10, num=5000):
    v = sato_tate(E, num)
    d, total_number_of_points = dist(v,b)
    return graph(d, b, total_number_of_points) + sin2()

g = graph_ellcurve(EllipticCurve('11a'),200,150000)
show(g)

Checking problem 5.3:

var('e1'); var('e2'); e3 = -e1 - e2; B = -e1*e2*e3; A = e1*e2+e2*e3+e1*e3

e1 e2

(e1-e2)^2*(e1-e3)^2*(e2-e3)^2 + 4*A^3 + 27*B^2

(e1 - e2)^2*(e1 + 2*e2)^2*(2*e1 + e2)^2 + 27*(e1 + e2)^2*e1^2*e2^2 - 4*((e1 + e2)*e1 + (e1 + e2)*e2 - e1*e2)^3

expand(_)

0

R=Integers(next_prime(12345678912));   a=R.random_element(); b=R.random_element(); c=R.random_element();
d=b^2-a^3-c*a; [a,b,c,d]; E=EllipticCurve([c,d]); E; E([a,b]).order().factor()

[7498018392, 7507292418, 12337758188, 9018389178] Elliptic Curve defined by y^2 = x^3 + 12337758188*x + 9018389178 over Ring of integers modulo 12345678923 7 * 19 * 41 * 566011

time EllipticCurve('4077a').congruence_number()

EllipticCurve('5077a').integral_points?

No object '.integral_points' currently defined.

E=EllipticCurve('11a')

E.integral_points?

File: /usr/local/sage/local/lib/python2.6/site-packages/sage/schemes/elliptic_curves/ell_rational_field.py

Type: <type ‘instancemethod’>

Definition: E.integral_points(mw_base=’auto’, both_signs=False, verbose=False)

Docstring:

Computes all integral points (up to sign) on this elliptic curve.

INPUT:

• mw_base - list of EllipticCurvePoint generating the Mordell-Weil group of E (default: ‘auto’ - calls E.gens())
• both_signs - True/False (default False): if True the output contains both P and -P, otherwise only one of each pair.
• verbose - True/False (default False): if True, some details of the computation are output

OUTPUT: A sorted list of all the integral points on E (up to sign unless both_signs is True)

Note

The complexity increases exponentially in the rank of curve E. The computation time (but not the output!) depends on the Mordell-Weil basis. If mw_base is given but is not a basis for the Mordell-Weil group (modulo torsion), integral points which are not in the subgroup generated by the given points will almost certainly not be listed.

EXAMPLES: A curve of rank 3 with no torsion points

sage: E=EllipticCurve([0,0,1,-7,6])
sage: P1=E.point((2,0)); P2=E.point((-1,3)); P3=E.point((4,6))
sage: a=E.integral_points([P1,P2,P3]); a
[(-3 : 0 : 1), (-2 : 3 : 1), (-1 : 3 : 1), (0 : 2 : 1), (1 : 0 : 1), (2 : 0 : 1), (3 : 3 : 1), (4 : 6 : 1), (8 : 21 : 1), (11 : 35 : 1), (14 : 51 : 1), (21 : 95 : 1), (37 : 224 : 1), (52 : 374 : 1), (93 : 896 : 1), (342 : 6324 : 1), (406 : 8180 : 1), (816 : 23309 : 1)]

sage: a = E.integral_points([P1,P2,P3], verbose=True)
Using mw_basis  [(2 : 0 : 1), (3 : -4 : 1), (8 : -22 : 1)]
e1,e2,e3:  -3.0124303725933... 1.0658205476962... 1.94660982489710
Minimal eigenvalue of height pairing matrix:  0.63792081458500...
x-coords of points on compact component with  -3 <=x<= 1
[-3, -2, -1, 0, 1]
x-coords of points on non-compact component with  2 <=x<= 6
[2, 3, 4]
starting search of remaining points using coefficient bound  5
x-coords of extra integral points:
[2, 3, 4, 8, 11, 14, 21, 37, 52, 93, 342, 406, 816]
Total number of integral points: 18

It is not necessary to specify mw_base; if it is not provided, then the Mordell-Weil basis must be computed, which may take much longer.

sage: E=EllipticCurve([0,0,1,-7,6])
sage: a=E.integral_points(both_signs=True); a
[(-3 : -1 : 1), (-3 : 0 : 1), (-2 : -4 : 1), (-2 : 3 : 1), (-1 : -4 : 1), (-1 : 3 : 1), (0 : -3 : 1), (0 : 2 : 1), (1 : -1 : 1), (1 : 0 : 1), (2 : -1 : 1), (2 : 0 : 1), (3 : -4 : 1), (3 : 3 : 1), (4 : -7 : 1), (4 : 6 : 1), (8 : -22 : 1), (8 : 21 : 1), (11 : -36 : 1), (11 : 35 : 1), (14 : -52 : 1), (14 : 51 : 1), (21 : -96 : 1), (21 : 95 : 1), (37 : -225 : 1), (37 : 224 : 1), (52 : -375 : 1), (52 : 374 : 1), (93 : -897 : 1), (93 : 896 : 1), (342 : -6325 : 1), (342 : 6324 : 1), (406 : -8181 : 1), (406 : 8180 : 1), (816 : -23310 : 1), (816 : 23309 : 1)]

An example with negative discriminant:

sage: EllipticCurve('900d1').integral_points()
[(-11 : 27 : 1), (-4 : 34 : 1), (4 : 18 : 1), (16 : 54 : 1)]

Another example with rank 5 and no torsion points

sage: E=EllipticCurve([-879984,319138704])
sage: P1=E.point((540,1188)); P2=E.point((576,1836))
sage: P3=E.point((468,3132)); P4=E.point((612,3132))
sage: P5=E.point((432,4428))
sage: a=E.integral_points([P1,P2,P3,P4,P5]); len(a)  # long time (400s!)
54

The bug reported on trac #4525 is now fixed:

sage: EllipticCurve('91b1').integral_points()
[(-1 : 3 : 1), (1 : 0 : 1), (3 : 4 : 1)]

sage: [len(e.integral_points(both_signs=False)) for e in cremona_curves([11..100])] # long time
[2, 0, 2, 3, 2, 1, 3, 0, 2, 4, 2, 4, 3, 0, 0, 1, 2, 1, 2, 0, 2, 1, 0, 1, 3, 3, 1, 1, 4, 2, 3, 2, 0, 0, 5, 3, 2, 2, 1, 1, 1, 0, 1, 3, 0, 1, 0, 1, 1, 3, 6, 1, 2, 2, 2, 0, 0, 2, 3, 1, 2, 2, 1, 1, 0, 3, 2, 1, 0, 1, 0, 1, 3, 3, 1, 1, 5, 1, 0, 1, 1, 0, 1, 2, 0, 2, 0, 1, 1, 3, 1, 2, 2, 4, 4, 2, 1, 0, 0, 5, 1, 0, 1, 2, 0, 2, 2, 0, 0, 0, 1, 0, 3, 1, 5, 1, 2, 4, 1, 0, 1, 0, 1, 0, 1, 0, 2, 2, 0, 0, 1, 0, 1, 1, 4, 1, 0, 1, 1, 0, 4, 2, 0, 1, 1, 2, 3, 1, 1, 1, 1, 6, 2, 1, 1, 0, 2, 0, 6, 2, 0, 4, 2, 2, 0, 0, 1, 2, 0, 2, 1, 0, 3, 1, 2, 1, 4, 6, 3, 2, 1, 0, 2, 2, 0, 0, 5, 4, 1, 0, 0, 1, 0, 2, 2, 0, 0, 2, 3, 1, 3, 1, 1, 0, 1, 0, 0, 1, 2, 2, 0, 2, 0, 0, 1, 2, 0, 0, 4, 1, 0, 1, 1, 0, 1, 2, 0, 1, 4, 3, 1, 2, 2, 1, 1, 1, 1, 6, 3, 3, 3, 3, 1, 1, 1, 1, 1, 0, 7, 3, 0, 1, 3, 2, 1, 0, 3, 2, 1, 0, 2, 2, 6, 0, 0, 6, 2, 2, 3, 3, 5, 5, 1, 0, 6, 1, 0, 3, 1, 1, 2, 3, 1, 2, 1, 1, 0, 1, 0, 1, 0, 5, 5, 2, 2, 0, 0, 1, 0, 0, 0, 0, 1, 1, 1, 1]

The bug reported at #4897 is now fixed:

sage: [P[0] for P in EllipticCurve([0,0,0,-468,2592]).integral_points()]
[-24, -18, -14, -6, -3, 4, 6, 18, 21, 24, 36, 46, 102, 168, 186, 381, 1476, 2034, 67246]

Note

This function uses the algorithm given in [Co1].

REFERENCES:

• [Co1] Cohen H., Number Theory Vol I: Tools and Diophantine Equations GTM 239, Springer 2007

AUTHORS:

• Michael Mardaus (2008-07)
• Tobias Nagel (2008-07)
• John Cremona (2008-07)

EllipticCurve('56b1').modular_degree()

4

1+1

2