Statistiques
| Révision :

root / pobysoPythonSage / src / pobyso.py @ 58

Historique | Voir | Annoter | Télécharger (36,49 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
"""
35
Create the equivalent to an enum for the Sollya function types.
36
"""
37
(SOLLYA_BASE_FUNC_ABS,
38
SOLLYA_BASE_FUNC_ACOS,
39
    SOLLYA_BASE_FUNC_ACOSH,
40
    SOLLYA_BASE_FUNC_ADD,
41
    SOLLYA_BASE_FUNC_ASIN,
42
    SOLLYA_BASE_FUNC_ASINH,
43
    SOLLYA_BASE_FUNC_ATAN,
44
    SOLLYA_BASE_FUNC_ATANH,
45
    SOLLYA_BASE_FUNC_CEIL,
46
    SOLLYA_BASE_FUNC_CONSTANT,
47
    SOLLYA_BASE_FUNC_COS,
48
    SOLLYA_BASE_FUNC_COSH,
49
    SOLLYA_BASE_FUNC_DIV,
50
    SOLLYA_BASE_FUNC_DOUBLE,
51
    SOLLYA_BASE_FUNC_DOUBLEDOUBLE,
52
    SOLLYA_BASE_FUNC_DOUBLEEXTENDED,
53
    SOLLYA_BASE_FUNC_ERF,
54
    SOLLYA_BASE_FUNC_ERFC,
55
    SOLLYA_BASE_FUNC_EXP,
56
    SOLLYA_BASE_FUNC_EXP_M1,
57
    SOLLYA_BASE_FUNC_FLOOR,
58
    SOLLYA_BASE_FUNC_FREE_VARIABLE,
59
    SOLLYA_BASE_FUNC_HALFPRECISION,
60
    SOLLYA_BASE_FUNC_LIBRARYCONSTANT,
61
    SOLLYA_BASE_FUNC_LIBRARYFUNCTION,
62
    SOLLYA_BASE_FUNC_LOG,
63
    SOLLYA_BASE_FUNC_LOG_10,
64
    SOLLYA_BASE_FUNC_LOG_1P,
65
    SOLLYA_BASE_FUNC_LOG_2,
66
    SOLLYA_BASE_FUNC_MUL,
67
    SOLLYA_BASE_FUNC_NEARESTINT,
68
    SOLLYA_BASE_FUNC_NEG,
69
    SOLLYA_BASE_FUNC_PI,
70
    SOLLYA_BASE_FUNC_POW,
71
    SOLLYA_BASE_FUNC_PROCEDUREFUNCTION,
72
    SOLLYA_BASE_FUNC_QUAD,
73
    SOLLYA_BASE_FUNC_SIN,
74
    SOLLYA_BASE_FUNC_SINGLE,
75
    SOLLYA_BASE_FUNC_SINH,
76
    SOLLYA_BASE_FUNC_SQRT,
77
    SOLLYA_BASE_FUNC_SUB,
78
    SOLLYA_BASE_FUNC_TAN,
79
    SOLLYA_BASE_FUNC_TANH,
80
SOLLYA_BASE_FUNC_TRIPLEDOUBLE) = map(int,xrange(44))
81
print "\nSuperficial pobyso check..."    
82
print "First constant - SOLLYA_BASE_FUNC_ABS: ", SOLLYA_BASE_FUNC_ABS
83
print "Last constant  - SOLLYA_BASE_FUNC_TRIPLEDOUBLE: ", SOLLYA_BASE_FUNC_TRIPLEDOUBLE
84

    
85
pobyso_max_arity = 9
86

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

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

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

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

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

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

    
127

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

    
514
def pobyso_get_subfunctions(expressionSo):
515
    """ Legacy function. See pobyso_get_subfunctions_so_sa. """
516
    return(pobyso_get_subfunctions_so_sa(expressionSo)) 
517

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

    
557
def pobyso_get_prec_so_sa():
558
    """
559
    Get the current default precision in Sollya.
560
    The return value is Sage/Python int.
561
    """
562
    precSo = sollya_lib_get_prec(None)
563
    precSa = c_int(0)
564
    sollya_lib_get_constant_as_int(byref(precSa), precSo)
565
    sollya_lib_clear_obj(precSo)
566
    return(int(precSa.value))
567

    
568
def pobyso_get_prec_of_constant(ctExpSo):
569
    """ Legacy function. See pobyso_get_prec_of_constant_so_sa. """
570
    return(pobyso_get_prec_of_constant_so_sa(ctExpSo))
571

    
572
def pobyso_get_prec_of_constant_so_sa(ctExpSo):
573
    prec = c_int(0)
574
    retc = sollya_lib_get_prec_of_constant(byref(prec), ctExpSo, None)
575
    if retc == 0:
576
        return(None)
577
    return(int(prec.value))
578

    
579
def pobyso_get_prec_of_range_so_sa(rangeSo):
580
    prec = c_int(0)
581
    retc = sollya_lib_get_prec_of_range(byref(prec), rangeSo, None)
582
    if retc == 0:
583
        return(None)
584
    return(int(prec.value))
585

    
586
def pobyso_infnorm_so_so(func, interval, file = None, intervalList = None):
587
    print "Do not use this function. User pobyso_supnorm_so_so instead."
588
    return(None)
589

    
590
def pobyso_lib_init():
591
    sollya_lib_init(None)
592
    
593
def pobyso_name_free_variable(freeVariableName):
594
    """ Legacy function. See pobyso_name_free_variable_sa_so. """
595
    pobyso_name_free_variable_sa_so(freeVariableName)
596

    
597
def pobyso_name_free_variable_sa_so(freeVariableName):
598
    sollya_lib_name_free_variable(freeVariableName)
599

    
600
def pobyso_parse_string(string):
601
    """ Legacy function. See pobyso_parse_string_sa_so. """
602
    return(pobyso_parse_string_sa_so(string))
603
 
604
def pobyso_parse_string_sa_so(string):
605
    return(sollya_lib_parse_string(string))
606

    
607
def pobyso_range(rnLowerBound, rnUpperBound):
608
    """ Legacy function. See pobyso_range_sa_so. """
609
    return(pobyso_range_sa_so(rnLowerBound, rnUpperBound)) 
610

    
611
def pobyso_range_sa_so(rnLowerBoundSa, rnUpperBoundSa):
612
    """
613
    Return a Sollya range from to 2 RealField elements.
614
    The Sollya range element has a sufficient precision to hold all
615
    the digits of the bounds.
616
    """
617
    # TODO: check the bounds.
618
    lbPrec = rnLowerBoundSa.parent().precision()
619
    ubPrec = rnLowerBoundSa.parent().precision()
620
    currentSollyaPrecSa = pobyso_get_prec_so_sa()
621
    maxPrecSa = max(lbPrec, ubPrec, currentSollyaPrecSa)
622
    # Change the current Sollya precision only if necessary.
623
    if maxPrecSa > currentSollyaPrecSa:
624
        currentPrecSo = sollya_lib_get_prec(None)
625
        newPrecSo = solly_lib_constant_from_uint64(maxPrecSa)
626
        sollya_lib_set_prec(newPrecSo)
627
    lowerBoundSo = sollya_lib_constant(get_rn_value(rnLowerBoundSa))
628
    upperBoundSo = sollya_lib_constant(get_rn_value(rnUpperBoundSa))
629
    rangeSo = sollya_lib_range(lowerBoundSo, upperBoundSo)
630
    currentPrecSo = sollya_lib_get_prec(None)
631
    if maxPrecSa > currentSollyaPrecSa:
632
        sollya_lib_set_prec(currentPrecSo)
633
        sollya_lib_clear_obj(currentPrecSo)
634
        sollya_lib_clear_obj(newPrecSo)
635
    sollya_lib_clear_obj(lowerBoundSo)
636
    sollya_lib_clear_obj(upperBoundSo)
637
    return(rangeSo)
638

    
639
def pobyso_range_so_sa(rangeSo, realIntervalField = None):
640
    if realIntervalField is None:
641
        precSa = pobyso_get_prec_of_range_so_sa(rangeSo)
642
        realIntervalField = RealIntervalField(precSa)
643
    intervalSa = \
644
        pobyso_get_interval_from_range_so_sa(rangeSo, realIntervalField) 
645
    return(intervalSa)
646

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

    
758
def pobyso_set_canonical_on():
759
    sollya_lib_set_canonical(sollya_lib_on())
760

    
761
def pobyso_set_prec(p):
762
    """ Legacy function. See pobyso_set_prec_sa_so. """
763
    return( pobyso_set_prec_sa_so(p))
764

    
765
def pobyso_set_prec_sa_so(p):
766
    a = c_int(p)
767
    precSo = c_void_p(sollya_lib_constant_from_int(a))
768
    sollya_lib_set_prec(precSo)
769

    
770
def pobyso_supnorm_so_so(polySo, funcSo, intervalSo, errorTypeSo, accuracySo):
771
    return(sollya_lib_supnorm(polySo, funcSo, intervalSo, errorTypeSo, \
772
                              accuracySo))
773

    
774
def pobyso_taylor_expansion_with_change_var_so_so(functionSo, degreeSo, rangeSo, \
775
                                                errorTypeSo, \
776
                                                sollyaPrecSo=None):
777
    """
778
    Compute the Taylor expansion with the variable change
779
    x -> (x-intervalCenter) included.
780
    """
781
    # No global change of the working precision.
782
    if not sollyaPrecSo is None:
783
        initialPrecSo = sollya_lib_get_prec(None)
784
        sollya_lib_set_prec(sollyaPrecSo)
785
    #
786
    intervalCenterSo = sollya_lib_mid(rangeSo)
787
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, \
788
                                         intervalCenterSo, \
789
                                         rangeSo, errorTypeSo, None)
