Statistiques
| Révision :

root / pobysoPythonSage / src / pobyso.py @ 154

Historique | Voir | Annoter | Télécharger (49,89 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
    # From Sage to Sollya bounds.
116
    lowerBoundSo = sollya_lib_constant(get_rn_value(rnLowerBoundSa),
117
                                       maxPrecSa)
118
    upperBoundSo = sollya_lib_constant(get_rn_value(rnUpperBoundSa),
119
                                       maxPrecSa)
120
    # From Sollya bounds to range.
121
    rangeSo = sollya_lib_range(lowerBoundSo, upperBoundSo)
122
    # Back to original precision.
123
    # Clean up
124
    sollya_lib_clear_obj(lowerBoundSo)
125
    sollya_lib_clear_obj(upperBoundSo)
126
    return rangeSo
127
# End pobyso_bounds_to_range_sa_so
128

    
129
def pobyso_build_function_sub_so_so(exp1So, exp2So):
130
    return(sollya_lib_build_function_sub(exp1So, exp2So))
131

    
132
def pobyso_change_var_in_function_so_so(funcSo, chvarExpSo):
133
    """
134
    Variable change in a function.
135
    """
136
    return(sollya_lib_evaluate(funcSo,chvarExpSo))
137
# End pobyso_change_var_in_function_so_so     
138

    
139
def pobyso_chebyshevform_so_so(functionSo, degreeSo, intervalSo):
140
    resultSo = sollya_lib_chebyshevform(functionSo, degreeSo, intervalSo)
141
    return(resultSo)
142
# End pobyso_chebyshevform_so_so.
143

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

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

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

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

    
252
def pobyso_constant_1():
253
    """
254
    Obvious.
255
    Legacy function. See pobyso_constant_so_so. 
256
    """
257
    return(pobyso_constant_1_sa_so())
258

    
259
def pobyso_constant_1_sa_so():
260
    """
261
    Obvious.
262
    """
263
    return(pobyso_constant_from_int_sa_so(1))
264

    
265
def pobyso_constant_from_int(anInt):
266
    """ Legacy function. See pobyso_constant_from_int_sa_so. """
267
    return(pobyso_constant_from_int_sa_so(anInt))
268

    
269
def pobyso_constant_from_int_sa_so(anInt):
270
    """
271
    Get a Sollya constant from a Sage int.
272
    """
273
    return(sollya_lib_constant_from_int(int(anInt)))
274

    
275
def pobyso_constant_from_int_so_sa(constSo):
276
    """
277
    Get a Sage int from a Sollya int constant.
278
    Usefull for precision or powers in polynomials.
279
    """
280
    constSa = c_int(0)
281
    sollya_lib_get_constant_as_int(byref(constSa), constSo)
282
    return(constSa.value)
283
# End pobyso_constant_from_int_so_sa
284

    
285
def pobyso_function_type_as_string(funcType):
286
    """ Legacy function. See pobyso_function_type_as_string_so_sa. """
287
    return(pobyso_function_type_as_string_so_sa(funcType))
288

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

    
385
def pobyso_get_constant(rnArgSa, constSo):
386
    """ Legacy function. See pobyso_get_constant_so_sa. """
387
    return(pobyso_get_constant_so_sa(rnArgSa, constSo))
388

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

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

    
437
def pobyso_get_free_variable_name():
438
    """ 
439
    Legacy function. See pobyso_get_free_variable_name_so_sa.
440
    """
441
    return(pobyso_get_free_variable_name_so_sa())
442

    
443
def pobyso_get_free_variable_name_so_sa():
444
    return(sollya_lib_get_free_variable_name())
445
    
446
def pobyso_get_function_arity(expressionSo):
447
    """ 
448
    Legacy function. See pobyso_get_function_arity_so_sa.
449
    """
450
    return(pobyso_get_function_arity_so_sa(expressionSo))
451

    
452
def pobyso_get_function_arity_so_sa(expressionSo):
453
    arity = c_int(0)
454
    sollya_lib_get_function_arity(byref(arity),expressionSo)
455
    return(int(arity.value))
456

    
457
def pobyso_get_head_function(expressionSo):
458
    """ 
459
    Legacy function. See pobyso_get_head_function_so_sa. 
460
    """
461
    return(pobyso_get_head_function_so_sa(expressionSo)) 
462

    
463
def pobyso_get_head_function_so_sa(expressionSo):
464
    functionType = c_int(0)
465
    sollya_lib_get_head_function(byref(functionType), expressionSo, None)
466
    return(int(functionType.value))
467

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

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

    
533
def pobyso_get_max_prec_of_exp(soExp):
534
    """ Legacy function. See pobyso_get_max_prec_of_exp_so_sa. """
535
    return(pobyso_get_max_prec_of_exp_so_sa(soExp))
536

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

    
574
def pobyso_get_min_prec_of_constant_so_sa(constExpSo):
575
    """
576
    Get the minimum precision necessary to represent the value of a Sollya
577
    constant.
578
    MPFR_MIN_PREC and powers of 2 are taken into account.
579
    We assume that constExpSo is a point
580
    """
581
    constExpAsRnSa = pobyso_get_constant_as_rn_so_sa(constExpSo)
582
    return(min_mpfr_size(get_rn_value(constExpAsRnSa)))
583

    
584
def pobyso_get_sage_exp_from_sollya_exp(sollyaExpSo, realField = RR):
585
    """ Legacy function. See pobyso_get_sage_exp_from_sollya_exp_so_sa. """
586
    return(pobyso_get_sage_exp_from_sollya_exp_so_sa(sollyaExpSo, \
587
                                                     realField = RR))
588

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

    
635
def pobyso_get_poly_sa_so(polySo, realFieldSa=None):
636
    """
637
    Create a Sollya polynomial from a Sage polynomial.
638
    """
639
    pass
640
# pobyso_get_poly_sa_so
641

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

    
662
    expressionSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, \
663
                                                             realFieldSa)
