Statistiques
| Révision :

root / pobysoPythonSage / src / pobyso.py @ 66

Historique | Voir | Annoter | Télécharger (37,03 ko)

1
"""
2
Actual functions to use in Sage
3
ST 2012-11-13
4

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

8
pobyso functions come in five flavors: 
9
- the _so_so (arguments and returned objects are pointers to Sollya objects, 
10
  includes the void function and the no arguments function that return a 
11
  pointer to a Sollya object);
12
- the _so_sa (argument are pointers to Sollya objects, returned objects are
13
  Sage/Python objects or, more generally, information is transfered from the
14
  Sollya world to Sage/Python world; e.g. functions without arguments that
15
  return a Sage/Python object);
16
- the _sa_so (arguments are Sage/Python objects, returned objects are 
17
  pointers to Sollya objects);
18
- the sa_sa (arguments and returned objects are all Sage/Python objects);
19
- a catch all flavor, without any suffix, (e. g. functions that have no argument 
20
  nor return value).
21
This classification is not always very strict. Conversion functions from Sollya
22
to Sage/Python are sometimes decorated with Sage/Python arguments to set
23
the precision. These functions remain in the so_sa category.
24
NOTES:
25
Reported errors in Eclipse come from the calls to
26
the Sollya library
27

28
ToDo (among other things): 
29
 -memory management.
30
"""
31
from ctypes import *
32
import re
33
from sage.symbolic.expression_conversions import polynomial
34
from sage.symbolic.expression_conversions import PolynomialConverter
35
"""
36
Create the equivalent to an enum for the Sollya function types.
37
"""
38
(SOLLYA_BASE_FUNC_ABS,
39
SOLLYA_BASE_FUNC_ACOS,
40
    SOLLYA_BASE_FUNC_ACOSH,
41
    SOLLYA_BASE_FUNC_ADD,
42
    SOLLYA_BASE_FUNC_ASIN,
43
    SOLLYA_BASE_FUNC_ASINH,
44
    SOLLYA_BASE_FUNC_ATAN,
45
    SOLLYA_BASE_FUNC_ATANH,
46
    SOLLYA_BASE_FUNC_CEIL,
47
    SOLLYA_BASE_FUNC_CONSTANT,
48
    SOLLYA_BASE_FUNC_COS,
49
    SOLLYA_BASE_FUNC_COSH,
50
    SOLLYA_BASE_FUNC_DIV,
51
    SOLLYA_BASE_FUNC_DOUBLE,
52
    SOLLYA_BASE_FUNC_DOUBLEDOUBLE,
53
    SOLLYA_BASE_FUNC_DOUBLEEXTENDED,
54
    SOLLYA_BASE_FUNC_ERF,
55
    SOLLYA_BASE_FUNC_ERFC,
56
    SOLLYA_BASE_FUNC_EXP,
57
    SOLLYA_BASE_FUNC_EXP_M1,
58
    SOLLYA_BASE_FUNC_FLOOR,
59
    SOLLYA_BASE_FUNC_FREE_VARIABLE,
60
    SOLLYA_BASE_FUNC_HALFPRECISION,
61
    SOLLYA_BASE_FUNC_LIBRARYCONSTANT,
62
    SOLLYA_BASE_FUNC_LIBRARYFUNCTION,
63
    SOLLYA_BASE_FUNC_LOG,
64
    SOLLYA_BASE_FUNC_LOG_10,
65
    SOLLYA_BASE_FUNC_LOG_1P,
66
    SOLLYA_BASE_FUNC_LOG_2,
67
    SOLLYA_BASE_FUNC_MUL,
68
    SOLLYA_BASE_FUNC_NEARESTINT,
69
    SOLLYA_BASE_FUNC_NEG,
70
    SOLLYA_BASE_FUNC_PI,
71
    SOLLYA_BASE_FUNC_POW,
72
    SOLLYA_BASE_FUNC_PROCEDUREFUNCTION,
73
    SOLLYA_BASE_FUNC_QUAD,
74
    SOLLYA_BASE_FUNC_SIN,
75
    SOLLYA_BASE_FUNC_SINGLE,
76
    SOLLYA_BASE_FUNC_SINH,
77
    SOLLYA_BASE_FUNC_SQRT,
78
    SOLLYA_BASE_FUNC_SUB,
79
    SOLLYA_BASE_FUNC_TAN,
80
    SOLLYA_BASE_FUNC_TANH,
81
SOLLYA_BASE_FUNC_TRIPLEDOUBLE) = map(int,xrange(44))
82
print "\nSuperficial pobyso check..."    
83
print "First constant - SOLLYA_BASE_FUNC_ABS: ", SOLLYA_BASE_FUNC_ABS
84
print "Last constant  - SOLLYA_BASE_FUNC_TRIPLEDOUBLE: ", SOLLYA_BASE_FUNC_TRIPLEDOUBLE
85

    
86
pobyso_max_arity = 9
87

    
88
def pobyso_autoprint(arg):
89
    sollya_lib_autoprint(arg,None)
90

    
91
def pobyso_autoprint_so_so(arg):
92
    sollya_lib_autoprint(arg,None)
93
    
94
def pobyso_absolute_so_so():
95
    return(sollya_lib_absolute(None))
96

    
97
def pobyso_build_function_sub_so_so(exp1So, exp2So):
98
    return(sollya_lib_build_function_sub(exp1So, exp2So))
99

    
100
def pobyso_cmp(rnArg, soCte):
101
    """
102
    Compare the MPFR value a RealNumber with that of a Sollya constant.
103
    
104
    Get the value of the Sollya constant into a RealNumber and compare
105
    using MPFR. Could be optimized by working directly with a mpfr_t
106
    for the intermediate number. 
107
    """
108
    precisionOfCte = c_int(0)
109
    # From the Sollya constant, create a local Sage RealNumber.
110
    sollya_lib_get_prec_of_constant(precisionOfCte, soCte) 
