Statistiques
| Révision :

root / pobysoPythonSage / src / pobyso.py @ 160

Historique | Voir | Annoter | Télécharger (51,34 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().replace('_SAGE_VAR_', ''))
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().replace('_SAGE_VAR_', '') + ')')
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_error_so():
286
    return sollya_lib_error(None)
287
# End pobyso_error().
288

    
289
def pobyso_function_type_as_string(funcType):
290
    """ Legacy function. See pobyso_function_type_as_string_so_sa. """
291
    return(pobyso_function_type_as_string_so_sa(funcType))
292

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

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

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

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

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

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

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

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

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

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

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

    
537
def pobyso_get_max_prec_of_exp(soExp):
538
    """ Legacy function. See pobyso_get_max_prec_of_exp_so_sa. """
539
    return(pobyso_get_max_prec_of_exp_so_sa(soExp))
540

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

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

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

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

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

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

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

    
678
def pobyso_get_subfunctions(expressionSo):
679
    """ Legacy function. See pobyso_get_subfunctions_so_sa. """
680
    return(pobyso_get_subfunctions_so_sa(expressionSo)) 
681

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

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

    
749
def pobyso_get_prec_of_constant(ctExpSo):
750
    """ Legacy function. See pobyso_get_prec_of_constant_so_sa. """
751
    return(pobyso_get_prec_of_constant_so_sa(ctExpSo))
752

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

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

    
768
def pobyso_guess_degree_sa_sa(functionSa, intervalSa, approxErrorSa, 
769
                              weightSa=None, degreeBoundSa=None):
770
    """
771
    Sa_sa variant of the solly_guessdegree function.
772
    Return 0 if something goes wrong.
773
    """
774
    functionAsStringSa = functionSa._assume_str().replace('_SAGE_VAR_', '')
775
    functionSo = pobyso_parse_string_sa_so(functionAsStringSa)
776
    if pobyso_is_error_so_sa(functionSo):
777
        sollya_lib_clear_obj(functionSo)
778
        return 0
779
    rangeSo = pobyso_interval_to_range_sa_so(intervalSa)
780
    # The approximation error is expected to be a floating point number.
781
    if pobyso_is_floating_point_number_sa_sa(approxErrorSa):
782
        approxErrorSo = pobyso_constant_sa_so(approxErrorSa)
783
    else:
784
        approxErrorSo = pobyso_constant_sa_so(RR(approxErrorSa))
785
    if not weightSa is None:
786
        weightAsStringSa = weightSa._assume_str().replace('_SAGE_VAR_', '')
787
        weightSo = pobyso_parse_string_sa_so(weightAsStringSa)
788
        if pobyso_is_error(weightSo):
789
            sollya_lib_clear_obj(functionSo)
790
            sollya_lib_clear_obj(rangeSo)
791
            sollya_lib_clear_obj(approxErrorSo)   
792
            sollya_lib_clear_obj(weightSo)
793
            return 0   
794
    else:
795
        weightSo = None
796
    if not degreeBoundSa is None:
797
        degreeBoundSo = pobyso_constant_from_int_sa_so(degreeBoundSa)
798
    else:
799
        degreeBoundSo = None
800
    guessedDegreeSa = pobyso_guess_degree_so_sa(functionSo,
801
                                              rangeSo,
802
                                              approxErrorSo,
803
                                              weightSo,
804
                                              degreeBoundSo)
805
    sollya_lib_clear_obj(functionSo)
806
    sollya_lib_clear_obj(rangeSo)
807
    sollya_lib_clear_obj(approxErrorSo)
808
    if not weightSo is None:
809
        sollya_lib_clear_obj(weightSo)
810
    if not degreeBoundSo is None:
811
        sollya_lib_clear_obj(degreeBoundSo)
812
    return guessedDegreeSa
813
# End poyso_guess_degree_sa_sa
814

    
815
def pobyso_guess_degree_so_sa(functionSo, rangeSo, errorSo, weightSo=None, \
816
                              degreeBoundSo=None):
817
    """
818
    Thin wrapper around the guessdegree function.
819
    Nevertheless, some precision control stuff has been appended.
820
    """
821
    # Deal with Sollya internal precision issues: if it is too small, 
822
    # compared with the error, increases it to about twice -log2(error).
