"""
@file pobyso.py
Actual functions to use in Sage
@author S.T. 
@date 2012-11-13

Command line syntax:
  use from Sage (via the "load" or the "attach" commands)

pobyso functions come in five flavors: 
- the _so_so (arguments and returned objects are pointers to Sollya objects, 
  includes the void function and the no arguments function that return a 
  pointer to a Sollya object);
- the _so_sa (argument are pointers to Sollya objects, returned objects are
  Sage/Python objects or, more generally, information is transfered from the
  Sollya world to Sage/Python world; e.g. functions without arguments that
  return a Sage/Python object);
- the _sa_so (arguments are Sage/Python objects, returned objects are 
  pointers to Sollya objects);
- the sa_sa (arguments and returned objects are all Sage/Python objects);
- a catch all flavor, without any suffix, (e. g. functions that have no 
  argument nor return value).

This classification is not always very strict. Conversion functions from Sollya
to Sage/Python are sometimes decorated with Sage/Python arguments to set
the precision. These functions remain in the so_sa category.

@note
Reported errors in Eclipse come from the calls to
the Sollya library

ToDo (among other things): 
 -memory management.
"""
from ctypes import *
import re
from sage.symbolic.expression_conversions import polynomial
from sage.symbolic.expression_conversions import PolynomialConverter
"""
Create the equivalent to an enum for the Sollya function types.
"""
(SOLLYA_BASE_FUNC_ABS,
SOLLYA_BASE_FUNC_ACOS,
    SOLLYA_BASE_FUNC_ACOSH,
    SOLLYA_BASE_FUNC_ADD,
    SOLLYA_BASE_FUNC_ASIN,
    SOLLYA_BASE_FUNC_ASINH,
    SOLLYA_BASE_FUNC_ATAN,
    SOLLYA_BASE_FUNC_ATANH,
    SOLLYA_BASE_FUNC_CEIL,
    SOLLYA_BASE_FUNC_CONSTANT,
    SOLLYA_BASE_FUNC_COS,
    SOLLYA_BASE_FUNC_COSH,
    SOLLYA_BASE_FUNC_DIV,
    SOLLYA_BASE_FUNC_DOUBLE,
    SOLLYA_BASE_FUNC_DOUBLEDOUBLE,
    SOLLYA_BASE_FUNC_DOUBLEEXTENDED,
    SOLLYA_BASE_FUNC_ERF,
    SOLLYA_BASE_FUNC_ERFC,
    SOLLYA_BASE_FUNC_EXP,
    SOLLYA_BASE_FUNC_EXP_M1,
    SOLLYA_BASE_FUNC_FLOOR,
    SOLLYA_BASE_FUNC_FREE_VARIABLE,
    SOLLYA_BASE_FUNC_HALFPRECISION,
    SOLLYA_BASE_FUNC_LIBRARYCONSTANT,
    SOLLYA_BASE_FUNC_LIBRARYFUNCTION,
    SOLLYA_BASE_FUNC_LOG,
    SOLLYA_BASE_FUNC_LOG_10,
    SOLLYA_BASE_FUNC_LOG_1P,
    SOLLYA_BASE_FUNC_LOG_2,
    SOLLYA_BASE_FUNC_MUL,
    SOLLYA_BASE_FUNC_NEARESTINT,
    SOLLYA_BASE_FUNC_NEG,
    SOLLYA_BASE_FUNC_PI,
    SOLLYA_BASE_FUNC_POW,
    SOLLYA_BASE_FUNC_PROCEDUREFUNCTION,
    SOLLYA_BASE_FUNC_QUAD,
    SOLLYA_BASE_FUNC_SIN,
    SOLLYA_BASE_FUNC_SINGLE,
    SOLLYA_BASE_FUNC_SINH,
    SOLLYA_BASE_FUNC_SQRT,
    SOLLYA_BASE_FUNC_SUB,
    SOLLYA_BASE_FUNC_TAN,
    SOLLYA_BASE_FUNC_TANH,
SOLLYA_BASE_FUNC_TRIPLEDOUBLE) = map(int,xrange(44))
sys.stderr.write("Superficial pobyso check...\n")
#print "First constant - SOLLYA_BASE_FUNC_ABS: ", SOLLYA_BASE_FUNC_ABS
#print "Last constant  - SOLLYA_BASE_FUNC_TRIPLEDOUBLE: ", SOLLYA_BASE_FUNC_TRIPLEDOUBLE

pobyso_max_arity = 9

def pobyso_absolute_so_so():
    return(sollya_lib_absolute(None))

def pobyso_autoprint(arg):
    sollya_lib_autoprint(arg, None)

def pobyso_autoprint_so_so(arg):
    sollya_lib_autoprint(arg,None)
    
def pobyso_bounds_to_range_sa_so(rnLowerBoundSa, rnUpperBoundSa, \
                                 precisionSa=None):
    """
    Return a Sollya range from to 2 RealField Sage elements.
    The Sollya range element has a sufficient precision to hold all
    the digits of the widest of the Sage bounds.
    """
    # Sanity check.
    if rnLowerBoundSa > rnUpperBoundSa:
        return None
    # Precision stuff.
    if precisionSa is None:
        # Check for the largest precision.
        lbPrecSa = rnLowerBoundSa.parent().precision()
        ubPrecSa = rnLowerBoundSa.parent().precision()
        maxPrecSa = max(lbPrecSa, ubPrecSa)
    else:
        maxPrecSa = precisionSa
    # From Sage to Sollya bounds.
#    lowerBoundSo = sollya_lib_constant(get_rn_value(rnLowerBoundSa),
#                                       maxPrecSa)
    lowerBoundSo = pobyso_constant_sa_so(rnLowerBoundSa,
                                         maxPrecSa)
    upperBoundSo = pobyso_constant_sa_so(rnUpperBoundSa,
                                       maxPrecSa)
    
    # From Sollya bounds to range.
    rangeSo = sollya_lib_range(lowerBoundSo, upperBoundSo)
    # Back to original precision.
    # Clean up
    sollya_lib_clear_obj(lowerBoundSo)
    sollya_lib_clear_obj(upperBoundSo)
    return rangeSo
# End pobyso_bounds_to_range_sa_so

def pobyso_build_end_elliptic_list_so_so(*args):
    """
    From Sollya argument objects, create a Sollya end elliptic list.
    Elements of the list are "eaten" (should not be cleared individually, 
    are cleared when the list is cleared).
    """
    if len(args) == 0:
        ## Called with an empty list produced "error".
        return sollya_lib_build_end_elliptic_list(None)
    index = 0
    ## One can not append elements to an elliptic list, prepend only is 
    #  permitted.
    for argument in reversed(args):
        if index == 0:
            listSo = sollya_lib_build_end_elliptic_list(argument, None)
        else:
            listSo = sollya_lib_prepend(argument, listSo)
        index += 1
    return listSo
        
# End pobyso_build_end_elliptic_list_so_so

def pobyso_build_function_sub_so_so(exp1So, exp2So):
    return sollya_lib_build_function_sub(exp1So, exp2So)

def pobyso_build_list_of_ints_sa_so(*args):
    """
    Build a Sollya list from Sage integral arguments. 
    """
    if len(args) == 0:
        return pobyso_build_list_so_so()
    ## Make a Sage list of integral constants.
    intsList = []
    for intElem in args:
        intsList.append(pobyso_constant_from_int_sa_so(intElem))
    return pobyso_build_list_so_so(*intsList) 
    
def pobyso_build_list_so_so(*args):
    """
    Make a Sollya list out of Sollya objects passed as arguments.
    If one wants to call it with a list argument, should prepend a "*"
    before the list variable name.
    Elements of the list are "eaten" (should not be cleared individually, 
    are cleared when the list is cleared).
    """
    if len(args) == 0:
        ## Called with an empty list produced "error".
        return sollya_lib_build_list(None)
    index = 0
    ## Append the Sollya elements one by one.
    for elementSo in args:
        if index == 0:
            listSo = sollya_lib_build_list(elementSo, None)
        else:
            listSo = sollya_lib_append(listSo, elementSo)
        index += 1
    return listSo
# End pobyso_build list_so_so
    

def pobyso_change_var_in_function_so_so(funcSo, chvarExpSo):
    """
    Variable change in a function.
    """
    return(sollya_lib_evaluate(funcSo,chvarExpSo))
# End pobyso_change_var_in_function_so_so     

def pobyso_chebyshevform_so_so(functionSo, degreeSo, intervalSo):
    resultSo = sollya_lib_chebyshevform(functionSo, degreeSo, intervalSo)
    return(resultSo)
# End pobyso_chebyshevform_so_so.

def pobyso_clear_obj(objSo):
    """
    Free a Sollya object's memory.
    Very thin wrapper around sollya_lib_clear_obj().
    """
    sollya_lib_clear_obj(objSo)
# End pobyso_clear_obj. 

def pobyso_clear_taylorform_sa_so(taylorFormSaSo):
    """
    This method is necessary to correctly clean up the memory from Taylor forms.
    These are made of a Sollya object, a Sollya object list, a Sollya object.
    For no clearly understood reason, sollya_lib_clear_object_list crashed 
    when applied to the object list.
    Here, we decompose it into Sage list of Sollya objects references and we
     clear them one at a time. 
    """
    sollya_lib_clear_obj(taylorFormSaSo[0])
    (coefficientsErrorsListSaSo, numElementsSa, isEndEllipticSa) = \
        pobyso_get_list_elements_so_so(taylorFormSaSo[1])
    for element in coefficientsErrorsListSaSo:
        sollya_lib_clear_obj(element)
    sollya_lib_clear_obj(taylorFormSaSo[1])
    sollya_lib_clear_obj(taylorFormSaSo[2])
# End pobyso_clear_taylorform_sa_so 

def pobyso_cmp(rnArgSa, cteSo):
    """
    Deprecated, use pobyso_cmp_sa_so_sa instead.
    """
    print "Deprecated, use pobyso_cmp_sa_so_sa instead."
    return pobyso_cmp_sa_so_sa(rnArgSa, cteSo)
# End  pobyso_cmp
    
def pobyso_cmp_sa_so_sa(rnArgSa, cteSo):
    """
    Compare the MPFR value a RealNumber with that of a Sollya constant.
    
    Get the value of the Sollya constant into a RealNumber and compare
    using MPFR. Could be optimized by working directly with a mpfr_t
    for the intermediate number. 
    """
    # Get the precision of the Sollya constant to build a Sage RealNumber 
    # with enough precision.to hold it.
    precisionOfCte = c_int(0)
    # From the Sollya constant, create a local Sage RealNumber.
    sollya_lib_get_prec_of_constant(precisionOfCte, cteSo) 
    #print "Precision of constant: ", precisionOfCte
    RRRR = RealField(precisionOfCte.value)
    rnLocalSa = RRRR(0)
    sollya_lib_get_constant(get_rn_value(rnLocalSa), cteSo)
    #
    ## Compare the Sage RealNumber version of the Sollya constant with rnArg
    #  through a direct comparison of underlying MPFR numbers.
    return cmp_rn_value(rnArgSa, rnLocal)
# End pobyso_smp_sa_so_sa

