Statistiques
| Révision :

root / pobysoPythonSage / src / pobyso.py @ 153

Historique | Voir | Annoter | Télécharger (47,5 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_absolute_so_so():
89
    return(sollya_lib_absolute(None))
90

    
91
def pobyso_autoprint(arg):
92
    sollya_lib_autoprint(arg,None)
93

    
94
def pobyso_autoprint_so_so(arg):
95
    sollya_lib_autoprint(arg,None)
96
    
97
def pobyso_bounds_to_range_sa_so(rnLowerBoundSa, rnUpperBoundSa, \
98
                                 precisionSa=None):
99
    """
100
    Return a Sollya range from to 2 RealField Sage elements.
101
    The Sollya range element has a sufficient precision to hold all
102
    the digits of the Sage bounds.
103
    """
104
    # Sanity check.
105
    if rnLowerBoundSa > rnUpperBoundSa:
106
        return None
107
    # Precision stuff.
108
    if precisionSa is None:
109
        # Check for the largest precision.
110
        lbPrecSa = rnLowerBoundSa.parent().precision()
111
        ubPrecSa = rnLowerBoundSa.parent().precision()
112
        maxPrecSa = max(lbPrecSa, ubPrecSa)
113
    else:
114
        maxPrecSa = precisionSa
115
    sollyaCurrentPrecSo = pobyso_get_prec_so()
116
    sollyaCurrentPrecSa = pobyso_constant_from_int_so_sa(sollyaCurrentPrecSo)
117
    # Change the current Sollya precision only if necessary.
118
    if maxPrecSa > sollyaCurrentPrecSa:
119
        pobyso_set_prec_sa_so(maxPrecSa)
120
    # From Sage to Sollya bounds.
121
    lowerBoundSo = sollya_lib_constant(get_rn_value(rnLowerBoundSa))
122
    upperBoundSo = sollya_lib_constant(get_rn_value(rnUpperBoundSa))
123
    # From Sollya bounds to range.
124
    rangeSo = sollya_lib_range(lowerBoundSo, upperBoundSo)
125
    # Back to original precision.
126
    if maxPrecSa > sollyaCurrentPrecSa:
127
        sollya_lib_set_prec(sollyaCurrentPrecSo)
128
    # Clean up
129
    sollya_lib_clear_obj(sollyaCurrentPrecSo)
130
    sollya_lib_clear_obj(lowerBoundSo)
131
    sollya_lib_clear_obj(upperBoundSo)
132
    return(rangeSo)
133
# End pobyso_bounds_to_range_sa_so
134

    
135
def pobyso_build_function_sub_so_so(exp1So, exp2So):
136
    return(sollya_lib_build_function_sub(exp1So, exp2So))
137

    
138
def pobyso_change_var_in_function_so_so(funcSo, chvarExpSo):
139
    """
140
    Variable change in a function.
141
    """
142
    return(sollya_lib_evaluate(funcSo,chvarExpSo))
143
# End pobyso_change_var_in_function_so_so     
144

    
145
def pobyso_chebyshevform_so_so(functionSo, degreeSo, intervalSo):
146
    resultSo = sollya_lib_chebyshevform(functionSo, degreeSo, intervalSo)
147
    return(resultSo)
148
# End pobyso_chebyshevform_so_so.
149

    
150
def pobyso_clear_taylorform_sa_so(taylorFormSaSo):
151
    """
152
    This method is necessary to correctly clean up the memory from Taylor forms.
153
    These are made of a Sollya object, a Sollya object list, a Sollya object.
154
    For no clearly understood reason, sollya_lib_clear_object_list crashed 
155
    when applied to the object list.
156
    Here, we decompose it into Sage list of Sollya objects references and we
157
     clear them one by one. 
158
    """
159
    sollya_lib_clear_obj(taylorFormSaSo[0])
160
    (coefficientsErrorsListSaSo, numElementsSa, isEndEllipticSa) = \
161
        pobyso_get_list_elements_so_so(taylorFormSaSo[1])
162
    for element in coefficientsErrorsListSaSo:
163
        sollya_lib_clear_obj(element)
164
    sollya_lib_clear_obj(taylorFormSaSo[1])
165
    sollya_lib_clear_obj(taylorFormSaSo[2])
166
# End pobyso_clear_taylorform_sa_so 
167

    
168
def pobyso_cmp(rnArgSa, cteSo):
169
    """
170
    Compare the MPFR value a RealNumber with that of a Sollya constant.
171
    
172
    Get the value of the Sollya constant into a RealNumber and compare
173
    using MPFR. Could be optimized by working directly with a mpfr_t
174
    for the intermediate number. 
175
    """
176
    # Get the precision of the Sollya constant to build a Sage RealNumber 
177
    # with enough precision.to hold it.
178
    precisionOfCte = c_int(0)
179
    # From the Sollya constant, create a local Sage RealNumber.
180
    sollya_lib_get_prec_of_constant(precisionOfCte, cteSo) 
181
    #print "Precision of constant: ", precisionOfCte
182
    RRRR = RealField(precisionOfCte.value)
183
    rnLocalSa = RRRR(0)
184
    sollya_lib_get_constant(get_rn_value(rnLocalSa), cteSo)
185
    #
186
    ## Compare the Sage RealNumber version of the Sollya constant with rnArg.
187
    return(cmp_rn_value(rnArgSa, rnLocal))
188
# End pobyso_smp
189

    
190
def pobyso_compute_pos_function_abs_val_bounds_sa_sa(funcSa, lowerBoundSa, \
191
                                                     upperBoundSa):
192
    """
193
    TODO: completely rework and test.
194
    """
195
    pobyso = pobyso_name_free_variable_sa_so(funcSa.variables()[0])
196
    funcSo = pobyso_parse_string(funcSa._assume_str())
197
    rangeSo = pobyso_range_sa_so(lowerBoundSa, upperBoundSa)
198
    infnormSo = pobyso_infnorm_so_so(funcSo,rangeSo)
199
    # Sollya return the infnorm as an interval.
200
    fMaxSa = pobyso_get_interval_from_range_so_sa(infnormSo)
201
    # Get the top bound and compute the binade top limit.
202
    fMaxUpperBoundSa = fMaxSa.upper()
203
    binadeTopLimitSa = 2**ceil(fMaxUpperBoundSa.log2())
204
    # Put up together the function to use to compute the lower bound.
205
    funcAuxSo = pobyso_parse_string(str(binadeTopLimitSa) +  \
206
                                    '-(' + f._assume_str() + ')')
207
    pobyso_autoprint(funcAuxSo)
208
    # Clear the Sollya range before a new call to infnorm and issue the call.
209
    sollya_lib_clear_obj(infnormSo)
210
    infnormSo = pobyso_infnorm_so_so(funcAuxSo,rangeSo)
211
    fMinSa = pobyso_get_interval_from_range_so_sa(infnormSo)
212
    sollya_lib_clear_obj(infnormSo)
213
    fMinLowerBoundSa = binadeTopLimitSa - fMinSa.lower()
214
    # Compute the maximum of the precisions of the different bounds.
215
    maxPrecSa = max([fMinLowerBoundSa.parent().precision(), \
216
                     fMaxUpperBoundSa.parent().precision()])
217
    # Create a RealIntervalField and create an interval with the "good" bounds.
218
    RRRI = RealIntervalField(maxPrecSa)
219
    imageIntervalSa = RRRI(fMinLowerBoundSa, fMaxUpperBoundSa)
220
    # Free the unneeded Sollya objects
221
    sollya_lib_clear_obj(funcSo)
222
    sollya_lib_clear_obj(funcAuxSo)
223
    sollya_lib_clear_obj(rangeSo)
224
    return(imageIntervalSa)
225
# End pobyso_compute_pos_function_abs_val_bounds_sa_sa
226

    
227
def pobyso_constant(rnArg):
228
    """ Legacy function. See pobyso_constant_sa_so. """
229
    return(pobyso_constant_sa_so(rnArg))
230
    
231
def pobyso_constant_sa_so(rnArgSa, precisionSa=None):
232
    """
233
    Create a Sollya constant from a Sage RealNumber.
234
    """
235
    # Precision stuff
236
    if precisionSa is None:
237
        precisionSa = rnArgSa.parent().precision()
238
    currentSollyaPrecisionSo = sollya_lib_get_prec()
239
    currentSollyaPrecisionSa = \
240
        pobyso_constant_from_int(currentSollyaPrecisionSo)
241
    # Sollya constant creation takes place here.
242
    if precisionSa > currentSollyaPrecisionSa:
243
        pobyso_set_prec_sa_so(precisionSa)
244
        constantSo = sollya_lib_constant(get_rn_value(rnArgSa))
245
        pobyso_set_prec_sa_so(currentSollyaPrecision)
246
    else:
247
        constantSo = sollya_lib_constant(get_rn_value(rnArgSa))
248
    sollya_lib_clear_obj(currentSollyaPrecisionSo)
249
    return(constantSo)
250
# End pobyso_constant_sa_so
251
     
252
def pobyso_constant_0_sa_so():
253
    """
254
    Obvious.
255
    """
256
    return(pobyso_constant_from_int_sa_so(0))
257

    
258
def pobyso_constant_1():
259
    """
260
    Obvious.
261
    Legacy function. See pobyso_constant_so_so. 
262
    """
263
    return(pobyso_constant_1_sa_so())
264

    
265
def pobyso_constant_1_sa_so():
266
    """
267
    Obvious.
268
    """
269
    return(pobyso_constant_from_int_sa_so(1))
270

    
271
def pobyso_constant_from_int(anInt):
272
    """ Legacy function. See pobyso_constant_from_int_sa_so. """
273
    return(pobyso_constant_from_int_sa_so(anInt))
274

    
275
def pobyso_constant_from_int_sa_so(anInt):
276
    """
277
    Get a Sollya constant from a Sage int.
278
    """
279
    return(sollya_lib_constant_from_int(int(anInt)))
280

    
281
def pobyso_constant_from_int_so_sa(constSo):
282
    """
283
    Get a Sage int from a Sollya int constant.
284
    Usefull for precision or powers in polynomials.
285
    """
286
    constSa = c_int(0)
287
    sollya_lib_get_constant_as_int(byref(constSa), constSo)
288
    return(constSa.value)
289
# End pobyso_constant_from_int_so_sa
290

    
291
def pobyso_function_type_as_string(funcType):
292
    """ Legacy function. See pobyso_function_type_as_string_so_sa. """
293
    return(pobyso_function_type_as_string_so_sa(funcType))
294

    
295
def pobyso_function_type_as_string_so_sa(funcType):
296
    """
297
    Numeric Sollya function codes -> Sage mathematical function names.
298
    Notice that pow -> ^ (a la Sage, not a la Python).
299
    """
300
    if funcType == SOLLYA_BASE_FUNC_ABS:
301
        return "abs"
302
    elif funcType == SOLLYA_BASE_FUNC_ACOS:
303
        return "arccos"
304
    elif funcType == SOLLYA_BASE_FUNC_ACOSH:
305
        return "arccosh"
306
    elif funcType == SOLLYA_BASE_FUNC_ADD:
307
        return "+"
308
    elif funcType == SOLLYA_BASE_FUNC_ASIN:
309
        return "arcsin"
310
    elif funcType == SOLLYA_BASE_FUNC_ASINH:
311
        return "arcsinh"
312
    elif funcType == SOLLYA_BASE_FUNC_ATAN:
313
        return "arctan"
314
    elif funcType == SOLLYA_BASE_FUNC_ATANH:
315
        return "arctanh"
316
    elif funcType == SOLLYA_BASE_FUNC_CEIL:
317
        return "ceil"
318
    elif funcType == SOLLYA_BASE_FUNC_CONSTANT:
319
        return "cte"
320
    elif funcType == SOLLYA_BASE_FUNC_COS:
321
        return "cos"
322
    elif funcType == SOLLYA_BASE_FUNC_COSH:
323
        return "cosh"
324
    elif funcType == SOLLYA_BASE_FUNC_DIV:
325
        return "/"
326
    elif funcType == SOLLYA_BASE_FUNC_DOUBLE:
327
        return "double"
328
    elif funcType == SOLLYA_BASE_FUNC_DOUBLEDOUBLE:
329
        return "doubleDouble"
330
    elif funcType == SOLLYA_BASE_FUNC_DOUBLEEXTENDED:
331
        return "doubleDxtended"
332
    elif funcType == SOLLYA_BASE_FUNC_ERF:
333
        return "erf"
334
    elif funcType == SOLLYA_BASE_FUNC_ERFC:
335
        return "erfc"
336
    elif funcType == SOLLYA_BASE_FUNC_EXP:
337
        return "exp"
338
    elif funcType == SOLLYA_BASE_FUNC_EXP_M1:
339
        return "expm1"
340
    elif funcType == SOLLYA_BASE_FUNC_FLOOR:
341
        return "floor"
342
    elif funcType == SOLLYA_BASE_FUNC_FREE_VARIABLE:
343
        return "freeVariable"
344
    elif funcType == SOLLYA_BASE_FUNC_HALFPRECISION:
345
        return "halfPrecision"
346
    elif funcType == SOLLYA_BASE_FUNC_LIBRARYCONSTANT:
347
        return "libraryConstant"
348
    elif funcType == SOLLYA_BASE_FUNC_LIBRARYFUNCTION:
349
        return "libraryFunction"
350
    elif funcType == SOLLYA_BASE_FUNC_LOG:
351
        return "log"
352
    elif funcType == SOLLYA_BASE_FUNC_LOG_10:
353
        return "log10"
354
    elif funcType == SOLLYA_BASE_FUNC_LOG_1P:
355
        return "log1p"
356
    elif funcType == SOLLYA_BASE_FUNC_LOG_2:
357
        return "log2"
358
    elif funcType == SOLLYA_BASE_FUNC_MUL:
359
        return "*"
360
    elif funcType == SOLLYA_BASE_FUNC_NEARESTINT:
361
        return "round"
362
    elif funcType == SOLLYA_BASE_FUNC_NEG:
363
        return "__neg__"
364
    elif funcType == SOLLYA_BASE_FUNC_PI:
365
        return "pi"
366
    elif funcType == SOLLYA_BASE_FUNC_POW:
367
        return "^"
368
    elif funcType == SOLLYA_BASE_FUNC_PROCEDUREFUNCTION:
369
        return "procedureFunction"
370
    elif funcType == SOLLYA_BASE_FUNC_QUAD:
371
        return "quad"
372
    elif funcType == SOLLYA_BASE_FUNC_SIN:
373
        return "sin"
374
    elif funcType == SOLLYA_BASE_FUNC_SINGLE:
375
        return "single"
376
    elif funcType == SOLLYA_BASE_FUNC_SINH:
377
        return "sinh"
378
    elif funcType == SOLLYA_BASE_FUNC_SQRT:
379
        return "sqrt"
380
    elif funcType == SOLLYA_BASE_FUNC_SUB:
381
        return "-"
382
    elif funcType == SOLLYA_BASE_FUNC_TAN:
383
        return "tan"
384
    elif funcType == SOLLYA_BASE_FUNC_TANH:
385
        return "tanh"
386
    elif funcType == SOLLYA_BASE_FUNC_TRIPLEDOUBLE:
387
        return "tripleDouble"
388
    else:
389
        return None
390

    
391
def pobyso_get_constant(rnArgSa, constSo):
392
    """ Legacy function. See pobyso_get_constant_so_sa. """
393
    return(pobyso_get_constant_so_sa(rnArgSa, constSo))
394

    
395
def pobyso_get_constant_so_sa(rnArgSa, constSo):
396
    """
397
    Set the value of rnArgSo to the value of constSo in MPFR_RNDN mode.
398
    rnArg must already exist and belong to some RealField.
399
    We assume that constSo points to a Sollya constant.
400
    """
401
    return(sollya_lib_get_constant(get_rn_value(rnArgSa), constSo))
402
    
403
def pobyso_get_constant_as_rn(ctExpSo):
404
    """ 
405
    Legacy function. See pobyso_get_constant_as_rn_so_sa. 
406
    """ 
407
    return(pobyso_get_constant_as_rn_so_sa(ctExpSo))
408
    
409
def pobyso_get_constant_as_rn_so_sa(constExpSo):
410
    """
411
    Get a Sollya constant as a Sage "real number".
412
    The precision of the floating-point number returned is that of the Sollya
413
    constant.
414
    """
415
    precisionSa  = pobyso_get_prec_of_constant_so_sa(constExpSo) 
416
    RRRR = RealField(precisionSa)
417
    rnSa = RRRR(0)
418
    sollya_lib_get_constant(get_rn_value(rnSa), constExpSo)
419
    return(rnSa)
420
# End pobyso_get_constant_as_rn_so_sa
421

    
422
def pobyso_get_constant_as_rn_with_rf(ctExp, realField):
423
    """ 
424
    Legacy function. See pobyso_get_constant_as_rn_with_rf_so_sa.
425
    """
426
    return(pobyso_get_constant_as_rn_with_rf_so_sa(ctExp, realField))
427
    
428
def pobyso_get_constant_as_rn_with_rf_so_sa(ctExpSo, realFieldSa = None):
429
    """
430
    Get a Sollya constant as a Sage "real number".
431
    If no real field is specified, the precision of the floating-point number 
432
    returned is that of the Sollya constant.
433
    Otherwise is is that of the real field. Hence rounding may happen.
434
    """
435
    if realFieldSa is None:
436
        sollyaPrecSa = pobyso_get_prec_of_constant_so_sa(ctExpSo)
437
        realFieldSa = RealField(sollyaPrecSa)
438
    rnSa = realFieldSa(0)
439
    sollya_lib_get_constant(get_rn_value(rnSa), ctExpSo)
440
    return(rnSa)
441
# End pobyso_get_constant_as_rn_with_rf_so_sa
442

    
443
def pobyso_get_free_variable_name():
444
    """ 
445
    Legacy function. See pobyso_get_free_variable_name_so_sa.
446
    """
447
    return(pobyso_get_free_variable_name_so_sa())
448

    
449
def pobyso_get_free_variable_name_so_sa():
450
    return(sollya_lib_get_free_variable_name())
451
    
452
def pobyso_get_function_arity(expressionSo):
453
    """ 
454
    Legacy function. See pobyso_get_function_arity_so_sa.
455
    """
456
    return(pobyso_get_function_arity_so_sa(expressionSo))
457

    
458
def pobyso_get_function_arity_so_sa(expressionSo):
459
    arity = c_int(0)
460
    sollya_lib_get_function_arity(byref(arity),expressionSo)
461
    return(int(arity.value))
462

    
463
def pobyso_get_head_function(expressionSo):
464
    """ 
465
    Legacy function. See pobyso_get_head_function_so_sa. 
466
    """
467
    return(pobyso_get_head_function_so_sa(expressionSo)) 
468

    
469
def pobyso_get_head_function_so_sa(expressionSo):
470
    functionType = c_int(0)
471
    sollya_lib_get_head_function(byref(functionType), expressionSo, None)
472
    return(int(functionType.value))
473

    
474
def pobyso_get_interval_from_range_so_sa(soRange, realIntervalFieldSa = None ):
475
    """
476
    Return the Sage interval corresponding to the Sollya range argument.
477
    If no reaIntervalField is passed as an argument, the interval bounds are not
478
    rounded: they are elements of RealIntervalField of the "right" precision
479
    to hold all the digits.
480
    """
481
    prec = c_int(0)
482
    if realIntervalFieldSa is None:
483
        retval = sollya_lib_get_prec_of_range(byref(prec), soRange, None)
484
        if retval == 0:
485
            return(None)
486
        realIntervalFieldSa = RealIntervalField(prec.value)
487
    intervalSa = realIntervalFieldSa(0,0)
488
    retval = \
489
        sollya_lib_get_interval_from_range(get_interval_value(intervalSa),\
490
                                           soRange)
491
    if retval == 0:
492
        return(None)
493
    return(intervalSa)
494
# End pobyso_get_interval_from_range_so_sa
495

    
496
def pobyso_get_list_elements(soObj):
497
    """ Legacy function. See pobyso_get_list_elements_so_so. """
498
    return(pobyso_get_list_elements_so_so(soObj))
499
 
500
def pobyso_get_list_elements_so_so(objectListSo):
501
    """
502
    Get the Sollya list elements as a Sage/Python array of Sollya objects.
503
    
504
    INPUT:
505
    - objectListSo: a Sollya list of Sollya objects.
506
    
507
    OUTPUT:
508
    - a Sage/Python tuple made of:
509
      - a Sage/Python list of Sollya objects,
510
      - a Sage/Python int holding the number of elements,
511
      - a Sage/Python int stating (!= 0) that the list is end-elliptic.
512
    NOTE::
513
        We recover the addresses of the Sollya object from the list of pointers
514
        returned by sollya_lib_get_list_elements. The list itself is freed.
515
    TODO::
516
        Figure out what to do with numElements since the number of elements
517
        can easily be recovered from the list itself. 
518
        Ditto for isEndElliptic.
519
    """
520
    listAddress = POINTER(c_longlong)()
521
    numElements = c_int(0)
522
    isEndElliptic = c_int(0)
523
    listAsSageList = []
524
    result = sollya_lib_get_list_elements(byref(listAddress),\
525
                                          byref(numElements),\
526
                                          byref(isEndElliptic),\
527
                                          objectListSo)
528
    if result == 0 :
529
        return None
530
    for i in xrange(0, numElements.value, 1):
531
       #listAsSageList.append(sollya_lib_copy_obj(listAddress[i]))
532
       listAsSageList.append(listAddress[i])
533
       # Clear each of the elements returned by Sollya.
534
       #sollya_lib_clear_obj(listAddress[i])
535
    # Free the list itself.   
536
    sollya_lib_free(listAddress)
537
    return(listAsSageList, numElements.value, isEndElliptic.value)
538

    
539
def pobyso_get_max_prec_of_exp(soExp):
540
    """ Legacy function. See pobyso_get_max_prec_of_exp_so_sa. """
541
    return(pobyso_get_max_prec_of_exp_so_sa(soExp))
542

    
543
def pobyso_get_max_prec_of_exp_so_sa(expSo):
544
    """
545
    Get the maximum precision used for the numbers in a Sollya expression.
546
    
547
    Arguments:
548
    soExp -- a Sollya expression pointer
549
    Return value:
550
    A Python integer
551
    TODO: 
552
    - error management;
553
    - correctly deal with numerical type such as DOUBLEEXTENDED.
554
    """
555
    maxPrecision = 0
556
    minConstPrec = 0
557
    currentConstPrec = 0
558
    operator = pobyso_get_head_function_so_sa(expSo)
559
    if (operator != SOLLYA_BASE_FUNC_CONSTANT) and \
560
    (operator != SOLLYA_BASE_FUNC_FREE_VARIABLE):
561
        (arity, subexpressions) = pobyso_get_subfunctions_so_sa(expSo)
562
        for i in xrange(arity):
563
            maxPrecisionCandidate = \
564
                pobyso_get_max_prec_of_exp_so_sa(subexpressions[i])
565
            if maxPrecisionCandidate > maxPrecision:
566
                maxPrecision = maxPrecisionCandidate
567
        return(maxPrecision)
568
    elif operator == SOLLYA_BASE_FUNC_CONSTANT:
569
        #minConstPrec = pobyso_get_min_prec_of_constant_so_sa(expSo)
570
        #currentConstPrec = pobyso_get_min_prec_of_constant_so_sa(soExp)
571
        #print minConstPrec, " - ", currentConstPrec 
572
        return(pobyso_get_min_prec_of_constant_so_sa(expSo))
573
    
574
    elif operator == SOLLYA_BASE_FUNC_FREE_VARIABLE:
575
        return(0)
576
    else:
577
        print "pobyso_get_max_prec_of_exp_so_sa: unexepected operator."
578
        return(0)
579

    
580
def pobyso_get_min_prec_of_constant_so_sa(constExpSo):
581
    """
582
    Get the minimum precision necessary to represent the value of a Sollya
583
    constant.
584
    MPFR_MIN_PREC and powers of 2 are taken into account.
585
    We assume that constExpSo is a point
586
    """
587
    constExpAsRnSa = pobyso_get_constant_as_rn_so_sa(constExpSo)
588
    return(min_mpfr_size(get_rn_value(constExpAsRnSa)))
589

    
590
def pobyso_get_sage_exp_from_sollya_exp(sollyaExpSo, realField = RR):
591
    """ Legacy function. See pobyso_get_sage_exp_from_sollya_exp_so_sa. """
592
    return(pobyso_get_sage_exp_from_sollya_exp_so_sa(sollyaExpSo, \
593
                                                     realField = RR))
594

    
595
def pobyso_get_sage_exp_from_sollya_exp_so_sa(sollyaExpSo, realFieldSa = RR):
596
    """
597
    Get a Sage expression from a Sollya expression. 
598
    Currently only tested with polynomials with floating-point coefficients.
599
    Notice that, in the returned polynomial, the exponents are RealNumbers.
600
    """
601
    #pobyso_autoprint(sollyaExp)
602
    operatorSa = pobyso_get_head_function_so_sa(sollyaExpSo)
603
    sollyaLibFreeVariableName = sollya_lib_get_free_variable_name()
604
    # Constants and the free variable are special cases.
605
    # All other operator are dealt with in the same way.
606
    if (operatorSa != SOLLYA_BASE_FUNC_CONSTANT) and \
607
       (operatorSa != SOLLYA_BASE_FUNC_FREE_VARIABLE):
608
        (aritySa, subexpressionsSa) = pobyso_get_subfunctions_so_sa(sollyaExpSo)
609
        if aritySa == 1:
610
            sageExpSa = eval(pobyso_function_type_as_string_so_sa(operatorSa) + \
611
            "(" + pobyso_get_sage_exp_from_sollya_exp_so_sa(subexpressionsSa[0], \
612
            realFieldSa) + ")")
613
        elif aritySa == 2:
614
            # We do not get through the preprocessor.
615
            # The "^" operator is then a special case.
616
            if operatorSa == SOLLYA_BASE_FUNC_POW:
617
                operatorAsStringSa = "**"
618
            else:
619
                operatorAsStringSa = \
620
                    pobyso_function_type_as_string_so_sa(operatorSa)
621
            sageExpSa = \
622
              eval("pobyso_get_sage_exp_from_sollya_exp_so_sa(subexpressionsSa[0], realFieldSa)"\
623
              + " " + operatorAsStringSa + " " + \
624
                   "pobyso_get_sage_exp_from_sollya_exp_so_sa(subexpressionsSa[1], realFieldSa)")
625
        # We do not know yet how to deal with arity >= 3 
626
        # (is there any in Sollya anyway?).
627
        else:
628
            sageExpSa = eval('None')
629
        return(sageExpSa)
630
    elif operatorSa == SOLLYA_BASE_FUNC_CONSTANT:
631
        #print "This is a constant"
632
        return pobyso_get_constant_as_rn_with_rf_so_sa(sollyaExpSo, realFieldSa)
633
    elif operatorSa == SOLLYA_BASE_FUNC_FREE_VARIABLE:
634
        #print "This is free variable"
635
        return(eval(sollyaLibFreeVariableName))
636
    else:
637
        print "Unexpected"
638
        return eval('None')
639
# End pobyso_get_sage_poly_from_sollya_poly
640

    
641
def pobyso_get_poly_sa_so(polySo, realFieldSa=None):
642
    """
643
    Create a Sollya polynomial from a Sage polynomial.
644
    """
645
    pass
646
# pobyso_get_poly_sa_so
647

    
648
def pobyso_get_poly_so_sa(polySo, realFieldSa=None):
649
    """
650
    Convert a Sollya polynomial into a Sage polynomial.
651
    We assume that the polynomial is in canonical form.
652
    If no realField is given, a RealField corresponding to the maximum 
653
    precision of the coefficients is internally computed.
654
    The real field is not returned but can be easily retrieved from 
655
    the polynomial itself.
656
    ALGORITHM:
657
    - (optional) compute the RealField of the coefficients;
658
    - convert the Sollya expression into a Sage expression;
659
    - convert the Sage expression into a Sage polynomial
660
    TODO: the canonical thing for the polynomial.
661
    """    
662
    if realFieldSa is None:
663
        expressionPrecSa = pobyso_get_max_prec_of_exp_so_sa(polySo)
664
        realFieldSa = RealField(expressionPrecSa)
665
    #print "Sollya expression before...",
666
    #pobyso_autoprint(polySo)
667

    
668
    expressionSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, \
669
                                                             realFieldSa)