823
    errorSa = pobyso_get_constant_as_rn_with_rf_so_sa(errorSo)
824
    log2ErrorSa = errorSa.log2()
825
    if log2ErrorSa < 0:
826
        neededPrecisionSa = int(2 * int(-log2ErrorSa) / 64) * 64
827
    else:
828
        neededPrecisionSa = int(2 * int(log2ErrorSa) / 64) * 64
829
    #print "Needed precision:", neededPrecisionSa
830
    currentPrecSa = pobyso_get_prec_so_sa()
831
    if neededPrecisionSa > currentPrecSa:
832
        currentPrecSo = pobyso_get_prec_so()
833
        pobyso_set_prec_sa_so(neededPrecisionSa)
834
    # weightSo and degreeBoundsSo are optional arguments.
835
    if weightSo is None:
836
        degreeRangeSo = sollya_lib_guessdegree(functionSo, rangeSo, errorSo, None)
837
    elif degreeBoundSo is None:
838
        degreeRangeSo =  sollya_lib_guessdegree(functionSo, rangeSo, \
839
                                                errorSo, weightSo, None)
840
    else:
841
        degreeRangeSo =  sollya_lib_guessdegree(functionSo, rangeSo, errorSo, \
842
                                                weightSo, degreeBoundSo, None)
843
    # Restore internal precision, if applicable.
844
    if neededPrecisionSa > currentPrecSa:
845
        pobyso_set_prec_so_so(currentPrecSo)
846
        sollya_lib_clear_obj(currentPrecSo)
847
    degreeIntervalSa = pobyso_range_to_interval_so_sa(degreeRangeSo)
848
    sollya_lib_clear_obj(degreeRangeSo)
849
    # When ok, both bounds match.
850
    # When the degree bound is too low, the upper bound is the degree
851
    # for which the error can be honored.
852
    # When it really goes wrong, the upper bound is infinity. 
853
    if degreeIntervalSa.lower() == degreeIntervalSa.upper():
854
        return int(degreeIntervalSa.lower())
855
    else:
856
        if degreeIntervalSa.upper().is_infinity():
857
            return None
858
        else:
859
            return int(degreeIntervalSa.upper())
860
    # End pobyso_guess_degree_so_sa
861

    
862
def pobyso_infnorm_so_so(func, interval, file = None, intervalList = None):
863
    print "Do not use this function. User pobyso_supnorm_so_so instead."
864
    return(None)
865

    
866
def pobyso_interval_to_range_sa_so(intervalSa, precisionSa=None):
867
    if precisionSa is None:
868
        precisionSa = intervalSa.parent().precision()
869
    intervalSo = pobyso_bounds_to_range_sa_so(intervalSa.lower(),\
870
                                              intervalSa.upper(),\
871
                                              precisionSa)
872
    return(intervalSo)
873
# End pobyso_interval_to_range_sa_so
874

    
875
def pobyso_is_error_so_sa(objSo):
876
    """
877
    Thin wrapper around the sollya_lib_obj_is_error() function.
878
    """
879
    if sollya_lib_obj_is_error(objSo) != 0:
880
        return True
881
    else:
882
        return False
883
# End pobyso_is_error-so_sa
884

    
885
def pobyso_is_floating_point_number_sa_sa(numberSa):
886
    """
887
    Check whether a Sage number is floating point
888
    """
889
    return numberSa.parent().__class__ is RR.__class__
890

    
891
def pobyso_lib_init():
892
    sollya_lib_init(None)
893

    
894
def pobyso_lib_close():
895
    sollya_lib_close(None)
896
    
897
def pobyso_name_free_variable(freeVariableNameSa):
898
    """ Legacy function. See pobyso_name_free_variable_sa_so. """
899
    pobyso_name_free_variable_sa_so(freeVariableNameSa)
900

    
901
def pobyso_name_free_variable_sa_so(freeVariableNameSa):
902
    """
903
    Set the free variable name in Sollya from a Sage string.
904
    """
905
    sollya_lib_name_free_variable(freeVariableNameSa)
906

    
907
def pobyso_parse_string(string):
908
    """ Legacy function. See pobyso_parse_string_sa_so. """
909
    return(pobyso_parse_string_sa_so(string))
910
 