664
    #print "...Sollya expression after.",
665
    #pobyso_autoprint(polySo)
666
    polyVariableSa = expressionSa.variables()[0]
667
    polyRingSa = realFieldSa[str(polyVariableSa)]
668
    #print polyRingSa
669
    # Do not use the polynomial(expressionSa, ring=polyRingSa) form!
670
    polynomialSa = polyRingSa(expressionSa)
671
    return(polynomialSa)
672
# End pobyso_get_sage_poly_from_sollya_poly
673

    
674
def pobyso_get_subfunctions(expressionSo):
675
    """ Legacy function. See pobyso_get_subfunctions_so_sa. """
676
    return(pobyso_get_subfunctions_so_sa(expressionSo)) 
677

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

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

    
745
def pobyso_get_prec_of_constant(ctExpSo):
746
    """ Legacy function. See pobyso_get_prec_of_constant_so_sa. """
747
    return(pobyso_get_prec_of_constant_so_sa(ctExpSo))
748

    
749
def pobyso_get_prec_of_constant_so_sa(ctExpSo):
750
    prec = c_int(0)
751
    retc = sollya_lib_get_prec_of_constant(byref(prec), ctExpSo, None)
752
    if retc == 0:
753
        return(None)
754
    return(int(prec.value))
755

    
756
def pobyso_get_prec_of_range_so_sa(rangeSo):
757
    prec = c_int(0)
758
    retc = sollya_lib_get_prec_of_range(byref(prec), rangeSo, None)
759
    if retc == 0:
760
        return(None)
761
    return(int(prec.value))
762
# End pobyso_get_prec_of_range_so_sa()
763

    
764
def pobyso_guess_degree_sa_sa(functionSa, intervalSa, errorSa, weightSa=None, \
765
                              degreeBoundSa=None):
766
    functionAsStringSa = functionSa._assume_str()
767
    functionSo = pobyso_parse_string_sa_so(functionAsStringSa)
768
    rangeSo = pobyso_interval_to_range_sa_so(intervalSa)
769
    errorSo = pobyso_constant_sa_so(errorSa)
770
    if not weightSa is None:
771
        weightAsStringSa = weightSa._assume_str()
772
        weightSo = pobyso_parse_string_sa_so(weightAsStringSa)
773
    else:
774
        weightSo = None
775
    if not degreeBoundSa is None:
776
        degreeBoundSo = pobyso_constant_from_int_sa_so(degreeBoundSa)
777
    else:
778
        degreeBoundSo = None