111
    #print "Precision of constant: ", precisionOfCte
112
    RRRR = RealField(precisionOfCte.value)
113
    rnLocal = RRRR(0)
114
    sollya_lib_get_constant(get_rn_value(rnLocal), soCte)
115
    #print "rnDummy: ", rnDummy
116
    # Compare the local Sage RealNumber with rnArg.
117
    return(cmp_rn_value(rnArg, rnLocal))
118

    
119
def pobyso_change_var_in_function_so_so(funcSo, chvarExpSo):
120
    return(sollya_lib_evaluate(funcSo,chvarExpSo))
121
    
122

    
123
def pobyso_chebyshevform_so_so(functionSo, degreeSo, intervalSo):
124
    resultSo = sollya_lib_chebyshevform(functionSo, degreeSo, intervalSo)
125
    return(resultSo)
126
        
127

    
128

    
129
def pobyso_compute_pos_function_abs_val_bounds_sa_sa(funcSa, lowerBoundSa, \
130
                                                     upperBoundSa):
131
    """
132
    TODO: set the variable name in Sollya.
133
    """
134
    funcSo = pobyso_parse_string(funcSa._assume_str())
135
    rangeSo = pobyso_range_sa_so(lowerBoundSa, upperBoundSa)
136
    infnormSo = pobyso_infnorm_so_so(funcSo,rangeSo)
137
    fMaxSa = pobyso_get_interval_from_range_so_sa(infnormSo)
138
    # Get the top bound and compute the binade top limit.
139
    fMaxUpperBoundSa = fMaxSa.upper()
140
    binadeTopLimitSa = 2**ceil(fMaxUpperBoundSa.log2())
141
    # Put up together the function to use to compute the lower bound.
142
    funcAuxSo = pobyso_parse_string(str(binadeTopLimitSa) +  \
143
                                    '-(' + f._assume_str() + ')')
144
    pobyso_autoprint(funcAuxSo)
145
    # Clear the Sollay range before a new call to infnorm and issue the call.
146
    sollya_lib_clear_obj(infnormSo)
147
    infnormSo = pobyso_infnorm_so_so(funcAuxSo,rangeSo)
148
    fMinSa = pobyso_get_interval_from_range_so_sa(infnormSo)
149
    sollya_lib_clear_obj(infnormSo)
150
    fMinLowerBoundSa = topBinadeLimit - fMinSa.lower()
151
    # Compute the maximum of the precisions of the different bounds.
152
    maxPrecSa = max([fMinLowerBoundSa.parent().precision(), \
153
                     fMaxUpperBoundSa.parent().precision()])
154
    # Create a RealIntervalField and create an interval with the "good" bounds.
155
    RRRI = RealIntervalField(maxPrecSa)
156
    imageIntervalSa = RRRI(fMinLowerBoundSa, fMaxUpperBoundSa)
157
    # Free the uneeded Sollya objects
158
    sollya_lib_clear_obj(funcSo)
159
    sollya_lib_clear_obj(funcAuxSo)
160
    sollya_lib_clear_obj(rangeSo)
161
    return(imageIntervalSa)
162
    # End pobyso_compute_function_abs_val_bounds_sa_sa
163

    
164
def pobyso_constant(rnArg):
165
    """ Legacy function. See pobyso_constant_sa_so. """
166
    return(pobyso_constant_sa_so(rnArg))
167
    
168
def pobyso_constant_sa_so(rnArg):
169
    """
170
    Create a Sollya constant from a RealNumber.
171
    """
172
    return(sollya_lib_constant(get_rn_value(rnArg)))
173
    
174
def pobyso_constant_0_sa_so():
175
    return(pobyso_constant_from_int_sa_so(0))
176

    
177
def pobyso_constant_1():
178
    """ Legacy function. See pobyso_constant_so_so. """
179
    return(pobyso_constant_1_sa_so())
180

    
181
def pobyso_constant_1_sa_so():
182
    return(pobyso_constant_from_int_sa_so(1))
183

    
184
def pobyso_constant_from_int(anInt):
185
    """ Legacy function. See pobyso_constant_from_int_sa_so. """
186
    return(pobyso_constant_from_int_sa_so(anInt))
187

    
188
def pobyso_constant_from_int_sa_so(anInt):
189
    return(sollya_lib_constant_from_int(int(anInt)))
190

    
191
def pobyso_function_type_as_string(funcType):
192
    """ Legacy function. See pobyso_function_type_as_string_so_sa. """
193
    return(pobyso_function_type_as_string_so_sa(funcType))
194

    
195
def pobyso_function_type_as_string_so_sa(funcType):
196
    """
197
    Numeric Sollya function codes -> Sage mathematical function names.
198
    Notice that pow -> ^ (a la Sage, not a la Python).
199
    """
200
    if funcType == SOLLYA_BASE_FUNC_ABS:
201
        return "abs"
202
    elif funcType == SOLLYA_BASE_FUNC_ACOS:
203
        return "arccos"
204
    elif funcType == SOLLYA_BASE_FUNC_ACOSH:
205
        return "arccosh"
206
    elif funcType == SOLLYA_BASE_FUNC_ADD:
207
        return "+"
208
    elif funcType == SOLLYA_BASE_FUNC_ASIN:
209
        return "arcsin"
210
    elif funcType == SOLLYA_BASE_FUNC_ASINH:
211
        return "arcsinh"
212
    elif funcType == SOLLYA_BASE_FUNC_ATAN:
213
        return "arctan"
214
    elif funcType == SOLLYA_BASE_FUNC_ATANH:
215
        return "arctanh"
216
    elif funcType == SOLLYA_BASE_FUNC_CEIL:
217
        return "ceil"
218
    elif funcType == SOLLYA_BASE_FUNC_CONSTANT:
219
        return "cte"
220
    elif funcType == SOLLYA_BASE_FUNC_COS:
221
        return "cos"
222
    elif funcType == SOLLYA_BASE_FUNC_COSH:
223
        return "cosh"