911
def pobyso_parse_string_sa_so(string):
912
    """
913
    Get the Sollya expression computed from a Sage string or
914
    a Sollya error object if parsing failed.
915
    """
916
    return(sollya_lib_parse_string(string))
917

    
918
def pobyso_range(rnLowerBound, rnUpperBound):
919
    """ Legacy function. See pobyso_range_sa_so. """
920
    return(pobyso_range_sa_so(rnLowerBound, rnUpperBound)) 
921

    
922

    
923
def pobyso_range_to_interval_so_sa(rangeSo, realIntervalFieldSa = None):
924
    """
925
    Get a Sage interval from a Sollya range.
926
    If no realIntervalField is given as a parameter, the Sage interval
927
    precision is that of the Sollya range.
928
    Otherwise, the precision is that of the realIntervalField. In this case
929
    rounding may happen.
930
    """
931
    if realIntervalFieldSa is None:
932
        precSa = pobyso_get_prec_of_range_so_sa(rangeSo)
933
        realIntervalFieldSa = RealIntervalField(precSa)
934
    intervalSa = \
935
        pobyso_get_interval_from_range_so_sa(rangeSo, realIntervalFieldSa) 
936
    return(intervalSa)
937

    
938
def pobyso_remez_canonical_sa_sa(func, \
939
                                 degree, \
940
                                 lowerBound, \
941
                                 upperBound, \
942
                                 weight = None, \
943
                                 quality = None):
944
    """
945
    All arguments are Sage/Python.
946
    The functions (func and weight) must be passed as expressions or strings.
947
    Otherwise the function fails. 
948
    The return value is a Sage polynomial.
949
    """
950
    var('zorglub')    # Dummy variable name for type check only. Type of 
951
    # zorglub is "symbolic expression".
952
    polySo = pobyso_remez_canonical_sa_so(func, \
953
                                 degree, \
954
                                 lowerBound, \
955
                                 upperBound, \
956
                                 weight, \
957
                                 quality)
958
    # String test
959
    if parent(func) == parent("string"):
960
        functionSa = eval(func)
961
    # Expression test.
962
    elif type(func) == type(zorglub):
963
        functionSa = func
964
    else:
965
        return None
966
    #
967
    maxPrecision = 0
968
    if polySo is None:
969
        return(None)
970
    maxPrecision = pobyso_get_max_prec_of_exp_so_sa(polySo)
971
    RRRRSa = RealField(maxPrecision)
972
    polynomialRingSa = RRRRSa[functionSa.variables()[0]]
973
    expSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, RRRRSa)
974
    polySa = polynomial(expSa, polynomialRingSa)
975
    sollya_lib_clear_obj(polySo)
976
    return(polySa)
977
# End pobyso_remez_canonical_sa_sa
978
    
979
def pobyso_remez_canonical(func, \
980
                           degree, \
981
                           lowerBound, \
982
                           upperBound, \
983
                           weight = "1", \
984
                           quality = None):
985
    """ Legacy function. See pobyso_remez_canonical_sa_so. """
986
    return(pobyso_remez_canonical_sa_so(func, \
987
                                        degree, \
988
                                        lowerBound, \
989
                                        upperBound, \
990
                                        weight, \
991
                                        quality))
992
def pobyso_remez_canonical_sa_so(func, \
993
                                 degree, \
994
                                 lowerBound, \
995
                                 upperBound, \
996
                                 weight = None, \
997
                                 quality = None):
998
    """
999
    All arguments are Sage/Python.
1000
    The functions (func and weight) must be passed as expressions or strings.
1001
    Otherwise the function fails. 
1002
    The return value is a pointer to a Sollya function.
1003
    """
1004
    var('zorglub')    # Dummy variable name for type check only. Type of
1005
    # zorglub is "symbolic expression".
1006
    currentVariableNameSa = None
1007
    # The func argument can be of different types (string, 
1008
    # symbolic expression...)
1009
    if parent(func) == parent("string"):
1010
        localFuncSa = eval(func)
1011
        if len(localFuncSa.variables()) > 0:
1012
            currentVariableNameSa = localFuncSa.variables()[0]
1013
            sollya_lib_name_free_variable(str(currentVariableNameSa))
1014
            functionSo = \
1015
              sollya_lib_parse_string(localFuncSa._assume_str().replace('_SAGE_VAR_', ''))
