Laboratoire de l'Informatique et du Parallélisme » Pobyso

#cimport cython

from libc.stdlib          cimport malloc,free

from sage.rings.real_mpfr cimport *

from sage.rings.real_mpfi cimport *

#

# Functions to manipulate the real RealNumber values.

#

# Import constants and declarations that are not exported

# from the Sage real_mpfr.pyx file.

# Notes:

# - notice how #define constants are imported.

# - before mpfr_min_prec() is declared, we need to "ctypedef" mpfr_prec_t

#   if we want to be consistent with MPFR notations. The actual type comes

#   from mpfr.h but, according to documentation, any type of same general

#   kind would have been OK.

cdef extern from "mpfr.h":

    ctypedef unsigned short mpfr_prec_t

    cdef mpfr_prec_t MPFR_PREC_MIN "MPFR_PREC_MIN"

    cdef mpfr_prec_t mpfr_min_prec(mpfr_t x)

cpdef int get_interval_value(RealIntervalFieldElement x):

    """

    Get the address of the value of a RealIntervalField Element.

    This address can be used to to set the value as well.

    """

    return(<int>&(x.value))

cpdef int get_rn_value(RealNumber x):

    """

    Get the address of the value of a RealNumber. This address can be

    used to to set the value as well.

    """

    return(<int>&(x.value))

cpdef set_interval_value(RealIntervalFieldElement x, int p):

    """

    Set the value of x from MPFI interval p.

    """

    # We use the "array trick" to workaround the pointer dereferencing

    # issue (in Cython there is no unary operator "*").

    #

    mpfi_set(x.value, (<mpfi_t *>p)[0])

cpdef set_rn_value(RealNumber x, int p, mp_rnd_t rnd=GMP_RNDN):

    """

    Set the value of x from MPFR number p.

    """

    # We use the "array trick" to workaround the pointer dereferencing

    # issue (in Cython there is no unary operator "*").

    mpfr_set(x.value, (<mpfr_t *>p)[0], rnd)

#

cpdef mpfr_prec_t min_rn_size(RealNumber x):

    """

    Compute the minimum number of bits necessary to represent

    the x Sage RealNumber's value.

    If x is a power of 2 (only 1 bit needed) return the MPFR

    MPFR_PREC_MIN value.

    """

    cdef mpfr_prec_t min_prec = mpfr_min_prec(x.value)

    if min_prec < MPFR_PREC_MIN:

        return(MPFR_PREC_MIN)

    else:

        return(min_prec)

#

cpdef min_mpfr_size(int p):

    """

    Compute the minimum number of bits necessary to represent

    the MPFR value pointed by the mpfr_t pointer p (hence the int type).

    If *p is a power of 2 (only 1 bit needed) return the MPFR

    MPFR_PREC_MIN value.

    """

    # We use the "array trick" to workaround the pointer dereferencing

    # issue (in Cython there is no unary operator "*").

    cdef mpfr_prec_t min_prec = mpfr_min_prec((<mpfr_t *>p)[0])

    if min_prec < MPFR_PREC_MIN:

        return(MPFR_PREC_MIN)

    else:

        return(min_prec)

cpdef free_mpfi_num(int mpfiNum):

    """

    Free a MPFI interval pointed by mpfiNum (hece the int type).

    """

    # We use the "array trick" to workaround the pointer dereferencing

    # issue (in Cython there is no unary operator "*").

    mpfi_clear((<mpfi_t*>mpfiNum)[0])

cpdef free_mpfr_num(int mpfrNum):

    """

    Free a MPFR number pointed by mpfrNum (hece the int type).

    """

    # We use the "array trick" to workaround the pointer dereferencing

    # issue (in Cython there is no unary operator "*").

    mpfr_clear((<mpfr_t*>mpfrNum)[0])

cpdef int new_interval_num(unsigned int prec):

    """

    Create a newly allocated MPFI interval to precision prec.

    """

    cdef mpfi_t *newValue = NULL

    print "Precision : ", prec

    newValue = <mpfi_t*>malloc(sizeof(mpfi_t))

    mpfi_init2(newValue[0], prec)

    print "Creation done!"

    return(<int>newValue)

cpdef int new_mpfr_num(unsigned int prec):

    cdef mpfr_t *newValue = NULL

    print "Precision : ", prec

    newValue = <mpfr_t*>malloc(sizeof(mpfr_t))

    mpfr_init2(newValue[0], prec)

    print "Creation done!"

    return(<int>newValue)

#

cpdef cmp_rn_value(RealNumber x, RealNumber y):

    """

    Compare two RN with the mpfr_cmp function

    """

    return(mpfr_cmp(x.value, y.value))

cpdef cmp_interval_value(RealIntervalFieldElement x,

                         RealIntervalFieldElement y):

    """

    Compare two intervals with the mpfi_cmp function

    """

    return(mpfi_cmp(x.value, y.value))

# Do not use this function!

# It changes the precision of the instance behind the curtain, but precision is

# class attribute. If the precision is increased, memory space is wasted and

# other problems may surface in subsequent operations.

# If the precision is decreased the semantics

# are broken since the expected precision (from the class attribute) is not

# the actual precision.

#

cpdef s_rn_value(RealNumber x, int p):

    mpfr_clear(x.value)

    mpfr_init2(x.value, mpfr_get_prec((<mpfr_t *>p)[0]))

    #x.__prec = mpfr_get_prec(x.value)

    printf("x.value prec: %d\n", mpfr_get_prec(x.value))

    mpfr_set(x.value, (<mpfr_t *>p)[0], GMP_RNDN)