Package sage :: Package schemes :: Package elliptic_curves :: Module formal_group :: Class EllipticCurveFormalGroup
[hide private]
[frames] | no frames]

Class EllipticCurveFormalGroup

source code

                      object --+    
                               |    
structure.sage_object.SageObject --+
                                   |
                                  EllipticCurveFormalGroup


The formal group associated to an elliptic curve.



Instance Methods [hide private]
 
__init__(self, E)
x.__init__(...) initializes x; see x.__class__.__doc__ for signature
source code
 
_repr_(self) source code
 
curve(self)
The elliptic curve this formal group is associated to.
source code
 
w(self, prec=20)
The formal group power series w.
source code
 
x(self, prec=20)
Return the formal series $x(t) = t/w(t)$ in terms of the local parameter $t = -x/y$ at infinity.
source code
 
y(self, prec=20)
Return the formal series $y(t) = -1/w(t)$ in terms of the local parameter $t = -x/y$ at infinity.
source code
 
differential(self, prec=20)
Returns the power series $f(t) = 1 + \cdots$ such that $f(t) dt$ is the usual invariant differential $dx/(2y + a_1 x + a_3)$.
source code
 
log(self, prec=20)
Returns the power series $f(t) = t + \cdots$ which is an isomorphism to the additive formal group.
source code
 
inverse(self, prec=20)
The formal group inverse law i(t), which satisfies F(t, i(t)) = 0.
source code
 
group_law(self, prec=10)
The formal group law.
source code
 
mult_by_n(self, n, prec=10)
The formal 'multiplication by n' endomorphism $[n]$.
source code
 
sigma(self, prec=10) source code

Inherited from structure.sage_object.SageObject: __hash__, __new__, __repr__, _axiom_, _axiom_init_, _gap_, _gap_init_, _gp_, _gp_init_, _interface_, _interface_init_, _interface_is_cached_, _kash_, _kash_init_, _macaulay2_, _macaulay2_init_, _magma_, _magma_init_, _maple_, _maple_init_, _mathematica_, _mathematica_init_, _maxima_, _maxima_init_, _octave_, _octave_init_, _pari_, _pari_init_, _r_init_, _sage_, _singular_, _singular_init_, category, db, dump, dumps, plot, rename, reset_name, save, version

Inherited from object: __delattr__, __getattribute__, __reduce__, __reduce_ex__, __setattr__, __str__

Properties [hide private]

Inherited from object: __class__

Method Details [hide private]

__init__(self, E)
(Constructor)

source code 
x.__init__(...) initializes x; see x.__class__.__doc__ for signature

Overrides: object.__init__
(inherited documentation)

curve(self)

source code 

The elliptic curve this formal group is associated to.

EXAMPLES:
    sage: E = EllipticCurve("37a")
    sage: F = E.formal_group()
    sage: F.curve()
    Elliptic Curve defined by y^2 + y = x^3 - x over Rational Field        

w(self, prec=20)

source code 

The formal group power series w.

INPUT:
    prec -- integer

OUTPUT:
    a power series with given precision

DETAILS:
    Return the formal power series
    $$
           w(t) = t^3 + a_1 t^4 + (a_2 + a_1^2) t^5 + \cdots
    $$
    to precision $O(t^prec)$ of Proposition IV.1.1 of [Silverman
    AEC1].  This is the formal expansion of $w = -1/y$ about the
    formal parameter $t = -x/y$ at $\\infty$.

    The result is cached, and a cached version is returned if
    possible.

WARNING:
    The resulting power series will have precision prec, but its
    parent PowerSeriesRing will have default precision 20 (or whatever
    the default default is).

ALGORITHM:
    Uses Newton's method to solve the elliptic curve equation
    at the origin. Complexity is roughly $O(M(n))$ where
    $n$ is the precision and $M(n)$ is the time required to multiply
    polynomials of length $n$ over the coefficient ring of $E$.

AUTHOR:
    -- David Harvey (2006-09-09): modified to use Newton's method
    instead of a recurrence formula.