1016
    # Expression test.
1017
    elif type(func) == type(zorglub):
1018
        # Until we are able to translate Sage expressions into Sollya 
1019
        # expressions : parse the string version.
1020
        if len(func.variables()) > 0:
1021
            currentVariableNameSa = func.variables()[0]
1022
            sollya_lib_name_free_variable(str(currentVariableNameSa))
1023
            functionSo = \
1024
              sollya_lib_parse_string(func._assume_str().replace('_SAGE_VAR_', ''))
1025
    else:
1026
        return(None)
1027
    if weight is None: # No weight given -> 1.
1028
        weightSo = pobyso_constant_1_sa_so()
1029
    elif parent(weight) == parent("string"): # Weight given as string: parse it.
1030
        weightSo = sollya_lib_parse_string(func)
1031
    elif type(weight) == type(zorglub): # Weight given as symbolice expression.
1032
        functionSo = \
1033
          sollya_lib_parse_string_sa_so(weight._assume_str().replace('_SAGE_VAR_', ''))
1034
    else:
1035
        return(None)
1036
    degreeSo = pobyso_constant_from_int(degree)
1037
    rangeSo = pobyso_bounds_to_range_sa_so(lowerBound, upperBound)
1038
    if not quality is None:
1039
        qualitySo= pobyso_constant_sa_so(quality)
1040
    else:
1041
        qualitySo = None
1042
        
1043
    remezPolySo = sollya_lib_remez(functionSo, \
1044
                                   degreeSo, \
1045
                                   rangeSo, \
1046
                                   weightSo, \
1047
                                   qualitySo, \
1048
                                   None)
1049
    sollya_lib_clear_obj(functionSo)
1050
    sollya_lib_clear_obj(degreeSo)
1051
    sollya_lib_clear_obj(rangeSo)
1052
    sollya_lib_clear_obj(weightSo)
1053
    if not qualitySo is None:
1054
        sollya_lib_clear_obj(qualitySo)
1055
    return(remezPolySo)
1056
# End pobyso_remez_canonical_sa_so
1057

    
1058
def pobyso_remez_canonical_so_so(funcSo, \
1059
                                 degreeSo, \
1060
                                 rangeSo, \
1061
                                 weightSo = pobyso_constant_1_sa_so(),\
1062
                                 qualitySo = None):
1063
    """
1064
    All arguments are pointers to Sollya objects.
1065
    The return value is a pointer to a Sollya function.
1066
    """
1067
    if not sollya_lib_obj_is_function(funcSo):
1068
        return(None)
1069
    return(sollya_lib_remez(funcSo, degreeSo, rangeSo, weightSo, qualitySo, None))
1070
    
1071
def pobyso_set_canonical_off():
1072
    sollya_lib_set_canonical(sollya_lib_off())
1073

    
1074
def pobyso_set_canonical_on():
1075
    sollya_lib_set_canonical(sollya_lib_on())
1076

    
1077
def pobyso_set_prec(p):
1078
    """ Legacy function. See pobyso_set_prec_sa_so. """
1079
    pobyso_set_prec_sa_so(p)
1080

    
1081
def pobyso_set_prec_sa_so(p):
1082
    a = c_int(p)
1083
    precSo = c_void_p(sollya_lib_constant_from_int(a))
1084
    sollya_lib_set_prec(precSo, None)
1085

    
1086
def pobyso_set_prec_so_so(newPrecSo):
1087
    sollya_lib_set_prec(newPrecSo, None)
1088

    
1089
def pobyso_supnorm_so_so(polySo, funcSo, intervalSo, errorTypeSo = None,\
1090
                         accuracySo = None):
1091
    """
1092
    Computes the supnorm of the approximation error between the given 
1093
    polynomial and function.
1094
    errorTypeSo defaults to "absolute".
1095
    accuracySo defaults to 2^(-40).
1096
    """
1097
    if errorTypeSo is None:
1098
        errorTypeSo = sollya_lib_absolute(None)
1099
        errorTypeIsNone = True
1100
    else:
1101
        errorTypeIsNone = False
1102
    #
1103
    if accuracySo is None:
1104
        # Notice the **!
1105
        accuracySo = pobyso_constant_sa_so(RR(2**(-40)))