def pobyso_compute_pos_function_abs_val_bounds_sa_sa(funcSa, lowerBoundSa, \
                                                     upperBoundSa):
    """
    TODO: completely rework and test.
    """
    pobyso = pobyso_name_free_variable_sa_so(funcSa.variables()[0])
    funcSo = pobyso_parse_string(funcSa._assume_str().replace('_SAGE_VAR_', ''))
    rangeSo = pobyso_range_sa_so(lowerBoundSa, upperBoundSa)
    infnormSo = pobyso_infnorm_so_so(funcSo,rangeSo)
    # Sollya return the infnorm as an interval.
    fMaxSa = pobyso_get_interval_from_range_so_sa(infnormSo)
    # Get the top bound and compute the binade top limit.
    fMaxUpperBoundSa = fMaxSa.upper()
    binadeTopLimitSa = 2**ceil(fMaxUpperBoundSa.log2())
    # Put up together the function to use to compute the lower bound.
    funcAuxSo = pobyso_parse_string(str(binadeTopLimitSa) +  \
                                    '-(' + f._assume_str().replace('_SAGE_VAR_', '') + ')')
    pobyso_autoprint(funcAuxSo)
    # Clear the Sollya range before a new call to infnorm and issue the call.
    sollya_lib_clear_obj(infnormSo)
    infnormSo = pobyso_infnorm_so_so(funcAuxSo,rangeSo)
    fMinSa = pobyso_get_interval_from_range_so_sa(infnormSo)
    sollya_lib_clear_obj(infnormSo)
    fMinLowerBoundSa = binadeTopLimitSa - fMinSa.lower()
    # Compute the maximum of the precisions of the different bounds.
    maxPrecSa = max([fMinLowerBoundSa.parent().precision(), \
                     fMaxUpperBoundSa.parent().precision()])
    # Create a RealIntervalField and create an interval with the "good" bounds.
    RRRI = RealIntervalField(maxPrecSa)
    imageIntervalSa = RRRI(fMinLowerBoundSa, fMaxUpperBoundSa)
    # Free the unneeded Sollya objects
    sollya_lib_clear_obj(funcSo)
    sollya_lib_clear_obj(funcAuxSo)
    sollya_lib_clear_obj(rangeSo)
    return(imageIntervalSa)
# End pobyso_compute_pos_function_abs_val_bounds_sa_sa

def pobyso_compute_precision_decay_ratio_function_sa_so():
    """
    Compute the precision decay ratio function for polynomial 
    coefficient progressive trucation.
    """
    functionText = """
    proc(deg, a, b, we, wq) 
    {
      k = we * (exp(x/a)-1) + wq * (b*x)^2 + (1-we-wq) * x;
      return k/k(d);
    };
    """
    return pobyso_parse_string_sa_so(functionText)
# End  pobyso_compute_precision_decay_ratio_function.


def pobyso_constant(rnArg):
    """ Legacy function. See pobyso_constant_sa_so. """
    return(pobyso_constant_sa_so(rnArg))
    
def pobyso_constant_sa_so(rnArgSa, precisionSa=None):
    """
    Create a Sollya constant from a Sage RealNumber.
    The sollya_lib_constant() function creates a constant
    with the same precision as the source.
    """
    ## Precision stuff. If one wants to change precisions,
    #  everything takes place in Sage. That only makes
    #  sense if one wants to reduce the precision.
    # TODO: revisit precision stuff with new technique.
    if not precisionSa is None:
        RRR = RealField(precisionSa)
        rnArgSa = RRR(rnArgSa)
    #print rnArgSa, rnArgSa.precision()
    # Sollya constant creation takes place here.
    return sollya_lib_constant(get_rn_value(rnArgSa))
# End pobyso_constant_sa_so
     
def pobyso_constant_0_sa_so():
    """
    Obvious.
    """
    return pobyso_constant_from_int_sa_so(0)

def pobyso_constant_1():
    """
    Obvious.
    Legacy function. See pobyso_constant_so_so. 
    """
    return pobyso_constant_1_sa_so()

def pobyso_constant_1_sa_so():
    """
    Obvious.
    """
    return(pobyso_constant_from_int_sa_so(1))

def pobyso_constant_from_int(anInt):
    """ Legacy function. See pobyso_constant_from_int_sa_so. """
    return pobyso_constant_from_int_sa_so(anInt)

def pobyso_constant_from_int_sa_so(anInt):
    """
    Get a Sollya constant from a Sage int.
    """
    return sollya_lib_constant_from_int64(long(anInt))

def pobyso_constant_from_int_so_sa(constSo):
    """
    Get a Sage int from a Sollya int constant.
    Usefull for precision or powers in polynomials.
    """
    constSa = c_long(0)
    sollya_lib_get_constant_as_int64(byref(constSa), constSo)
    return constSa.value
# End pobyso_constant_from_int_so_sa

def pobyso_constant_from_mpq_sa_so(rationalSa):
    """
    Make a Sollya constant from Sage rational.
    The Sollya constant is an unevaluated expression.
    Hence no precision argument is needed.
    It is better to leave this way since Sollya has its own
    optimized evaluation mecanism that tries very hard to
    return exact values or at least faithful ones.
    """
    ratExprSo = \
        sollya_lib_constant_from_mpq(sgmp_get_rational_value(rationalSa))
    return ratExprSo
# End pobyso_constant_from_mpq_sa_so.

def pobyso_constant_sollya_prec_sa_so(rnArgSa):
    """
    Create a Sollya constant from a Sage RealNumber at the
    current precision in Sollya.
    """
    currentSollyaPrecSa = pobyso_get_prec_so_sa()
    return pobyso_constant_sa_so(rnArgSa, currentSollyaPrecSa)
# End pobyso_constant_sollya_prec_sa_so

def pobyso_end_elliptic_list_so_sa_so(objectsListSo, intCountSa):
    """
    Create a Sollya end elliptic list made of the objectListSo[0] to
     objectsListSo[intCountSa-1] objects.
    """     
    return sollya_lib_end_elliptic_list(objectSo, int(intCountSa))

def pobyso_error_so():
    return sollya_lib_error(None)
# End pobyso_error().

def pobyso_evaluate_so_so(funcSo, argumentSo):
    """
    Evaluates funcSo for arguemntSo through sollya_lib_evaluate().
    """
    return sollya_lib_evaluate(funcSo, argumentSo)
# End pobyso_evaluate_so_so.
#
def pobyso_diff_so_so(funcSo):
    """
    Very thin wrapper around sollya_lib_diff.
    """
    ## TODO: add a check to make sure funcSo is a functional expression.
    return sollya_lib_diff(funcSo)

def pobyso_dirty_find_zeros_so_so(funcSo, rangeSo):
    """
    Thin wrapper over sollya_lib_dirtyfindzeros()
    """
    return sollya_lib_dirtyfindzeros(funcSo, rangeSo)
# End pobys_dirty_find_zeros

def pobyso_dirty_inf_norm_so_so(funcSo, rangeSo, preSo=None):
    """
    Thin wrapper around sollya_dirtyinfnorm().
    """
    # TODO: manage the precision change.
    
    return sollya_lib_dirtyinfnorm(funcSo, rangeSo)
# End pobyso_dirty_inf_norm_so_so

def pobyso_float_list_so_sa(listSo):
    """
    Return a Sollya list of floating-point numbers as a Sage list of 
    floating-point numbers.
    TODO: add a test to make sure that each element of the list is a constant.
    """
    listSa   = []
    ## The function returns none if the list is empty or an error has happened.
    retVal = pobyso_get_list_elements_so_so(listSo)
    if retVal is None:
        return listSa
    ## Just in case the interface is changed and an empty list is returned 
    #  instead of None.
    elif len(retVal) == 0:
        return listSa
    else:
        ## Remember pobyso_get_list_elements_so_so returns more information
        #  than just the elements of the list (# elements, is_elliptic)
        listSaSo, numElements, isEndElliptic = retVal
    ## Return an empty list.
    if numElements == 0:
        return listSa
    ## Search first for the maximum precision of the elements
    maxPrecSa = 0
    for floatSo in listSaSo:
        #pobyso_autoprint(floatSo)
        curPrecSa =  pobyso_get_prec_of_constant_so_sa(floatSo)
        if curPrecSa > maxPrecSa:
            maxPrecSa = curPrecSa
    ##
    RF = RealField(maxPrecSa)
    ##
    for floatSo in listSaSo:
        listSa.append(pobyso_get_constant_as_rn_with_rf_so_sa(floatSo))
    return listSa
# End pobyso_float_list_so_sa

def pobyso_float_poly_sa_so(polySa, precSa = None):
    """
    Create a Sollya polynomial from a Sage RealField polynomial.
    """
    ## TODO: filter arguments.
    ## Precision. If a precision is given, convert the polynomial
    #  into the right polynomial field. If not convert it straight
    #  to Sollya.
    sollyaPrecChanged = False 
    (initialSollyaPrecSo, initialSollyaPrecSa) = pobyso_get_prec_so_so_sa()
    if precSa is None:
        precSa = polySa.parent().base_ring().precision()
    if (precSa > initialSollyaPrecSa):
        assert precSa >= 2, "Precision change <2 requested"
        if precSa <= 2:
            print inspect.stack()[0][3], ": precision change <= 2 requested"
        precSo = pobyso_constant_from_int(precSa)
        pobyso_set_prec_so_so(precSo)
        sollya_lib_clear_obj(precSo)
        sollyaPrecChanged = True    
    ## Get exponents and coefficients.
    exponentsSa     = polySa.exponents()
    coefficientsSa  = polySa.coefficients()
    ## Build the polynomial.
    polySo = None
    for coefficientSa, exponentSa in zip(coefficientsSa, exponentsSa):
        #print coefficientSa.n(prec=precSa), exponentSa
        coefficientSo = \
            pobyso_constant_sa_so(coefficientSa)
        #pobyso_autoprint(coefficientSo)
        exponentSo = \
            pobyso_constant_from_int_sa_so(exponentSa)
        #pobyso_autoprint(exponentSo)
        monomialSo = sollya_lib_build_function_pow(
                       sollya_lib_build_function_free_variable(),
                       exponentSo)
        polyTermSo = sollya_lib_build_function_mul(coefficientSo,
                                                       monomialSo)
        if polySo is None:
            polySo = polyTermSo
        else:
            polySo = sollya_lib_build_function_add(polySo, polyTermSo)
    if sollyaPrecChanged:
        pobyso_set_prec_so_so(initialSollyaPrecSo)
        sollya_lib_clear_obj(initialSollyaPrecSo)
    return polySo
# End pobyso_float_poly_sa_so    

def pobyso_float_poly_so_sa(polySo, realFieldSa=None):
    """
    Convert a Sollya polynomial into a Sage floating-point polynomial.
    If no realField is given, a RealField corresponding to the maximum 
    precision of the coefficients is internally computed.
    The real field is not returned but can be easily retrieved from 
    the polynomial itself.
    ALGORITHM:
    - (optional) compute the RealField of the coefficients;
    - convert the Sollya expression into a Sage expression;
    - convert the Sage expression into a Sage polynomial
    """    
    if realFieldSa is None:
        expressionPrecSa = pobyso_get_max_prec_of_exp_so_sa(polySo)
        #print "Maximum precision of Sollya polynomial coefficients:", expressionPrecSa
        if expressionPrecSa < 2 or expressionPrecSa > 2147483391:
            print "Maximum degree of expression:", expressionPrecSa
        realFieldSa      = RealField(expressionPrecSa)
    #print "Sollya expression before...",
    #pobyso_autoprint(polySo)

    expressionSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo,
                                                             realFieldSa)
    #print "...Sollya expression after."
    #pobyso_autoprint(polySo)
    polyVariableSa = expressionSa.variables()[0]
    polyRingSa     = realFieldSa[str(polyVariableSa)]
    #print polyRingSa
    # Do not use the polynomial(expressionSa, ring=polyRingSa) form!
    polynomialSa = polyRingSa(expressionSa)
    polyCoeffsListSa = polynomialSa.coefficients()
    #for coeff in polyCoeffsListSa:
    #    print coeff.abs().n()
    return polynomialSa
# End pobyso_float_poly_so_sa

def pobyso_free_variable():
    """
    Ultra thin wrapper around the sollya_lib_function_build_free_variable function.
    """
    return sollya_lib_build_function_free_variable()
    