779
    guessedDegreeSa = pobyso_guess_degree_so_sa(functionSo,
780
                                              rangeSo,
781
                                              errorSo,
782
                                              weightSo,
783
                                              degreeBoundSo)
784
    sollya_lib_clear_obj(functionSo)
785
    sollya_lib_clear_obj(rangeSo)
786
    sollya_lib_clear_obj(errorSo)
787
    if not weightSo is None:
788
        sollya_lib_clear_obj(weightSo)
789
    if not degreeBoundSo is None:
790
        sollya_lib_clear_obj(degreeBoundSo)
791
    return guessedDegreeSa
792
# End poyso_guess_degree_sa_sa
793

    
794
def pobyso_guess_degree_so_sa(functionSo, rangeSo, errorSo, weightSo=None, \
795
                              degreeBoundSo=None):
796
    """
797
    Thin wrapper around the guessdegree function.
798
    Nevertheless, some precision control stuff has been appended.
799
    """
800
    # Deal with Sollya internal precision issues: if it is too small, 
801
    # compared with the error, increases it to about twice -log2(error).
802
    errorSa = pobyso_get_constant_as_rn_with_rf_so_sa(errorSo)
803
    log2ErrorSa = errorSa.log2()
804
    if log2ErrorSa < 0:
805
        neededPrecisionSa = int(2 * int(-log2ErrorSa) / 64) * 64
806
    else:
807
        neededPrecisionSa = int(2 * int(log2ErrorSa) / 64) * 64
808
    #print "Needed precision:", neededPrecisionSa
809
    currentPrecSa = pobyso_get_prec_so_sa()
810
    if neededPrecisionSa > currentPrecSa:
811
        currentPrecSo = pobyso_get_prec_so()
812
        pobyso_set_prec_sa_so(neededPrecisionSa)
813
    # weightSo and degreeBoundsSo are optional arguments.
814
    if weightSo is None:
815
        degreeRangeSo = sollya_lib_guessdegree(functionSo, rangeSo, errorSo, None)
816
    elif degreeBoundSo is None:
817
        degreeRangeSo =  sollya_lib_guessdegree(functionSo, rangeSo, \
818
                                                errorSo, weightSo, None)
819
    else:
820
        degreeRangeSo =  sollya_lib_guessdegree(functionSo, rangeSo, errorSo, \
821
                                                weightSo, degreeBoundSo, None)
822
    # Restore internal precision, if applicable.
823
    if neededPrecisionSa > currentPrecSa:
824
        pobyso_set_prec_so_so(currentPrecSo)
825
        sollya_lib_clear_obj(currentPrecSo)
826
    degreeIntervalSa = pobyso_range_to_interval_so_sa(degreeRangeSo)
827
    sollya_lib_clear_obj(degreeRangeSo)
828
    # When ok, both bounds match.
829
    # When the degree bound is too low, the upper bound is the degree
830
    # for which the error can be honored.
831
    # When it really goes wrong, the upper bound is infinity. 
832
    if degreeIntervalSa.lower() == degreeIntervalSa.upper():
833
        return int(degreeIntervalSa.lower())
834
    else:
835
        if degreeIntervalSa.upper().is_infinity():
836
            return None
837
        else:
838
            return int(degreeIntervalSa.upper())
839
    # End pobyso_guess_degree_so_sa
840

    
841
def pobyso_infnorm_so_so(func, interval, file = None, intervalList = None):
842
    print "Do not use this function. User pobyso_supnorm_so_so instead."
843
    return(None)
844

    
845
def pobyso_interval_to_range_sa_so(intervalSa, precisionSa=None):
846
    if precisionSa is None:
847
        precisionSa = intervalSa.parent().precision()
848
    intervalSo = pobyso_bounds_to_range_sa_so(intervalSa.lower(),\
849
                                              intervalSa.upper(),\
850
                                              precisionSa)
851
    return(intervalSo)
852
# End pobyso_interval_to_range_sa_so
853

    
854
def pobyso_lib_init():
855
    sollya_lib_init(None)
856

    
857
def pobyso_lib_close():
858
    sollya_lib_close(None)
859
    
860
def pobyso_name_free_variable(freeVariableNameSa):
861
    """ Legacy function. See pobyso_name_free_variable_sa_so. """
862
    pobyso_name_free_variable_sa_so(freeVariableNameSa)