1106
        accuracyIsNone = True
1107
    else:
1108
        accuracyIsNone = False
1109
    pobyso_autoprint(accuracySo)
1110
    resultSo = \
1111
        sollya_lib_supnorm(polySo, funcSo, intervalSo, errorTypeSo, \
1112
                              accuracySo)
1113
    if errorTypeIsNone:
1114
        sollya_lib_clear_obj(errorTypeSo)
1115
    if accuracyIsNone:
1116
        sollya_lib_clear_obj(accuracySo)
1117
    return resultSo
1118
# End pobyso_supnorm_so_so
1119

    
1120
def pobyso_taylor_expansion_with_change_var_so_so(functionSo, degreeSo, \
1121
                                                  rangeSo, \
1122
                                                  errorTypeSo=None, \
1123
                                                  sollyaPrecSo=None):
1124
    """
1125
    Compute the Taylor expansion with the variable change
1126
    x -> (x-intervalCenter) included.
1127
    """
1128
    # No global change of the working precision.
1129
    if not sollyaPrecSo is None:
1130
        initialPrecSo = sollya_lib_get_prec(None)
1131
        sollya_lib_set_prec(sollyaPrecSo)
1132
    #
1133
    # Error type stuff: default to absolute.
1134
    if errorTypeSo is None:
1135
        errorTypeIsNone = True
1136
        errorTypeSo = sollya_lib_absolute(None)
1137
    else:
1138
        errorTypeIsNone = False
1139
    intervalCenterSo = sollya_lib_mid(rangeSo)
1140
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, \
1141
                                         intervalCenterSo, \
1142
                                         rangeSo, errorTypeSo, None)
1143
    # taylorFormListSaSo is a Python list of Sollya objects references that 
1144
    # are copies of the elements of taylorFormSo.
1145
    # pobyso_get_list_elements_so_so clears taylorFormSo.
1146
    (taylorFormListSo, numElements, isEndElliptic) = \
1147
        pobyso_get_list_elements_so_so(taylorFormSo)
1148
    polySo = taylorFormListSo[0]
1149
    errorRangeSo = taylorFormListSo[2]
1150
    maxErrorSo = sollya_lib_sup(errorRangeSo)
1151
    changeVarExpSo = sollya_lib_build_function_sub(\
1152
                       sollya_lib_build_function_free_variable(),\
1153
                       sollya_lib_copy_obj(intervalCenterSo))
1154
    polyVarChangedSo = sollya_lib_evaluate(polySo, changeVarExpSo) 
1155
    sollya_lib_clear_obj(changeVarExpSo)
1156
    # If changed, reset the Sollya working precision.
1157
    if not sollyaPrecSo is None:
1158
        sollya_lib_set_prec(initialPrecSo)
1159
        sollya_lib_clear_obj(initialPrecSo)
1160
    if errorTypeIsNone:
1161
        sollya_lib_clear_obj(errorTypeSo)
1162
    sollya_lib_clear_obj(taylorFormSo)
1163
    # Do not clear maxErrorSo.
1164
    return((polyVarChangedSo, intervalCenterSo, maxErrorSo))
1165
# end pobyso_taylor_expansion_with_change_var_so_so
1166

    
1167
def pobyso_taylor_expansion_no_change_var_so_so(functionSo, degreeSo, rangeSo,
1168
                                                errorTypeSo=None, 
1169
                                                sollyaPrecSo=None):
1170
    """
1171
    Compute the Taylor expansion without the variable change
1172
    x -> x-intervalCenter.
1173
    """
1174
    # No global change of the working precision.
1175
    if not sollyaPrecSo is None:
1176
        initialPrecSo = sollya_lib_get_prec(None)
1177
        sollya_lib_set_prec(sollyaPrecSo)
1178
    # Error type stuff: default to absolute.
1179
    if errorTypeSo is None:
1180
        errorTypeIsNone = True
1181
        errorTypeSo = sollya_lib_absolute(None)
1182
    else:
1183
        errorTypeIsNone = False
1184
    intervalCenterSo = sollya_lib_mid(rangeSo, None)
1185
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo,
1186
                                         intervalCenterSo,
1187
                                         rangeSo, errorTypeSo, None)
1188
    # taylorFormListSaSo is a Python list of Sollya objects references that 