def pobyso_function_type_as_string(funcType):
    """ Legacy function. See pobyso_function_type_as_string_so_sa. """
    return(pobyso_function_type_as_string_so_sa(funcType))

def pobyso_function_type_as_string_so_sa(funcType):
    """
    Numeric Sollya function codes -> Sage mathematical function names.
    Notice that pow -> ^ (a la Sage, not a la Python).
    """
    if funcType == SOLLYA_BASE_FUNC_ABS:
        return "abs"
    elif funcType == SOLLYA_BASE_FUNC_ACOS:
        return "arccos"
    elif funcType == SOLLYA_BASE_FUNC_ACOSH:
        return "arccosh"
    elif funcType == SOLLYA_BASE_FUNC_ADD:
        return "+"
    elif funcType == SOLLYA_BASE_FUNC_ASIN:
        return "arcsin"
    elif funcType == SOLLYA_BASE_FUNC_ASINH:
        return "arcsinh"
    elif funcType == SOLLYA_BASE_FUNC_ATAN:
        return "arctan"
    elif funcType == SOLLYA_BASE_FUNC_ATANH:
        return "arctanh"
    elif funcType == SOLLYA_BASE_FUNC_CEIL:
        return "ceil"
    elif funcType == SOLLYA_BASE_FUNC_CONSTANT:
        return "cte"
    elif funcType == SOLLYA_BASE_FUNC_COS:
        return "cos"
    elif funcType == SOLLYA_BASE_FUNC_COSH:
        return "cosh"
    elif funcType == SOLLYA_BASE_FUNC_DIV:
        return "/"
    elif funcType == SOLLYA_BASE_FUNC_DOUBLE:
        return "double"
    elif funcType == SOLLYA_BASE_FUNC_DOUBLEDOUBLE:
        return "doubleDouble"
    elif funcType == SOLLYA_BASE_FUNC_DOUBLEEXTENDED:
        return "doubleDxtended"
    elif funcType == SOLLYA_BASE_FUNC_ERF:
        return "erf"
    elif funcType == SOLLYA_BASE_FUNC_ERFC:
        return "erfc"
    elif funcType == SOLLYA_BASE_FUNC_EXP:
        return "exp"
    elif funcType == SOLLYA_BASE_FUNC_EXP_M1:
        return "expm1"
    elif funcType == SOLLYA_BASE_FUNC_FLOOR:
        return "floor"
    elif funcType == SOLLYA_BASE_FUNC_FREE_VARIABLE:
        return "freeVariable"
    elif funcType == SOLLYA_BASE_FUNC_HALFPRECISION:
        return "halfPrecision"
    elif funcType == SOLLYA_BASE_FUNC_LIBRARYCONSTANT:
        return "libraryConstant"
    elif funcType == SOLLYA_BASE_FUNC_LIBRARYFUNCTION:
        return "libraryFunction"
    elif funcType == SOLLYA_BASE_FUNC_LOG:
        return "log"
    elif funcType == SOLLYA_BASE_FUNC_LOG_10:
        return "log10"
    elif funcType == SOLLYA_BASE_FUNC_LOG_1P:
        return "log1p"
    elif funcType == SOLLYA_BASE_FUNC_LOG_2:
        return "log2"
    elif funcType == SOLLYA_BASE_FUNC_MUL:
        return "*"
    elif funcType == SOLLYA_BASE_FUNC_NEARESTINT:
        return "round"
    elif funcType == SOLLYA_BASE_FUNC_NEG:
        return "__neg__"
    elif funcType == SOLLYA_BASE_FUNC_PI:
        return "pi"
    elif funcType == SOLLYA_BASE_FUNC_POW:
        return "^"
    elif funcType == SOLLYA_BASE_FUNC_PROCEDUREFUNCTION:
        return "procedureFunction"
    elif funcType == SOLLYA_BASE_FUNC_QUAD:
        return "quad"
    elif funcType == SOLLYA_BASE_FUNC_SIN:
        return "sin"
    elif funcType == SOLLYA_BASE_FUNC_SINGLE:
        return "single"
    elif funcType == SOLLYA_BASE_FUNC_SINH:
        return "sinh"
    elif funcType == SOLLYA_BASE_FUNC_SQRT:
        return "sqrt"
    elif funcType == SOLLYA_BASE_FUNC_SUB:
        return "-"
    elif funcType == SOLLYA_BASE_FUNC_TAN:
        return "tan"
    elif funcType == SOLLYA_BASE_FUNC_TANH:
        return "tanh"
    elif funcType == SOLLYA_BASE_FUNC_TRIPLEDOUBLE:
        return "tripleDouble"
    else:
        return None

def pobyso_get_constant(rnArgSa, constSo):
    """ Legacy function. See pobyso_get_constant_so_sa. """
    return pobyso_get_constant_so_sa(rnArgSa, constSo)
# End pobyso_get_constant

def pobyso_get_constant_so_sa(rnArgSa, constSo):
    """
    Set the value of rnArgSo to the value of constSo in MPFR_RNDN mode.
    rnArg must already exist and belong to some RealField.
    We assume that constSo points to a Sollya constant.
    """
    outcome = sollya_lib_get_constant(get_rn_value(rnArgSa), constSo)
    if outcome == 0: # Failure because constSo is not a constant expression.
        return None
    else:
        return outcome
# End  pobyso_get_constant_so_sa
   
def pobyso_get_constant_as_rn(ctExpSo):
    """ 
    Legacy function. See pobyso_get_constant_as_rn_so_sa. 
    """ 
    return(pobyso_get_constant_as_rn_so_sa(ctExpSo))
    
def pobyso_get_constant_as_rn_so_sa(constExpSo):
    """
    Get a Sollya constant as a Sage "real number".
    The precision of the floating-point number returned is that of the Sollya
    constant.
    """
    #print "Before computing precision of variable..."
    #pobyso_autoprint(constExpSo)
    precisionSa  = pobyso_get_prec_of_constant_so_sa(constExpSo)
    #print "precisionSa:", precisionSa
    ## If the expression can not be exactly converted, None is returned.
    #  In this case opt for the Sollya current expression.
    if precisionSa is None:
        precisionSa = pobyso_get_prec_so_sa()
    RRRR = RealField(precisionSa)
    rnSa = RRRR(0)
    outcome = sollya_lib_get_constant(get_rn_value(rnSa), constExpSo)
    if outcome == 0:
        return None
    else:
        return rnSa
# End pobyso_get_constant_as_rn_so_sa

def pobyso_get_constant_as_rn_with_rf(ctExp, realField):
    """ 
    Legacy function. See pobyso_get_constant_as_rn_with_rf_so_sa.
    """
    return pobyso_get_constant_as_rn_with_rf_so_sa(ctExp, realField)
# End pobyso_get_constant_as_rn_with_rf
    
def pobyso_get_constant_as_rn_with_rf_so_sa(ctExpSo, realFieldSa = None):
    """
    Get a Sollya constant as a Sage "real number".
    If no real field is specified, the precision of the floating-point number 
    returned is that of the Sollya constant.
    Otherwise is is that of the real field. Hence rounding may happen.
    """
    if realFieldSa is None:
        return pobyso_get_constant_as_rn_so_sa(ctExpSo)
    rnSa = realFieldSa(0)
    outcome = sollya_lib_get_constant(get_rn_value(rnSa), ctExpSo)
    if outcome == 0:
        return None
    else:
        return rnSa
# End pobyso_get_constant_as_rn_with_rf_so_sa

def pobyso_get_free_variable_name():
    """ 
    Legacy function. See pobyso_get_free_variable_name_so_sa.
    """
    return(pobyso_get_free_variable_name_so_sa())

def pobyso_get_free_variable_name_so_sa():
    return sollya_lib_get_free_variable_name()
    
def pobyso_get_function_arity(expressionSo):
    """ 
    Legacy function. See pobyso_get_function_arity_so_sa.
    """
    return(pobyso_get_function_arity_so_sa(expressionSo))

def pobyso_get_function_arity_so_sa(expressionSo):
    arity = c_int(0)
    sollya_lib_get_function_arity(byref(arity),expressionSo)
    return int(arity.value)

def pobyso_get_head_function(expressionSo):
    """ 
    Legacy function. See pobyso_get_head_function_so_sa. 
    """
    return(pobyso_get_head_function_so_sa(expressionSo)) 

def pobyso_get_head_function_so_sa(expressionSo):
    functionType = c_int(0)
    sollya_lib_get_head_function(byref(functionType), expressionSo)
    return int(functionType.value)

def pobyso_get_interval_from_range_so_sa(soRange, realIntervalFieldSa = None ):
    """
    Return the Sage interval corresponding to the Sollya range argument.
    If no reaIntervalField is passed as an argument, the interval bounds are not
    rounded: they are elements of RealIntervalField of the "right" precision
    to hold all the digits.
    """
    prec = c_int(0)
    if realIntervalFieldSa is None:
        retval = sollya_lib_get_prec_of_range(byref(prec), soRange, None)
        if retval == 0:
            return None
        realIntervalFieldSa = RealIntervalField(prec.value)
    intervalSa = realIntervalFieldSa(0,0)
    retval = \
        sollya_lib_get_interval_from_range(get_interval_value(intervalSa),\
                                           soRange)
    if retval == 0:
        return None
    return intervalSa
# End pobyso_get_interval_from_range_so_sa

def pobyso_get_list_elements(soObj):
    """ Legacy function. See pobyso_get_list_elements_so_so. """
    return pobyso_get_list_elements_so_so(soObj)
 
def pobyso_get_list_elements_so_so(objectListSo):
    """
    Get the Sollya list elements as a Sage/Python array of Sollya objects.
    
    INPUT:
    - objectListSo: a Sollya list of Sollya objects.
    
    OUTPUT:
    - a Sage/Python tuple made of:
      - a Sage/Python list of Sollya objects,
      - a Sage/Python int holding the number of elements,
      - a Sage/Python int stating (!= 0) that the list is end-elliptic.
    NOTE::
        We recover the addresses of the Sollya object from the list of pointers
        returned by sollya_lib_get_list_elements. The list itself is freed.
    TODO::
        Figure out what to do with numElements since the number of elements
        can easily be recovered from the list itself. 
        Ditto for isEndElliptic.
    """
    listAddress = POINTER(c_longlong)()
    numElements = c_int(0)
    isEndElliptic = c_int(0)
    listAsSageList = []
    result = sollya_lib_get_list_elements(byref(listAddress),\
                                          byref(numElements),\
                                          byref(isEndElliptic),\
                                          objectListSo)
    if result == 0 :
        return None
    for i in xrange(0, numElements.value, 1):
       #listAsSageList.append(sollya_lib_copy_obj(listAddress[i]))
       listAsSageList.append(listAddress[i])
       # Clear each of the elements returned by Sollya.
       #sollya_lib_clear_obj(listAddress[i])
    # Free the list itself.   
    sollya_lib_free(listAddress)
    return (listAsSageList, numElements.value, isEndElliptic.value)

def pobyso_get_max_prec_of_exp(soExp):
    """ Legacy function. See pobyso_get_max_prec_of_exp_so_sa. """
    return pobyso_get_max_prec_of_exp_so_sa(soExp)