670
    #print "...Sollya expression after.",
671
    #pobyso_autoprint(polySo)
672
    polyVariableSa = expressionSa.variables()[0]
673
    polyRingSa = realFieldSa[str(polyVariableSa)]
674
    #print polyRingSa
675
    # Do not use the polynomial(expressionSa, ring=polyRingSa) form!
676
    polynomialSa = polyRingSa(expressionSa)
677
    return(polynomialSa)
678
# End pobyso_get_sage_poly_from_sollya_poly
679

    
680
def pobyso_get_subfunctions(expressionSo):
681
    """ Legacy function. See pobyso_get_subfunctions_so_sa. """
682
    return(pobyso_get_subfunctions_so_sa(expressionSo)) 
683

    
684
def pobyso_get_subfunctions_so_sa(expressionSo):
685
    """
686
    Get the subfunctions of an expression.
687
    Return the number of subfunctions and the list of subfunctions addresses.
688
    S.T.: Could not figure out another way than that ugly list of declarations
689
    to recover the addresses of the subfunctions.
690
    We limit ourselves to arity 8 functions. 
691
    """
692
    subf0 = c_int(0)
693
    subf1 = c_int(0)
694
    subf2 = c_int(0)
695
    subf3 = c_int(0)
696
    subf4 = c_int(0)