1189
    # are copies of the elements of taylorFormSo.
1190
    # pobyso_get_list_elements_so_so clears taylorFormSo.
1191
    (taylorFormListSaSo, numElementsSa, isEndEllipticSa) = \
1192
        pobyso_get_list_elements_so_so(taylorFormSo)
1193
    polySo = sollya_lib_copy_obj(taylorFormListSaSo[0])
1194
    #print "Num elements:", numElementsSa
1195
    sollya_lib_clear_obj(taylorFormSo)
1196
    #polySo = taylorFormListSaSo[0]
1197
    #errorRangeSo = sollya_lib_copy_obj(taylorFormListSaSo[2])
1198
    errorRangeSo = taylorFormListSaSo[2]
1199
    # No copy_obj needed here: a new object is created.
1200
    maxErrorSo = sollya_lib_sup(errorRangeSo)
1201
    # If changed, reset the Sollya working precision.
1202
    if not sollyaPrecSo is None:
1203
        sollya_lib_set_prec(initialPrecSo)
1204
        sollya_lib_clear_obj(initialPrecSo)
1205
    if errorTypeIsNone:
1206
        sollya_lib_clear_obj(errorTypeSo)
1207
    pobyso_clear_taylorform_sa_so(taylorFormListSaSo)
1208
    return((polySo, intervalCenterSo, maxErrorSo))
1209
# end pobyso_taylor_expansion_no_change_var_so_so
1210

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

    
1215
def pobyso_taylor_so_so(functionSo, degreeSo, pointSo):
1216
    return(sollya_lib_taylor(functionSo, degreeSo, pointSo))
1217
    
1218
def pobyso_taylorform(function, degree, point = None, 
1219
                      interval = None, errorType=None):
1220
    """ Legacy function. See pobyso_taylorform_sa_sa;"""
1221
    
1222
def pobyso_taylorform_sa_sa(functionSa, \
1223
                            degreeSa, \
1224
                            pointSa, \
1225
                            intervalSa=None, \
1226
                            errorTypeSa=None, \
1227
                            precisionSa=None):
1228
    """
1229
    Compute the Taylor form of 'degreeSa' for 'functionSa' at 'pointSa' 
1230
    for 'intervalSa' with 'errorTypeSa' (a string) using 'precisionSa'. 
1231
    point: must be a Real or a Real interval.
1232
    return the Taylor form as an array
1233
    TODO: take care of the interval and of the point when it is an interval;
1234
          when errorType is not None;
1235
          take care of the other elements of the Taylor form (coefficients 
1236
          errors and delta.
1237
    """
1238
    # Absolute as the default error.
1239
    if errorTypeSa is None:
1240
        errorTypeSo = sollya_lib_absolute()
1241
    elif errorTypeSa == "relative":
1242
        errorTypeSo = sollya_lib_relative()
1243
    elif errortypeSa == "absolute":
1244
        errorTypeSo = sollya_lib_absolute()
1245
    else:
1246
        # No clean up needed.
1247
        return None
1248
    # Global precision stuff
1249
    precisionChangedSa = False
1250
    currentSollyaPrecSo = pobyso_get_prec_so()
1251
    currentSollyaPrecSa = pobyso_constant_from_int_so_sa(currentSollyaPrecSo)
1252
    if not precisionSa is None:
1253
        if precisionSa > currentSollyaPrecSa:
1254
            pobyso_set_prec_sa_so(precisionSa)
1255
            precisionChangedSa = True
1256
            
1257
    if len(functionSa.variables()) > 0:
1258
        varSa = functionSa.variables()[0]
1259
        pobyso_name_free_variable_sa_so(str(varSa))
1260
    # In any case (point or interval) the parent of pointSa has a precision
1261
    # method.
1262
    pointPrecSa = pointSa.parent().precision()
1263
    if precisionSa > pointPrecSa:
1264
        pointPrecSa = precisionSa
1265
    # In any case (point or interval) pointSa has a base_ring() method.
1266
    pointBaseRingString = str(pointSa.base_ring())
1267
    if re.search('Interval', pointBaseRingString) is None: # Point
1268
        pointSo = pobyso_constant_sa_so(pointSa, pointPrecSa)
1269
    else: # Interval.
1270
        pointSo = pobyso_interval_to_range_sa_so(pointSa, pointPrecSa)