def pobyso_get_max_prec_of_exp_so_sa(expSo):
    """
    Get the maximum precision used for the numbers in a Sollya expression.
    
    Arguments:
    soExp -- a Sollya expression pointer
    Return value:
    A Python integer
    TODO: 
    - error management;
    - correctly deal with numerical type such as DOUBLEEXTENDED.
    """
    if expSo is None:
        print inspect.stack()[0][3], ": expSo is None."
        return 0
    maxPrecision = 0
    minConstPrec = 0
    currentConstPrec = 0
    #pobyso_autoprint(expSo)
    operator = pobyso_get_head_function_so_sa(expSo)
    if (operator != SOLLYA_BASE_FUNC_CONSTANT) and \
    (operator != SOLLYA_BASE_FUNC_FREE_VARIABLE):
        (arity, subexpressions) = pobyso_get_subfunctions_so_sa(expSo)
        for i in xrange(arity):
            maxPrecisionCandidate = \
                pobyso_get_max_prec_of_exp_so_sa(subexpressions[i])
            if maxPrecisionCandidate > maxPrecision:
                maxPrecision = maxPrecisionCandidate
        return maxPrecision
    elif operator == SOLLYA_BASE_FUNC_CONSTANT:
        #minConstPrec = pobyso_get_min_prec_of_constant_so_sa(expSo)
        #currentConstPrec = pobyso_get_min_prec_of_constant_so_sa(soExp)
        #print minConstPrec, " - ", currentConstPrec 
        return pobyso_get_min_prec_of_constant_so_sa(expSo)
    
    elif operator == SOLLYA_BASE_FUNC_FREE_VARIABLE:
        return 0
    else:
        print "pobyso_get_max_prec_of_exp_so_sa: unexepected operator."
        return 0

def pobyso_get_min_prec_of_constant_so_sa(constExpSo):
    """
    Get the minimum precision necessary to represent the value of a Sollya
    constant.
    MPFR_MIN_PREC and powers of 2 are taken into account.
    We assume that constExpSo is a pointer to a Sollay constant expression.
    """
    constExpAsRnSa = pobyso_get_constant_as_rn_so_sa(constExpSo)
    return(min_mpfr_size(get_rn_value(constExpAsRnSa)))

def pobyso_get_poly_so_sa(polySo, realFieldSa=None):
    """
    Convert a Sollya polynomial into a Sage polynomial.
    Legacy function. Use pobyso_float_poly_so_sa() instead.
    """    
    return pobyso_float_poly_so_sa(polySo,realFieldSa)
# End pobyso_get_poly_so_sa

def pobyso_get_prec():
    """ Legacy function. See pobyso_get_prec_so_sa(). """
    return pobyso_get_prec_so_sa()

def pobyso_get_prec_so():
    """
    Get the current default precision in Sollya.
    The return value is a Sollya object.
    Usefull when modifying the precision back and forth by avoiding
    extra conversions.
    """
    return sollya_lib_get_prec(None)
    
def pobyso_get_prec_so_sa():
    """
    Get the current default precision in Sollya.
    The return value is Sage/Python int.
    """
    precSo = sollya_lib_get_prec()
    precSa = pobyso_constant_from_int_so_sa(precSo)
    sollya_lib_clear_obj(precSo)
    return precSa
# End pobyso_get_prec_so_sa.

def pobyso_get_prec_so_so_sa():
    """
    Return the current precision both as a Sollya object and a
    Sage integer as hybrid tuple.
    To avoid multiple calls for precision manipulations. 
    """
    precSo = sollya_lib_get_prec()
    precSa = pobyso_constant_from_int_so_sa(precSo)
    return (precSo, int(precSa))
    
def pobyso_get_prec_of_constant(ctExpSo):
    """ Legacy function. See pobyso_get_prec_of_constant_so_sa. """
    return pobyso_get_prec_of_constant_so_sa(ctExpSo)

def pobyso_get_prec_of_constant_so_sa(ctExpSo):
    """
    Tries to find a precision to represent ctExpSo without rounding.
    If not possible, returns None.
    """
    #print "Entering pobyso_get_prec_of_constant_so_sa..."
    prec = c_int(0)
    retc = sollya_lib_get_prec_of_constant(byref(prec), ctExpSo, None)
    if retc == 0:
        #print "pobyso_get_prec_of_constant_so_sa failed."
        return None
    #print "...exiting pobyso_get_prec_of_constant_so_sa."
    return int(prec.value)

def pobyso_get_prec_of_range_so_sa(rangeSo):
    """
    Returns the number of bits elements of a range are coded with.
    """
    prec = c_int(0)
    retc = sollya_lib_get_prec_of_range(byref(prec), rangeSo, None)
    if retc == 0:
        return(None)
    return int(prec.value)
# End pobyso_get_prec_of_range_so_sa()

def pobyso_get_sage_exp_from_sollya_exp(sollyaExpSo, realField = RR):
    """ Legacy function. See pobyso_get_sage_exp_from_sollya_exp_so_sa. """
    return pobyso_get_sage_exp_from_sollya_exp_so_sa(sollyaExpSo, 
                                                     realField = RR)

def pobyso_get_sage_exp_from_sollya_exp_so_sa(sollyaExpSo, realFieldSa = RR):
    """
    Get a Sage expression from a Sollya expression. 
    Currently only tested with polynomials with floating-point coefficients.
    Notice that, in the returned polynomial, the exponents are RealNumbers.
    """
    #pobyso_autoprint(sollyaExp)
    operatorSa = pobyso_get_head_function_so_sa(sollyaExpSo)
    sollyaLibFreeVariableName = sollya_lib_get_free_variable_name()
    ## Get rid of the "_"'s in "_x_", if any.
    sollyaLibFreeVariableName = re.sub('_', '', sollyaLibFreeVariableName)
    # Constants and the free variable are special cases.
    # All other operator are dealt with in the same way.
    if (operatorSa != SOLLYA_BASE_FUNC_CONSTANT) and \
       (operatorSa != SOLLYA_BASE_FUNC_FREE_VARIABLE):
        (aritySa, subexpressionsSa) = pobyso_get_subfunctions_so_sa(sollyaExpSo)
        if aritySa == 1:
            sageExpSa = eval(pobyso_function_type_as_string_so_sa(operatorSa) + \
            "(" + pobyso_get_sage_exp_from_sollya_exp_so_sa(subexpressionsSa[0], \
            realFieldSa) + ")")
        elif aritySa == 2:
            # We do not get through the preprocessor.
            # The "^" operator is then a special case.
            if operatorSa == SOLLYA_BASE_FUNC_POW:
                operatorAsStringSa = "**"
            else:
                operatorAsStringSa = \
                    pobyso_function_type_as_string_so_sa(operatorSa)
            sageExpSa = \
              eval("pobyso_get_sage_exp_from_sollya_exp_so_sa(subexpressionsSa[0], realFieldSa)"\
              + " " + operatorAsStringSa + " " + \
                   "pobyso_get_sage_exp_from_sollya_exp_so_sa(subexpressionsSa[1], realFieldSa)")
        # We do not know yet how to deal with arity >= 3 
        # (is there any in Sollya anyway?).
        else:
            sageExpSa = eval('None')
        return sageExpSa
    elif operatorSa == SOLLYA_BASE_FUNC_CONSTANT:
        #print "This is a constant"
        return pobyso_get_constant_as_rn_with_rf_so_sa(sollyaExpSo, realFieldSa)
    elif operatorSa == SOLLYA_BASE_FUNC_FREE_VARIABLE:
        #print "This is the free variable"
        return eval(sollyaLibFreeVariableName)
    else:
        print "Unexpected"
        return eval('None')
# End pobyso_get_sage_exp_from_sollya_exp_so_sa


def pobyso_get_subfunctions(expressionSo):
    """ Legacy function. See pobyso_get_subfunctions_so_sa. """
    return pobyso_get_subfunctions_so_sa(expressionSo) 
# End pobyso_get_subfunctions.
 
def pobyso_get_subfunctions_so_sa(expressionSo):
    """
    Get the subfunctions of an expression.
    Return the number of subfunctions and the list of subfunctions addresses.
    S.T.: Could not figure out another way than that ugly list of declarations
    to recover the addresses of the subfunctions.
    We limit ourselves to arity 8 functions. 
    """
    subf0 = c_int(0)
    subf1 = c_int(0)
    subf2 = c_int(0)
    subf3 = c_int(0)
    subf4 = c_int(0)
    subf5 = c_int(0)
    subf6 = c_int(0)
    subf7 = c_int(0)
    subf8 = c_int(0)
    arity = c_int(0)
    nullPtr = POINTER(c_int)()
    sollya_lib_get_subfunctions(expressionSo, byref(arity), \
      byref(subf0), byref(subf1), byref(subf2), byref(subf3), \
      byref(subf4), byref(subf5),\
      byref(subf6), byref(subf7), byref(subf8), nullPtr, None) 
#    byref(cast(subfunctions[0], POINTER(c_int))), \
#    byref(cast(subfunctions[0], POINTER(c_int))), \
#    byref(cast(subfunctions[2], POINTER(c_int))), \
#    byref(cast(subfunctions[3], POINTER(c_int))), \
#    byref(cast(subfunctions[4], POINTER(c_int))), \
#    byref(cast(subfunctions[5], POINTER(c_int))), \
#    byref(cast(subfunctions[6], POINTER(c_int))), \
#    byref(cast(subfunctions[7], POINTER(c_int))), \
#    byref(cast(subfunctions[8], POINTER(c_int))), nullPtr)
    subfunctions = [subf0, subf1, subf2, subf3, subf4, subf5, subf6, subf7, \
                    subf8]
    subs = []
    if arity.value > pobyso_max_arity:
        return(0,[])
    for i in xrange(arity.value):
        subs.append(int(subfunctions[i].value))
        #print subs[i]
    return (int(arity.value), subs)
# End pobyso_get_subfunctions_so_sa
    
def pobyso_guess_degree_sa_sa(functionSa, intervalSa, approxErrorSa, 
                              weightSa=None, degreeBoundSa=None):
    """
    Sa_sa variant of the solly_guessdegree function.
    Return 0 if something goes wrong.
    """
    functionAsStringSa = functionSa._assume_str().replace('_SAGE_VAR_', '')
    functionSo = pobyso_parse_string_sa_so(functionAsStringSa)
    if pobyso_is_error_so_sa(functionSo):
        sollya_lib_clear_obj(functionSo)
        return 0
    rangeSo = pobyso_interval_to_range_sa_so(intervalSa)
    # The approximation error is expected to be a floating point number.
    if pobyso_is_floating_point_number_sa_sa(approxErrorSa):
        approxErrorSo = pobyso_constant_sa_so(approxErrorSa)
    else:
        approxErrorSo = pobyso_constant_sa_so(RR(approxErrorSa))
    if not weightSa is None:
        weightAsStringSa = weightSa._assume_str().replace('_SAGE_VAR_', '')
        weightSo = pobyso_parse_string_sa_so(weightAsStringSa)
        if pobyso_is_error_so_sa(weightSo):
            sollya_lib_clear_obj(functionSo)
            sollya_lib_clear_obj(rangeSo)
            sollya_lib_clear_obj(approxErrorSo)   
            sollya_lib_clear_obj(weightSo)
            return 0   
    else:
        weightSo = None
    if not degreeBoundSa is None:
        degreeBoundSo = pobyso_constant_from_int_sa_so(degreeBoundSa)
    else:
        degreeBoundSo = None
    guessedDegreeSa = pobyso_guess_degree_so_sa(functionSo,
                                                rangeSo,
                                                approxErrorSo,
                                                weightSo,
                                                degreeBoundSo)
    sollya_lib_clear_obj(functionSo)
    sollya_lib_clear_obj(rangeSo)
    sollya_lib_clear_obj(approxErrorSo)
    if not weightSo is None:
        sollya_lib_clear_obj(weightSo)
    if not degreeBoundSo is None:
        sollya_lib_clear_obj(degreeBoundSo)
    return guessedDegreeSa
# End poyso_guess_degree_sa_sa