224
    elif funcType == SOLLYA_BASE_FUNC_DIV:
225
        return "/"
226
    elif funcType == SOLLYA_BASE_FUNC_DOUBLE:
227
        return "double"
228
    elif funcType == SOLLYA_BASE_FUNC_DOUBLEDOUBLE:
229
        return "doubleDouble"
230
    elif funcType == SOLLYA_BASE_FUNC_DOUBLEEXTENDED:
231
        return "doubleDxtended"
232
    elif funcType == SOLLYA_BASE_FUNC_ERF:
233
        return "erf"
234
    elif funcType == SOLLYA_BASE_FUNC_ERFC:
235
        return "erfc"
236
    elif funcType == SOLLYA_BASE_FUNC_EXP:
237
        return "exp"
238
    elif funcType == SOLLYA_BASE_FUNC_EXP_M1:
239
        return "expm1"
240
    elif funcType == SOLLYA_BASE_FUNC_FLOOR:
241
        return "floor"
242
    elif funcType == SOLLYA_BASE_FUNC_FREE_VARIABLE:
243
        return "freeVariable"
244
    elif funcType == SOLLYA_BASE_FUNC_HALFPRECISION:
245
        return "halfPrecision"
246
    elif funcType == SOLLYA_BASE_FUNC_LIBRARYCONSTANT:
247
        return "libraryConstant"
248
    elif funcType == SOLLYA_BASE_FUNC_LIBRARYFUNCTION:
249
        return "libraryFunction"
250
    elif funcType == SOLLYA_BASE_FUNC_LOG:
251
        return "log"
252
    elif funcType == SOLLYA_BASE_FUNC_LOG_10:
253
        return "log10"
254
    elif funcType == SOLLYA_BASE_FUNC_LOG_1P:
255
        return "log1p"
256
    elif funcType == SOLLYA_BASE_FUNC_LOG_2:
257
        return "log2"
258
    elif funcType == SOLLYA_BASE_FUNC_MUL:
259
        return "*"
260
    elif funcType == SOLLYA_BASE_FUNC_NEARESTINT:
261
        return "round"
262
    elif funcType == SOLLYA_BASE_FUNC_NEG:
263
        return "__neg__"
264
    elif funcType == SOLLYA_BASE_FUNC_PI:
265
        return "pi"
266
    elif funcType == SOLLYA_BASE_FUNC_POW:
267
        return "^"
268
    elif funcType == SOLLYA_BASE_FUNC_PROCEDUREFUNCTION:
269
        return "procedureFunction"
270
    elif funcType == SOLLYA_BASE_FUNC_QUAD:
271
        return "quad"
272
    elif funcType == SOLLYA_BASE_FUNC_SIN:
273
        return "sin"
274
    elif funcType == SOLLYA_BASE_FUNC_SINGLE:
275
        return "single"
276
    elif funcType == SOLLYA_BASE_FUNC_SINH:
277
        return "sinh"
278
    elif funcType == SOLLYA_BASE_FUNC_SQRT:
279
        return "sqrt"
280
    elif funcType == SOLLYA_BASE_FUNC_SUB:
281
        return "-"
282
    elif funcType == SOLLYA_BASE_FUNC_TAN:
283
        return "tan"
284
    elif funcType == SOLLYA_BASE_FUNC_TANH:
285
        return "tanh"
286
    elif funcType == SOLLYA_BASE_FUNC_TRIPLEDOUBLE:
287
        return "tripleDouble"
288
    else:
289
        return None
290

    
291
def pobyso_get_constant(rnArg, soConst):
292
    """ Legacy function. See pobyso_get_constant_so_sa. """
293
    return(pobyso_get_constant_so_sa(rnArg, soConst))
294

    
295
def pobyso_get_constant_so_sa(rnArg, soConst):
296
    """
297
    Set the value of rnArg to the value of soConst in MPFR_RNDN mode.
298
    rnArg must already exist and belong to some RealField.
299
    We assume that soConst points to a Sollya constant.
300
    """
301
    return(sollya_lib_get_constant(get_rn_value(rnArg), soConst))
302
    
303
def pobyso_get_constant_as_rn(ctExpSo):
304
    """ Legacy function. See pobyso_get_constant_as_rn_so_sa. """ 
305
    return(pobyso_get_constant_as_rn_so_sa(ctExpSo))
306
    
307
def pobyso_get_constant_as_rn_so_sa(constExpSo):
308
    precisionSa  = pobyso_get_prec_of_constant_so_sa(constExpSo) 
309
    RRRR = RealField(precisionSa)
310
    rnSa = RRRR(0)
311
    sollya_lib_get_constant(get_rn_value(rnSa), constExpSo)
312
    return(rnSa)
313

    
314
def pobyso_get_constant_as_rn_with_rf(ctExp, realField):
315
    """ Legacy function. See pobyso_get_constant_as_rn_with_rf_so_sa."""
316
    return(pobyso_get_constant_as_rn_with_rf_so_sa(ctExp, realField))
317
    
318
def pobyso_get_constant_as_rn_with_rf_so_sa(ctExpSo, realFieldSa = None):
319
    if realFieldSa is None:
320
        sollyaPrecSa = pobyso_get_prec_so_sa()
321
        realFieldSa = RealField(sollyaPrecSa)
322
    rnSa = realFieldSa(0)
323
    sollya_lib_get_constant(get_rn_value(rnSa), ctExpSo)
324
    return(rnSa)
325

    
326
def pobyso_get_free_variable_name():
327
    """ Legacy function. See pobyso_get_free_variable_name_so_sa."""
328
    return(pobyso_get_free_variable_name_so_sa())
329

    
330
def pobyso_get_free_variable_name_so_sa():
331
    return(sollya_lib_get_free_variable_name())
332
    
333
def pobyso_get_function_arity(expressionSo):
334
    """ Legacy function. See pobyso_get_function_arity_so_sa."""
335
    return(pobyso_get_function_arity_so_sa(expressionSo))