697
    subf5 = c_int(0)
698
    subf6 = c_int(0)
699
    subf7 = c_int(0)
700
    subf8 = c_int(0)
701
    arity = c_int(0)
702
    nullPtr = POINTER(c_int)()
703
    sollya_lib_get_subfunctions(expressionSo, byref(arity), \
704
      byref(subf0), byref(subf1), byref(subf2), byref(subf3), \
705
      byref(subf4), byref(subf5),\
706
      byref(subf6), byref(subf7), byref(subf8), nullPtr, None) 
707
#    byref(cast(subfunctions[0], POINTER(c_int))), \
708
#    byref(cast(subfunctions[0], POINTER(c_int))), \
709
#    byref(cast(subfunctions[2], POINTER(c_int))), \
710
#    byref(cast(subfunctions[3], POINTER(c_int))), \
711
#    byref(cast(subfunctions[4], POINTER(c_int))), \
712
#    byref(cast(subfunctions[5], POINTER(c_int))), \
713
#    byref(cast(subfunctions[6], POINTER(c_int))), \
714
#    byref(cast(subfunctions[7], POINTER(c_int))), \
715
#    byref(cast(subfunctions[8], POINTER(c_int))), nullPtr)
716
    subfunctions = [subf0, subf1, subf2, subf3, subf4, subf5, subf6, subf7, \
717
                    subf8]