def pobyso_guess_degree_so_sa(functionSo, rangeSo, errorSo, weightSo=None, \
                              degreeBoundSo=None):
    """
    Thin wrapper around the guessdegree function.
    Nevertheless, some precision control stuff has been appended.
    """
    # Deal with Sollya internal precision issues: if it is too small, 
    # compared with the error, increases it to about twice -log2(error).
    errorSa = pobyso_get_constant_as_rn_with_rf_so_sa(errorSo)
    log2ErrorSa = errorSa.log2()
    if log2ErrorSa < 0:
        neededPrecisionSa = int(2 * int(-log2ErrorSa) / 64) * 64
    else:
        neededPrecisionSa = int(2 * int(log2ErrorSa) / 64) * 64
    #print "Needed precision:", neededPrecisionSa
    sollyaPrecisionHasChanged = False
    (initialPrecSo, initialPrecSa) = pobyso_get_prec_so_so_sa()
    if neededPrecisionSa > initialPrecSa:
        if neededPrecisionSa <= 2:
            print inspect.stack()[0][3], ": precision change <= 2 requested."
        pobyso_set_prec_sa_so(neededPrecisionSa)
        sollyaPrecisionHasChanged = True
    #print "Guessing degree..."
    # weightSo and degreeBoundsSo are optional arguments.
    # As declared, sollya_lib_guessdegree must take 5 arguments.
    if weightSo is None:
        degreeRangeSo = sollya_lib_guessdegree(functionSo, rangeSo, errorSo, 
                                               0, 0, None)
    elif degreeBoundSo is None:
        degreeRangeSo =  sollya_lib_guessdegree(functionSo, rangeSo, \
                                                errorSo, weightSo, 0, None)
    else:
        degreeRangeSo =  sollya_lib_guessdegree(functionSo, rangeSo, errorSo, \
                                                weightSo, degreeBoundSo, None)
    #print "...degree guess done."
    # Restore internal precision, if applicable.
    if sollyaPrecisionHasChanged:
        pobyso_set_prec_so_so(initialPrecSo)
        sollya_lib_clear_obj(initialPrecSo)
    degreeIntervalSa = pobyso_range_to_interval_so_sa(degreeRangeSo)
    sollya_lib_clear_obj(degreeRangeSo)
    # When ok, both bounds match.
    # When the degree bound is too low, the upper bound is the degree
    # for which the error can be honored.
    # When it really goes wrong, the upper bound is infinity. 
    if degreeIntervalSa.lower() == degreeIntervalSa.upper():
        return int(degreeIntervalSa.lower())
    else:
        if degreeIntervalSa.upper().is_infinity():
            return None
        else:
            return int(degreeIntervalSa.upper())
    # End pobyso_guess_degree_so_sa

def pobyso_inf_so_so(intervalSo):
    """
    Very thin wrapper around sollya_lib_inf().
    """
    return sollya_lib_inf(intervalSo)
# End pobyso_inf_so_so.
    
def pobyso_infnorm_so_so(func, interval, file = None, intervalList = None):
    print "Do not use this function. User pobyso_supnorm_so_so instead."
    return None

def pobyso_interval_to_range_sa_so(intervalSa, precisionSa=None):
    if precisionSa is None:
        precisionSa = intervalSa.parent().precision()
    intervalSo = pobyso_bounds_to_range_sa_so(intervalSa.lower(),\
                                              intervalSa.upper(),\
                                              precisionSa)
    return intervalSo
# End pobyso_interval_to_range_sa_so

def pobyso_is_error_so_sa(objSo):
    """
    Thin wrapper around the sollya_lib_obj_is_error() function.
    """
    if sollya_lib_obj_is_error(objSo) != 0:
        return True
    else:
        return False
# End pobyso_is_error-so_sa

def pobyso_is_floating_point_number_sa_sa(numberSa):
    """
    Check whether a Sage number is floating point.
    Exception stuff added because numbers other than
    floating-point ones do not have the is_real() attribute.
    """
    try:
        return numberSa.is_real()
    except AttributeError:
        return False
# End pobyso_is_floating_piont_number_sa_sa

def pobyso_is_range_so_sa(rangeCandidateSo):
    """
    Thin wrapper over sollya_lib_is_range.
    """
    return sollya_lib_obj_is_range(rangeCandidateSo) != 0
    
# End pobyso_is_range_so_sa


def pobyso_lib_init():
    sollya_lib_init(None)

def pobyso_lib_close():
    sollya_lib_close(None)

def pobyso_name_free_variable(freeVariableNameSa):
    """ Legacy function. See pobyso_name_free_variable_sa_so. """
    pobyso_name_free_variable_sa_so(freeVariableNameSa)

def pobyso_name_free_variable_sa_so(freeVariableNameSa):
    """
    Set the free variable name in Sollya from a Sage string.
    """
    sollya_lib_name_free_variable(freeVariableNameSa)

def pobyso_parse_string(string):
    """ Legacy function. See pobyso_parse_string_sa_so. """
    return pobyso_parse_string_sa_so(string)
 
def pobyso_parse_string_sa_so(string):
    """
    Get the Sollya expression computed from a Sage string or
    a Sollya error object if parsing failed.
    """
    return sollya_lib_parse_string(string)

def pobyso_precision_so_sa(ctExpSo):
    """
    Computes the necessary precision to represent a number.
    If x is not zero, it can be uniquely written as x = m · 2e 
    where m is an odd integer and e is an integer. 
    precision(x) returns the number of bits necessary to write m 
    in binary (i.e. ceil(log2(m))).
    """
    #TODO: take care of the special case: 0, @NaN@, @Inf@
    precisionSo = sollya_lib_precision(ctExpSo)
    precisionSa = pobyso_constant_from_int_so_sa(precisionSo)
    sollya_lib_clear_obj(precisionSo)
    return precisionSa
# End pobyso_precision_so_sa

def pobyso_polynomial_coefficients_progressive_round_so_so(polySo,
                                                           funcSo,
                                                           icSo,
                                                           intervalSo,
                                                           itpSo,
                                                           ftpSo,
                                                           maxPrecSo,
                                                           maxErrSo,
                                                           debug=False):
    if debug:
        print "Input arguments:"
        pobyso_autoprint(polySo)
        pobyso_autoprint(funcSo)
        pobyso_autoprint(icSo)
        pobyso_autoprint(intervalSo)
        pobyso_autoprint(itpSo)
        pobyso_autoprint(ftpSo)
        pobyso_autoprint(maxPrecSo)
        pobyso_autoprint(maxErrSo)
        print "________________"
    
    ## Higher order function see: 
    #  http://effbot.org/pyfaq/how-do-you-make-a-higher-order-function-in-python.htm
    def precision_decay_ratio_function(degreeSa):
        def outer(x):
            def inner(x):
                we = 3/8
                wq = 2/8
                a  = 2.2
                b  = 2 
                return we*(exp(x/a)-1) +  wq*((b*x)**2) + (1-we-wq)*x
            return  inner(x)/inner(degreeSa)
        return outer
    
    #
    degreeSa        = pobyso_polynomial_degree_so_sa(polySo)
    ratio           = precision_decay_ratio_function(degreeSa)
    itpSa           = pobyso_constant_from_int_so_sa(itpSo)
    ftpSa           = pobyso_constant_from_int_so_sa(ftpSo)
    maxPrecSa       = pobyso_constant_from_int_so_sa(maxPrecSo)
    maxErrSa        = pobyso_get_constant_as_rn_so_sa(maxErrSo)
    if debug:
        print "degreeSa:", degreeSa
        print "ratio:", ratio
        print "itpsSa:", itpSa
        print "ftpSa:", ftpSa
        print "maxPrecSa:", maxPrecSa
        print "maxErrSa:", maxErrSa
    lastResPolySo   = None
    lastInfNormSo   = None
    #print "About to enter the while loop..."
    while True: 
        resPolySo   = pobyso_constant_0_sa_so()
        pDeltaSa    = ftpSa - itpSa
        for indexSa in reversed(xrange(0,degreeSa+1)):
            #print "Index:", indexSa
            indexSo = pobyso_constant_from_int_sa_so(indexSa)
            #pobyso_autoprint(indexSo)
            #print ratio(indexSa)
            ctpSa = floor(ftpSa - (pDeltaSa * ratio(indexSa)))
            ctpSo = pobyso_constant_from_int_sa_so(ctpSa)
            if debug:
                print "Index:", indexSa, " - Target precision:", 
                pobyso_autoprint(ctpSo)
            cmonSo  = \
                sollya_lib_build_function_mul(sollya_lib_coeff(polySo, indexSo),
                                      sollya_lib_build_function_pow( \
                                          sollya_lib_build_function_free_variable(), \
                                          indexSo))
            #pobyso_autoprint(cmonSo)
            cmonrSo = pobyso_round_coefficients_single_so_so(cmonSo, ctpSo)
            sollya_lib_clear_obj(cmonSo)
            #pobyso_autoprint(cmonrSo)
            resPolySo = sollya_lib_build_function_add(resPolySo,
                                                      cmonrSo)
            #pobyso_autoprint(resPolySo)
        # End for index
        freeVarSo     = sollya_lib_build_function_free_variable()
        changeVarSo   = sollya_lib_sub(freeVarSo, icSo)
        resPolyCvSo   = sollya_lib_evaluate(resPolySo, changeVarSo)
        errFuncSo = sollya_lib_build_function_sub(sollya_lib_copy_obj(funcSo), 
                                                  resPolyCvSo)
        infNormSo = sollya_lib_dirtyinfnorm(errFuncSo, intervalSo)
        cerrSa    = pobyso_get_constant_as_rn_so_sa(infNormSo)
        if debug:
            print "Infnorm (Sollya):", 
            pobyso_autoprint(infNormSo)
        sollya_lib_clear_obj(errFuncSo)
        #print "Infnorm  (Sage):", cerrSa
        if (cerrSa > maxErrSa):
            if debug:
                print "Error is too large."
            if lastResPolySo is None:
                if debug:
                    print "Enlarging prec."
                ntpSa = floor(ftpSa + ftpSa/50)
                ## Can't enlarge (numerical)
                if ntpSa == ftpSa:
                    sollya_lib_clear_obj(resPolySo)
                    return None
                ## Can't enlarge (not enough precision left)
                if ntpSa > maxPrecSa:
                    sollya_lib_clear_obj(resPolySo)
                    return None
                ftpSa = ntpSa
                continue
            ## One enlargement took place.
            else:
                if debug:
                    print "Exit with the last before last polynomial."
                    print "Precision of highest degree monomial:", itpSa
                    print "Precision of constant term          :", ftpSa
                sollya_lib_clear_obj(resPolySo)
                sollya_lib_clear_obj(infNormSo)
                return (lastResPolySo, lastInfNormSo)
        # cerrSa <= maxErrSa: scrap more bits, possibly.
        else:
            if debug:
                print "Error is too small"
            if cerrSa <= (maxErrSa/2):
                if debug: 
                    print "Shrinking prec."
                ntpSa = floor(ftpSa - ftpSa/50)
                ## Can't shrink (numerical)
                if ntpSa == ftpSa:
                    if not lastResPolySo is None:
                        sollya_lib_clear_obj(lastResPolySo)
                    if not lastInfNormSo is None:
                        sollya_lib_clear_obj(lastInfNormSo)
                    if debug:
                        print "Exit because can't shrink anymore (numerically)."
                        print "Precision of highest degree monomial:", itpSa
                        print "Precision of constant term          :", ftpSa
                    return (resPolySo, infNormSo)
                ## Can't shrink (not enough precision left)
                if ntpSa <= itpSa:
                    if not lastResPolySo is None:
                        sollya_lib_clear_obj(lastResPolySo)
                    if not lastInfNormSo is None:
                        sollya_lib_clear_obj(lastInfNormSo)
                        print "Exit because can't shrink anymore (no bits left)."
                        print "Precision of highest degree monomial:", itpSa
                        print "Precision of constant term          :", ftpSa
                    return (resPolySo, infNormSo)
                ftpSa = ntpSa
                if not lastResPolySo is None:
                    sollya_lib_clear_obj(lastResPolySo)
                if not lastInfNormSo is None:
                    sollya_lib_clear_obj(lastInfNormSo)
                lastResPolySo = resPolySo
                lastInfNormSo = infNormSo
                continue
            else: # Error is not that small, just return
                if not lastResPolySo is None:
                    sollya_lib_clear_obj(lastResPolySo)
                if not lastInfNormSo is None:
                    sollya_lib_clear_obj(lastInfNormSo)
                if debug:
                    print "Exit normally."
                    print "Precision of highest degree monomial:", itpSa
                    print "Precision of constant term          :", ftpSa
                return (resPolySo, infNormSo)                
    # End wile True
    return None