336

    
337
def pobyso_get_function_arity_so_sa(expressionSo):
338
    arity = c_int(0)
339
    sollya_lib_get_function_arity(byref(arity),expressionSo)
340
    return(int(arity.value))
341

    
342
def pobyso_get_head_function(expressionSo):
343
    """ Legacy function. See pobyso_get_head_function_so_sa. """
344
    return(pobyso_get_head_function_so_sa(expressionSo)) 
345

    
346
def pobyso_get_head_function_so_sa(expressionSo):
347
    functionType = c_int(0)
348
    sollya_lib_get_head_function(byref(functionType), expressionSo, None)
349
    return(int(functionType.value))
350

    
351
def pobyso_get_interval_from_range_so_sa(soRange, realIntervalFieldSa = None ):
352
    """
353
    Return the Sage interval corresponding to the Sollya range argument.
354
    If no reaInterval lField is passed as argument, the interval bounds are not
355
    rounded: they are elements of RealIntervalField of the "right" precision
356
    to hold all the digits.
357
    """
358
    prec = c_int(0)
359
    if realIntervalFieldSa is None:
360
        retval = sollya_lib_get_prec_of_range(byref(prec), soRange, None)
361
        if retval == 0:
362
            return(None)
363
        realIntervalFieldSa = RealIntervalField(prec.value)
364
    intervalSa = realIntervalFieldSa(0,0)
365
    retval = \
366
        sollya_lib_get_interval_from_range(get_interval_value(intervalSa),\
367
                                           soRange)
368
    if retval == 0:
369
        return(None)
370
    return(intervalSa)
371
# End pobyso_get_interval_from_range_so_sa
372

    
373
def pobyso_get_list_elements(soObj):
374
    """ Legacy function. See pobyso_get_list_elements_so_so. """
375
    return(pobyso_get_list_elements_so_so(soObj))
376
 
377
def pobyso_get_list_elements_so_so(soObj):
378
    """
379
    Get the list elements as a Sage/Python array of Sollya objects.
380
    The other data returned are also Sage/Python objects.
381
    """
382
    listAddress = POINTER(c_longlong)()
383
    numElements = c_int(0)
384
    isEndElliptic = c_int(0)
385
    listAsList = []
386
    result = sollya_lib_get_list_elements(byref(listAddress),\
387
                                          byref(numElements),\
388
                                          byref(isEndElliptic),\
389
                                          soObj)
390
    if result == 0 :
391
        return None
392
    for i in xrange(0, numElements.value, 1):
393
        listAsList.append(listAddress[i])
394
    return(listAsList, numElements.value, isEndElliptic.value)
395

    
396
def pobyso_get_max_prec_of_exp(soExp):
397
    """ Legacy function. See pobyso_get_max_prec_of_exp_so_sa. """
398
    return(pobyso_get_max_prec_of_exp_so_sa(soExp))
399

    
400
def pobyso_get_max_prec_of_exp_so_sa(soExp):
401
    """
402
    Get the maximum precision used for the numbers in a Sollya expression.
403
    
404
    Arguments:
405
    soExp -- a Sollya expression pointer
406
    Return value:
407
    A Python integer
408
    TODO: 
409
    - error management;
410
    - correctly deal with numerical type such as DOUBLEEXTENDED.
411
    """
412
    maxPrecision = 0
413
    minConstPrec = 0
414
    currentConstPrec = 0
415
    operator = pobyso_get_head_function_so_sa(soExp)
416
    if (operator != SOLLYA_BASE_FUNC_CONSTANT) and \
417
    (operator != SOLLYA_BASE_FUNC_FREE_VARIABLE):
418
        (arity, subexpressions) = pobyso_get_subfunctions_so_sa(soExp)
419
        for i in xrange(arity):
420
            maxPrecisionCandidate = \
421
                pobyso_get_max_prec_of_exp_so_sa(subexpressions[i])
422
            if maxPrecisionCandidate > maxPrecision:
423
                maxPrecision = maxPrecisionCandidate
424
        return(maxPrecision)
425
    elif operator == SOLLYA_BASE_FUNC_CONSTANT:
426
        minConstPrec = pobyso_get_min_prec_of_constant_so_sa(soExp)
427
        #currentConstPrec = pobyso_get_min_prec_of_constant_so_sa(soExp)
428
        #print minConstPrec, " - ", currentConstPrec 
429
        return(pobyso_get_min_prec_of_constant_so_sa(soExp))
430
    
431
    elif operator == SOLLYA_BASE_FUNC_FREE_VARIABLE:
432
        return(0)
433
    else:
434
        print "pobyso_get_max_prec_of_exp_so_sa: unexepected operator."
435
        return(0)
436

    
437
def pobyso_get_min_prec_of_constant_so_sa(soConstExp):
438
    """
439
    Get the minimum precision necessary to represent the value of a Sollya
440
    constant.
441
    MPFR_MIN_PREC and powers of 2 are taken into account.
442
    We assume that soCteExp is a point
443
    """
444
    constExpAsRn = pobyso_get_constant_as_rn_so_sa(soConstExp)
445
    return(min_mpfr_size(get_rn_value(constExpAsRn)))
446

    
447
def pobyso_get_sage_exp_from_sollya_exp(sollyaExp, realField = RR):
448
    """ Legacy function. See pobyso_get_sage_exp_from_sollya_exp_so_sa. """
449
    return(pobyso_get_sage_exp_from_sollya_exp_so_sa(sollyaExp, realField = RR))
450

    
451
def pobyso_get_sage_exp_from_sollya_exp_so_sa(sollyaExp, realField = RR):
452
    """
453
    Get a Sage expression from a Sollya expression. 
454
    Currently only tested with polynomials with floating-point coefficients.
455
    Notice that, in the returned polynomial, the exponents are RealNumbers.
456
    """
457
    #pobyso_autoprint(sollyaExp)
458
    operator = pobyso_get_head_function_so_sa(sollyaExp)