718
    subs = []
719
    if arity.value > pobyso_max_arity:
720
        return(0,[])
721
    for i in xrange(arity.value):
722
        subs.append(int(subfunctions[i].value))
723
        #print subs[i]
724
    return(int(arity.value), subs)
725
    
726
def pobyso_get_prec():
727
    """ Legacy function. See pobyso_get_prec_so_sa(). """
728
    return(pobyso_get_prec_so_sa())
729

    
730
def pobyso_get_prec_so():
731
    """
732
    Get the current default precision in Sollya.
733
    The return value is a Sollya object.
734
    Usefull when modifying the precision back and forth by avoiding
735
    extra conversions.
736
    """
737
    return(sollya_lib_get_prec(None))
738
    
739
def pobyso_get_prec_so_sa():
740
    """
741
    Get the current default precision in Sollya.
742
    The return value is Sage/Python int.
743
    """
744
    precSo = sollya_lib_get_prec(None)
745
    precSa = c_int(0)
746
    sollya_lib_get_constant_as_int(byref(precSa), precSo)
747
    sollya_lib_clear_obj(precSo)
748
    return(int(precSa.value))
749
# End pobyso_get_prec_so_sa.
750

    
751
def pobyso_get_prec_of_constant(ctExpSo):
752
    """ Legacy function. See pobyso_get_prec_of_constant_so_sa. """
