source: sage/algebras/quaternion_algebra.py @ 4069:ce1f01cf88d7

Revision 4069:ce1f01cf88d7, 14.0 KB checked in by William Stein <wstein@…>, 6 years ago (diff)

merge
Line 
1"""
2Quaternion algebras
3
4AUTHOR: David Kohel, 2005-09
5
6TESTS:
7    sage: A = QuaternionAlgebra(QQ, -1,-1, names=list('ijk'))
8    sage: A == loads(dumps(A))
9    True
10    sage: i, j, k = A.gens()
11    sage: i == loads(dumps(i))
12    True
13
14"""
15
16#*****************************************************************************
17#  Copyright (C) 2005 David Kohel <kohel@maths.usyd.edu>
18#
19#  Distributed under the terms of the GNU General Public License (GPL)
20#
21#    This code is distributed in the hope that it will be useful,
22#    but WITHOUT ANY WARRANTY; without even the implied warranty
23#    of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
24#
25#  See the GNU General Public License for more details; the full text
26#  is available at:
27#
28#                  http://www.gnu.org/licenses/
29#*****************************************************************************
30
31from sage.misc.misc import mul
32from sage.rings.arith import kronecker, GCD, hilbert_symbol
33from sage.rings.integer import Integer
34from sage.rings.all import (IntegerRing, RationalField, PolynomialRing, is_Field)
35from sage.modules.free_module import FreeModule
36from sage.modules.free_module import VectorSpace
37from sage.matrix.matrix_space import MatrixSpace
38from sage.algebras.free_algebra import FreeAlgebra
39from sage.algebras.free_algebra_quotient import FreeAlgebraQuotient
40from sage.algebras.algebra_element import AlgebraElement
41from sage.algebras.free_algebra_quotient_element import FreeAlgebraQuotientElement
42from sage.algebras.quaternion_algebra_element import QuaternionAlgebraElement, QuaternionAlgebraElement_fast
43
44import weakref
45
46def sign(x):
47    if x < 0:
48        return -1
49    elif x > 0:
50        return 1
51    else:
52        return 0
53
54def fundamental_discriminant(D):
55    D = Integer(D)
56    D = D.square_free_part()
57    if D%4 in (0,1):
58        return D
59    return 4*D
60
61def ramified_primes(a,b):
62    a = Integer(a); b = Integer(b)
63    if a.is_square() or b.is_square() or (a+b).is_square():
64        return [ ]
65    a = a.square_free_part()
66    b = b.square_free_part()
67    c = Integer(GCD(a,b))
68    if c != 1:
69        p = c.factor()[0][0]
70        ram_prms = ramified_primes(p,-(b//p))
71        for p in ramified_primes(a//p,b):
72            if p in ram_prms:
73                ram_prms.remove(p)
74            else:
75                ram_prms.append(p)
76        ram_prms.sort()
77        return ram_prms
78    ram_prms = [ ]
79    S1 = [ p[0] for p in abs(a).factor() ]
80    for p in S1:
81        if p == 2 and b%4 == 3:
82            if kronecker(a+b,p) == -1:
83                ram_prms.append(p)
84        elif kronecker(b,p) == -1:
85            ram_prms.append(p)
86    S2 = [ p[0] for p in abs(b).factor() ]
87    for q in S2:
88        if q == 2 and a%4 == 3:
89            if kronecker(a+b,q) == -1:
90                ram_prms.append(q)
91        elif kronecker(a,q) == -1:
92            ram_prms.append(q)
93    if not 2 in ram_prms and a%4 == 3 and b%4 == 3:
94        ram_prms.append(2)
95    ram_prms.sort()
96    return ram_prms
97
98def ramified_primes_from_discs(D1,D2,T):
99    M = Integer(GCD([D1,D2,T]))
100    D3 = (T**2 - D1*D2)//4
101    facs = D3.factor()
102    D1 = fundamental_discriminant(D1)
103    D2 = fundamental_discriminant(D2)
104    D3 = fundamental_discriminant(D3)
105    ram_prms = []
106    for pow in facs:
107        p = pow[0]
108        if pow[1]%2 == 1: 
109            chi = (kronecker(D,p) for D in (D1,D2,D3))
110            if -1 in chi:
111                ram_prms.append(p)
112            elif not 1 in chi and hilbert_symbol(D1,D3,p) == -1:
113                ram_prms.append(p)
114        elif D1%p == 0 and D2%p == 0:
115            chi = (kronecker(D,p) for D in (D1,D2,D3))
116            if hilbert_symbol(D1,D3,p) == -1:
117                ram_prms.append(p)
118    return ram_prms
119
120_cache = {}
121
122def QuaternionAlgebra(K, a, b, names=['i','j','k'], denom=1):
123    """
124    Return the quaternion algebra over $K$ generated by $i$, $j$, and $k$
125    such that $i^2 = a$, $j^2 = b$, and $ij=-ji=k$.
126
127    INPUT:
128        K -- field
129        a -- element of K
130        b -- element of K
131        names -- list of three strings
132        denom -- (optional, default 1)
133
134    EXAMPLES:
135        sage: A.<i,j,k> = QuaternionAlgebra(QQ, -1,-1)
136        sage: i^2
137        -1
138        sage: j^2
139        -1
140        sage: i*j
141        k
142        sage: j*i
143        -k
144        sage: (i + j + k)^2
145        -3
146        sage: A.ramified_primes()
147        [2]
148    """
149    if not is_Field(K):
150        raise ValueError, "Base ring K (= %s) must be a field."%K
151
152    if K(2) == 0:
153        raise ValueError, "Base ring K (= %s) must be a field of characteristic different from 2."%K
154
155    try: 
156        a = K(a)
157    except TypeError:
158        raise ValueError, "Arguments a = %s and b = %s must coerce into K (= %s)."%(a,b,K)
159    try: 
160        b = K(b)
161    except TypeError:
162        raise ValueError, "Arguments a = %s and b = %s must coerce into K (= %s)."%(a,b,K)
163
164    if a == 0 or b == 0:
165        raise ValueError, "Arguments a = %s and b = %s must be nonzero."%(a,b)
166
167    if denom == 0:
168        raise ValueError, "Argument denom (= %s) must be a nonzero."%denom
169
170    if isinstance(names, list):
171        names = tuple(names)
172
173    key = (K, a, b, names, denom) 
174    global _cache
175    if _cache.has_key(key):
176        A = _cache[key]()
177        if not A is None:
178            return A
179
180    if K is RationalField():
181        prms = ramified_primes(a,b)
182        H = QuaternionAlgebra_generic(K, [2,0,0,0], prms)
183    else:
184        H = QuaternionAlgebra_generic(K, [2,0,0,0])
185    A = FreeAlgebra(K,3, names=names)
186    F = A.monoid()
187    mons = [ F(1) ] + [ F.gen(i) for i in range(3) ] 
188    M = MatrixSpace(K,4)
189    m = denom
190    c = a/m
191    d = b/m
192    mats = [
193        M([0,1,0,0, a,0,0,0, 0,0,0,-m, 0,0,-c,0]),
194        M([0,0,1,0, 0,0,0,m, b,0,0,0, 0,d,0,0]),
195        M([0,0,0,1, 0,0,c,0, 0,-d,0,0, -c*d,0,0,0]) ]
196    FreeAlgebraQuotient.__init__(H, A, mons=mons, mats=mats, names=names)
197    _cache[key] = weakref.ref(H)
198    return H
199
200def QuaternionAlgebraWithInnerProduct(K, norms, traces, names):
201    """
202    """
203    (n1,n2,n3) = norms
204    (t1,t2,t3,t12,t13,t23) = traces
205    T = t1*t23 + t2*t13 - t3*t12
206    N = (t1*t2 - t12)*t13*t23 \
207        + t23*(t23 - t2*t3)*n1 \
208        + t13*(t13 - t1*t3)*n2 \
209        + t12*(t12 - t1*t2)*n3 \
210        + t1**2*n2*n3 + t2**2*n1*n3 + t3**2*n1*n2 - 4*n1*n2*n3
211    if True and K(2).is_unit(): 
212        D = T**2 - 4*N
213        try:
214            S = D.square_root()
215        except ArithmeticError:
216            raise ValueError, "Invalid inner product input."
217        assert bool
218        x = (T + S)/2
219    else:
220        # In characteristic 2 we can't diagonalize the quadratic
221        # form so we solve for its roots.
222        X = PolynomialRing(K).gen()
223        try:
224            x = (X**2 - T*X + N).roots()[0][0]
225        except IndexError:
226            raise ValueError, "Invalid inner product input."
227    assert (x**2 - T*x + N) == 0
228    M = MatrixSpace(K,5)
229    m = M([ [    2,  t1,  t2, t3, t1*t2-t12 ],
230            [   t1, 2*n1,  t12, t13, t2*n1  ],
231            [   t2,  t12, 2*n2, t23, t1*n2  ],
232            [   t3,  t13,  t23, 2*n3,  x ], 
233            [ t1*t2-t12, t2*n1, t1*n2, x, 2*n1*n2] ])
234    v = m.kernel().gen(0)
235    r4 = -1/v[4]
236    V = VectorSpace(K,4)
237    vij = V([ r4*v[i] for i in range(4) ])
238    (s0,s1,s2,s3) = vij.list()
239    r3 = 1/s3
240    vik = r3 * (V([ s1*n1, -(s0+s1*t1), -n1, 0 ]) + (t1-s2)*vij)
241    vkj = r3 * (V([ s2*n2, -n2, -(s0+s2*t2), 0 ]) + (t2-s1)*vij)
242    vji = V([ -t12, t2, t1, 0 ]) - vij
243    vki = V([ -t13, t3, 0, t1 ]) - vik
244    vjk = V([ -t23, 0, t3, t2 ]) - vkj
245    H = QuaternionAlgebra_generic(K, [2,t1,t2,t3])
246    A = FreeAlgebra(K,3, names=names)
247    F = A.monoid()
248    mons = [ F(1) ] + [ F.gen(i) for i in range(3) ] 
249    M = MatrixSpace(K,4)
250    mi = M([ [0,1,0,0], [-n1,t1,0,0], vji.list(), vki.list() ])
251    mj = M([ [0,0,1,0], vij.list(), [-n2,0,t2,0], vkj.list() ])
252    mk = M([ [0,0,0,1], vik.list(), vjk.list(), [-n3,0,0,t3] ])
253    mats = [mi,mj,mk]
254    FreeAlgebraQuotient.__init__(H, A, mons=mons, mats=mats, names=names)
255    return H
256
257def QuaternionAlgebraWithGramMatrix(K, gram, names):
258    """
259    """
260    if not isinstance(gram,Matrix) or not gram.is_symmetric:
261        raise AttributeError, "Argument gram (= %s) must be a symmetric matrix"%gram
262    (q0,t01,t02,t03,_,q1,t12,t13,_,_,q3,t23,_,_,_,q4) = gram.list()
263    n1 = q1/2
264    n2 = q2/2
265    n3 = q3/2
266    return QuaternionAlgebraWithInnerProduct(K,norms,traces,names=names)
267
268def QuaternionAlgebraWithDiscriminants(D1, D2, T, names=['i','j','k'], M=2):
269    r"""
270    Return the quaternion algebra over the rationals generated by $i$,
271    $j$, and $k = (ij - ji)/M$ where $\Z[i]$, $\Z[j]$, and $\Z[k]$ are
272    quadratic suborders of discriminants $D_1$, $D_2$, and $D_3 = (D_1
273    D_2 - T^2)/M^2$, respectively.  The traces of $i$ and $j$ are
274    chosen in $\{0,1\}$.
275
276    The integers $D_1$, $D_2$ and $T$ must all be even or all odd, and
277    $D_1$, $D_2$ and $D_3$ must each be the discriminant of some
278    quadratic order, i.e. nonsquare integers = 0, 1 (mod 4).
279
280    INPUT:
281        D1 -- Integer
282        D2 -- Integer
283        T  -- Integer
284        M -- Integer (default: 2)
285       
286    OUTPUT:
287        A quaternion algebra.
288
289    EXAMPLES:
290        sage: A = QuaternionAlgebraWithDiscriminants(-7,-47,1, names=['i','j','k'])
291        sage: print A
292        Quaternion algebra with generators (i, j, k) over Rational Field
293        sage: i, j, k = A.gens()
294        sage: i**2
295        -2 + i
296        sage: j**2
297        -12 + j
298        sage: k**2
299        -24 + k
300        sage: i.minimal_polynomial('x')
301        x^2 - x + 2
302        sage: j.minimal_polynomial('x')
303        x^2 - x + 12
304    """
305    QQ = RationalField()
306    ZZ = IntegerRing()
307    if M != 2:
308        raise ValueError, "Not implemented: Argument denom (= %s) must be 2."%M
309    # adapt types of input arguments
310    D1 = ZZ(D1); D2 = ZZ(D2); T = ZZ(T)
311    if (D1*D2)%M**2-(T**2)%M**2 != 0:
312        raise ArithmeticError, \
313              "Each of (D1 D2 - T^2) = %s must be 0 mod M^2 = %s"%(D3,M**2)
314    # Check that the D_i are 0 or 1 mod 4 and non-square:
315    if D1.is_square() or D2.is_square():
316        raise ArithmeticError, "D1 (=%s) and D2 (=%s) must be nonsquare."%(D1,D2)
317    t1 = D1%4; t2 = D2%4
318    if t1 in [2,3] or t2 in [2,3]:
319        raise ArithmeticError, "D1 (=%s) and D2 (=%s) must be in {0,1} mod 4."%(D1,D2)
320    n1 = (t1-D1)//4; n2 = (t2-D2)//4
321    t12 = (t1*t2 - T)//2
322    A = FreeAlgebra(RationalField(),3,names=names)
323    i, j, k = A.monoid().gens()
324    # Right matrix action on algebra:
325    MQ = MatrixSpace(QQ,4)
326    mi = MQ([0,1,0,0, -n1,t1,0,0, -t12,t2,t1,-1, -n1*t2,t1*t2-t12,n1,0])
327    mj = MQ([0,0,1,0, 0,0,0,1, -n2,0,t2,0, 0,-n2,0,t2])
328    # N.B. mk = mi*mj
329    t3 = t1*t2-t12
330    mk = MQ([0,0,0,1, 0,0,-n1,t1,  -n2*t1,n2,t3,0, -n1*n2,0,0,t3])
331    mats = [ mi, mj, mk ] 
332    prms = ramified_primes_from_discs(D1,D2,T)
333    H = QuaternionAlgebra_generic(QQ, [2,t1,t2,t3], prms)
334    A = FreeAlgebra(QQ,3, names=names)
335    F = A.monoid()
336    mons = [ F(1) ] + [ F.gen(i) for i in range(3) ] 
337    FreeAlgebraQuotient.__init__(H, A, mons=mons, mats=mats, names=names)
338    return H
339
340class QuaternionAlgebra_generic(FreeAlgebraQuotient):
341
342    def __init__(self, K, basis_traces = None, ramified_primes=None):
343        """
344        """
345        if ramified_primes != None:
346            self._ramified_primes = ramified_primes
347        if basis_traces != None:
348            self._basis_traces = basis_traces
349       
350    def __call__(self, x):
351        if hasattr(self, 'discriminants'):
352            if isinstance(x, QuaternionAlgebraElement_fast) and x.parent() is self: 
353                return x
354            return QuaternionAlgebraElement_fast(self,x)
355        else:
356            if isinstance(x, QuaternionAlgebraElement) and x.parent() is self: 
357                return x
358            return QuaternionAlgebraElement(self,x)
359       
360    def __repr__(self):
361        return "Quaternion algebra with generators %s over %s"%(
362            self.gens(), self.base_ring())
363
364    def gen(self,i):
365        """
366        The i-th generator of the quaternion algebra.
367        """
368        if i < 0 or not i < 4: 
369            raise IndexError, \
370                "Argument i (= %s) must be between 0 and %s."%(i, 3)
371        K = self.base_ring()
372        F = self.free_algebra().monoid()
373        return QuaternionAlgebraElement(self,{F.gen(i):K(1)})
374
375    def basis(self):
376        return (self(1), self([0,1,0,0]), self([0,0,1,0]), self([0,0,0,1]))
377       
378    def discriminant(self):
379        return self.gram_matrix().determinant()
380
381    def gram_matrix(self):
382        """
383        The Gram matrix of the inner product determined by the norm.
384        """
385        try: 
386            M = self.__vector_space.inner_product_matrix()
387        except: 
388            K = self.base_ring()
389            M = MatrixSpace(K,4)(0)
390            B = self.basis()
391            for i in range(4):
392                x = B[i]
393                M[i,i] = 2*(x.reduced_norm())
394                for j in range(i+1,4):
395                    y = B[j]
396                    c = (x * y.conjugate()).reduced_trace()
397                    M[i,j] = c
398                    M[j,i] = c
399            self.__vector_space = VectorSpace(K,4,inner_product_matrix = M)
400        return M
401       
402    def inner_product_matrix(self):
403        return self.gram_matrix()
404
405    def is_commutative(self):
406        """
407        EXAMPLES:
408            sage: Q.<i,j,k> = QuaternionAlgebra(QQ, -3,-7)
409            sage: Q.is_commutative()
410            False
411        """
412        return False
413   
414    def ramified_primes(self):
415        try:
416            return self._ramified_primes
417        except:
418            raise AttributeError, "Ramified primes have not been computed."
419
420    def random_element(self):
421        K = self.base_ring()
422        return self([ K.random_element() for _ in range(4) ])
423
424    def vector_space(self):
425        try:
426            V = self.__vector_space
427        except:
428            M = self.gram_matrix() # induces construction of V
429            V = self.__vector_space
430        return V
431
432
433def QuaternionAlgebra_fast(K, a, b, names="ijk"):
434    H = QuaternionAlgebra(K, a, b, names)
435    H.discriminants = (a,b)
436    return H
437
438class QuaternionAlgebra_faster(QuaternionAlgebra_generic):
439   
440    def __call__(self, x):
441        if isinstance(x, QuaternionAlgebraElement_fast) and x.parent() is self: 
442            return x
443        return QuaternionAlgebraElement_fast(self,x)
Note: See TracBrowser for help on using the repository browser.