459
    # Constants and the free variable are special cases.
460
    # All other operator are dealt with in the same way.
461
    if (operator != SOLLYA_BASE_FUNC_CONSTANT) and \
462
       (operator != SOLLYA_BASE_FUNC_FREE_VARIABLE):
463
        (arity, subexpressions) = pobyso_get_subfunctions_so_sa(sollyaExp)
464
        if arity == 1:
465
            sageExp = eval(pobyso_function_type_as_string_so_sa(operator) + \
466
            "(" + pobyso_get_sage_exp_from_sollya_exp_so_sa(subexpressions[0], \
467
            realField) + ")")
468
        elif arity == 2:
469
            # We do not get through the preprocessor.
470
            # The "^" operator is then a special case.
471
            if operator == SOLLYA_BASE_FUNC_POW:
472
                operatorAsString = "**"
473
            else:
474
                operatorAsString = \
475
                    pobyso_function_type_as_string_so_sa(operator)
476
            sageExp = \
477
              eval("pobyso_get_sage_exp_from_sollya_exp_so_sa(subexpressions[0], realField)"\
478
              + " " + operatorAsString + " " + \
479
                   "pobyso_get_sage_exp_from_sollya_exp_so_sa(subexpressions[1], realField)")
480
        # We do not know yet how to deal with arity >= 3 
481
        # (is there any in Sollya anyway?).
482
        else:
483
            sageExp = eval('None')
484
        return(sageExp)
485
    elif operator == SOLLYA_BASE_FUNC_CONSTANT:
486
        #print "This is a constant"
487
        return pobyso_get_constant_as_rn_with_rf_so_sa(sollyaExp, realField)
488
    elif operator == SOLLYA_BASE_FUNC_FREE_VARIABLE:
489
        #print "This is free variable"
490
        return(eval(sollya_lib_get_free_variable_name()))
491
    else:
492
        print "Unexpected"
493
        return eval('None')
494
# End pobyso_get_sage_poly_from_sollya_poly
495
def pobyso_get_poly_so_sa(polySo, realFieldSa=None):
496
    """
497
    Convert a Sollya polynomial into a Sage polynomial.
498
    We assume that the polynomial is in canonical form.
499
    If no realField is given, a RealField corresponding to the maximum precision 
500
    of the coefficients is internally computed.
501
    It is not returned but can be easily retrieved from the polynomial itself.
502
    Main steps:
503
    - (optional) compute the RealField of the coefficients;
504
    - convert the Sollya expression into a Sage expression;
505
    - convert the Sage expression into a Sage polynomial
506
    TODO: the canonical thing for the polynomial.
507
    """    
508
    if realFieldSa is None:
509
        expressionPrecSa = pobyso_get_max_prec_of_exp_so_sa(polySo)
510
        realFieldSa = RealField(expressionPrecSa)
511
    #print "Sollya expression before...",
512
    #pobyso_autoprint(polySo)
513

    
514
    expressionSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, \
515
                                                             realFieldSa)
516
    #print "...Sollya expression after.",
517
    #pobyso_autoprint(polySo)
518
    polyVariableSa = expressionSa.variables()[0]
519
    polyRingSa = realFieldSa[str(polyVariableSa)]
520
    print polyRingSa
521
    # Do not use the polynomial(expressionSa, ring=polyRingSa) form!
522
    polynomialSa = polyRingSa(expressionSa)
523
    return(polynomialSa)
524
# End pobyso_get_sage_poly_from_sollya_poly
525

    
526
def pobyso_get_subfunctions(expressionSo):
527
    """ Legacy function. See pobyso_get_subfunctions_so_sa. """
528
    return(pobyso_get_subfunctions_so_sa(expressionSo)) 
529

    
530
def pobyso_get_subfunctions_so_sa(expressionSo):
531
    """
532
    Get the subfunctions of an expression.
533
    Return the number of subfunctions and the list of subfunctions addresses.
534
    S.T.: Could not figure out another way than that ugly list of declarations
535
    to recover the addresses of the subfunctions. 
536
    """
537
    subf0 = c_int(0)
538
    subf1 = c_int(0)
539
    subf2 = c_int(0)
540
    subf3 = c_int(0)
541
    subf4 = c_int(0)
542
    subf5 = c_int(0)
543
    subf6 = c_int(0)
544
    subf7 = c_int(0)
545
    subf8 = c_int(0)
546
    arity = c_int(0)
547
    nullPtr = POINTER(c_int)()
548
    sollya_lib_get_subfunctions(expressionSo, byref(arity), \
549
    byref(subf0), byref(subf1), byref(subf2), byref(subf3), byref(subf4), byref(subf5),\
550
     byref(subf6), byref(subf7), byref(subf8), nullPtr, None) 
551
#    byref(cast(subfunctions[0], POINTER(c_int))), byref(cast(subfunctions[0], POINTER(c_int))), \
552
#    byref(cast(subfunctions[2], POINTER(c_int))), byref(cast(subfunctions[3], POINTER(c_int))), \
553
#    byref(cast(subfunctions[4], POINTER(c_int))), byref(cast(subfunctions[5], POINTER(c_int))), \
554
#    byref(cast(subfunctions[6], POINTER(c_int))), byref(cast(subfunctions[7], POINTER(c_int))), \
555
#    byref(cast(subfunctions[8], POINTER(c_int))), nullPtr)
556
    subfunctions = [subf0, subf1, subf2, subf3, subf4, subf5, subf6, subf7, subf8]
557
    subs = []
558
    if arity.value > pobyso_max_arity:
559
        return(0,[])
560
    for i in xrange(arity.value):
561
        subs.append(int(subfunctions[i].value))
562
        #print subs[i]
563
    return(int(arity.value), subs)
564
    
565
def pobyso_get_prec():
566
    """ Legacy function. See pobyso_get_prec_so_sa(). """
567
    return(pobyso_get_prec_so_sa())
