Statistiques
| Révision :

root / pobysoPythonSage / src / pobyso.py @ 123

Historique | Voir | Annoter | Télécharger (46,76 ko)

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

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

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

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

    
86
pobyso_max_arity = 9
87

    
88
def pobyso_absolute_so_so():
89
    return(sollya_lib_absolute(None))
90

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

    
768
def pobyso_infnorm_so_so(func, interval, file = None, intervalList = None):
769
    print "Do not use this function. User pobyso_supnorm_so_so instead."
770
    return(None)
771

    
772
def pobyso_interval_to_range_sa_so(intervalSa, precisionSa=None):
773
    if precisionSa is None:
774
        precisionSa = intervalSa.parent().precision()
775
    intervalSo = pobyso_bounds_to_range_sa_so(intervalSa.lower(),\
776
                                              intervalSa.upper(),\
777
                                              precisionSa)
778
    return(intervalSo)
779
# End pobyso_interval_to_range_sa_so
780

    
781
def pobyso_lib_init():
782
    sollya_lib_init(None)
783

    
784
def pobyso_lib_close():
785
    sollya_lib_close(None)
786
    
787
def pobyso_name_free_variable(freeVariableNameSa):
788
    """ Legacy function. See pobyso_name_free_variable_sa_so. """
789
    pobyso_name_free_variable_sa_so(freeVariableNameSa)
790

    
791
def pobyso_name_free_variable_sa_so(freeVariableNameSa):
792
    """
793
    Set the free variable name in Sollya from a Sage string.
794
    """
795
    sollya_lib_name_free_variable(freeVariableNameSa)
796

    
797
def pobyso_parse_string(string):
798
    """ Legacy function. See pobyso_parse_string_sa_so. """
799
    return(pobyso_parse_string_sa_so(string))
800
 
801
def pobyso_parse_string_sa_so(string):
802
    """
803
    Get the Sollya expression computed from a Sage string.
804
    """
805
    return(sollya_lib_parse_string(string))
806

    
807
def pobyso_range(rnLowerBound, rnUpperBound):
808
    """ Legacy function. See pobyso_range_sa_so. """
809
    return(pobyso_range_sa_so(rnLowerBound, rnUpperBound)) 
810

    
811

    
812
def pobyso_range_to_interval_so_sa(rangeSo, realIntervalFieldSa = None):
813
    """
814
    Get a Sage interval from a Sollya range.
815
    If no realIntervalField is given as a parameter, the Sage interval
816
    precision is that of the Sollya range.
817
    Otherwise, the precision is that of the realIntervalField. In this case
818
    rounding may happen.
819
    """
820
    if realIntervalFieldSa is None:
821
        precSa = pobyso_get_prec_of_range_so_sa(rangeSo)
822
        realIntervalFieldSa = RealIntervalField(precSa)
823
    intervalSa = \
824
        pobyso_get_interval_from_range_so_sa(rangeSo, realIntervalFieldSa) 
825
    return(intervalSa)
826

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

    
944
def pobyso_remez_canonical_so_so(funcSo, \
945
                                 degreeSo, \
946
                                 rangeSo, \
947
                                 weightSo = pobyso_constant_1_sa_so(),\
948
                                 qualitySo = None):
949
    """
950
    All arguments are pointers to Sollya objects.
951
    The return value is a pointer to a Sollya function.
952
    """
953
    if not sollya_lib_obj_is_function(funcSo):
954
        return(None)
955
    return(sollya_lib_remez(funcSo, degreeSo, rangeSo, weightSo, qualitySo, None))
956
    
957
def pobyso_set_canonical_off():
958
    sollya_lib_set_canonical(sollya_lib_off())
959

    
960
def pobyso_set_canonical_on():
961
    sollya_lib_set_canonical(sollya_lib_on())
962

    
963
def pobyso_set_prec(p):
964
    """ Legacy function. See pobyso_set_prec_sa_so. """
965
    pobyso_set_prec_sa_so(p)
966

    
967
def pobyso_set_prec_sa_so(p):
968
    a = c_int(p)
969
    precSo = c_void_p(sollya_lib_constant_from_int(a))
970
    sollya_lib_set_prec(precSo, None)
971

    
972
def pobyso_set_prec_so_so(newPrecSo):
973
    sollya_lib_set_prec(newPrecSo, None)
974

    
975
def pobyso_supnorm_so_so(polySo, funcSo, intervalSo, errorTypeSo = None,\
976
                         accuracySo = None):