790
    (taylorFormListSo, numElements, isEndElliptic) = \
791
        pobyso_get_list_elements_so_so(taylorFormSo)
792
    polySo = taylorFormListSo[0]
793
    errorRangeSo = taylorFormListSo[2]
794
    maxErrorSo = sollya_lib_sup(errorRangeSo)
795
    changeVarExpSo = sollya_lib_build_function_sub(\
796
                       sollya_lib_build_function_free_variable(),\
797
                       sollya_lib_copy_obj(intervalCenterSo))
798
    polyVarChangedSo = sollya_lib_evaluate(polySo, changeVarExpSo) 
799
    # If changed, reset the Sollya working precision.
800
    if not sollyaPrecSo is None:
801
        sollya_lib_set_prec(initialPrecSo)
802
        sollya_lib_clear_obj(initailPrecSo)
803
    return([polyVarChangedSo, intervalCenterSo, maxErrorSo])
804
# end pobyso_taylor_expansion_with_change_var_so_so
805

    
806
def pobyso_taylor_expansion_no_change_var_so_so(functionSo, degreeSo, rangeSo, \
807
                                                errorTypeSo, \
808
                                                sollyaPrecSo=None):
809
    """
810
    Compute the Taylor expansion without the variable change
811
    x -> x-intervalCenter.
812
    """