568

    
569
def pobyso_get_prec_so_sa():
570
    """
571
    Get the current default precision in Sollya.
572
    The return value is Sage/Python int.
573
    """
574
    precSo = sollya_lib_get_prec(None)
575
    precSa = c_int(0)
576
    sollya_lib_get_constant_as_int(byref(precSa), precSo)
577
    sollya_lib_clear_obj(precSo)
578
    return(int(precSa.value))
579

    
580
def pobyso_get_prec_of_constant(ctExpSo):
581
    """ Legacy function. See pobyso_get_prec_of_constant_so_sa. """
582
    return(pobyso_get_prec_of_constant_so_sa(ctExpSo))
583

    
584
def pobyso_get_prec_of_constant_so_sa(ctExpSo):
585
    prec = c_int(0)
586
    retc = sollya_lib_get_prec_of_constant(byref(prec), ctExpSo, None)
587
    if retc == 0:
588
        return(None)
589
    return(int(prec.value))
590

    
591
def pobyso_get_prec_of_range_so_sa(rangeSo):
592
    prec = c_int(0)
593
    retc = sollya_lib_get_prec_of_range(byref(prec), rangeSo, None)
594
    if retc == 0:
595
        return(None)
596
    return(int(prec.value))
597

    
598
def pobyso_infnorm_so_so(func, interval, file = None, intervalList = None):
599
    print "Do not use this function. User pobyso_supnorm_so_so instead."
600
    return(None)
601

    
602
def pobyso_lib_init():
603
    sollya_lib_init(None)
604
    
605
def pobyso_name_free_variable(freeVariableName):
606
    """ Legacy function. See pobyso_name_free_variable_sa_so. """
607
    pobyso_name_free_variable_sa_so(freeVariableName)
608

    
609
def pobyso_name_free_variable_sa_so(freeVariableName):
610
    sollya_lib_name_free_variable(freeVariableName)
611

    
612
def pobyso_parse_string(string):
613
    """ Legacy function. See pobyso_parse_string_sa_so. """
614
    return(pobyso_parse_string_sa_so(string))
615
 
616
def pobyso_parse_string_sa_so(string):
617
    return(sollya_lib_parse_string(string))
618

    
619
def pobyso_range(rnLowerBound, rnUpperBound):
620
    """ Legacy function. See pobyso_range_sa_so. """
621
    return(pobyso_range_sa_so(rnLowerBound, rnUpperBound)) 
622

    
623
def pobyso_bounds_to_range_sa_so(rnLowerBoundSa, rnUpperBoundSa):
624
    """
625
    Return a Sollya range from to 2 RealField elements.
626
    The Sollya range element has a sufficient precision to hold all
627
    the digits of the bounds.
628
    """
629
    if rnLowerBoundSa > rnUpperBoundSa:
630
        return None
631
    lbPrec = rnLowerBoundSa.parent().precision()
632
    ubPrec = rnLowerBoundSa.parent().precision()
633
    currentSollyaPrecSa = pobyso_get_prec_so_sa()
634
    maxPrecSa = max(lbPrec, ubPrec, currentSollyaPrecSa)
635
    # Change the current Sollya precision only if necessary.
636
    if maxPrecSa > currentSollyaPrecSa:
637
        currentPrecSo = sollya_lib_get_prec(None)
638
        newPrecSo = solly_lib_constant_from_uint64(maxPrecSa)
639
        sollya_lib_set_prec(newPrecSo)
640
    lowerBoundSo = sollya_lib_constant(get_rn_value(rnLowerBoundSa))
641
    upperBoundSo = sollya_lib_constant(get_rn_value(rnUpperBoundSa))
642
    rangeSo = sollya_lib_range(lowerBoundSo, upperBoundSo)
643
    currentPrecSo = sollya_lib_get_prec(None)
644
    if maxPrecSa > currentSollyaPrecSa:
645
        sollya_lib_set_prec(currentPrecSo)
646
        sollya_lib_clear_obj(currentPrecSo)
647
        sollya_lib_clear_obj(newPrecSo)
648
    sollya_lib_clear_obj(lowerBoundSo)
649
    sollya_lib_clear_obj(upperBoundSo)
650
    return(rangeSo)
651

    
652
def pobyso_range_to_interval_so_sa(rangeSo, realIntervalField = None):
653
    if realIntervalField is None:
654
        precSa = pobyso_get_prec_of_range_so_sa(rangeSo)
655
        realIntervalField = RealIntervalField(precSa)
656
    intervalSa = \
657
        pobyso_get_interval_from_range_so_sa(rangeSo, realIntervalField) 
658
    return(intervalSa)
659

    
660
def pobyso_remez_canonical_sa_sa(func, \
661
                                 degree, \
662
                                 lowerBound, \
663
                                 upperBound, \
664
                                 weight = None, \
665
                                 quality = None):
666
    """
667
    All arguments are Sage/Python.
668
    The functions (func and weight) must be passed as expressions or strings.
669
    Otherwise the function fails. 
670
    The return value is a pointer is a Sage polynomial.
671
    """
672
    var('zorglub')    # Dummy variable name for type check only.
673
    polySo = pobyso_remez_canonical_sa_so(func, \
674
                                 degree, \
675
                                 lowerBound, \
676
                                 upperBound, \
677
                                 weight = None, \
678
                                 quality = None)
679
    if parent(func) == parent("string"):
680
        functionSa = eval(func)
681
    # Expression test.
682
    elif type(func) == type(zorglub):
683
        functionSa = func
684
    maxPrecision = 0
685
    if polySo is None:
686
        return(None)
687
    maxPrecision = pobyso_get_max_prec_of_exp_so_sa(polySo)
688
    RRRR = RealField(maxPrecision)
689
    polynomialRing = RRRR[functionSa.variables()[0]]
690
    expSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, RRRR)
691
    polySa = polynomial(expSa, polynomialRing)
692
    return(polySa)
693
    