863

    
864
def pobyso_name_free_variable_sa_so(freeVariableNameSa):
865
    """
866
    Set the free variable name in Sollya from a Sage string.
867
    """
868
    sollya_lib_name_free_variable(freeVariableNameSa)
869

    
870
def pobyso_parse_string(string):
871
    """ Legacy function. See pobyso_parse_string_sa_so. """
872
    return(pobyso_parse_string_sa_so(string))
873
 
874
def pobyso_parse_string_sa_so(string):
875
    """
876
    Get the Sollya expression computed from a Sage string.
877
    """
878
    return(sollya_lib_parse_string(string))
879

    
880
def pobyso_range(rnLowerBound, rnUpperBound):
881
    """ Legacy function. See pobyso_range_sa_so. """
882
    return(pobyso_range_sa_so(rnLowerBound, rnUpperBound)) 
883

    
884

    
885
def pobyso_range_to_interval_so_sa(rangeSo, realIntervalFieldSa = None):
886
    """
887
    Get a Sage interval from a Sollya range.
888
    If no realIntervalField is given as a parameter, the Sage interval
889
    precision is that of the Sollya range.
890
    Otherwise, the precision is that of the realIntervalField. In this case
891
    rounding may happen.
892
    """
893
    if realIntervalFieldSa is None:
894
        precSa = pobyso_get_prec_of_range_so_sa(rangeSo)
895
        realIntervalFieldSa = RealIntervalField(precSa)
896
    intervalSa = \
897
        pobyso_get_interval_from_range_so_sa(rangeSo, realIntervalFieldSa) 
898
    return(intervalSa)
899

    
900
def pobyso_remez_canonical_sa_sa(func, \
901
                                 degree, \
902
                                 lowerBound, \
903
                                 upperBound, \
904
                                 weight = None, \
905
                                 quality = None):
906
    """
907
    All arguments are Sage/Python.
908
    The functions (func and weight) must be passed as expressions or strings.
909
    Otherwise the function fails. 
910
    The return value is a Sage polynomial.
911
    """
912
    var('zorglub')    # Dummy variable name for type check only. Type of 
913
    # zorglub is "symbolic expression".
914
    polySo = pobyso_remez_canonical_sa_so(func, \
915
                                 degree, \
916
                                 lowerBound, \
917
                                 upperBound, \
918
                                 weight, \
919
                                 quality)
920
    # String test
921
    if parent(func) == parent("string"):
922
        functionSa = eval(func)
923
    # Expression test.
924
    elif type(func) == type(zorglub):
925
        functionSa = func
926
    else:
927
        return None
928
    #
929
    maxPrecision = 0
930
    if polySo is None:
931
        return(None)
932
    maxPrecision = pobyso_get_max_prec_of_exp_so_sa(polySo)
933
    RRRRSa = RealField(maxPrecision)
934
    polynomialRingSa = RRRRSa[functionSa.variables()[0]]
935
    expSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, RRRRSa)
936
    polySa = polynomial(expSa, polynomialRingSa)
937
    sollya_lib_clear_obj(polySo)
938
    return(polySa)
939
# End pobyso_remez_canonical_sa_sa
940
    
941
def pobyso_remez_canonical(func, \
942
                           degree, \
943
                           lowerBound, \
944
                           upperBound, \
945
                           weight = "1", \
946
                           quality = None):
947
    """ Legacy function. See pobyso_remez_canonical_sa_so. """
948
    return(pobyso_remez_canonical_sa_so(func, \
949
                                        degree, \
950
                                        lowerBound, \
951
                                        upperBound, \
952
                                        weight, \
953
                                        quality))
954
def pobyso_remez_canonical_sa_so(func, \
955
                                 degree, \
956
                                 lowerBound, \
957
                                 upperBound, \
958
                                 weight = None, \
959
                                 quality = None):
960
    """
961
    All arguments are Sage/Python.
962
    The functions (func and weight) must be passed as expressions or strings.
963
    Otherwise the function fails. 
964
    The return value is a pointer to a Sollya function.
965
    """
966
    var('zorglub')    # Dummy variable name for type check only. Type of
967
    # zorglub is "symbolic expression".
968
    currentVariableNameSa = None
969
    # The func argument can be of different types (string, 
970
    # symbolic expression...)
971
    if parent(func) == parent("string"):