# End pobyso_polynomial_coefficients_progressive_truncate_so_so.

def pobyso_polynomial_degree_so_sa(polySo):
    """
    Return the degree of a Sollya polynomial as a Sage int.
    """
    degreeSo = sollya_lib_degree(polySo)
    return pobyso_constant_from_int_so_sa(degreeSo)
# End pobyso_polynomial_degree_so_sa

def pobyso_polynomial_degree_so_so(polySo):
    """
    Thin wrapper around lib_sollya_degree().
    """
    return sollya_lib_degree(polySo)
# End pobyso_polynomial_degree_so_so
                                                                  
def pobyso_range(rnLowerBound, rnUpperBound):
    """ Legacy function. See pobyso_range_sa_so. """
    return pobyso_range_sa_so(rnLowerBound, rnUpperBound) 

def pobyso_range_from_bounds_sa_so(rnLowerBound, rnUpperBound, precSa = None):
    """
    Create a Sollya range from 2 Sage real numbers as bounds
    """
    # TODO check precision stuff.
    sollyaPrecChanged = False
    (initialSollyaPrecSo, initialSollyaPrecSa) = \
        pobyso_get_prec_so_so_sa()
    if precSa is None:
        precSa = max(rnLowerBound.parent().prec(), rnUpperBound.parent().prec())
    if precSa > initialSollyaPrecSa:
        if precSa <= 2:
            print inspect.stack()[0][3], ": precision change <= 2 requested."
        pobyso_set_prec_sa_so(precSa)
        sollyaPrecChanged = True
    rangeSo = sollya_lib_range_from_bounds(get_rn_value(rnLowerBound),
                                           get_rn_value(rnUpperBound))
    if sollyaPrecChanged:
        pobyso_set_prec_so_so(initialSollyaPrecSo)
        sollya_lib_clear_obj(initialSollyaPrecSo)
    return rangeSo
# End pobyso_range_from_bounds_sa_so

def pobyso_range_max_abs_so_so(rangeSo):
    """
    Return, as Sollya constant, the maximum absolute value of a Sollay range.
    """
    lowerBoundSo = sollya_lib_inf(rangeSo)
    upperBoundSo = sollya_lib_sup(rangeSo)
    #
    lowerBoundSo = sollya_lib_build_function_abs(lowerBoundSo)
    upperBoundSo = sollya_lib_build_function_abs(upperBoundSo)
    #pobyso_autoprint(lowerBoundSo)
    #pobyso_autoprint(upperBoundSo)
    #
    maxAbsSo = sollya_lib_max(lowerBoundSo, upperBoundSo, None)
    #sollya_lib_clear_obj(lowerBoundSo)
    #sollya_lib_clear_obj(upperBoundSo)
    return maxAbsSo
# End pobyso_range_max_abs_so_so
 
def pobyso_range_to_interval_so_sa(rangeSo, realIntervalFieldSa = None):
    """
    Get a Sage interval from a Sollya range.
    If no realIntervalField is given as a parameter, the Sage interval
    precision is that of the Sollya range.
    Otherwise, the precision is that of the realIntervalField. In this case
    rounding may happen.
    """
    if realIntervalFieldSa is None:
        precSa = pobyso_get_prec_of_range_so_sa(rangeSo)
        realIntervalFieldSa = RealIntervalField(precSa)
    intervalSa = \
        pobyso_get_interval_from_range_so_sa(rangeSo, realIntervalFieldSa) 
    return intervalSa
# End pobyso_range_to_interval_so_sa
#
def pobyso_relative_so_so():
    """
    Very thin wrapper around the sollya_lib_relative function.
    """
    return sollya_lib_relative()
# End pobyso_relative_so_so
#
def pobyso_rat_poly_sa_so(polySa, precSa = None):
    """
    Create a Sollya polynomial from a Sage rational polynomial.
    """
    ## TODO: filter arguments.
    ## Precision. If no precision is given, use the current precision
    #  of Sollya.
    if precSa is None:
        precSa =  pobyso_get_prec_so_sa()
    #print "Precision:",  precSa
    RRR = RealField(precSa)
    ## Create a Sage polynomial in the "right" precision.
    P_RRR = RRR[polySa.variables()[0]]
    polyFloatSa = P_RRR(polySa)
    ## Make sure no precision is provided: pobyso_float_poly_sa_so will
    #  recover it all by itself and not make an extra conversion.
    return pobyso_float_poly_sa_so(polyFloatSa)
    
# End pobyso_rat_poly_sa_so    
    
def pobyso_remez_canonical_sa_sa(func, \
                                 degree, \
                                 lowerBound, \
                                 upperBound, \
                                 weight = None, \
                                 quality = None):
    """
    All arguments are Sage/Python.
    The functions (func and weight) must be passed as expressions or strings.
    Otherwise the function fails. 
    The return value is a Sage polynomial.
    """
    var('zorglub')    # Dummy variable name for type check only. Type of 
    # zorglub is "symbolic expression".
    polySo = pobyso_remez_canonical_sa_so(func, \
                                 degree, \
                                 lowerBound, \
                                 upperBound, \
                                 weight, \
                                 quality)
    # String test
    if parent(func) == parent("string"):
        functionSa = eval(func)
    # Expression test.
    elif type(func) == type(zorglub):
        functionSa = func
    else:
        return None
    #
    maxPrecision = 0
    if polySo is None:
        return(None)
    maxPrecision = pobyso_get_max_prec_of_exp_so_sa(polySo)
    RRRRSa = RealField(maxPrecision)
    polynomialRingSa = RRRRSa[functionSa.variables()[0]]
    expSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, RRRRSa)
    polySa = polynomial(expSa, polynomialRingSa)
    sollya_lib_clear_obj(polySo)
    return(polySa)
# End pobyso_remez_canonical_sa_sa
    
def pobyso_remez_canonical(func, \
                           degree, \
                           lowerBound, \
                           upperBound, \
                           weight = "1", \
                           quality = None):
    """ Legacy function. See pobyso_remez_canonical_sa_so. """
    return(pobyso_remez_canonical_sa_so(func, \
                                        degree, \
                                        lowerBound, \
                                        upperBound, \
                                        weight, \
                                        quality))
# End pobyso_remez_canonical.

def pobyso_remez_canonical_sa_so(func, \
                                 degree, \
                                 lowerBound, \
                                 upperBound, \
                                 weight = None, \
                                 quality = None):
    """
    All arguments are Sage/Python.
    The functions (func and weight) must be passed as expressions or strings.
    Otherwise the function fails. 
    The return value is a pointer to a Sollya function.
    lowerBound and upperBound mus be reals.
    """
    var('zorglub')    # Dummy variable name for type check only. Type of
    # zorglub is "symbolic expression".
    currentVariableNameSa = None
    # The func argument can be of different types (string, 
    # symbolic expression...)
    if parent(func) == parent("string"):
        localFuncSa = sage_eval(func,globals())
        if len(localFuncSa.variables()) > 0:
            currentVariableNameSa = localFuncSa.variables()[0]
            sollya_lib_name_free_variable(str(currentVariableNameSa))
            functionSo = \
              sollya_lib_parse_string(localFuncSa._assume_str().replace('_SAGE_VAR_', ''))
    # Expression test.
    elif type(func) == type(zorglub):
        # Until we are able to translate Sage expressions into Sollya 
        # expressions : parse the string version.
        if len(func.variables()) > 0:
            currentVariableNameSa = func.variables()[0]
            sollya_lib_name_free_variable(str(currentVariableNameSa))
            functionSo = \
              sollya_lib_parse_string(func._assume_str().replace('_SAGE_VAR_', ''))
    else:
        return(None)
    if weight is None: # No weight given -> 1.
        weightSo = pobyso_constant_1_sa_so()
    elif parent(weight) == parent("string"): # Weight given as string: parse it.
        weightSo = sollya_lib_parse_string(func)
    elif type(weight) == type(zorglub): # Weight given as symbolice expression.
        functionSo = \
          sollya_lib_parse_string_sa_so(weight._assume_str().replace('_SAGE_VAR_', ''))
    else:
        return(None)
    degreeSo = pobyso_constant_from_int(degree)
    rangeSo = pobyso_bounds_to_range_sa_so(lowerBound, upperBound)
    if not quality is None:
        qualitySo= pobyso_constant_sa_so(quality)
    else:
        qualitySo = None
        
    remezPolySo = sollya_lib_remez(functionSo, \
                                   degreeSo, \
                                   rangeSo, \
                                   weightSo, \
                                   qualitySo, \
                                   None)
    sollya_lib_clear_obj(functionSo)
    sollya_lib_clear_obj(degreeSo)
    sollya_lib_clear_obj(rangeSo)
    sollya_lib_clear_obj(weightSo)
    if not qualitySo is None:
        sollya_lib_clear_obj(qualitySo)
    return(remezPolySo)
# End pobyso_remez_canonical_sa_so

def pobyso_remez_canonical_so_so(funcSo, \
                                 degreeSo, \
                                 rangeSo, \
                                 weightSo = pobyso_constant_1_sa_so(),\
                                 qualitySo = None):
    """
    All arguments are pointers to Sollya objects.
    The return value is a pointer to a Sollya function.
    """
    if not sollya_lib_obj_is_function(funcSo):
        return(None)
    return(sollya_lib_remez(funcSo, degreeSo, rangeSo, weightSo, qualitySo, None))
# End pobyso_remez_canonical_so_so.

def pobyso_remez_exponents_list_sa_so(func,     \
                                 exponentsList, \
                                 lowerBound,    \
                                 upperBound,    \
                                 weight = None, \
                                 quality = None):
    """
    All arguments are Sage/Python.
    The functions (func and weight) must be passed as expressions or strings.
    Otherwise the function fails. 
    The return value is a pointer to a Sollya function.
    lowerBound and upperBound mus be reals.
    """
    var('zorglub')    # Dummy variable name for type check only. Type of
    # zorglub is "symbolic expression".
    currentVariableNameSa = None
    # The func argument can be of different types (string, 
    # symbolic expression...)
    if parent(func) == parent("string"):
        localFuncSa = sage_eval(func,globals())
        if len(localFuncSa.variables()) > 0:
            currentVariableNameSa = localFuncSa.variables()[0]
            sollya_lib_name_free_variable(str(currentVariableNameSa))
            functionSo = \
              sollya_lib_parse_string(localFuncSa._assume_str().replace('_SAGE_VAR_', ''))
    # Expression test.
    elif type(func) == type(zorglub):
        # Until we are able to translate Sage expressions into Sollya 
        # expressions : parse the string version.
        if len(func.variables()) > 0:
            currentVariableNameSa = func.variables()[0]
            sollya_lib_name_free_variable(str(currentVariableNameSa))
            functionSo = \
              sollya_lib_parse_string(func._assume_str().replace('_SAGE_VAR_', ''))
    else:
        return(None)
    ## Deal with the weight, much in the same way as with the function.
    if weight is None: # No weight given -> 1.
        weightSo = pobyso_constant_1_sa_so()
    elif parent(weight) == parent("string"): # Weight given as string: parse it.
        weightSo = sollya_lib_parse_string(func)
    elif type(weight) == type(zorglub): # Weight given as symbolice expression.
        functionSo = \
          sollya_lib_parse_string(weight._assume_str().replace('_SAGE_VAR_', '',100))
    else:
        return(None)
    rangeSo = pobyso_bounds_to_range_sa_so(lowerBound, upperBound)
    if not quality is None:
        qualitySo= pobyso_constant_sa_so(quality)
    else:
        qualitySo = None
    #
    ## Tranform the Sage list of exponents into a Sollya list.
    exponentsListSo = pobyso_build_list_of_ints_sa_so(*exponentsList)
    remezPolySo = sollya_lib_remez(functionSo, \
                                   exponentsListSo, \
                                   rangeSo, \
                                   weightSo, \
                                   qualitySo, \
                                   None)
    sollya_lib_clear_obj(functionSo)
    sollya_lib_clear_obj(exponentsListSo)
    sollya_lib_clear_obj(rangeSo)
    sollya_lib_clear_obj(weightSo)
    if not qualitySo is None:
        sollya_lib_clear_obj(qualitySo)
    return(remezPolySo)
