Ticket #7577: mpolynomial_rr_libsingular.patch

File mpolynomial_rr_libsingular.patch, 15.7 KB (added by malb, 12 years ago)
  • module_list.py

    # HG changeset patch
    # User Martin Albrecht <malb@informatik.uni-bremen.de>
    # Date 1259708951 0
    # Node ID 043ba443bd743f212d704018423a96d8d1682cd3
    # Parent  75513ee47886beb483985250de3ff630524f7e1b
    moving multivariate polynomials over RR to libSingular
    
    diff -r 75513ee47886 -r 043ba443bd74 module_list.py
    a b  
    448448   
    449449    Extension('sage.libs.singular.singular',
    450450              sources = ['sage/libs/singular/singular.pyx'],
    451               libraries = ['m', 'readline', 'singular', 'givaro', 'gmpxx', 'gmp'],
     451              libraries = ['m', 'readline', 'singular', 'givaro', 'gmpxx', 'gmp', 'mpfr'],
    452452              language="c++",
    453453              include_dirs = [SAGE_ROOT +'/local/include/singular'],
    454454              depends = [SAGE_ROOT + "/local/include/libsingular.h"]),
  • sage/libs/mpfr.pxd

    diff -r 75513ee47886 -r 043ba443bd74 sage/libs/mpfr.pxd
    a b  
    11
    22from sage.libs.gmp.all cimport mp_exp_t, mp_size_t, gmp_randstate_t, mpz_t, mpq_t
    33
     4from sage.libs.gmp.types cimport mpf_t
     5
    46cdef extern from "mpfr.h":
    57    ctypedef struct __mpfr_struct:
    68        pass
     
    4244    # int mpfr_set_decimal64 (mpfr_t rop, _Decimal64 op, mp_rnd_t rnd)
    4345    int mpfr_set_z (mpfr_t rop, mpz_t op, mp_rnd_t rnd)
    4446    int mpfr_set_q (mpfr_t rop, mpq_t op, mp_rnd_t rnd)
    45     # int mpfr_set_f (mpfr_t rop, mpf_t op, mp_rnd_t rnd)
     47    int mpfr_set_f (mpfr_t rop, mpf_t op, mp_rnd_t rnd)
    4648    int mpfr_set_ui_2exp (mpfr_t rop, unsigned long int op, mp_exp_t e, mp_rnd_t rnd)
    4749    int mpfr_set_si_2exp (mpfr_t rop, long int op, mp_exp_t e, mp_rnd_t rnd)
    4850    # int mpfr_set_uj_2exp (mpfr_t rop, uintmax_t op, intmax_t e, mp_rnd_t rnd)
     
    7678    # uintmax_t mpfr_get_uj (mpfr_t op, mp_rnd_t rnd)
    7779    mp_exp_t mpfr_get_z_exp (mpz_t rop, mpfr_t op)
    7880    void mpfr_get_z (mpz_t rop, mpfr_t op, mp_rnd_t rnd)
    79     # int mpfr_get_f (mpf_t rop, mpfr_t op, mp_rnd_t rnd)
     81    int mpfr_get_f (mpf_t rop, mpfr_t op, mp_rnd_t rnd)
    8082    char * mpfr_get_str (char *str, mp_exp_t *expptr, int b, size_t n, mpfr_t op, mp_rnd_t rnd)
    8183    void mpfr_free_str (char *str)
    8284    bint mpfr_fits_ulong_p (mpfr_t op, mp_rnd_t rnd)
  • sage/libs/singular/ring.pyx

    diff -r 75513ee47886 -r 043ba443bd74 sage/libs/singular/ring.pyx
    a b  
    3030from sage.rings.number_field.number_field_base cimport NumberField
    3131from sage.rings.rational_field import RationalField
    3232from sage.rings.ring import FiniteField as FiniteField_generic
     33from sage.rings.real_mpfr import is_RealField
    3334
    3435from sage.rings.polynomial.term_order import TermOrder
    3536from sage.rings.polynomial.multi_polynomial_libsingular cimport MPolynomial_libsingular, MPolynomialRing_libsingular
     
    209210            ringflaga = <__mpz_struct*>omAlloc(sizeof(__mpz_struct))
    210211            mpz_init_set_ui(ringflaga, characteristic)
    211212            ringflagb = 1
     213    elif is_RealField(base_ring):
     214        characteristic = -1
    212215    else:
    213216        raise NotImplementedError("Base ring is not supported.")
    214217
     
    233236
    234237        _ring.minpoly=<number*>nmp
    235238
     239    if characteristic == -1:
     240        _ring.float_len = base_ring.prec()/3.5
     241        _ring.float_len2 = 0
     242
    236243    nblcks = len(order.blocks)
    237244    offset = 0
    238245
  • sage/libs/singular/singular-cdefs.pxi

    diff -r 75513ee47886 -r 043ba443bd74 sage/libs/singular/singular-cdefs.pxi
    a b  
    146146        unsigned long ringflagb
    147147        int pCompIndex # index of components
    148148
     149        short float_len
     150        short float_len2
     151
    149152        n_Procs_s*    cf
    150153        int ref
    151154
     
    960963cdef extern from "syz.h":
    961964    ctypedef struct syStrategy "ssyStrategy":
    962965        short references
     966
     967
     968cdef extern from "mpr_complex.h":
     969    ctypedef struct gmp_float "gmp_float":
     970        mpf_t* (*mpfp)()
     971
     972    gmp_float * new_gmp_float_from_mpf "new gmp_float"(mpf_t val)
  • sage/libs/singular/singular.pyx

    diff -r 75513ee47886 -r 043ba443bd74 sage/libs/singular/singular.pyx
    a b  
    4747from sage.structure.parent_base cimport ParentWithBase
    4848from sage.rings.polynomial.multi_polynomial_libsingular cimport MPolynomial_libsingular
    4949
     50from sage.libs.mpfr cimport mpfr_set_f, mpfr_set_f, mpfr_init2, mpfr_t, mpfr_rnd_t, mpfr_get_f
     51from sage.libs.gmp.types cimport mpf_t
     52from sage.libs.gmp.mpf cimport mpf_init2
     53from sage.libs.singular.decl cimport gmp_float, new_gmp_float_from_mpf
     54from sage.rings.real_mpfr cimport RealField, RealNumber
     55
     56cdef RealNumber si2sa_RR(number *n, ring *_ring, RealField base):
     57    """
     58    TESTS::
     59
     60        sage: P.<x,y,z> = RR[]
     61        sage: P(0.3).lc()
     62        0.300000000000000
     63        sage: P(1).lc()
     64        1.00000000000000
     65        sage: P(0).lc()
     66        0.000000000000000
     67        sage: P(-0.3).lc()
     68        -0.300000000000000
     69        sage: type(P(3.0).lc())
     70        <type 'sage.rings.real_mpfr.RealNumber'>
     71    """
     72    cdef RealNumber ret = PY_NEW(RealNumber)
     73    mpfr_init2(ret.value, base.prec())
     74    ret._parent = base
     75    cdef gmp_float *_n = <gmp_float*>(n)
     76    cdef mpf_t *__n = <mpf_t *>_n.mpfp()
     77    mpfr_set_f(ret.value,__n[0], base.rnd)
     78    return ret
     79
    5080cdef Rational si2sa_QQ(number *n, ring *_ring):
    5181    """
    5282    TESTS::
     
    315345
    316346    return base(_ring.cf.nInt(n))
    317347
     348cdef number * sa2si_RR(RealNumber n, ring *_ring):
     349    """
     350    TESTS::
     351
     352        sage: P.<x,y,z> = RR[]
     353        sage: P(0) + P(0.5) - 0.5
     354        0
     355        sage: P(0.5) + 3.0/5.0 - 3.0/5.0
     356        0.5
     357        sage: P(2.0/3.0) + 1.0/4.0 - 1/4.0
     358        0.666666666666667
     359        sage: P(12345678901234567890/23.0) + 5/2 - 5/2
     360        0.536768647879764e + 18
     361    """
     362    cdef mpf_t _r
     363    mpf_init2(_r, (<RealField>n._parent).prec())
     364    mpfr_get_f(_r, n.value, (<RealField>n._parent).rnd)
     365    cdef gmp_float *_n = new_gmp_float_from_mpf(_r)
     366    return <number*>_n
     367   
    318368cdef number *sa2si_QQ(Rational r, ring *_ring):
    319369    """
    320370    TESTS::
     
    563613            return base(nInt(n))
    564614        return si2sa_ZZmod(n, _ring, base)
    565615
     616    elif PY_TYPE_CHECK(base, RealField):
     617        return si2sa_RR(n, _ring, base)
     618
    566619    else:
    567620        raise ValueError, "cannot convert from SINGULAR number"
    568621
     
    592645        if _ring.ringtype == 0:
    593646            return n_Init(int(elem),_ring)
    594647        return sa2si_ZZmod(elem, _ring)
     648    elif PY_TYPE_CHECK(elem._parent, RealField):
     649        return sa2si_RR(elem, _ring)
     650   
    595651    else:
    596652        raise ValueError, "cannot convert to SINGULAR number"
    597653
  • sage/rings/polynomial/multi_polynomial.pyx

    diff -r 75513ee47886 -r 043ba443bd74 sage/rings/polynomial/multi_polynomial.pyx
    a b  
    2525        TESTS::
    2626
    2727            sage: type(RR['x,y'])
    28             <class 'sage.rings.polynomial.multi_polynomial_ring.MPolynomialRing_polydict_domain'>
     28            <type 'sage.rings.polynomial.multi_polynomial_libsingular.MPolynomialRing_libsingular'>
    2929            sage: type(RR['x, y'](0))
    30             <class 'sage.rings.polynomial.multi_polynomial_element.MPolynomial_polydict'>
     30            <type 'sage.rings.polynomial.multi_polynomial_libsingular.MPolynomial_libsingular'>
    3131
    3232            sage: int(RR['x,y'](0)) # indirect doctest
    3333            0
  • sage/rings/polynomial/multi_polynomial_element.py

    diff -r 75513ee47886 -r 043ba443bd74 sage/rings/polynomial/multi_polynomial_element.py
    a b  
    700700            sage: R.<x,y> = RR[]
    701701            sage: f=x*y+5
    702702            sage: c=f.coefficient({x:0,y:0}); c
    703             5.00000000000000
     703            5
    704704            sage: parent(c)
    705705            Multivariate Polynomial Ring in x, y over Real Field with 53 bits of precision
    706706       
  • sage/rings/polynomial/multi_polynomial_ideal.py

    diff -r 75513ee47886 -r 043ba443bd74 sage/rings/polynomial/multi_polynomial_ideal.py
    a b  
    18861886        There are two real intersections::
    18871887       
    18881888            sage: I.variety(ring=RR)
    1889             [{y: 0.361103080528647, x: 2.76929235423863},
    1890              {y: 1.00000000000000, x: 1.00000000000000}]
     1889            [{y: 1.00000000000000, x: 1.00000000000000},
     1890             {y: 0.361103080528647, x: 2.76929235423863}]
    18911891            sage: I.variety(ring=AA)
    18921892            [{x: 2.769292354238632?, y: 0.3611030805286474?},
    18931893             {x: 1, y: 1}]
  • sage/rings/polynomial/multi_polynomial_libsingular.pyx

    diff -r 75513ee47886 -r 043ba443bd74 sage/rings/polynomial/multi_polynomial_libsingular.pyx
    a b  
    199199from sage.structure.element cimport RingElement
    200200from sage.structure.element cimport ModuleElement
    201201from sage.structure.element cimport Element
    202 from sage.structure.element cimport CommutativeRingElement
     202from sage.structure.element cimport RingElement
    203203
    204204from sage.structure.factorization import Factorization
    205205from sage.structure.sequence import Sequence
     
    426426            else:
    427427                raise TypeError("incompatable parents.")
    428428           
    429         elif PY_TYPE_CHECK(element, CommutativeRingElement):
     429        elif PY_TYPE_CHECK(element, RingElement):
    430430            # base ring elements
    431431            if  <Parent>element.parent() is base_ring:
    432432                # shortcut for GF(p)
     
    10771077            # singular converts to bits from base_10 in mpr_complex.cc by:
    10781078            #  size_t bits = 1 + (size_t) ((float)digits * 3.5);
    10791079            precision = base_ring.precision()
    1080             digits = sage.rings.arith.ceil((2*precision - 2)/7.0)
     1080            digits = sage.functions.other.ceil((2*precision - 2)/7.0)
    10811081            self.__singular = singular.ring("(real,%d,0)"%digits, _vars, order=order)
    10821082
    10831083        elif is_ComplexField(base_ring):
    10841084            # singular converts to bits from base_10 in mpr_complex.cc by:
    10851085            #  size_t bits = 1 + (size_t) ((float)digits * 3.5);
    10861086            precision = base_ring.precision()
    1087             digits = sage.rings.arith.ceil((2*precision - 2)/7.0)
     1087            digits = sage.functions.other.ceil((2*precision - 2)/7.0)
    10881088            self.__singular = singular.ring("(complex,%d,0,I)"%digits, _vars,  order=order)
    10891089
    10901090        elif base_ring.is_prime_field():
  • sage/rings/real_mpfr.pxd

    diff -r 75513ee47886 -r 043ba443bd74 sage/rings/real_mpfr.pxd
    a b  
    44from sage.libs.mpfr cimport *
    55
    66include '../ext/cdefs.pxi'
    7 include '../libs/pari/decl.pxi'
    87
    98cimport sage.rings.ring
    109import  sage.rings.ring
     
    2726    cdef char init
    2827    cdef RealNumber _new(self)
    2928    cdef _set(self, x, int base)
    30     cdef _set_from_GEN_REAL(self, GEN g)
    3129    cdef RealNumber abs(RealNumber self)
  • sage/rings/real_mpfr.pyx

    diff -r 75513ee47886 -r 043ba443bd74 sage/rings/real_mpfr.pyx
    a b  
    116116include '../ext/interrupt.pxi'
    117117include "../ext/stdsage.pxi"
    118118include "../ext/random.pxi"
     119include '../libs/pari/decl.pxi'
    119120
    120121cimport sage.rings.ring
    121122import  sage.rings.ring
     
    249250        RealField_cache[prec, sci_not, rnd] = R = RealField(prec=prec, sci_not=sci_not, rnd=rnd)
    250251        return R
    251252
     253
    252254cdef class RealField(sage.rings.ring.Field):
    253255    """
    254256    An approximation to the field of real numbers using floating point
     
    858860
    859861cdef class RealLiteral(RealNumber)
    860862
     863cdef _set_from_GEN_REAL(RealNumber self, GEN g):
     864    """
     865    EXAMPLES::
     866
     867        sage: rt2 = sqrt(pari('2.0'))
     868        sage: rt2
     869        1.41421356237310
     870        sage: rt2.python()
     871        1.414213562373095048801688724              # 32-bit
     872        1.4142135623730950488016887242096980786    # 64-bit
     873        sage: rt2.python().prec()
     874        96                                         # 32-bit
     875        128                                        # 64-bit
     876        sage: pari(rt2.python()) == rt2
     877        True
     878        sage: for i in xrange(1, 1000):
     879        ...       assert(sqrt(pari(i)) == pari(sqrt(pari(i)).python()))
     880        sage: (-3.1415)._pari_().python()
     881        -3.14150000000000000
     882    """
     883    cdef int sgn
     884    sgn = signe(g)
     885
     886    if sgn == 0:
     887        mpfr_set_ui(self.value, 0, GMP_RNDN)
     888        return
     889
     890    cdef int wordsize
     891
     892    if sage.misc.misc.is_64_bit:
     893        wordsize = 64
     894    else:
     895        wordsize = 32
     896
     897    cdef mpz_t mantissa
     898    mpz_init(mantissa)
     899
     900    mpz_import(mantissa, lg(g) - 2, 1, wordsize/8, 0, 0, &g[2])
     901
     902    cdef int exponent
     903    exponent = expo(g)
     904
     905    # Round to nearest for best results when setting a low-precision
     906    # MPFR from a high-precision GEN
     907    mpfr_set_z(self.value, mantissa, GMP_RNDN)
     908    mpfr_mul_2si(self.value, self.value, exponent - wordsize * (lg(g) - 2) + 1, GMP_RNDN)
     909
     910    if sgn < 0:
     911        mpfr_neg(self.value, self.value, GMP_RNDN)
     912
     913    mpz_clear(mantissa)
     914
    861915cdef class RealNumber(sage.structure.element.RingElement):
    862916    """
    863917    A floating point approximation to a real number using any specified
     
    9861040            mpfr_set_q(self.value, (<Rational>x).value, parent.rnd)
    9871041        elif PY_TYPE_CHECK(x, gen) and typ((<gen>x).g) == t_REAL:
    9881042            _gen = x
    989             self._set_from_GEN_REAL(_gen.g)
     1043            _set_from_GEN_REAL(self, _gen.g)
    9901044        elif isinstance(x, int):
    9911045            mpfr_set_si(self.value, x, parent.rnd)
    9921046        elif isinstance(x, long):
     
    10081062                else:
    10091063                    raise TypeError, "Unable to convert x (='%s') to real number."%s
    10101064               
    1011     cdef _set_from_GEN_REAL(self, GEN g):
    1012         """
    1013         EXAMPLES::
    1014        
    1015             sage: rt2 = sqrt(pari('2.0'))
    1016             sage: rt2
    1017             1.41421356237310
    1018             sage: rt2.python()
    1019             1.414213562373095048801688724              # 32-bit
    1020             1.4142135623730950488016887242096980786    # 64-bit
    1021             sage: rt2.python().prec()
    1022             96                                         # 32-bit
    1023             128                                        # 64-bit
    1024             sage: pari(rt2.python()) == rt2
    1025             True
    1026             sage: for i in xrange(1, 1000):
    1027             ...       assert(sqrt(pari(i)) == pari(sqrt(pari(i)).python()))
    1028             sage: (-3.1415)._pari_().python()
    1029             -3.14150000000000000
    1030         """
    1031         cdef int sgn
    1032         sgn = signe(g)
    1033 
    1034         if sgn == 0:
    1035             mpfr_set_ui(self.value, 0, GMP_RNDN)
    1036             return
    1037 
    1038         cdef int wordsize
    1039 
    1040         if sage.misc.misc.is_64_bit:
    1041             wordsize = 64
    1042         else:
    1043             wordsize = 32
    1044 
    1045         cdef mpz_t mantissa
    1046         mpz_init(mantissa)
    1047 
    1048         mpz_import(mantissa, lg(g) - 2, 1, wordsize/8, 0, 0, &g[2])
    1049 
    1050         cdef int exponent
    1051         exponent = expo(g)
    1052 
    1053         # Round to nearest for best results when setting a low-precision
    1054         # MPFR from a high-precision GEN
    1055         mpfr_set_z(self.value, mantissa, GMP_RNDN)
    1056         mpfr_mul_2si(self.value, self.value, exponent - wordsize * (lg(g) - 2) + 1, GMP_RNDN)
    1057 
    1058         if sgn < 0:
    1059             mpfr_neg(self.value, self.value, GMP_RNDN)
    1060 
    1061         mpz_clear(mantissa)
    1062 
    10631065    def __reduce__(self):
    10641066        """
    10651067        EXAMPLES::