694
def pobyso_remez_canonical(func, \
695
                           degree, \
696
                           lowerBound, \
697
                           upperBound, \
698
                           weight = "1", \
699
                           quality = None):
700
    """ Legacy function. See pobyso_remez_canonical_sa_so. """
701
    return(pobyso_remez_canonical_sa_so(func, \
702
                                        degree, \
703
                                        lowerBound, \
704
                                        upperBound, \
705
                                        weight, \
706
                                        quality))
707
def pobyso_remez_canonical_sa_so(func, \
708
                                 degree, \
709
                                 lowerBound, \
710
                                 upperBound, \
711
                                 weight = None, \
712
                                 quality = None):
713
    """
714
    All arguments are Sage/Python.
715
    The functions (func and weight) must be passed as expressions or strings.
716
    Otherwise the function fails. 
717
    The return value is a pointer to a Sollya function.
718
    """
719
    var('zorglub')    # Dummy variable name for type check only.
720
    currentVariableName = None
721
    # The func argument can be of different types (string, 
722
    # symbolic expression...)
723
    if parent(func) == parent("string"):
724
        functionSo = sollya_lib_parse_string(func)
725
    # Expression test.
726
    elif type(func) == type(zorglub):
727
        # Until we are able to translate Sage expressions into Sollya 
728
        # expressions : parse the string version.
729
        currentVariableName = func.variables()[0]
730
        sollya_lib_name_free_variable(str(currentVariableName))
731
        functionSo = sollya_lib_parse_string(func._assume_str())
732
    else:
733
        return(None)
734
    if weight is None:
735
        weightSo = pobyso_constant_1_sa_so()
736
    elif parent(weight) == parent("string"):
737
        weightSo = sollya_lib_parse_string(func)
738
    elif type(weight) == type(zorglub): 
739
        functionSo = sollya_lib_parse_string_sa_so(weight._assume_str())
740
    else:
741
        return(None)
742
    degreeSo = pobyso_constant_from_int(degree)
743
    rangeSo = pobyso_range_sa_so(lowerBound, upperBound)
744
    if not quality is None:
745
        qualitySo= pobyso_constant_sa_so(quality)
746
    else:
747
        qualitySo = None
748
    return(sollya_lib_remez(functionSo, \
749
                            degreeSo, \
750
                            rangeSo, \
751
                            weightSo, \
752
                            qualitySo, \
753
                            None))
754
    
755
def pobyso_remez_canonical_so_so(funcSo, \
756
                                 degreeSo, \
757
                                 rangeSo, \
758
                                 weightSo = pobyso_constant_1_sa_so(),\
759
                                 qualitySo = None):
760
    """
761
    All arguments are pointers to Sollya objects.
762
    The return value is a pointer to a Sollya function.
763
    """
764
    if not sollya_lib_obj_is_function(funcSo):
765
        return(None)
766
    return(sollya_lib_remez(funcSo, degreeSo, rangeSo, weightSo, qualitySo, None))
767
    
768
def pobyso_set_canonical_off():
769
    sollya_lib_set_canonical(sollya_lib_off())
770

    
771
def pobyso_set_canonical_on():
772
    sollya_lib_set_canonical(sollya_lib_on())
773

    
774
def pobyso_set_prec(p):
775
    """ Legacy function. See pobyso_set_prec_sa_so. """
776
    return( pobyso_set_prec_sa_so(p))
777

    
778
def pobyso_set_prec_sa_so(p):
779
    a = c_int(p)
780
    precSo = c_void_p(sollya_lib_constant_from_int(a))
781
    sollya_lib_set_prec(precSo)
782

    
783
def pobyso_supnorm_so_so(polySo, funcSo, intervalSo, errorTypeSo, accuracySo):
784
    return(sollya_lib_supnorm(polySo, funcSo, intervalSo, errorTypeSo, \
785
                              accuracySo))
786

    
787
def pobyso_taylor_expansion_with_change_var_so_so(functionSo, degreeSo, rangeSo, \
788
                                                errorTypeSo, \
789
                                                sollyaPrecSo=None):
790
    """
791
    Compute the Taylor expansion with the variable change
792
    x -> (x-intervalCenter) included.
793
    """
794
    # No global change of the working precision.
795
    if not sollyaPrecSo is None:
796
        initialPrecSo = sollya_lib_get_prec(None)
797
        sollya_lib_set_prec(sollyaPrecSo)
798
    #
799
    intervalCenterSo = sollya_lib_mid(rangeSo)
800
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, \
801
                                         intervalCenterSo, \
802
                                         rangeSo, errorTypeSo, None)
803
    (taylorFormListSo, numElements, isEndElliptic) = \
804
        pobyso_get_list_elements_so_so(taylorFormSo)
805
    polySo = taylorFormListSo[0]
806
    errorRangeSo = taylorFormListSo[2]
807
    maxErrorSo = sollya_lib_sup(errorRangeSo)
808
    changeVarExpSo = sollya_lib_build_function_sub(\
809
                       sollya_lib_build_function_free_variable(),\
810
                       sollya_lib_copy_obj(intervalCenterSo))
811
    polyVarChangedSo = sollya_lib_evaluate(polySo, changeVarExpSo) 
812
    # If changed, reset the Sollya working precision.
813
    if not sollyaPrecSo is None:
814
        sollya_lib_set_prec(initialPrecSo)
815
        sollya_lib_clear_obj(initailPrecSo)
816
    return((polyVarChangedSo, intervalCenterSo, maxErrorSo))
817
# end pobyso_taylor_expansion_with_change_var_so_so
818

    
819
def pobyso_taylor_expansion_no_change_var_so_so(functionSo, degreeSo, rangeSo, \
820
                                                errorTypeSo, \
821
                                                sollyaPrecSo=None):
822
    """
823
    Compute the Taylor expansion without the variable change
824
    x -> x-intervalCenter.
825
    """
826
    # No global change of the working precision.
827
    if not sollyaPrecSo is None:
828
        initialPrecSo = sollya_lib_get_prec(None)