753
    return(pobyso_get_prec_of_constant_so_sa(ctExpSo))
754

    
755
def pobyso_get_prec_of_constant_so_sa(ctExpSo):
756
    prec = c_int(0)
757
    retc = sollya_lib_get_prec_of_constant(byref(prec), ctExpSo, None)
758
    if retc == 0:
759
        return(None)
760
    return(int(prec.value))
761

    
762
def pobyso_get_prec_of_range_so_sa(rangeSo):
763
    prec = c_int(0)
764
    retc = sollya_lib_get_prec_of_range(byref(prec), rangeSo, None)
765
    if retc == 0:
766
        return(None)
767
    return(int(prec.value))
768
# End pobyso_get_prec_of_range_so_sa()
769

    
770
def pobyso_guess_degree_so_sa(functionSo, rangeSo, errorSo, weightSo=None, \
771
                              degreeBoundsSo=None):
772
    # weightSo and degreeBoundsSo are optional arguments.
773
    if weightSo is None:
774
        degreeRangeSo = sollya_lib_guessdegree(functionSo, rangeSo, errorSo)
775
    elif degreeBounds is None:
776
        degreeRangeSo =  sollya_lib_guessdegree(functionSo, rangeSo, \
777
                                                errorSo, weightSo)
778
    else:
779
        degreeRangeSo =  sollya_lib_guessdegree(functionSo, rangeSo, errorSo, \
780
                                                weightSo, degreeBoundsSo)
781
    pobyso_range_to_interval_so_sa(degreeRangeSo)
782
# End pobyso_guess_degree_so_sa
783

    
784
def pobyso_infnorm_so_so(func, interval, file = None, intervalList = None):
785
    print "Do not use this function. User pobyso_supnorm_so_so instead."
786
    return(None)
787

    
788
def pobyso_interval_to_range_sa_so(intervalSa, precisionSa=None):
789
    if precisionSa is None:
790
        precisionSa = intervalSa.parent().precision()
791
    intervalSo = pobyso_bounds_to_range_sa_so(intervalSa.lower(),\
792
                                              intervalSa.upper(),\
793
                                              precisionSa)
794
    return(intervalSo)
795
# End pobyso_interval_to_range_sa_so
796

    
797
def pobyso_lib_init():
798
    sollya_lib_init(None)
799

    
800
def pobyso_lib_close():
801
    sollya_lib_close(None)
802
    