EXAMPLES:
    sage: e = EllipticCurve([0, 0, 1, -1, 0])
    sage: e.formal_group().w(10)
     t^3 + t^6 - t^7 + 2*t^9 + O(t^10)

  Check that caching works:
    sage: e = EllipticCurve([3, 2, -4, -2, 5])
    sage: e.formal_group().w(20)
     t^3 + 3*t^4 + 11*t^5 + 35*t^6 + 101*t^7 + 237*t^8 + 312*t^9 - 949*t^10 - 10389*t^11 - 57087*t^12 - 244092*t^13 - 865333*t^14 - 2455206*t^15 - 4366196*t^16 + 6136610*t^17 + 109938783*t^18 + 688672497*t^19 + O(t^20)
    sage: e.formal_group().w(7)
     t^3 + 3*t^4 + 11*t^5 + 35*t^6 + O(t^7)
    sage: e.formal_group().w(35)
     t^3 + 3*t^4 + 11*t^5 + 35*t^6 + 101*t^7 + 237*t^8 + 312*t^9 - 949*t^10 - 10389*t^11 - 57087*t^12 - 244092*t^13 - 865333*t^14 - 2455206*t^15 - 4366196*t^16 + 6136610*t^17 + 109938783*t^18 + 688672497*t^19 + 3219525807*t^20 + 12337076504*t^21 + 38106669615*t^22 + 79452618700*t^23 - 33430470002*t^24 - 1522228110356*t^25 - 10561222329021*t^26 - 52449326572178*t^27 - 211701726058446*t^28 - 693522772940043*t^29 - 1613471639599050*t^30 - 421817906421378*t^31 + 23651687753515182*t^32 + 181817896829144595*t^33 + 950887648021211163*t^34 + O(t^35)

x(self, prec=20)

source code 

Return the formal series $x(t) = t/w(t)$ in terms of the local
parameter $t = -x/y$ at infinity.

INPUT:
    prec -- integer

OUTPUT:
    a laurent series with given precision

DETAILS:
    Return the formal series
    $$
           x(t) = t^{-2} - a_1 t^{-1} - a_2 - a_3 t - \cdots
    $$
    to precision $O(t^prec)$ of page 113 of [Silverman
    AEC1].

WARNING:
    The resulting series will have precision prec, but its
    parent PowerSeriesRing will have default precision 20 (or whatever
    the default default is).

EXAMPLES:
    sage: EllipticCurve([0, 0, 1, -1, 0]).formal_group().x(10)
     t^-2 - t + t^2 - t^4 + 2*t^5 - t^6 - 2*t^7 + 6*t^8 - 6*t^9 + O(t^10)

y(self, prec=20)

source code 

Return the formal series $y(t) = -1/w(t)$ in terms of the local
parameter $t = -x/y$ at infinity.

INPUT:
    prec -- integer

OUTPUT:
    a laurent series with given precision

DETAILS:
    Return the formal series
    $$
           y(t) = - t^{-3} + a_1 t^{-2} + a_2 t + a_3 + \cdots
    $$
    to precision $O(t^prec)$ of page 113 of [Silverman
    AEC1].

    The result is cached, and a cached version is returned if
    possible.

WARNING:
    The resulting series will have precision prec, but its
    parent PowerSeriesRing will have default precision 20 (or whatever
    the default default is).

EXAMPLES:
    sage: EllipticCurve([0, 0, 1, -1, 0]).formal_group().y(10)
     -t^-3 + 1 - t + t^3 - 2*t^4 + t^5 + 2*t^6 - 6*t^7 + 6*t^8 + 3*t^9 + O(t^10)

differential(self, prec=20)

source code 

Returns the power series $f(t) = 1 + \cdots$ such that $f(t) dt$ is
the usual invariant differential $dx/(2y + a_1 x + a_3)$.

INPUT:
   prec -- nonnegative integer, answer will be returned $O(t^{\var{prec}})$

OUTPUT:
    a power series with given precision

DETAILS:
    Return the formal series
    $$
           f(t) = 1 + a_1 t + ({a_1}^2 + a_2) t^2 + \cdots
    $$
    to precision $O(t^prec)$ of page 113 of [Silverman
    AEC1].

    The result is cached, and a cached version is returned if
    possible.

WARNING:
    The resulting series will have precision prec, but its
    parent PowerSeriesRing will have default precision 20 (or whatever
    the default default is).

EXAMPLES:
   sage: EllipticCurve([-1, 1/4]).formal_group().differential(15)
    1 - 2*t^4 + 3/4*t^6 + 6*t^8 - 5*t^10 - 305/16*t^12 + 105/4*t^14 + O(t^15)
   sage: EllipticCurve(Integers(53), [-1, 1/4]).formal_group().differential(15)
    1 + 51*t^4 + 14*t^6 + 6*t^8 + 48*t^10 + 24*t^12 + 13*t^14 + O(t^15)

AUTHOR:
   -- David Harvey (2006-09-10): factored out of log

log(self, prec=20)

source code 

Returns the power series $f(t) = t + \cdots$ which is an isomorphism
to the additive formal group.

Generally this only makes sense in characteristic zero, although the
terms before $t^p$ may work in characteristic $p$.

INPUT:
   prec -- nonnegative integer

OUTPUT:
    a power series with given precision

EXAMPLES:
   sage: EllipticCurve([-1, 1/4]).formal_group().log(15)
    t - 2/5*t^5 + 3/28*t^7 + 2/3*t^9 - 5/11*t^11 - 305/208*t^13 + O(t^15)