# End pobyso_remez_exponentsList_sa_so
#
def pobyso_round_coefficients_so_so(polySo, truncFormatListSo):
    """
    A wrapper around the "classical" sollya_lib_roundcoefficients: a Sollya
    polynomial and a Sollya list are given as arguments.
    """
    return sollya_lib_roundcoefficients(polySo, truncFormatListSo)

def pobyso_round_coefficients_progressive_so_so(polySo, 
                                                funcSo,
                                                precSo,
                                                intervalSo,
                                                icSo,
                                                currentApproxErrorSo,
                                                approxAccurSo,
                                                debug=False):
    if debug:
        print "Input arguments:"
        print "Polynomial: ", ; pobyso_autoprint(polySo)
        print "Function: ", ; pobyso_autoprint(funcSo)
        print "Internal precision: ", ; pobyso_autoprint(precSo)
        print "Interval: ", ; pobyso_autoprint(intervalSo)
        print "Current approximation error: ", ; pobyso_autoprint(currentApproxErrorSo)
        print "Requested approxiation error: ", ; pobyso_autoprint(approxAccurSo)
        print "________________"
    approxAccurSa        = pobyso_get_constant_as_rn_so_sa(approxAccurSo)
    currentApproxErrorSa = pobyso_get_constant_as_rn_so_sa(currentApproxErrorSo)
    ## If the current approximation error is too close to the target, there is
    #  no possible gain.
    if currentApproxErrorSa >= approxAccurSa / 2:
        #### Do not return the initial argument but copies: caller may free 
        #    as inutile after call.
        return (sollya_lib_copy_obj(polySo),
                sollya_lib_copy_obj(currentApproxErrorSo))
    degreeSa             = pobyso_polynomial_degree_so_sa(polySo)
    intervalSa           = pobyso_range_to_interval_so_sa(intervalSo)
    
    if debug:
        print "degreeSa              :", degreeSa
        print "intervalSa            :", intervalSa.str(style='brackets')
        print "currentApproxErrorSa  :", currentApproxErrorSa 
        print "approxAccurSa         :", approxAccurSa 
    ### Start with a 0 value expression.
    radiusSa = intervalSa.absolute_diameter() / 2
    if debug:
        print "log2(radius):", RR(radiusSa).log2()
    iterIndex = 0
    while True: 
        resPolySo = pobyso_constant_0_sa_so()
        roundedPolyApproxAccurSa = approxAccurSa / 2
        currentRadiusPowerSa = 1 
        for degree in xrange(0,degreeSa + 1):
            #### At round 0, use the agressive formula. At round 1, run the
            #    proved formula.
            if iterIndex == 0:
                roundingPowerSa = \
                    floor(((currentRadiusPowerSa/roundedPolyApproxAccurSa)*(degree+1)).log2())
            else:
                roundingPowerSa = \
                    floor(((currentRadiusPowerSa/roundedPolyApproxAccurSa)*(degreeSa+1)).log2())
            ## Under extreme conditions the above formulas can evaluate under 2, which is the
            #  minimal precision of an MPFR number.
            if roundingPowerSa < 2:
                roundingPowerSa = 2
            if debug:
                print "roundedPolyApproxAccurSa", roundedPolyApproxAccurSa
                print "currentRadiusPowerSa", currentRadiusPowerSa
                print "Current rounding exponent:", roundingPowerSa
            currentRadiusPowerSa *= radiusSa
            index1So = pobyso_constant_from_int_sa_so(degree)
            index2So = pobyso_constant_from_int_sa_so(degree)
            ### Create a monomial with:
            #   - the coefficient in the initial monomial at the current degrree;
            #   - the current exponent;
            #   - the free variable.
            cmonSo  = \
                sollya_lib_build_function_mul(sollya_lib_coeff(polySo, index1So),
                                              sollya_lib_build_function_pow( \
                                                sollya_lib_build_function_free_variable(), \
                                                index2So))
            roundingPowerSo = pobyso_constant_from_int_sa_so(roundingPowerSa)
            cmonrSo = pobyso_round_coefficients_single_so_so(cmonSo, roundingPowerSo)
            sollya_lib_clear_obj(cmonSo)
            ### Add to the result polynomial.
            resPolySo = sollya_lib_build_function_add(resPolySo,
                                                      cmonrSo)
        # End for.
        ### Check the new polynomial.
        freeVarSo     = sollya_lib_build_function_free_variable()
        changeVarSo   = sollya_lib_sub(freeVarSo, icSo)
        resPolyCvSo   = sollya_lib_evaluate(resPolySo, changeVarSo)
        errFuncSo     = sollya_lib_build_function_sub(sollya_lib_copy_obj(funcSo), 
                                                      resPolyCvSo)
        infNormSo = sollya_lib_dirtyinfnorm(errFuncSo, intervalSo)
        ### This also clears resPolyCvSo.
        sollya_lib_clear_obj(errFuncSo)
        cerrSa    = pobyso_get_constant_as_rn_so_sa(infNormSo)
        if debug:
            print "Error of the new polynomial:", cerrSa
        ### If at round 1, return the initial polynomial error. This should
        #   never happen since the rounding algorithm is proved. But some 
        #   circumstances may break it (e.g. internal precision of tools).
        if cerrSa > approxAccurSa:
            if iterIndex > 0: # Round 1 and beyond.
                sollya_lib_clear_obj(resPolySo)
                sollya_lib_clear_obj(infNormSo)
                #### Do not return the arguments but copies: the caller may free
                #    free the former as inutile after call.
                return (sollya_lib_copy_obj(polySo), 
                        sollya_lib_copy_obj(currentApproxErrorSo))
            else: # Round 0, got round 1
                sollya_lib_clear_obj(resPolySo)
                sollya_lib_clear_obj(infNormSo)
                iterIndex += 1
                continue
        ### If get here it is because cerrSa <= approxAccurSa
        ### Approximation error of the new polynomial is acceptable.
        return (resPolySo, infNormSo)
    # End while True
# End pobyso_round_coefficients_progressive_so_so
    
def pobyso_round_coefficients_single_so_so(polySo, commonPrecSo):
    """
    Create a rounded coefficients polynomial from polynomial argument to
    the number of bits in size argument.
    All coefficients are set to the same precision.
    """
    ## TODO: check arguments.
    endEllipListSo = pobyso_build_end_elliptic_list_so_so(commonPrecSo)
    polySo = sollya_lib_roundcoefficients(polySo, endEllipListSo, None)
    sollya_lib_clear_obj(endEllipListSo)
    #sollya_lib_clear_obj(endEllipListSo)
    return polySo
    
# End pobyso_round_coefficients_single_so_so

def pobyso_set_canonical_off():
    sollya_lib_set_canonical(sollya_lib_off())

def pobyso_set_canonical_on():
    sollya_lib_set_canonical(sollya_lib_on())

def pobyso_set_prec(p):
    """ Legacy function. See pobyso_set_prec_sa_so. """
    pobyso_set_prec_sa_so(p)

def pobyso_set_prec_sa_so(p):
    #a = c_int(p)
    #precSo = c_void_p(sollya_lib_constant_from_int(a))
    #precSo = sollya_lib_constant_from_int(a)
    precSo = pobyso_constant_from_int_sa_so(p)
    sollya_lib_set_prec(precSo)
    sollya_lib_clear_obj(precSo)
# End pobyso_set_prec_sa_so.

def pobyso_set_prec_so_so(newPrecSo):
    sollya_lib_set_prec(newPrecSo)
# End pobyso_set_prec_so_so.

def pobyso_inf_so_so(intervalSo):
    """
    Very thin wrapper around sollya_lib_inf().
    """
    return sollya_lib_inf(intervalSo)
# End pobyso_inf_so_so.
#   
def pobyso_supnorm_so_sa(polySo, funcSo, intervalSo, errorTypeSo = None,\
                         accuracySo = None, realFieldSa = None):
    """
    Computes the supremum norm from Solly input arguments and returns a
    Sage floating-point number whose precision is set by the only Sage argument.
    
    The returned value is the maximum of the absolute values of the range
    elements  returned by the Sollya supnorm functions
    """
    supNormRangeSo = pobyso_supnorm_so_so(polySo, 
                                          funcSo, 
                                          intervalSo, 
                                          errorTypeSo,
                                          accuracySo)
    supNormSo = pobyso_range_max_abs_so_so(supNormRangeSo)
    sollya_lib_clear_obj(supNormRangeSo)
    #pobyso_autoprint(supNormSo)
    supNormSa = pobyso_get_constant_as_rn_with_rf_so_sa(supNormSo, realFieldSa)
    sollya_lib_clear_obj(supNormSo)
    return supNormSa 
# End pobyso_supnorm_so_sa.
#
def pobyso_supnorm_so_so(polySo, funcSo, intervalSo, errorTypeSo = None,\
                         accuracySo = None):
    """
    Computes the supnorm of the approximation error between the given 
    polynomial and function. Attention: returns a range!
    errorTypeSo defaults to "absolute".
    accuracySo defaults to 2^(-40).
    """
    if errorTypeSo is None:
        errorTypeSo = sollya_lib_absolute(None)
        errorTypeIsNone = True
    else:
        errorTypeIsNone = False
    #
    if accuracySo is None:
        # Notice the **: we are in Pythonland here!
        accuracySo = pobyso_constant_sa_so(RR(2**(-40)))
        accuracyIsNone = True
    else:
        accuracyIsNone = False
    #pobyso_autoprint(accuracySo)
    resultSo = \
        sollya_lib_supnorm(polySo, funcSo, intervalSo, errorTypeSo, \
                              accuracySo)
    if errorTypeIsNone:
        sollya_lib_clear_obj(errorTypeSo)
    if accuracyIsNone:
        sollya_lib_clear_obj(accuracySo)
    return resultSo