972
        localFuncSa = eval(func)
973
        if len(localFuncSa.variables()) > 0:
974
            currentVariableNameSa = localFuncSa.variables()[0]
975
            sollya_lib_name_free_variable(str(currentVariableNameSa))
976
            functionSo = sollya_lib_parse_string(localFuncSa._assume_str())
977
    # Expression test.
978
    elif type(func) == type(zorglub):
979
        # Until we are able to translate Sage expressions into Sollya 
980
        # expressions : parse the string version.
981
        if len(func.variables()) > 0:
982
            currentVariableNameSa = func.variables()[0]
983
            sollya_lib_name_free_variable(str(currentVariableNameSa))
984
            functionSo = sollya_lib_parse_string(func._assume_str())
985
    else:
986
        return(None)
987
    if weight is None: # No weight given -> 1.
988
        weightSo = pobyso_constant_1_sa_so()
989
    elif parent(weight) == parent("string"): # Weight given as string: parse it.
990
        weightSo = sollya_lib_parse_string(func)
991
    elif type(weight) == type(zorglub): # Weight given as symbolice expression.
992
        functionSo = sollya_lib_parse_string_sa_so(weight._assume_str())
993
    else:
994
        return(None)
995
    degreeSo = pobyso_constant_from_int(degree)
996
    rangeSo = pobyso_bounds_to_range_sa_so(lowerBound, upperBound)
997
    if not quality is None:
998
        qualitySo= pobyso_constant_sa_so(quality)
999
    else:
1000
        qualitySo = None
1001
        
1002
    remezPolySo = sollya_lib_remez(functionSo, \
1003
                                   degreeSo, \
1004
                                   rangeSo, \
1005
                                   weightSo, \
1006
                                   qualitySo, \
1007
                                   None)
1008
    sollya_lib_clear_obj(functionSo)
1009
    sollya_lib_clear_obj(degreeSo)
1010
    sollya_lib_clear_obj(rangeSo)
1011
    sollya_lib_clear_obj(weightSo)
1012
    if not qualitySo is None:
1013
        sollya_lib_clear_obj(qualitySo)
1014
    return(remezPolySo)
1015
# End pobyso_remez_canonical_sa_so
1016

    
1017
def pobyso_remez_canonical_so_so(funcSo, \
1018
                                 degreeSo, \
1019
                                 rangeSo, \
1020
                                 weightSo = pobyso_constant_1_sa_so(),\
1021
                                 qualitySo = None):
1022
    """
1023
    All arguments are pointers to Sollya objects.
1024
    The return value is a pointer to a Sollya function.
1025
    """
1026
    if not sollya_lib_obj_is_function(funcSo):
1027
        return(None)
1028
    return(sollya_lib_remez(funcSo, degreeSo, rangeSo, weightSo, qualitySo, None))
1029
    
1030
def pobyso_set_canonical_off():
1031
    sollya_lib_set_canonical(sollya_lib_off())
1032

    
1033
def pobyso_set_canonical_on():
1034
    sollya_lib_set_canonical(sollya_lib_on())
1035

    
1036
def pobyso_set_prec(p):
1037
    """ Legacy function. See pobyso_set_prec_sa_so. """
1038
    pobyso_set_prec_sa_so(p)
1039

    
1040
def pobyso_set_prec_sa_so(p):
1041
    a = c_int(p)
1042
    precSo = c_void_p(sollya_lib_constant_from_int(a))
1043
    sollya_lib_set_prec(precSo, None)
1044

    
1045
def pobyso_set_prec_so_so(newPrecSo):
1046
    sollya_lib_set_prec(newPrecSo, None)
1047

    
1048
def pobyso_supnorm_so_so(polySo, funcSo, intervalSo, errorTypeSo = None,\
1049
                         accuracySo = None):
1050
    """
1051
    Computes the supnorm of the approximation error between the given 
1052
    polynomial and function.
1053
    errorTypeSo defaults to "absolute".
1054
    accuracySo defaults to 2^(-40).
1055
    """
1056
    if errorTypeSo is None:
1057
        errorTypeSo = sollya_lib_absolute(None)
1058
        errorTypeIsNone = True
1059
    else:
1060
        errorTypeIsNone = False
1061
    #
1062
    if accuracySo is None:
1063
        # Notice the **!