803
def pobyso_name_free_variable(freeVariableNameSa):
804
    """ Legacy function. See pobyso_name_free_variable_sa_so. """
805
    pobyso_name_free_variable_sa_so(freeVariableNameSa)
806

    
807
def pobyso_name_free_variable_sa_so(freeVariableNameSa):
808
    """
809
    Set the free variable name in Sollya from a Sage string.
810
    """
811
    sollya_lib_name_free_variable(freeVariableNameSa)
812

    
813
def pobyso_parse_string(string):
814
    """ Legacy function. See pobyso_parse_string_sa_so. """
815
    return(pobyso_parse_string_sa_so(string))
816
 
817
def pobyso_parse_string_sa_so(string):
818
    """
819
    Get the Sollya expression computed from a Sage string.
820
    """
821
    return(sollya_lib_parse_string(string))
822

    
823
def pobyso_range(rnLowerBound, rnUpperBound):
824
    """ Legacy function. See pobyso_range_sa_so. """
825
    return(pobyso_range_sa_so(rnLowerBound, rnUpperBound)) 
826

    
827

    
828
def pobyso_range_to_interval_so_sa(rangeSo, realIntervalFieldSa = None):
829
    """
830
    Get a Sage interval from a Sollya range.
831
    If no realIntervalField is given as a parameter, the Sage interval
832
    precision is that of the Sollya range.
833
    Otherwise, the precision is that of the realIntervalField. In this case
834
    rounding may happen.
835
    """
836
    if realIntervalFieldSa is None:
837
        precSa = pobyso_get_prec_of_range_so_sa(rangeSo)
838
        realIntervalFieldSa = RealIntervalField(precSa)
839
    intervalSa = \
840
        pobyso_get_interval_from_range_so_sa(rangeSo, realIntervalFieldSa) 
841
    return(intervalSa)
842

    
843
def pobyso_remez_canonical_sa_sa(func, \
844
                                 degree, \
845
                                 lowerBound, \
846
                                 upperBound, \
847
                                 weight = None, \
848
                                 quality = None):
849
    """
850
    All arguments are Sage/Python.
851
    The functions (func and weight) must be passed as expressions or strings.
852
    Otherwise the function fails. 
853
    The return value is a Sage polynomial.
854
    """
855
    var('zorglub')    # Dummy variable name for type check only. Type of 
856
    # zorglub is "symbolic expression".
857
    polySo = pobyso_remez_canonical_sa_so(func, \
858
                                 degree, \
859
                                 lowerBound, \
860
                                 upperBound, \
861
                                 weight, \
862
                                 quality)
863
    # String test
864
    if parent(func) == parent("string"):
865
        functionSa = eval(func)
866
    # Expression test.
867
    elif type(func) == type(zorglub):
868
        functionSa = func
869
    else:
870
        return None
871
    #
872
    maxPrecision = 0
873
    if polySo is None:
874
        return(None)
875
    maxPrecision = pobyso_get_max_prec_of_exp_so_sa(polySo)
876
    RRRRSa = RealField(maxPrecision)
877
    polynomialRingSa = RRRRSa[functionSa.variables()[0]]
878
    expSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, RRRRSa)
879
    polySa = polynomial(expSa, polynomialRingSa)
880
    sollya_lib_clear_obj(polySo)
881
    return(polySa)
882
# End pobyso_remez_canonical_sa_sa
883
    
884
def pobyso_remez_canonical(func, \
885
                           degree, \
886
                           lowerBound, \
887
                           upperBound, \
888
                           weight = "1", \
889
                           quality = None):
890
    """ Legacy function. See pobyso_remez_canonical_sa_so. """
891
    return(pobyso_remez_canonical_sa_so(func, \
892
                                        degree, \
893
                                        lowerBound, \
894
                                        upperBound, \
895
                                        weight, \
896
                                        quality))
897
def pobyso_remez_canonical_sa_so(func, \
898
                                 degree, \
899
                                 lowerBound, \
900
                                 upperBound, \
901
                                 weight = None, \
902
                                 quality = None):
903
    """
904
    All arguments are Sage/Python.
905
    The functions (func and weight) must be passed as expressions or strings.
906
    Otherwise the function fails. 
907
    The return value is a pointer to a Sollya function.
908
    """
909
    var('zorglub')    # Dummy variable name for type check only. Type of
910
    # zorglub is "symbolic expression".
911
    currentVariableNameSa = None
912
    # The func argument can be of different types (string, 
913
    # symbolic expression...)
914
    if parent(func) == parent("string"):
915
        localFuncSa = eval(func)
916
        if len(localFuncSa.variables()) > 0:
917
            currentVariableNameSa = localFuncSa.variables()[0]
918
            sollya_lib_name_free_variable(str(currentVariableNameSa))
919
            functionSo = sollya_lib_parse_string(localFuncSa._assume_str())
920
    # Expression test.
921
    elif type(func) == type(zorglub):
922
        # Until we are able to translate Sage expressions into Sollya 
923
        # expressions : parse the string version.
924
        if len(func.variables()) > 0:
925
            currentVariableNameSa = func.variables()[0]
926
            sollya_lib_name_free_variable(str(currentVariableNameSa))
927
            functionSo = sollya_lib_parse_string(func._assume_str())
928
    else:
929
        return(None)
930
    if weight is None: # No weight given -> 1.
931
        weightSo = pobyso_constant_1_sa_so()
932
    elif parent(weight) == parent("string"): # Weight given as string: parse it.
933
        weightSo = sollya_lib_parse_string(func)
934
    elif type(weight) == type(zorglub): # Weight given as symbolice expression.
935
        functionSo = sollya_lib_parse_string_sa_so(weight._assume_str())
936
    else:
937
        return(None)
938
    degreeSo = pobyso_constant_from_int(degree)
939
    rangeSo = pobyso_bounds_to_range_sa_so(lowerBound, upperBound)
940
    if not quality is None:
941
        qualitySo= pobyso_constant_sa_so(quality)
942
    else:
943
        qualitySo = None
944
        
945
    remezPolySo = sollya_lib_remez(functionSo, \
946
                                   degreeSo, \
947
                                   rangeSo, \
948
                                   weightSo, \
949
                                   qualitySo, \
950
                                   None)
951
    sollya_lib_clear_obj(functionSo)
952
    sollya_lib_clear_obj(degreeSo)
953
    sollya_lib_clear_obj(rangeSo)
954
    sollya_lib_clear_obj(weightSo)
955
    if not qualitySo is None:
956
        sollya_lib_clear_obj(qualitySo)
957
    return(remezPolySo)
958
# End pobyso_remez_canonical_sa_so
959

    
960
def pobyso_remez_canonical_so_so(funcSo, \
961
                                 degreeSo, \
962
                                 rangeSo, \
963
                                 weightSo = pobyso_constant_1_sa_so(),\
964
                                 qualitySo = None):
965
    """
966
    All arguments are pointers to Sollya objects.
967
    The return value is a pointer to a Sollya function.
968
    """
969
    if not sollya_lib_obj_is_function(funcSo):
970
        return(None)