# End pobyso_supnorm_so_so
#
def pobyso_taylor_expansion_no_change_var_so_so(functionSo, 
                                                degreeSo, 
                                                rangeSo,
                                                errorTypeSo=None, 
                                                sollyaPrecSo=None):
    """
    Compute the Taylor expansion without the variable change
    x -> x-intervalCenter.
    """
    # Change internal Sollya precision, if needed.
    (initialSollyaPrecSo, initialSollyaPrecSa) = pobyso_get_prec_so_so_sa()
    sollyaPrecChanged = False
    if sollyaPrecSo is None:
        pass
    else:
        sollya_lib_set_prec(sollyaPrecSo)
        sollyaPrecChanged = True
    # Error type stuff: default to absolute.
    if errorTypeSo is None:
        errorTypeIsNone = True
        errorTypeSo = sollya_lib_absolute(None)
    else:
        errorTypeIsNone = False
    intervalCenterSo = sollya_lib_mid(rangeSo, None)
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo,
                                         intervalCenterSo,
                                         rangeSo, errorTypeSo, None)
    # taylorFormListSaSo is a Python list of Sollya objects references that 
    # are copies of the elements of taylorFormSo.
    # pobyso_get_list_elements_so_so clears taylorFormSo.
    (taylorFormListSaSo, numElementsSa, isEndEllipticSa) = \
        pobyso_get_list_elements_so_so(taylorFormSo)
    polySo = sollya_lib_copy_obj(taylorFormListSaSo[0])
    #print "Num elements:", numElementsSa
    sollya_lib_clear_obj(taylorFormSo)
    #polySo = taylorFormListSaSo[0]
    #errorRangeSo = sollya_lib_copy_obj(taylorFormListSaSo[2])
    errorRangeSo = taylorFormListSaSo[2]
    # No copy_obj needed here: a new objects are created.
    maxErrorSo    = sollya_lib_sup(errorRangeSo)
    minErrorSo    = sollya_lib_inf(errorRangeSo)
    absMaxErrorSo = sollya_lib_abs(maxErrorSo)
    absMinErrorSo = sollya_lib_abs(minErrorSo)
    sollya_lib_clear_obj(maxErrorSo)
    sollya_lib_clear_obj(minErrorSo)
    absMaxErrorSa = pobyso_get_constant_as_rn_so_sa(absMaxErrorSo)
    absMinErrorSa = pobyso_get_constant_as_rn_so_sa(absMinErrorSo)
    # If changed, reset the Sollya working precision.
    if sollyaPrecChanged:
        sollya_lib_set_prec(initialSollyaPrecSo)
    sollya_lib_clear_obj(initialSollyaPrecSo)
    if errorTypeIsNone:
        sollya_lib_clear_obj(errorTypeSo)
    pobyso_clear_taylorform_sa_so(taylorFormListSaSo)
    if absMaxErrorSa > absMinErrorSa:
        sollya_lib_clear_obj(absMinErrorSo)
        return (polySo, intervalCenterSo, absMaxErrorSo)
    else:
        sollya_lib_clear_obj(absMaxErrorSo)
        return (polySo, intervalCenterSo, absMinErrorSo)
# end pobyso_taylor_expansion_no_change_var_so_so

def pobyso_taylor_expansion_with_change_var_so_so(functionSo, degreeSo, \
                                                  rangeSo, \
                                                  errorTypeSo=None, \
                                                  sollyaPrecSo=None):
    """
    Compute the Taylor expansion with the variable change
    x -> (x-intervalCenter) included.
    """
    # Change Sollya internal precision, if need.
    sollyaPrecChanged = False
    (initialSollyaPrecSo, initialSollyaPrecSa) = pobyso_get_prec_so_sos_sa()
    if sollyaPrecSo is None:
        pass
    else:
        sollya_lib_set_prec(sollyaPrecSo)
        sollyaPrecChanged = True
    #
    # Error type stuff: default to absolute.
    if errorTypeSo is None:
        errorTypeIsNone = True
        errorTypeSo = sollya_lib_absolute(None)
    else:
        errorTypeIsNone = False
    intervalCenterSo = sollya_lib_mid(rangeSo)
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, \
                                         intervalCenterSo, \
                                         rangeSo, errorTypeSo, None)
    # taylorFormListSaSo is a Python list of Sollya objects references that 
    # are copies of the elements of taylorFormSo.
    # pobyso_get_list_elements_so_so clears taylorFormSo.
    (taylorFormListSo, numElements, isEndElliptic) = \
        pobyso_get_list_elements_so_so(taylorFormSo)
    polySo = taylorFormListSo[0]
    errorRangeSo = taylorFormListSo[2]
    maxErrorSo    = sollya_lib_sup(errorRangeSo)
    minErrorSo    = sollya_lib_inf(errorRangeSo)
    absMaxErrorSo = sollya_lib_abs(maxErrorSo)
    absMinErrorSo = sollya_lib_abs(minErrorSo)
    sollya_lib_clear_obj(maxErrorSo)
    sollya_lib_clear_obj(minErrorSo)
    absMaxErrorSa = pobyso_get_constant_as_rn_so_sa(absMaxErrorSo)
    absMinErrorSa = pobyso_get_constant_as_rn_so_sa(absMinErrorSo)
    changeVarExpSo = sollya_lib_build_function_sub(\
                       sollya_lib_build_function_free_variable(),\
                       sollya_lib_copy_obj(intervalCenterSo))
    polyVarChangedSo = sollya_lib_evaluate(polySo, changeVarExpSo)
    sollya_lib_clear_obj(polySo) 
    sollya_lib_clear_obj(changeVarExpSo)
    # If changed, reset the Sollya working precision.
    if sollyaPrecChanged:
        sollya_lib_set_prec(initialSollyaPrecSo)
    sollya_lib_clear_obj(initialSollyaPrecSo)
    if errorTypeIsNone:
        sollya_lib_clear_obj(errorTypeSo)
    sollya_lib_clear_obj(taylorFormSo)
    # Do not clear maxErrorSo.
    if absMaxErrorSa > absMinErrorSa:
        sollya_lib_clear_obj(absMinErrorSo)
        return((polyVarChangedSo, intervalCenterSo, absMaxErrorSo))
    else:
        sollya_lib_clear_obj(absMaxErrorSo)
        return((polyVarChangedSo, intervalCenterSo, absMinErrorSo))
# end pobyso_taylor_expansion_with_change_var_so_so

def pobyso_taylor(function, degree, point):
    """ Legacy function. See pobysoTaylor_so_so. """
    return(pobyso_taylor_so_so(function, degree, point))

def pobyso_taylor_so_so(functionSo, degreeSo, pointSo):
    return(sollya_lib_taylor(functionSo, degreeSo, pointSo))
    
def pobyso_taylorform(function, degree, point = None, 
                      interval = None, errorType=None):
    """ Legacy function. See pobyso_taylorform_sa_sa;"""
    
def pobyso_taylorform_sa_sa(functionSa, \
                            degreeSa, \
                            pointSa, \
                            intervalSa=None, \
                            errorTypeSa=None, \
                            precisionSa=None):
    """
    Compute the Taylor form of 'degreeSa' for 'functionSa' at 'pointSa' 
    for 'intervalSa' with 'errorTypeSa' (a string) using 'precisionSa'. 
    point: must be a Real or a Real interval.
    return the Taylor form as an array
    TODO: take care of the interval and of the point when it is an interval;
          when errorType is not None;
          take care of the other elements of the Taylor form (coefficients 
          errors and delta.
    """
    # Absolute as the default error.
    if errorTypeSa is None:
        errorTypeSo = sollya_lib_absolute()
    elif errorTypeSa == "relative":
        errorTypeSo = sollya_lib_relative()
    elif errortypeSa == "absolute":
        errorTypeSo = sollya_lib_absolute()
    else:
        # No clean up needed.
        return None
    # Global precision stuff
    sollyaPrecisionChangedSa = False
    (initialSollyaPrecSo, initialSollyaPrecSa) = pobyso_get_prec_so_so_sa()
    if precisionSa is None:
        precSa = initialSollyaPrecSa
    else:
        if precSa > initialSollyaPrecSa:
            if precSa <= 2:
                print inspect.stack()[0][3], ":precision change <= 2 requested."
            pobyso_set_prec_sa_so(precSa)
            sollyaPrecisionChangedSa = True
    #        
    if len(functionSa.variables()) > 0:
        varSa = functionSa.variables()[0]
        pobyso_name_free_variable_sa_so(str(varSa))
    # In any case (point or interval) the parent of pointSa has a precision
    # method.
    pointPrecSa = pointSa.parent().precision()
    if precSa > pointPrecSa:
        pointPrecSa = precSa
    # In any case (point or interval) pointSa has a base_ring() method.
    pointBaseRingString = str(pointSa.base_ring())
    if re.search('Interval', pointBaseRingString) is None: # Point
        pointSo = pobyso_constant_sa_so(pointSa, pointPrecSa)
    else: # Interval.
        pointSo = pobyso_interval_to_range_sa_so(pointSa, pointPrecSa)
    # Sollyafy the function.
    functionSo = pobyso_parse_string_sa_so(functionSa._assume_str().replace('_SAGE_VAR_', ''))
    if sollya_lib_obj_is_error(functionSo):
        print "pobyso_tailorform: function string can't be parsed!"
        return None
    # Sollyafy the degree
    degreeSo = sollya_lib_constant_from_int(int(degreeSa))
    # Sollyafy the point
    # Call Sollya
    taylorFormSo = \
        sollya_lib_taylorform(functionSo, degreeSo, pointSo, errorTypeSo,\
                                         None)
    sollya_lib_clear_obj(functionSo)
    sollya_lib_clear_obj(degreeSo)
    sollya_lib_clear_obj(pointSo)
    sollya_lib_clear_obj(errorTypeSo)
    (tfsAsList, numElements, isEndElliptic) = \
            pobyso_get_list_elements_so_so(taylorFormSo)
    polySo = tfsAsList[0]
    maxPrecision = pobyso_get_max_prec_of_exp_so_sa(polySo)
    polyRealField = RealField(maxPrecision)
    expSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, polyRealField)
    if sollyaPrecisionChangedSa:
        sollya_lib_set_prec(initialSollyaPrecSo)
    sollya_lib_clear_obj(initialSollyaPrecSo)
    polynomialRing = polyRealField[str(varSa)]
    polySa = polynomial(expSa, polynomialRing)
    taylorFormSa = [polySa]
    # Final clean-up.
    sollya_lib_clear_obj(taylorFormSo)
    return(taylorFormSa)
# End pobyso_taylor_form_sa_sa

def pobyso_taylorform_so_so(functionSo, degreeSo, pointSo, intervalSo=None, \
                            errorTypeSo=None):
    createdErrorType = False
    if errorTypeSo is None:
        errorTypeSo = sollya_lib_absolute()
        createdErrorType = True
    else:
        #TODO: deal with the other case.
        pass
    if intervalSo is None:
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
                                         errorTypeSo, None)
    else:
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
                                         intervalSo, errorTypeSo, None)
    if createdErrorType:
        sollya_lib_clear_obj(errorTypeSo)
    return resultSo
        

def pobyso_univar_polynomial_print_reverse(polySa):
    """ Legacy function. See pobyso_univar_polynomial_print_reverse_sa_sa. """
    return(pobyso_univar_polynomial_print_reverse_sa_sa(polySa))

def pobyso_univar_polynomial_print_reverse_sa_sa(polySa):
    """
    Return the string representation of a univariate polynomial with
    monomials ordered in the x^0..x^n order of the monomials.
    Remember: Sage
    """
    polynomialRing = polySa.base_ring()
    # A very expensive solution:
    # -create a fake multivariate polynomial field with only one variable,
    #   specifying a negative lexicographical order;
    mpolynomialRing = PolynomialRing(polynomialRing.base(), \
                                     polynomialRing.variable_name(), \
                                     1, order='neglex')
    # - convert the univariate argument polynomial into a multivariate
    #   version;
    p = mpolynomialRing(polySa)
    # - return the string representation of the converted form.
    # There is no simple str() method defined for p's class.
    return(p.__str__())
#
#print pobyso_get_prec()  
pobyso_set_prec(165)
#print pobyso_get_prec()  
#a=100
#print type(a)
#id(a)
#print "Max arity: ", pobyso_max_arity
#print "Function tripleDouble (43) as a string: ", pobyso_function_type_as_string(43)
#print "Function None (44) as a string: ", pobyso_function_type_as_string(44)
sys.stderr.write("\t...Pobyso check done.\n")