1064
        accuracySo = pobyso_constant_sa_so(RR(2**(-40)))
1065
        accuracyIsNone = True
1066
    else:
1067
        accuracyIsNone = False
1068
    pobyso_autoprint(accuracySo)
1069
    resultSo = \
1070
        sollya_lib_supnorm(polySo, funcSo, intervalSo, errorTypeSo, \
1071
                              accuracySo)
1072
    if errorTypeIsNone:
1073
        sollya_lib_clear_obj(errorTypeSo)
1074
    if accuracyIsNone:
1075
        sollya_lib_clear_obj(accuracySo)
1076
    return resultSo
1077
# End pobyso_supnorm_so_so
1078

    
1079
def pobyso_taylor_expansion_with_change_var_so_so(functionSo, degreeSo, \
1080
                                                  rangeSo, \
1081
                                                  errorTypeSo=None, \
1082
                                                  sollyaPrecSo=None):
1083
    """
1084
    Compute the Taylor expansion with the variable change
1085
    x -> (x-intervalCenter) included.
1086
    """
1087
    # No global change of the working precision.
1088
    if not sollyaPrecSo is None:
1089
        initialPrecSo = sollya_lib_get_prec(None)
1090
        sollya_lib_set_prec(sollyaPrecSo)
1091
    #
1092
    # Error type stuff: default to absolute.
1093
    if errorTypeSo is None:
1094
        errorTypeIsNone = True
1095
        errorTypeSo = sollya_lib_absolute(None)
1096
    else:
1097
        errorTypeIsNone = False
1098
    intervalCenterSo = sollya_lib_mid(rangeSo)
1099
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, \
1100
                                         intervalCenterSo, \
1101
                                         rangeSo, errorTypeSo, None)
1102
    # taylorFormListSaSo is a Python list of Sollya objects references that 
1103
    # are copies of the elements of taylorFormSo.
1104
    # pobyso_get_list_elements_so_so clears taylorFormSo.
1105
    (taylorFormListSo, numElements, isEndElliptic) = \
1106
        pobyso_get_list_elements_so_so(taylorFormSo)
1107
    polySo = taylorFormListSo[0]
1108
    errorRangeSo = taylorFormListSo[2]
1109
    maxErrorSo = sollya_lib_sup(errorRangeSo)
1110
    changeVarExpSo = sollya_lib_build_function_sub(\
1111
                       sollya_lib_build_function_free_variable(),\
1112
                       sollya_lib_copy_obj(intervalCenterSo))
1113
    polyVarChangedSo = sollya_lib_evaluate(polySo, changeVarExpSo) 
1114
    sollya_lib_clear_obj(changeVarExpSo)
1115
    # If changed, reset the Sollya working precision.
1116
    if not sollyaPrecSo is None:
1117
        sollya_lib_set_prec(initialPrecSo)
1118
        sollya_lib_clear_obj(initialPrecSo)
1119
    if errorTypeIsNone:
1120
        sollya_lib_clear_obj(errorTypeSo)
1121
    sollya_lib_clear_obj(taylorFormSo)
1122
    # Do not clear maxErrorSo.
1123
    return((polyVarChangedSo, intervalCenterSo, maxErrorSo))
1124
# end pobyso_taylor_expansion_with_change_var_so_so
1125

    
1126
def pobyso_taylor_expansion_no_change_var_so_so(functionSo, degreeSo, rangeSo,
1127
                                                errorTypeSo=None, 
1128
                                                sollyaPrecSo=None):
1129
    """
1130
    Compute the Taylor expansion without the variable change
1131
    x -> x-intervalCenter.
1132
    """
1133
    # No global change of the working precision.
1134
    if not sollyaPrecSo is None:
1135
        initialPrecSo = sollya_lib_get_prec(None)
1136
        sollya_lib_set_prec(sollyaPrecSo)
1137
    # Error type stuff: default to absolute.
1138
    if errorTypeSo is None:
1139
        errorTypeIsNone = True
1140
        errorTypeSo = sollya_lib_absolute(None)
1141
    else:
1142
        errorTypeIsNone = False
1143
    intervalCenterSo = sollya_lib_mid(rangeSo, None)
1144
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo,
1145
                                         intervalCenterSo,
1146
                                         rangeSo, errorTypeSo, None)
