Ticket #8938: trac_8938.patch

File trac_8938.patch, 11.8 KB (added by Francis Clarke, 12 years ago)
  • sage/libs/singular/polynomial.pyx

    # HG changeset patch
    # User Francis Clarke <F.Clarke@Swansea.ac.uk>
    # Date 1273438278 -3600
    # Node ID 8e87962ff463035428df2bfcff94b088d85eb770
    # Parent  e2ccb846f2962cbe254f534ececfd0fbc9ff5045
    #8938: fix latex for multivariate polynomials and elements of QQbar
    
    diff -r e2ccb846f296 -r 8e87962ff463 sage/libs/singular/polynomial.pyx
    a b  
    378378    - ``p`` - a Singular polynomial
    379379    - ``r`` - a Singular ring
    380380
    381     EXAMPLE::
     381    EXAMPLES::
    382382
    383383        sage: P.<x,y,z> = QQ[]
    384384        sage: latex(x) # indirect doctest
    385385        x
    386386        sage: latex(10*x^2 + 1/2*y)
    387387        10 x^{2} + \frac{1}{2} y
     388        sage: C5.<z> = CyclotomicField(5)
     389        sage: P.<s, t> = C5[]
     390        sage: latex(s + z^3*s + z*s*t + 2 + z) # indirect doctest
     391        z s t + \left(z^{3} + 1\right) s + z + 2
     392
    388393    """
    389     poly = ""
     394    poly = " "
    390395    cdef long e,j
    391396    cdef int n = r.N
     397    atomic_repr = base.is_atomic_repr()
    392398    while p:
    393         sign_switch = False
    394399
    395400        # First determine the multinomial:
    396401        multi = ""
     
    401406            if e > 1:
    402407                multi += "^{%d}"%e
    403408        multi = multi.lstrip().rstrip()
     409        if len(multi) > 0:
     410            multi = "|" + multi
    404411
    405         # Next determine coefficient of multinomial
     412        # Now determine the coefficient of the multinomial
    406413        c =  si2sa( p_GetCoeff(p, r), r, base)
    407         if len(multi) == 0:
    408             multi = latex(c)
    409         elif c != 1:
    410             if  c == -1:
    411                 if len(poly) > 0:
    412                     sign_switch = True
    413                 else:
    414                     multi = "- %s"%(multi)
    415             else:
    416                 multi = "%s %s"%(latex(c),multi)
     414        x = y = latex(c)
     415        if len(poly) > 1:
     416            poly += " + "
     417        if y.find("-") == 0:
     418            y = y[1:]
     419        if not atomic_repr and len(multi) > 0 and (y.find("+") != -1 or y.find("-") != -1):
     420            x = "\\left(%s\\right)" % x
    417421
    418         # Now add on coefficiented multinomials
    419         if len(poly) > 0:
    420             if sign_switch:
    421                 poly = poly + " - "
    422             else:
    423                 poly = poly + " + "
    424         poly = poly + multi
     422        #  and add on the term
     423        poly += "%s %s"%(x, multi)
    425424
    426425        p = pNext(p)
    427426
     427    poly = poly.replace(" + -", " - ")
     428    poly = re.sub(" 1(\.0+)? \|", " ", poly)
     429    poly = re.sub(" -1(\.0+)? \|", " -", poly)
     430    poly = poly.replace("|", "")
    428431    poly = poly.lstrip().rstrip()
    429     poly = poly.replace("+ -","- ")
    430432
    431433    if len(poly) == 0:
    432434        return "0"
  • sage/rings/polynomial/multi_polynomial_libsingular.pyx

    diff -r e2ccb846f296 -r 8e87962ff463 sage/rings/polynomial/multi_polynomial_libsingular.pyx
    a b  
    21382138            sage: P.<x,y,z> = QQ[]
    21392139            sage: f = - 1*x^2*y - 25/27 * y^3 - z^2
    21402140            sage: latex(f)
    2141             - x^{2} y - \frac{25}{27} y^{3} - z^{2}
     2141            -x^{2} y - \frac{25}{27} y^{3} - z^{2}
     2142            sage: C3.<w> = CyclotomicField(3)
     2143            sage: P.<s, t> = C3[]
     2144            sage: latex(s + w*s - w*s*t + 2*w)
     2145            -w s t + \left(w + 1\right) s + 2 w
    21422146        """
    21432147        cdef ring *_ring = (<MPolynomialRing_libsingular>self._parent)._ring
    21442148        gens = self.parent().latex_variable_names()
  • sage/rings/polynomial/polydict.pyx

    diff -r e2ccb846f296 -r 8e87962ff463 sage/rings/polynomial/polydict.pyx
    a b  
    4545include "../../ext/stdsage.pxi"
    4646include '../../ext/cdefs.pxi'
    4747
     48import re
    4849import copy
    4950from sage.structure.element import generic_power
    5051from sage.misc.misc import cputime
     
    388389            f[var] += i
    389390            H[ETuple(f)] = val
    390391        return PolyDict(H, zero=self.__zero, force_etuples=False)
    391            
    392     def latex(PolyDict self, vars, atomic_exponents=True, atomic_coefficients=True, cmpfn = None):
     392               
     393    def latex(PolyDict self, vars, atomic_coefficients=True, cmpfn = None):
    393394        """
    394395        Return a nice polynomial latex representation of this PolyDict, where
    395396        the vars are substituted in.
     
    397398        INPUT:
    398399       
    399400        - ``vars`` -- list
    400         - ``atomic_exponents`` -- bool (default: True)
    401401        - ``atomic_coefficients`` -- bool (default: True)
    402        
     402        - ``cmpfn`` -- comparison function for ordering monomials
     403
    403404        EXAMPLES::
    404405
    405406            sage: from sage.rings.polynomial.polydict import PolyDict       
    406407            sage: f = PolyDict({(2,3):2, (1,2):3, (2,1):4})
    407408            sage: f.latex(['a','WW'])
    408409            '2 a^{2}WW^{3} + 4 a^{2}WW + 3 aWW^{2}'
    409 
    410         When ``atomic_exponents`` is False, the exponents are surrounded in
    411         parenthesis, since ^ has such high precedence::
    412 
    413             # I've removed fractional exponent support in ETuple when moving to a sparse C integer array
    414             #sage: f = PolyDict({(2/3,3,5):2, (1,2,1):3, (2,1,1):4}, force_int_exponents=False)       
    415             #sage: f.latex(['a','b','c'], atomic_exponents=False)
    416             #'4 a^{2}bc + 3 ab^{2}c + 2 a^{2/3}b^{3}c^{5}'
    417410        """
    418         n = len(vars)
    419         poly = ""
     411        poly = " "
    420412        E = self.__repn.keys()
    421         if cmpfn:
    422             E.sort(cmp = cmpfn, reverse=True)
    423         else:
    424             E.sort(reverse=True)
    425         try:
    426             pos_one = self.__zero.parent()(1)
    427             neg_one = -pos_one
    428         except AttributeError:
    429             # probably self.__zero is not a ring element
    430             pos_one = 1
    431             neg_one = -1
     413        E.sort(cmp = cmpfn, reverse=True)
    432414        for e in E:
    433415            c = self.__repn[e]
    434             if c != self.__zero:
    435                 sign_switch = False
     416            x = y = latex(c)
     417            if x != '0':
    436418                # First determine the multinomial:
    437419                multi = ""
    438420                for j in e.nonzero_positions(sort=True):
    439                     multi = multi + vars[j]
     421                    multi += vars[j]
    440422                    if e[j] != 1:
    441                         multi =  multi +"^{%s}"%e[j]
    442                 # Next determine coefficient of multinomial
    443                 if len(multi) == 0:
    444                     multi = latex(c)
    445                 elif c == neg_one:
    446                     # handle -1 specially because it's a pain
    447                     if len(poly) > 0:
    448                         sign_switch = True
    449                     else:
    450                         multi = "-%s"%(multi)
    451                 elif c != pos_one:
    452                     if not atomic_coefficients:
    453                         c = latex(c)
    454                         if c.find("+") != -1 or c.find("-") != -1 or c.find(" ") != -1:
    455                             c = "(%s)"%c
    456                     multi = "%s %s"%(c,multi)
     423                        multi += "^{%s}" % e[j]
     424                multi = multi.lstrip().rstrip()
     425                if len(multi) > 0:
     426                    multi = "|" + multi
    457427
    458                 # Now add on coefficiented multinomials
    459                 if len(poly) > 0:
    460                     if sign_switch:
    461                         poly = poly + " - "
    462                     else:
    463                         poly = poly + " + "
    464                 poly = poly + multi
     428                # Now determine the coefficient of the multinomial
     429                if len(poly) > 1:
     430                    poly += " + "
     431                if y.find("-") == 0:
     432                    y = y[1:]
     433                if not atomic_coefficients and len(multi) > 0 and (y.find("+") != -1 or y.find("-") != -1):
     434                    x = "\\left(%s\\right)" % x
     435
     436                #  and add on the term
     437                poly += "%s %s"%(x, multi)
     438
     439        poly = poly.replace(" + -", " - ")
     440        poly = re.sub(" 1(\.0+)? \|", " ", poly)
     441        poly = re.sub(" -1(\.0+)? \|", " -", poly)
     442        poly = poly.replace("|", "")
    465443        poly = poly.lstrip().rstrip()
    466         poly = poly.replace("+ -","- ")
     444
    467445        if len(poly) == 0:
    468446            return "0"
    469447        return poly
    470    
    471448
    472449    def poly_repr(PolyDict self, vars, atomic_exponents=True, atomic_coefficients=True, cmpfn = None):
    473450        """
  • sage/rings/qqbar.py

    diff -r e2ccb846f296 -r 8e87962ff463 sage/rings/qqbar.py
    a b  
    495495import infinity
    496496from sage.misc.functional import cyclotomic_polynomial
    497497
     498from sage.misc.latex import latex
     499
    498500CC = ComplexField()
    499501CIF = ComplexIntervalField()
    500502
     
    22082210        else:
    22092211            return repr(RIF(self._value))
    22102212
     2213    def _latex_(self):
     2214        """
     2215        Returns the latex representation of this number.
     2216
     2217        EXAMPLES::
     2218
     2219            sage: latex(AA(22/7))
     2220            \frac{22}{7}
     2221            sage: latex(QQbar(I))
     2222            i
     2223            sage: latex(-QQbar.zeta(4) + 5)
     2224            -i + 5
     2225            sage: latex(QQbar.zeta(17))
     2226            0.9324722294043558? + 0.3612416661871530?i
     2227            sage: latex(AA(19).sqrt())
     2228            4.358898943540674?
     2229        """
     2230        if self._descr.is_rational():
     2231            return latex(self._descr._value)
     2232        if isinstance(self._descr, ANRootOfUnity) and self._descr._angle == QQ_1_4:
     2233            x = latex(self._descr._scale)
     2234            if x == '1' or x == '-1':
     2235                x = x[:-1]
     2236            return '%si' % x
     2237        if isinstance(self._descr, ANExtensionElement) and self._descr._generator is QQbar_I_generator:
     2238            s = latex(self._descr._value)
     2239            return s.replace('I', 'i')
     2240        if self.parent() is QQbar:
     2241            return latex(CIF(self._value))
     2242        else:
     2243            return latex(RIF(self._value))
     2244
    22112245    def _sage_input_(self, sib, coerce):
    22122246        r"""
    22132247        Produce an expression which will reproduce this value when evaluated.
  • sage/schemes/generic/affine_space.py

    diff -r e2ccb846f296 -r 8e87962ff463 sage/schemes/generic/affine_space.py
    a b  
    351351
    352352            sage: A.<x, y> = AffineSpace(2, ZZ)
    353353            sage: A._latex_generic_point([y-x^2])
    354             '\\left(- x^{2} + y\\right)'
     354            '\\left(-x^{2} + y\\right)'
    355355            sage: A._latex_generic_point()
    356356            '\\left(x, y\\right)'
    357357        """
  • sage/schemes/generic/projective_space.py

    diff -r e2ccb846f296 -r 8e87962ff463 sage/schemes/generic/projective_space.py
    a b  
    414414
    415415            sage: P.<x, y, z> = ProjectiveSpace(2, ZZ)
    416416            sage: P._latex_generic_point([z*y-x^2])
    417             '\\left(- x^{2} + y z\\right)'
     417            '\\left(-x^{2} + y z\\right)'
    418418            sage: P._latex_generic_point()
    419419            '\\left(x : y : z\\right)'
    420420        """
  • sage/schemes/plane_curves/curve.py

    diff -r e2ccb846f296 -r 8e87962ff463 sage/schemes/plane_curves/curve.py
    a b  
    2323            sage: x,y,z = PolynomialRing(QQ, 3, names='x,y,z').gens()
    2424            sage: C = Curve(y^2*z - x^3 - 17*x*z^2 + y*z^2)
    2525            sage: latex(C)
    26             - x^{3} + y^{2} z - 17 x z^{2} + y z^{2}
     26            -x^{3} + y^{2} z - 17 x z^{2} + y z^{2}
    2727
    2828            sage: A2 = AffineSpace(2, QQ, names=['x','y'])
    2929            sage: x, y = A2.coordinate_ring().gens()
    3030            sage: C = Curve(y^2 - x^3 - 17*x + y)
    3131            sage: latex(C)
    32             - x^{3} + y^{2} - 17 x + y
     32            -x^{3} + y^{2} - 17 x + y
    3333        """
    3434        return latex(self.defining_polynomial())
    3535