977
    """
978
    Computes the supnorm of the approximation error between the given 
979
    polynomial and function.
980
    errorTypeSo defaults to "absolute".
981
    accuracySo defaults to 2^(-40).
982
    """
983
    if errorTypeSo is None:
984
        errorTypeSo = sollya_lib_absolute(None)
985
        errorTypeIsNone = True
986
    else:
987
        errorTypeIsNone = False
988
    #
989
    if accuracySo is None:
990
        # Notice the **!
991
        accuracySo = pobyso_constant_sa_so(RR(2**(-40)))
992
        accuracyIsNone = True
993
    else:
994
        accuracyIsNone = False
995
    pobyso_autoprint(accuracySo)
996
    resultSo = \
997
        sollya_lib_supnorm(polySo, funcSo, intervalSo, errorTypeSo, \
998
                              accuracySo)
999
    if errorTypeIsNone:
1000
        sollya_lib_clear_obj(errorTypeSo)
1001
    if accuracyIsNone:
1002
        sollya_lib_clear_obj(accuracySo)
1003
    return resultSo
1004
# End pobyso_supnorm_so_so
1005

    
1006
def pobyso_taylor_expansion_with_change_var_so_so(functionSo, degreeSo, \
1007
                                                  rangeSo, \
1008
                                                  errorTypeSo=None, \
1009
                                                  sollyaPrecSo=None):
1010
    """
1011
    Compute the Taylor expansion with the variable change
1012
    x -> (x-intervalCenter) included.
1013
    """
1014
    # No global change of the working precision.
1015
    if not sollyaPrecSo is None:
1016
        initialPrecSo = sollya_lib_get_prec(None)
1017
        sollya_lib_set_prec(sollyaPrecSo)
1018
    #
1019
    # Error type stuff: default to absolute.
1020
    if errorTypeSo is None:
1021
        errorTypeIsNone = True
1022
        errorTypeSo = sollya_lib_absolute(None)
1023
    else:
1024
        errorTypeIsNone = False
1025
    intervalCenterSo = sollya_lib_mid(rangeSo)
1026
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, \
1027
                                         intervalCenterSo, \
1028
                                         rangeSo, errorTypeSo, None)
1029
    # taylorFormListSaSo is a Python list of Sollya objects references that 
1030
    # are copies of the elements of taylorFormSo.
1031
    # pobyso_get_list_elements_so_so clears taylorFormSo.
1032
    (taylorFormListSo, numElements, isEndElliptic) = \
1033
        pobyso_get_list_elements_so_so(taylorFormSo)
1034
    polySo = taylorFormListSo[0]
1035
    errorRangeSo = taylorFormListSo[2]
1036
    maxErrorSo = sollya_lib_sup(errorRangeSo)
1037
    changeVarExpSo = sollya_lib_build_function_sub(\
1038
                       sollya_lib_build_function_free_variable(),\
1039
                       sollya_lib_copy_obj(intervalCenterSo))
1040
    polyVarChangedSo = sollya_lib_evaluate(polySo, changeVarExpSo) 
1041
    sollya_lib_clear_obj(changeVarExpSo)
1042
    # If changed, reset the Sollya working precision.
1043
    if not sollyaPrecSo is None:
1044
        sollya_lib_set_prec(initialPrecSo)
1045
        sollya_lib_clear_obj(initialPrecSo)
1046
    if errorTypeIsNone:
1047
        sollya_lib_clear_obj(errorTypeSo)
1048
    sollya_lib_clear_obj(taylorFormSo)
1049
    # Do not clear maxErrorSo.
1050
    return((polyVarChangedSo, intervalCenterSo, maxErrorSo))
1051
# end pobyso_taylor_expansion_with_change_var_so_so
1052

    
1053
def pobyso_taylor_expansion_no_change_var_so_so(functionSo, degreeSo, rangeSo,
1054
                                                errorTypeSo=None, 
1055
                                                sollyaPrecSo=None):
1056
    """
1057
    Compute the Taylor expansion without the variable change
1058
    x -> x-intervalCenter.
1059
    """
1060
    # No global change of the working precision.
1061
    if not sollyaPrecSo is None:
1062
        initialPrecSo = sollya_lib_get_prec(None)
1063
        sollya_lib_set_prec(sollyaPrecSo)
1064
    # Error type stuff: default to absolute.
1065
    if errorTypeSo is None:
1066
        errorTypeIsNone = True