1147
    # taylorFormListSaSo is a Python list of Sollya objects references that 
1148
    # are copies of the elements of taylorFormSo.
1149
    # pobyso_get_list_elements_so_so clears taylorFormSo.
1150
    (taylorFormListSaSo, numElementsSa, isEndEllipticSa) = \
1151
        pobyso_get_list_elements_so_so(taylorFormSo)
1152
    polySo = sollya_lib_copy_obj(taylorFormListSaSo[0])
1153
    #print "Num elements:", numElementsSa
1154
    sollya_lib_clear_obj(taylorFormSo)
1155
    #polySo = taylorFormListSaSo[0]
1156
    #errorRangeSo = sollya_lib_copy_obj(taylorFormListSaSo[2])
1157
    errorRangeSo = taylorFormListSaSo[2]
1158
    # No copy_obj needed here: a new object is created.
1159
    maxErrorSo = sollya_lib_sup(errorRangeSo)
1160
    # If changed, reset the Sollya working precision.
1161
    if not sollyaPrecSo is None:
1162
        sollya_lib_set_prec(initialPrecSo)
1163
        sollya_lib_clear_obj(initialPrecSo)
1164
    if errorTypeIsNone:
1165
        sollya_lib_clear_obj(errorTypeSo)
1166
    pobyso_clear_taylorform_sa_so(taylorFormListSaSo)
1167
    return((polySo, intervalCenterSo, maxErrorSo))
1168
# end pobyso_taylor_expansion_no_change_var_so_so
1169

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

    
1174
def pobyso_taylor_so_so(functionSo, degreeSo, pointSo):
1175
    return(sollya_lib_taylor(functionSo, degreeSo, pointSo))
1176
    
1177
def pobyso_taylorform(function, degree, point = None, 
1178
                      interval = None, errorType=None):
1179
    """ Legacy function. See pobyso_taylorform_sa_sa;"""
1180
    
1181
def pobyso_taylorform_sa_sa(functionSa, \
1182
                            degreeSa, \
1183
                            pointSa, \
1184
                            intervalSa=None, \
1185
                            errorTypeSa=None, \
1186
                            precisionSa=None):
1187
    """
1188
    Compute the Taylor form of 'degreeSa' for 'functionSa' at 'pointSa' 
1189
    for 'intervalSa' with 'errorTypeSa' (a string) using 'precisionSa'. 
1190
    point: must be a Real or a Real interval.
1191
    return the Taylor form as an array
1192
    TODO: take care of the interval and of the point when it is an interval;
1193
          when errorType is not None;
1194
          take care of the other elements of the Taylor form (coefficients 
1195
          errors and delta.
1196
    """
1197
    # Absolute as the default error.
1198
    if errorTypeSa is None:
1199
        errorTypeSo = sollya_lib_absolute()
1200
    elif errorTypeSa == "relative":
1201
        errorTypeSo = sollya_lib_relative()
1202
    elif errortypeSa == "absolute":
1203
        errorTypeSo = sollya_lib_absolute()
1204
    else:
1205
        # No clean up needed.
1206
        return None
1207
    # Global precision stuff
1208
    precisionChangedSa = False
1209
    currentSollyaPrecSo = pobyso_get_prec_so()
1210
    currentSollyaPrecSa = pobyso_constant_from_int_so_sa(currentSollyaPrecSo)
1211
    if not precisionSa is None:
1212
        if precisionSa > currentSollyaPrecSa:
1213
            pobyso_set_prec_sa_so(precisionSa)
1214
            precisionChangedSa = True
1215
            
1216
    if len(functionSa.variables()) > 0:
1217
        varSa = functionSa.variables()[0]
1218
        pobyso_name_free_variable_sa_so(str(varSa))
1219
    # In any case (point or interval) the parent of pointSa has a precision
1220
    # method.
1221
    pointPrecSa = pointSa.parent().precision()
1222
    if precisionSa > pointPrecSa:
1223
        pointPrecSa = precisionSa
1224
    # In any case (point or interval) pointSa has a base_ring() method.
1225
    pointBaseRingString = str(pointSa.base_ring())
1226
    if re.search('Interval', pointBaseRingString) is None: # Point
1227
        pointSo = pobyso_constant_sa_so(pointSa, pointPrecSa)
1228
    else: # Interval.