AUTHOR:
   -- David Harvey (2006-09-10): rewrote to use differential

inverse(self, prec=20)

source code 

The formal group inverse law i(t), which satisfies F(t, i(t)) = 0.

INPUT:
    prec -- integer

OUTPUT:
    a power series with given precision

DETAILS:
    Return the formal power series
    $$
           i(t) = - t + a_1 t^2 + \cdots
    $$
    to precision $O(t^prec)$ of page 114 of [Silverman AEC1].

    The result is cached, and a cached version is returned if
    possible.

WARNING:
    The resulting power series will have precision prec, but its
    parent PowerSeriesRing will have default precision 20 (or whatever
    the default default is).

EXAMPLES:
    sage: e = EllipticCurve([1, 2])
    sage: F = e.formal_group().group_law(5)
    sage: i = e.formal_group().inverse(5)
    sage: Fx = F.base_base_extend(i.parent())
    sage: Fx (i) (i.parent().gen())
    O(t^5)

group_law(self, prec=10)

source code 

The formal group law.

INPUT:
    prec -- integer

OUTPUT:
    a power series with given precision in ZZ[[ ZZ[['t1']], 't2']]

DETAILS:
    Return the formal power series
    $$
           F(t1, t2) = t1 + t2 - a1 t1 t2 - \cdots
    $$
    to precision $O(t^prec)$ of page 115 of [Silverman AEC1].

    The result is cached, and a cached version is returned if
    possible.

WARNING:
    The resulting power series will have precision prec, but its
    parent PowerSeriesRing will have default precision 20 (or whatever
    the default default is).

AUTHOR:
    -- Nick Alexander: minor fixes, docstring

EXAMPLES:
    sage: e = EllipticCurve([1, 2])
    sage: F = e.formal_group().group_law(5); F
     t1 + O(t1^5) + (1 - 2*t1^4 + O(t1^5))*t2 + (-4*t1^3 + O(t1^5))*t2^2 + (-4*t1^2 - 30*t1^4 + O(t1^5))*t2^3 + (-2*t1 - 30*t1^3 + O(t1^5))*t2^4 + O(t2^5)
    sage: i = e.formal_group().inverse(5)
    sage: Fx = F.base_base_extend(i.parent())
    sage: Fx (i.parent().gen()) (i)
     O(t^5)

    Let's ensure caching with changed precision is working:
    sage: e.formal_group().group_law(4)
     t1 + O(t1^4) + (1 + O(t1^4))*t2 + (-4*t1^3 + O(t1^4))*t2^2 + (-4*t1^2 + O(t1^4))*t2^3 + O(t2^4)

mult_by_n(self, n, prec=10)

source code 

The formal 'multiplication by n' endomorphism $[n]$.

INPUT:
    prec -- integer

OUTPUT:
    a power series with given precision

DETAILS:
    Return the formal power series
    $$
           [n](t) = n t + \cdots
    $$
    to precision $O(t^prec)$ of Proposition 2.3 of [Silverman AEC1].

WARNING:
    The resulting power series will have precision prec, but its
    parent PowerSeriesRing will have default precision 20 (or whatever
    the default default is).

AUTHOR:
    -- Nick Alexander: minor fixes, docstring
    -- David Harvey (2007-03): faster algorithm for char 0 field case

EXAMPLES:
    sage: e = EllipticCurve([1, 2, 3, 4, 6])
    sage: e.formal_group().mult_by_n(0, 5)
     O(t^5)
    sage: e.formal_group().mult_by_n(1, 5)
     t + O(t^5)

We verify an identity of low degree:

    sage: none = e.formal_group().mult_by_n(-1, 5)
    sage: two = e.formal_group().mult_by_n(2, 5)
    sage: ntwo = e.formal_group().mult_by_n(-2, 5)
    sage: ntwo - none(two)
     O(t^5)
    sage: ntwo - two(none)
     O(t^5)

It's quite fast:
    sage: E = EllipticCurve("37a"); F = E.formal_group()
    sage: F.mult_by_n(100, 20)
    100*t - 49999950*t^4 + 3999999960*t^5 + 14285614285800*t^7 - 2999989920000150*t^8 + 133333325333333400*t^9 - 3571378571674999800*t^10 + 1402585362624965454000*t^11 - 146666057066712847999500*t^12 + 5336978000014213190385000*t^13 - 519472790950932256570002000*t^14 + 93851927683683567270392002800*t^15 - 6673787211563812368630730325175*t^16 + 320129060335050875009191524993000*t^17 - 45670288869783478472872833214986000*t^18 + 5302464956134111125466184947310391600*t^19 + O(t^20)