971
    return(sollya_lib_remez(funcSo, degreeSo, rangeSo, weightSo, qualitySo, None))
972
    
973
def pobyso_set_canonical_off():
974
    sollya_lib_set_canonical(sollya_lib_off())
975

    
976
def pobyso_set_canonical_on():
977
    sollya_lib_set_canonical(sollya_lib_on())
978

    
979
def pobyso_set_prec(p):
980
    """ Legacy function. See pobyso_set_prec_sa_so. """
981
    pobyso_set_prec_sa_so(p)
982

    
983
def pobyso_set_prec_sa_so(p):
984
    a = c_int(p)
985
    precSo = c_void_p(sollya_lib_constant_from_int(a))
986
    sollya_lib_set_prec(precSo, None)
987

    
988
def pobyso_set_prec_so_so(newPrecSo):
989
    sollya_lib_set_prec(newPrecSo, None)
990

    
991
def pobyso_supnorm_so_so(polySo, funcSo, intervalSo, errorTypeSo = None,\
992
                         accuracySo = None):
993
    """
994
    Computes the supnorm of the approximation error between the given 
995
    polynomial and function.
996
    errorTypeSo defaults to "absolute".
997
    accuracySo defaults to 2^(-40).
998
    """
999
    if errorTypeSo is None:
1000
        errorTypeSo = sollya_lib_absolute(None)
1001
        errorTypeIsNone = True
1002
    else:
1003
        errorTypeIsNone = False
1004
    #
1005
    if accuracySo is None:
1006
        # Notice the **!
1007
        accuracySo = pobyso_constant_sa_so(RR(2**(-40)))
1008
        accuracyIsNone = True
1009
    else:
1010
        accuracyIsNone = False
1011
    pobyso_autoprint(accuracySo)
1012
    resultSo = \
1013
        sollya_lib_supnorm(polySo, funcSo, intervalSo, errorTypeSo, \
1014
                              accuracySo)
1015
    if errorTypeIsNone:
1016
        sollya_lib_clear_obj(errorTypeSo)
1017
    if accuracyIsNone:
1018
        sollya_lib_clear_obj(accuracySo)
1019
    return resultSo
1020
# End pobyso_supnorm_so_so
1021

    
1022
def pobyso_taylor_expansion_with_change_var_so_so(functionSo, degreeSo, \
1023
                                                  rangeSo, \
1024
                                                  errorTypeSo=None, \
1025
                                                  sollyaPrecSo=None):
1026
    """
1027
    Compute the Taylor expansion with the variable change
1028
    x -> (x-intervalCenter) included.
1029
    """
1030
    # No global change of the working precision.
1031
    if not sollyaPrecSo is None:
1032
        initialPrecSo = sollya_lib_get_prec(None)
1033
        sollya_lib_set_prec(sollyaPrecSo)
1034
    #
1035
    # Error type stuff: default to absolute.
1036
    if errorTypeSo is None:
1037
        errorTypeIsNone = True
1038
        errorTypeSo = sollya_lib_absolute(None)
1039
    else:
1040
        errorTypeIsNone = False
1041
    intervalCenterSo = sollya_lib_mid(rangeSo)
1042
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, \
1043
                                         intervalCenterSo, \
1044
                                         rangeSo, errorTypeSo, None)
1045
    # taylorFormListSaSo is a Python list of Sollya objects references that 
1046
    # are copies of the elements of taylorFormSo.
1047
    # pobyso_get_list_elements_so_so clears taylorFormSo.
1048
    (taylorFormListSo, numElements, isEndElliptic) = \
1049
        pobyso_get_list_elements_so_so(taylorFormSo)
1050
    polySo = taylorFormListSo[0]
1051
    errorRangeSo = taylorFormListSo[2]
1052
    maxErrorSo = sollya_lib_sup(errorRangeSo)
1053
    changeVarExpSo = sollya_lib_build_function_sub(\
1054
                       sollya_lib_build_function_free_variable(),\
1055
                       sollya_lib_copy_obj(intervalCenterSo))
1056
    polyVarChangedSo = sollya_lib_evaluate(polySo, changeVarExpSo) 
1057
    sollya_lib_clear_obj(changeVarExpSo)
1058
    # If changed, reset the Sollya working precision.
1059
    if not sollyaPrecSo is None:
1060
        sollya_lib_set_prec(initialPrecSo)
1061
        sollya_lib_clear_obj(initialPrecSo)
1062
    if errorTypeIsNone:
1063
        sollya_lib_clear_obj(errorTypeSo)
1064
    sollya_lib_clear_obj(taylorFormSo)
1065
    # Do not clear maxErrorSo.
1066
    return((polyVarChangedSo, intervalCenterSo, maxErrorSo))
1067
# end pobyso_taylor_expansion_with_change_var_so_so
1068

    
1069
def pobyso_taylor_expansion_no_change_var_so_so(functionSo, degreeSo, rangeSo,
1070
                                                errorTypeSo=None, 
1071
                                                sollyaPrecSo=None):
1072
    """
1073
    Compute the Taylor expansion without the variable change
1074
    x -> x-intervalCenter.
1075
    """
1076
    # No global change of the working precision.
1077
    if not sollyaPrecSo is None:
1078
        initialPrecSo = sollya_lib_get_prec(None)
1079
        sollya_lib_set_prec(sollyaPrecSo)
1080
    # Error type stuff: default to absolute.
1081
    if errorTypeSo is None:
1082
        errorTypeIsNone = True
1083
        errorTypeSo = sollya_lib_absolute(None)
1084
    else:
1085
        errorTypeIsNone = False
1086
    intervalCenterSo = sollya_lib_mid(rangeSo, None)
1087
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo,
1088
                                         intervalCenterSo,
1089
                                         rangeSo, errorTypeSo, None)
1090
    # taylorFormListSaSo is a Python list of Sollya objects references that 
1091
    # are copies of the elements of taylorFormSo.
1092
    # pobyso_get_list_elements_so_so clears taylorFormSo.
1093
    (taylorFormListSaSo, numElementsSa, isEndEllipticSa) = \
1094
        pobyso_get_list_elements_so_so(taylorFormSo)
1095
    polySo = sollya_lib_copy_obj(taylorFormListSaSo[0])
1096
    #print "Num elements:", numElementsSa
1097
    sollya_lib_clear_obj(taylorFormSo)
1098
    #polySo = taylorFormListSaSo[0]
1099
    #errorRangeSo = sollya_lib_copy_obj(taylorFormListSaSo[2])
1100
    errorRangeSo = taylorFormListSaSo[2]
1101
    # No copy_obj needed here: a new object is created.
1102
    maxErrorSo = sollya_lib_sup(errorRangeSo)
1103
    # If changed, reset the Sollya working precision.
1104
    if not sollyaPrecSo is None:
1105
        sollya_lib_set_prec(initialPrecSo)
1106
        sollya_lib_clear_obj(initialPrecSo)
1107
    if errorTypeIsNone:
1108
        sollya_lib_clear_obj(errorTypeSo)
1109
    pobyso_clear_taylorform_sa_so(taylorFormListSaSo)
1110
    return((polySo, intervalCenterSo, maxErrorSo))