1271
    # Sollyafy the function.
1272
    functionSo = pobyso_parse_string_sa_so(functionSa._assume_str().replace('_SAGE_VAR_', ''))
1273
    if sollya_lib_obj_is_error(functionSo):
1274
        print "pobyso_tailorform: function string can't be parsed!"
1275
        return None
1276
    # Sollyafy the degree
1277
    degreeSo = sollya_lib_constant_from_int(int(degreeSa))
1278
    # Sollyafy the point
1279
    # Call Sollya
1280
    taylorFormSo = \
1281
        sollya_lib_taylorform(functionSo, degreeSo, pointSo, errorTypeSo,\
1282
                                         None)
1283
    sollya_lib_clear_obj(functionSo)
1284
    sollya_lib_clear_obj(degreeSo)
1285
    sollya_lib_clear_obj(pointSo)
1286
    sollya_lib_clear_obj(errorTypeSo)
1287
    (tfsAsList, numElements, isEndElliptic) = \
1288
            pobyso_get_list_elements_so_so(taylorFormSo)
1289
    polySo = tfsAsList[0]
1290
    maxPrecision = pobyso_get_max_prec_of_exp_so_sa(polySo)
1291
    polyRealField = RealField(maxPrecision)
1292
    expSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, polyRealField)
1293
    if precisionChangedSa:
1294
        sollya_lib_set_prec(currentSollyaPrecSo)
1295
        sollya_lib_clear_obj(currentSollyaPrecSo)
1296
    polynomialRing = polyRealField[str(varSa)]
1297
    polySa = polynomial(expSa, polynomialRing)
1298
    taylorFormSa = [polySa]
1299
    # Final clean-up.
1300
    sollya_lib_clear_obj(taylorFormSo)
1301
    return(taylorFormSa)
1302
# End pobyso_taylor_form_sa_sa
1303

    
1304
def pobyso_taylorform_so_so(functionSo, degreeSo, pointSo, intervalSo=None, \
1305
                            errorTypeSo=None):
1306
    createdErrorType = False
1307
    if errorTypeSo is None:
1308
        errorTypeSo = sollya_lib_absolute()
1309
        createdErrorType = True
1310
    else:
1311
        #TODO: deal with the other case.
1312
        pass
1313
    if intervalSo is None:
1314
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
1315
                                         errorTypeSo, None)
1316
    else:
1317
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
1318
                                         intervalSo, errorTypeSo, None)
1319
    if createdErrorType:
1320
        sollya_lib_clear_obj(errorTypeSo)
1321
    return(resultSo)
1322
        
1323

    
1324
def pobyso_univar_polynomial_print_reverse(polySa):
1325
    """ Legacy function. See pobyso_univar_polynomial_print_reverse_sa_sa. """
1326
    return(pobyso_univar_polynomial_print_reverse_sa_sa(polySa))
1327

    
1328
def pobyso_univar_polynomial_print_reverse_sa_sa(polySa):
1329
    """
1330
    Return the string representation of a univariate polynomial with
1331
    monomials ordered in the x^0..x^n order of the monomials.
1332
    Remember: Sage
1333
    """
1334
    polynomialRing = polySa.base_ring()
1335
    # A very expensive solution:
1336
    # -create a fake multivariate polynomial field with only one variable,
1337
    #   specifying a negative lexicographical order;
1338
    mpolynomialRing = PolynomialRing(polynomialRing.base(), \
1339
                                     polynomialRing.variable_name(), \
1340
                                     1, order='neglex')
1341
    # - convert the univariate argument polynomial into a multivariate
1342
    #   version;
1343
    p = mpolynomialRing(polySa)
1344
    # - return the string representation of the converted form.
1345
    # There is no simple str() method defined for p's class.
1346
    return(p.__str__())
1347
#
1348
print pobyso_get_prec()  
1349
pobyso_set_prec(165)
1350
print pobyso_get_prec()  
1351
a=100
1352
print type(a)
1353
id(a)
1354
print "Max arity: ", pobyso_max_arity
1355
print "Function tripleDouble (43) as a string: ", pobyso_function_type_as_string(43)
1356
print "Function None (44) as a string: ", pobyso_function_type_as_string(44)
1357
print "...Pobyso check done"