1229
        pointSo = pobyso_interval_to_range_sa_so(pointSa, pointPrecSa)
1230
    # Sollyafy the function.
1231
    functionSo = pobyso_parse_string_sa_so(functionSa._assume_str())
1232
    if sollya_lib_obj_is_error(functionSo):
1233
        print "pobyso_tailorform: function string can't be parsed!"
1234
        return None
1235
    # Sollyafy the degree
1236
    degreeSo = sollya_lib_constant_from_int(int(degreeSa))
1237
    # Sollyafy the point
1238
    # Call Sollya
1239
    taylorFormSo = \
1240
        sollya_lib_taylorform(functionSo, degreeSo, pointSo, errorTypeSo,\
1241
                                         None)
1242
    sollya_lib_clear_obj(functionSo)
1243
    sollya_lib_clear_obj(degreeSo)
1244
    sollya_lib_clear_obj(pointSo)
1245
    sollya_lib_clear_obj(errorTypeSo)
1246
    (tfsAsList, numElements, isEndElliptic) = \
1247
            pobyso_get_list_elements_so_so(taylorFormSo)
1248
    polySo = tfsAsList[0]
1249
    maxPrecision = pobyso_get_max_prec_of_exp_so_sa(polySo)
1250
    polyRealField = RealField(maxPrecision)
1251
    expSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, polyRealField)
1252
    if precisionChangedSa:
1253
        sollya_lib_set_prec(currentSollyaPrecSo)
1254
        sollya_lib_clear_obj(currentSollyaPrecSo)
1255
    polynomialRing = polyRealField[str(varSa)]
1256
    polySa = polynomial(expSa, polynomialRing)
1257
    taylorFormSa = [polySa]
1258
    # Final clean-up.
1259
    sollya_lib_clear_obj(taylorFormSo)
1260
    return(taylorFormSa)
1261
# End pobyso_taylor_form_sa_sa
1262

    
1263
def pobyso_taylorform_so_so(functionSo, degreeSo, pointSo, intervalSo=None, \
1264
                            errorTypeSo=None):
1265
    createdErrorType = False
1266
    if errorTypeSo is None:
1267
        errorTypeSo = sollya_lib_absolute()
1268
        createdErrorType = True
1269
    else:
1270
        #TODO: deal with the other case.
1271
        pass
1272
    if intervalSo is None:
1273
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
1274
                                         errorTypeSo, None)
1275
    else:
1276
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
1277
                                         intervalSo, errorTypeSo, None)
1278
    if createdErrorType:
1279
        sollya_lib_clear_obj(errorTypeSo)
1280
    return(resultSo)
1281
        
1282

    
1283
def pobyso_univar_polynomial_print_reverse(polySa):
1284
    """ Legacy function. See pobyso_univar_polynomial_print_reverse_sa_sa. """
1285
    return(pobyso_univar_polynomial_print_reverse_sa_sa(polySa))
1286

    
1287
def pobyso_univar_polynomial_print_reverse_sa_sa(polySa):
1288
    """
1289
    Return the string representation of a univariate polynomial with
1290
    monomials ordered in the x^0..x^n order of the monomials.
1291
    Remember: Sage
1292
    """
1293
    polynomialRing = polySa.base_ring()
1294
    # A very expensive solution:
1295
    # -create a fake multivariate polynomial field with only one variable,
1296
    #   specifying a negative lexicographical order;
1297
    mpolynomialRing = PolynomialRing(polynomialRing.base(), \
1298
                                     polynomialRing.variable_name(), \
1299
                                     1, order='neglex')
1300
    # - convert the univariate argument polynomial into a multivariate
1301
    #   version;
1302
    p = mpolynomialRing(polySa)
1303
    # - return the string representation of the converted form.
1304
    # There is no simple str() method defined for p's class.
1305
    return(p.__str__())
1306
#
1307
print pobyso_get_prec()  
1308
pobyso_set_prec(165)
1309
print pobyso_get_prec()  
1310
a=100
1311
print type(a)
1312
id(a)
1313
print "Max arity: ", pobyso_max_arity
1314
print "Function tripleDouble (43) as a string: ", pobyso_function_type_as_string(43)
1315
print "Function None (44) as a string: ", pobyso_function_type_as_string(44)
1316
print "...Pobyso check done"