1067
        errorTypeSo = sollya_lib_absolute(None)
1068
    else:
1069
        errorTypeIsNone = False
1070
    intervalCenterSo = sollya_lib_mid(rangeSo, None)
1071
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo,
1072
                                         intervalCenterSo,
1073
                                         rangeSo, errorTypeSo, None)
1074
    # taylorFormListSaSo is a Python list of Sollya objects references that 
1075
    # are copies of the elements of taylorFormSo.
1076
    # pobyso_get_list_elements_so_so clears taylorFormSo.
1077
    (taylorFormListSaSo, numElementsSa, isEndEllipticSa) = \
1078
        pobyso_get_list_elements_so_so(taylorFormSo)
1079
    polySo = sollya_lib_copy_obj(taylorFormListSaSo[0])
1080
    #print "Num elements:", numElementsSa
1081
    sollya_lib_clear_obj(taylorFormSo)
1082
    #polySo = taylorFormListSaSo[0]
1083
    #errorRangeSo = sollya_lib_copy_obj(taylorFormListSaSo[2])
1084
    errorRangeSo = taylorFormListSaSo[2]
1085
    # No copy_obj needed here: a new object is created.
1086
    maxErrorSo = sollya_lib_sup(errorRangeSo)
1087
    # If changed, reset the Sollya working precision.
1088
    if not sollyaPrecSo is None:
1089
        sollya_lib_set_prec(initialPrecSo)
1090
        sollya_lib_clear_obj(initialPrecSo)
1091
    if errorTypeIsNone:
1092
        sollya_lib_clear_obj(errorTypeSo)
1093
    pobyso_clear_taylorform_sa_so(taylorFormListSaSo)
1094
    return((polySo, intervalCenterSo, maxErrorSo))
1095
# end pobyso_taylor_expansion_no_change_var_so_so
1096

    
1097
def pobyso_taylor(function, degree, point):
1098
    """ Legacy function. See pobysoTaylor_so_so. """
1099
    return(pobyso_taylor_so_so(function, degree, point))
1100

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

    
1190
def pobyso_taylorform_so_so(functionSo, degreeSo, pointSo, intervalSo=None, \
1191
                            errorTypeSo=None):
1192
    createdErrorType = False
1193
    if errorTypeSo is None:
1194
        errorTypeSo = sollya_lib_absolute()
1195
        createdErrorType = True
1196
    else:
1197
        #TODO: deal with the other case.
1198
        pass
1199
    if intervalSo is None:
1200
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
1201
                                         errorTypeSo, None)
1202
    else:
1203
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
1204
                                         intervalSo, errorTypeSo, None)
1205
    if createdErrorType:
1206
        sollya_lib_clear_obj(errorTypeSo)
1207
    return(resultSo)
1208
        
1209

    
1210
def pobyso_univar_polynomial_print_reverse(polySa):
1211
    """ Legacy function. See pobyso_univar_polynomial_print_reverse_sa_sa. """
1212
    return(pobyso_univar_polynomial_print_reverse_sa_sa(polySa))
1213

    
1214
def pobyso_univar_polynomial_print_reverse_sa_sa(polySa):
1215
    """
1216
    Return the string representation of a univariate polynomial with
1217
    monomials ordered in the x^0..x^n order of the monomials.
1218
    Remember: Sage
1219
    """
1220
    polynomialRing = polySa.base_ring()
1221
    # A very expensive solution:
1222
    # -create a fake multivariate polynomial field with only one variable,
1223
    #   specifying a negative lexicographical order;
1224
    mpolynomialRing = PolynomialRing(polynomialRing.base(), \
1225
                                     polynomialRing.variable_name(), \
1226
                                     1, order='neglex')
1227
    # - convert the univariate argument polynomial into a multivariate
1228
    #   version;
1229
    p = mpolynomialRing(polySa)
1230
    # - return the string representation of the converted form.
1231
    # There is no simple str() method defined for p's class.
1232
    return(p.__str__())
1233
#
1234
print pobyso_get_prec()  
1235
pobyso_set_prec(165)
1236
print pobyso_get_prec()  
1237
a=100
1238
print type(a)
1239
id(a)
1240
print "Max arity: ", pobyso_max_arity
1241
print "Function tripleDouble (43) as a string: ", pobyso_function_type_as_string(43)
1242
print "Function None (44) as a string: ", pobyso_function_type_as_string(44)
1243
print "...Pobyso check done"