813
    # No global change of the working precision.
814
    if not sollyaPrecSo is None:
815
        initialPrecSo = sollya_lib_get_prec(None)
816
        sollya_lib_set_prec(sollyaPrecSo)
817
    #
818
    intervalCenterSo = sollya_lib_mid(rangeSo)
819
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, \
820
                                         intervalCenterSo, \
821
                                         rangeSo, errorTypeSo, None)
822
    (taylorFormListSo, numElements, isEndElliptic) = \
823
        pobyso_get_list_elements_so_so(taylorFormSo)
824
    polySo = taylorFormListSo[0]
825
    errorRangeSo = taylorFormListSo[2]
826
    maxErrorSo = sollya_lib_sup(errorRangeSo)
827
    # If changed, reset the Sollya working precision.
828
    if not sollyaPrecSo is None:
829
        sollya_lib_set_prec(initialPrecSo)
830
        sollya_lib_clear_obj(initailPrecSo)
831
    return([polySo, intervalCenterSo, maxErrorSo])
832
# end pobyso_taylor_expansion_no_change_var_so_so
833

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

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

    
902
def pobyso_taylorform_so_so(functionSo, degreeSo, pointSo, intervalSo=None, \
903
                            errorTypeSo=None):
904
    createdErrorType = False
905
    if errorTypeSo is None:
906
        errorTypeSo = sollya_lib_absolute()
907
        createdErrorType = True
908
    else:
909
        #TODO: deal with the other case.
910
        pass
911
    if intervalSo is None:
912
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
913
                                         errorTypeSo, None)
914
    else:
915
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
916
                                         intervalSo, errorTypeSo, None)
917
    if createdErrorType:
918
        sollya_lib_clear_obj(errorTypeSo)
919
    return(resultSo)
920
        
921

    
922
def pobyso_univar_polynomial_print_reverse(polySa):
923
    """ Legacy function. See pobyso_univar_polynomial_print_reverse_sa_sa. """
924
    return(pobyso_univar_polynomial_print_reverse_sa_sa(polySa))
925

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