1111
# end pobyso_taylor_expansion_no_change_var_so_so
1112

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

    
1117
def pobyso_taylor_so_so(functionSo, degreeSo, pointSo):
1118
    return(sollya_lib_taylor(functionSo, degreeSo, pointSo))
1119
    
1120
def pobyso_taylorform(function, degree, point = None, 
1121
                      interval = None, errorType=None):
1122
    """ Legacy function. See pobyso_taylorform_sa_sa;"""
1123
    
1124
def pobyso_taylorform_sa_sa(functionSa, \
1125
                            degreeSa, \
1126
                            pointSa, \
1127
                            intervalSa=None, \
1128
                            errorTypeSa=None, \
1129
                            precisionSa=None):
1130
    """
1131
    Compute the Taylor form of 'degreeSa' for 'functionSa' at 'pointSa' 
1132
    for 'intervalSa' with 'errorTypeSa' (a string) using 'precisionSa'. 
1133
    point: must be a Real or a Real interval.
1134
    return the Taylor form as an array
1135
    TODO: take care of the interval and of the point when it is an interval;
1136
          when errorType is not None;
1137
          take care of the other elements of the Taylor form (coefficients 
1138
          errors and delta.
1139
    """
1140
    # Absolute as the default error.
1141
    if errorTypeSa is None:
1142
        errorTypeSo = sollya_lib_absolute()
1143
    elif errorTypeSa == "relative":
1144
        errorTypeSo = sollya_lib_relative()
1145
    elif errortypeSa == "absolute":
1146
        errorTypeSo = sollya_lib_absolute()
1147
    else:
1148
        # No clean up needed.
1149
        return None
1150
    # Global precision stuff
1151
    precisionChangedSa = False
1152
    currentSollyaPrecSo = pobyso_get_prec_so()
1153
    currentSollyaPrecSa = pobyso_constant_from_int_so_sa(currentSollyaPrecSo)
1154
    if not precisionSa is None:
1155
        if precisionSa > currentSollyaPrecSa:
1156
            pobyso_set_prec_sa_so(precisionSa)
1157
            precisionChangedSa = True
1158
            
1159
    if len(functionSa.variables()) > 0:
1160
        varSa = functionSa.variables()[0]
1161
        pobyso_name_free_variable_sa_so(str(varSa))
1162
    # In any case (point or interval) the parent of pointSa has a precision
1163
    # method.
1164
    pointPrecSa = pointSa.parent().precision()
1165
    if precisionSa > pointPrecSa:
1166
        pointPrecSa = precisionSa
1167
    # In any case (point or interval) pointSa has a base_ring() method.
1168
    pointBaseRingString = str(pointSa.base_ring())
1169
    if re.search('Interval', pointBaseRingString) is None: # Point
1170
        pointSo = pobyso_constant_sa_so(pointSa, pointPrecSa)
1171
    else: # Interval.
1172
        pointSo = pobyso_interval_to_range_sa_so(pointSa, pointPrecSa)
1173
    # Sollyafy the function.
1174
    functionSo = pobyso_parse_string_sa_so(functionSa._assume_str())
1175
    if sollya_lib_obj_is_error(functionSo):
1176
        print "pobyso_tailorform: function string can't be parsed!"
1177
        return None
1178
    # Sollyafy the degree
1179
    degreeSo = sollya_lib_constant_from_int(int(degreeSa))
1180
    # Sollyafy the point
1181
    # Call Sollya
1182
    taylorFormSo = \
1183
        sollya_lib_taylorform(functionSo, degreeSo, pointSo, errorTypeSo,\
1184
                                         None)
1185
    sollya_lib_clear_obj(functionSo)
1186
    sollya_lib_clear_obj(degreeSo)
1187
    sollya_lib_clear_obj(pointSo)
1188
    sollya_lib_clear_obj(errorTypeSo)
1189
    (tfsAsList, numElements, isEndElliptic) = \
1190
            pobyso_get_list_elements_so_so(taylorFormSo)
1191
    polySo = tfsAsList[0]
1192
    maxPrecision = pobyso_get_max_prec_of_exp_so_sa(polySo)
1193
    polyRealField = RealField(maxPrecision)
1194
    expSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, polyRealField)
1195
    if precisionChangedSa:
1196
        sollya_lib_set_prec(currentSollyaPrecSo)
1197
        sollya_lib_clear_obj(currentSollyaPrecSo)
1198
    polynomialRing = polyRealField[str(varSa)]
1199
    polySa = polynomial(expSa, polynomialRing)
1200
    taylorFormSa = [polySa]
1201
    # Final clean-up.
1202
    sollya_lib_clear_obj(taylorFormSo)
1203
    return(taylorFormSa)
1204
# End pobyso_taylor_form_sa_sa
1205

    
1206
def pobyso_taylorform_so_so(functionSo, degreeSo, pointSo, intervalSo=None, \
1207
                            errorTypeSo=None):
1208
    createdErrorType = False
1209
    if errorTypeSo is None:
1210
        errorTypeSo = sollya_lib_absolute()
1211
        createdErrorType = True
1212
    else:
1213
        #TODO: deal with the other case.
1214
        pass
1215
    if intervalSo is None:
1216
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
1217
                                         errorTypeSo, None)
1218
    else:
1219
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
1220
                                         intervalSo, errorTypeSo, None)
1221
    if createdErrorType:
1222
        sollya_lib_clear_obj(errorTypeSo)
1223
    return(resultSo)
1224
        
1225

    
1226
def pobyso_univar_polynomial_print_reverse(polySa):
1227
    """ Legacy function. See pobyso_univar_polynomial_print_reverse_sa_sa. """
1228
    return(pobyso_univar_polynomial_print_reverse_sa_sa(polySa))
1229

    
1230
def pobyso_univar_polynomial_print_reverse_sa_sa(polySa):
1231
    """
1232
    Return the string representation of a univariate polynomial with
1233
    monomials ordered in the x^0..x^n order of the monomials.
1234
    Remember: Sage
1235
    """
1236
    polynomialRing = polySa.base_ring()
1237
    # A very expensive solution:
1238
    # -create a fake multivariate polynomial field with only one variable,
1239
    #   specifying a negative lexicographical order;
1240
    mpolynomialRing = PolynomialRing(polynomialRing.base(), \
1241
                                     polynomialRing.variable_name(), \
1242
                                     1, order='neglex')
1243
    # - convert the univariate argument polynomial into a multivariate
1244
    #   version;
1245
    p = mpolynomialRing(polySa)
1246
    # - return the string representation of the converted form.
1247
    # There is no simple str() method defined for p's class.
1248
    return(p.__str__())
1249
#
1250
print pobyso_get_prec()  
1251
pobyso_set_prec(165)
1252
print pobyso_get_prec()  
1253
a=100
1254
print type(a)
1255
id(a)
1256
print "Max arity: ", pobyso_max_arity
1257
print "Function tripleDouble (43) as a string: ", pobyso_function_type_as_string(43)
1258
print "Function None (44) as a string: ", pobyso_function_type_as_string(44)
1259
print "...Pobyso check done"