829
        sollya_lib_set_prec(sollyaPrecSo)
830
    #
831
    intervalCenterSo = sollya_lib_mid(rangeSo)
832
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, \
833
                                         intervalCenterSo, \
834
                                         rangeSo, errorTypeSo, None)
835
    (taylorFormListSo, numElements, isEndElliptic) = \
836
        pobyso_get_list_elements_so_so(taylorFormSo)
837
    polySo = taylorFormListSo[0]
838
    errorRangeSo = taylorFormListSo[2]
839
    maxErrorSo = sollya_lib_sup(errorRangeSo)
840
    # If changed, reset the Sollya working precision.
841
    if not sollyaPrecSo is None:
842
        sollya_lib_set_prec(initialPrecSo)
843
        sollya_lib_clear_obj(initialPrecSo)
844
    return((polySo, intervalCenterSo, maxErrorSo))
845
# end pobyso_taylor_expansion_no_change_var_so_so
846

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

    
851
def pobyso_taylor_so_so(functionSo, degreeSo, pointSo):
852
    return(sollya_lib_taylor(functionSo, degreeSo, pointSo))
853
    
854
def pobyso_taylorform(function, degree, point = None, interval = None, errorType=None):
855
    """ Legacy function. See ;"""
856
    
857
def pobyso_taylorform_sa_sa(functionSa, \
858
                            degree, \
859
                            point, \
860
                            precision, \
861
                            interval=None, \
862
                            errorType=None):
863
    """
864
    Compute the Taylor form of 'degree' for 'functionSa' at 'point' 
865
    for 'interval' with 'errorType'. 
866
    point: must be a Real or a Real interval.
867
    return the Taylor form as an array
868
    TODO: take care of the interval and of point when it is an interval;
869
          when errorType is not None;
870
          take care of the other elements of the Taylor form (coefficients errors and
871
          delta.
872
    """
873
    # Absolute as the default error.
874
    if errorType is None:
875
        errorTypeSo = sollya_lib_absolute()
876
    else:
877
        #TODO: deal with the other case.
878
        pass
879
    varSa = functionSa.variables()[0]
880
    pointBaseRingString = str(point.base_ring())
881
    if not re.search('Real', pointBaseRingString):
882
        return None
883
    # Call Sollya but first "sollyafy" the arguments.
884
    pobyso_name_free_variable_sa_so(str(varSa))
885
    #pobyso_set_prec_sa_so(300)
886
    # Sollyafy the function.
887
    functionSo = pobyso_parse_string_sa_so(functionSa._assume_str())
888
    if sollya_lib_obj_is_error(functionSo):
889
        print "pobyso_tailorform: function string can't be parsed!"
890
        return None
891
    # Sollyafy the degree
892
    degreeSo = sollya_lib_constant_from_int(int(degree))
893
    # Sollyafy the point
894
    if not re.search('Interval', pointBaseRingString):
895
        pointSo  = pobyso_constant_sa_so(point)
896
    else:
897
        # TODO: deal with the interval case.
898
        pass
899
    # Call Sollya
900
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, errorTypeSo,\
901
                                         None)
902
    (tfsAsList, numElements, isEndElliptic) = \
903
            pobyso_get_list_elements_so_so(taylorFormSo)
904
    polySo = tfsAsList[0]
905
    maxPrecision = pobyso_get_max_prec_of_exp_so_sa(polySo)
906
    polyRealField = RealField(maxPrecision)
907
    expSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, polyRealField)
908
    sollya_lib_close()
909
    polynomialRing = polyRealField[str(varSa)]
910
    polySa = polynomial(expSa, polynomialRing)
911
    taylorFormSa = [polySa]
912
    return(taylorFormSa)
913
# End pobyso_taylor_form_sa_sa
914

    
915
def pobyso_taylorform_so_so(functionSo, degreeSo, pointSo, intervalSo=None, \
916
                            errorTypeSo=None):
917
    createdErrorType = False
918
    if errorTypeSo is None:
919
        errorTypeSo = sollya_lib_absolute()
920
        createdErrorType = True
921
    else:
922
        #TODO: deal with the other case.
923
        pass
924
    if intervalSo is None:
925
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
926
                                         errorTypeSo, None)
927
    else:
928
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
929
                                         intervalSo, errorTypeSo, None)
930
    if createdErrorType:
931
        sollya_lib_clear_obj(errorTypeSo)
932
    return(resultSo)
933
        
934

    
935
def pobyso_univar_polynomial_print_reverse(polySa):
936
    """ Legacy function. See pobyso_univar_polynomial_print_reverse_sa_sa. """
937
    return(pobyso_univar_polynomial_print_reverse_sa_sa(polySa))
938

    
939
def pobyso_univar_polynomial_print_reverse_sa_sa(polySa):
940
    """
941
    Return the string representation of a univariate polynomial with
942
    monomials ordered in the x^0..x^n order of the monomials.
943
    Remember: Sage
944
    """
945
    polynomialRing = polySa.base_ring()
946
    # A very expensive solution:
947
    # -create a fake multivariate polynomial field with only one variable,
948
    #   specifying a negative lexicographical order;
949
    mpolynomialRing = PolynomialRing(polynomialRing.base(), \
950
                                     polynomialRing.variable_name(), \
951
                                     1, order='neglex')
952
    # - convert the univariate argument polynomial into a multivariate
953
    #   version;
954
    p = mpolynomialRing(polySa)
955
    # - return the string representation of the converted form.
956
    # There is no simple str() method defined for p's class.
957
    return(p.__str__())
958
#
959
print pobyso_get_prec()  
960
pobyso_set_prec(165)
961
print pobyso_get_prec()  
962
a=100
963
print type(a)
964
id(a)
965
print "Max arity: ", pobyso_max_arity
966
print "Function tripleDouble (43) as a string: ", pobyso_function_type_as_string(43)
967
print "Function None (44) as a string: ", pobyso_function_type_as_string(44)
968
print "...Pobyso check done"