Statistiques
| Révision :

root / pobysoPythonSage / src / pobyso.py @ 56

Historique | Voir | Annoter | Télécharger (33,59 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, includes
10
  the void function and the no arguments function that return a pointer to a Sollya
11
  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 Sollya
14
  world to Sage/Python world);
15
- the _sa_so (arguments are Sage/Python objects, returned objects are 
16
  pointers to Sollya objects);
17
- the sa_sa (arguments and returned objects are all Sage/Python objects);
18
- a catch all flavor, without any suffix, (e. g. functions that have no argument 
19
  nor return value).
20
NOTES:
21
Reported errors in Eclipse come from the calls to
22
the Sollya library
23

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

    
81
pobyso_max_arity = 9
82

    
83
def pobyso_autoprint(arg):
84
    sollya_lib_autoprint(arg,None)
85

    
86
def pobyso_autoprint_so_so(arg):
87
    sollya_lib_autoprint(arg,None)
88
    
89
def pobyso_absolute_so_so():
90
    return(sollya_lib_absolute(None))
91

    
92
def pobyso_build_function_sub_so_so(exp1So, exp2So):
93
    return(sollya_lib_build_function_sub(exp1So, exp2So))
94

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

    
114
def pobyso_change_var_in_function_so_so(funcSo, chvarExpSo):
115
    return(sollya_lib_evaluate(funcSo,chvarExpSo))
116
    
117

    
118
def pobyso_chebyshevform_so_so(functionSo, degreeSo, intervalSo):
119
    resultSo = sollya_lib_chebyshevform(functionSo, degreeSo, intervalSo)
120
    return(resultSo)
121
        
122

    
123

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

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

    
172
def pobyso_constant_1():
173
    """ Legacy function. See pobyso_constant_so_so. """
174
    return(pobyso_constant_1_sa_so())
175

    
176
def pobyso_constant_1_sa_so():
177
    return(pobyso_constant_from_int_sa_so(1))
178

    
179
def pobyso_constant_from_int(anInt):
180
    """ Legacy function. See pobyso_constant_from_int_sa_so. """
181
    return(pobyso_constant_from_int_sa_so(anInt))
182

    
183
def pobyso_constant_from_int_sa_so(anInt):
184
    return(sollya_lib_constant_from_int(int(anInt)))
185

    
186
def pobyso_function_type_as_string(funcType):
187
    """ Legacy function. See pobyso_function_type_as_string_so_sa. """
188
    return(pobyso_function_type_as_string_so_sa(funcType))
189

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

    
286
def pobyso_get_constant(rnArg, soConst):
287
    """ Legacy function. See pobyso_get_constant_so_sa. """
288
    return(pobyso_get_constant_so_sa(rnArg, soConst))
289

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

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

    
321
def pobyso_get_free_variable_name():
322
    """ Legacy function. See pobyso_get_free_variable_name_so_sa."""
323
    return(pobyso_get_free_variable_name_so_sa())
324

    
325
def pobyso_get_free_variable_name_so_sa():
326
    return(sollya_lib_get_free_variable_name())
327
    
328
def pobyso_get_function_arity(expressionSo):
329
    """ Legacy function. See pobyso_get_function_arity_so_sa."""
330
    return(pobyso_get_function_arity_so_sa(expressionSo))
331

    
332
def pobyso_get_function_arity_so_sa(expressionSo):
333
    arity = c_int(0)
334
    sollya_lib_get_function_arity(byref(arity),expressionSo)
335
    return(int(arity.value))
336

    
337
def pobyso_get_head_function(expressionSo):
338
    """ Legacy function. See pobyso_get_head_function_so_sa. """
339
    return(pobyso_get_head_function_so_sa(expressionSo)) 
340

    
341
def pobyso_get_head_function_so_sa(expressionSo):
342
    functionType = c_int(0)
343
    sollya_lib_get_head_function(byref(functionType), expressionSo, None)
344
    return(int(functionType.value))
345

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

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

    
391
def pobyso_get_max_prec_of_exp(soExp):
392
    """ Legacy function. See pobyso_get_max_prec_of_exp_so_sa. """
393
    return(pobyso_get_max_prec_of_exp_so_sa(soExp))
394

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

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

    
442
def pobyso_get_sage_exp_from_sollya_exp(sollyaExp, realField = RR):
443
    """ Legacy function. See pobyso_get_sage_exp_from_sollya_exp_so_sa. """
444
    return(pobyso_get_sage_exp_from_sollya_exp_so_sa(sollyaExp, realField = RR))
445

    
446
def pobyso_get_sage_exp_from_sollya_exp_so_sa(sollyaExp, realField = RR):
447
    """
448
    Get a Sage expression from a Sollya expression. 
449
    Currently only tested with polynomials with floating-point coefficients.
450
    Notice that, in the returned polynomial, the exponents are RealNumbers.
451
    """
452
    #pobyso_autoprint(sollyaExp)
453
    operator = pobyso_get_head_function_so_sa(sollyaExp)
454
    # Constants and the free variable are special cases.
455
    # All other operator are dealt with in the same way.
456
    if (operator != SOLLYA_BASE_FUNC_CONSTANT) and \
457
       (operator != SOLLYA_BASE_FUNC_FREE_VARIABLE):
458
        (arity, subexpressions) = pobyso_get_subfunctions_so_sa(sollyaExp)
459
        if arity == 1:
460
            sageExp = eval(pobyso_function_type_as_string_so_sa(operator) + \
461
            "(" + pobyso_get_sage_exp_from_sollya_exp_so_sa(subexpressions[0], \
462
            realField) + ")")
463
        elif arity == 2:
464
            if operator == SOLLYA_BASE_FUNC_POW:
465
                operatorAsString = "**"
466
            else:
467
                operatorAsString = \
468
                    pobyso_function_type_as_string_so_sa(operator)
469
            sageExp = \
470
              eval("pobyso_get_sage_exp_from_sollya_exp_so_sa(subexpressions[0], realField)"\
471
              + " " + operatorAsString + " " + \
472
                   "pobyso_get_sage_exp_from_sollya_exp_so_sa(subexpressions[1], realField)")
473
        # We do not know yet how to deal with arity > 3 (is there any in Sollya anyway?).
474
        else:
475
            sageExp = eval('None')
476
        return(sageExp)
477
    elif operator == SOLLYA_BASE_FUNC_CONSTANT:
478
        #print "This is a constant"
479
        return pobyso_get_constant_as_rn_with_rf_so_sa(sollyaExp, realField)
480
    elif operator == SOLLYA_BASE_FUNC_FREE_VARIABLE:
481
        #print "This is free variable"
482
        return(eval(sollya_lib_get_free_variable_name()))
483
    else:
484
        print "Unexpected"
485
        return eval('None')
486
# End pobyso_get_sage_poly_from_sollya_poly
487
    
488
def pobyso_get_subfunctions(expressionSo):
489
    """ Legacy function. See pobyso_get_subfunctions_so_sa. """
490
    return(pobyso_get_subfunctions_so_sa(expressionSo)) 
491

    
492
def pobyso_get_subfunctions_so_sa(expressionSo):
493
    """
494
    Get the subfunctions of an expression.
495
    Return the number of subfunctions and the list of subfunctions addresses.
496
    S.T.: Could not figure out another way than that ugly list of declarations
497
    to recover the addresses of the subfunctions. 
498
    """
499
    subf0 = c_int(0)
500
    subf1 = c_int(0)
501
    subf2 = c_int(0)
502
    subf3 = c_int(0)
503
    subf4 = c_int(0)
504
    subf5 = c_int(0)
505
    subf6 = c_int(0)
506
    subf7 = c_int(0)
507
    subf8 = c_int(0)
508
    arity = c_int(0)
509
    nullPtr = POINTER(c_int)()
510
    sollya_lib_get_subfunctions(expressionSo, byref(arity), \
511
    byref(subf0), byref(subf1), byref(subf2), byref(subf3), byref(subf4), byref(subf5),\
512
     byref(subf6), byref(subf7), byref(subf8), nullPtr, None) 
513
#    byref(cast(subfunctions[0], POINTER(c_int))), byref(cast(subfunctions[0], POINTER(c_int))), \
514
#    byref(cast(subfunctions[2], POINTER(c_int))), byref(cast(subfunctions[3], POINTER(c_int))), \
515
#    byref(cast(subfunctions[4], POINTER(c_int))), byref(cast(subfunctions[5], POINTER(c_int))), \
516
#    byref(cast(subfunctions[6], POINTER(c_int))), byref(cast(subfunctions[7], POINTER(c_int))), \
517
#    byref(cast(subfunctions[8], POINTER(c_int))), nullPtr)
518
    subfunctions = [subf0, subf1, subf2, subf3, subf4, subf5, subf6, subf7, subf8]
519
    subs = []
520
    if arity.value > pobyso_max_arity:
521
        return(0,[])
522
    for i in xrange(arity.value):
523
        subs.append(int(subfunctions[i].value))
524
        #print subs[i]
525
    return(int(arity.value), subs)
526
    
527
def pobyso_get_prec():
528
    """ Legacy function. See pobyso_get_prec_so_sa(). """
529
    return(pobyso_get_prec_so_sa())
530

    
531
def pobyso_get_prec_so_sa():
532
    """
533
    Get the current default precision in Sollya.
534
    The return value is Sage/Python int.
535
    """
536
    precSo = sollya_lib_get_prec(None)
537
    precSa = c_int(0)
538
    sollya_lib_get_constant_as_int(byref(precSa), precSo)
539
    sollya_lib_clear_obj(precSo)
540
    return(int(precSa.value))
541

    
542
def pobyso_get_prec_of_constant(ctExpSo):
543
    """ Legacy function. See pobyso_get_prec_of_constant_so_sa. """
544
    return(pobyso_get_prec_of_constant_so_sa(ctExpSo))
545

    
546
def pobyso_get_prec_of_constant_so_sa(ctExpSo):
547
    prec = c_int(0)
548
    retc = sollya_lib_get_prec_of_constant(byref(prec), ctExpSo, None)
549
    if retc == 0:
550
        return(None)
551
    return(int(prec.value))
552

    
553
def pobyso_get_prec_of_range_so_sa(rangeSo):
554
    prec = c_int(0)
555
    retc = sollya_lib_get_prec_of_range(byref(prec), rangeSo, None)
556
    if retc == 0:
557
        return(None)
558
    return(int(prec.value))
559

    
560
def pobyso_infnorm_so_so(func, interval, file = None, intervalList = None):
561
    print "Do not use this function. User pobyso_supnorm_so_so instead."
562
    return(None)
563

    
564
def pobyso_lib_init():
565
    sollya_lib_init(None)
566
    
567
def pobyso_name_free_variable(freeVariableName):
568
    """ Legacy function. See pobyso_name_free_variable_sa_so. """
569
    pobyso_name_free_variable_sa_so(freeVariableName)
570

    
571
def pobyso_name_free_variable_sa_so(freeVariableName):
572
    sollya_lib_name_free_variable(freeVariableName)
573

    
574
def pobyso_parse_string(string):
575
    """ Legacy function. See pobyso_parse_string_sa_so. """
576
    return(pobyso_parse_string_sa_so(string))
577
 
578
def pobyso_parse_string_sa_so(string):
579
    return(sollya_lib_parse_string(string))
580

    
581
def pobyso_range(rnLowerBound, rnUpperBound):
582
    """ Legacy function. See pobyso_range_sa_so. """
583
    return(pobyso_range_sa_so(rnLowerBound, rnUpperBound)) 
584

    
585
def pobyso_range_sa_so(rnLowerBoundSa, rnUpperBoundSa):
586
    """
587
    Return a Sollya range from to 2 RealField elements.
588
    The Sollya range element has a sufficient precision to hold all
589
    the digits of the bounds.
590
    """
591
    # TODO: check the bounds.
592
    lbPrec = rnLowerBoundSa.parent().precision()
593
    ubPrec = rnLowerBoundSa.parent().precision()
594
    currentSollyaPrecSa = pobyso_get_prec_so_sa()
595
    maxPrecSa = max(lbPrec, ubPrec, currentSollyaPrecSa)
596
    # Change the current Sollya precision only if necessary.
597
    if maxPrecSa > currentSollyaPrecSa:
598
        currentPrecSo = sollya_lib_get_prec(None)
599
        newPrecSo = solly_lib_constant_from_uint64(maxPrecSa)
600
        sollya_lib_set_prec(newPrecSo)
601
    lowerBoundSo = sollya_lib_constant(get_rn_value(rnLowerBoundSa))
602
    upperBoundSo = sollya_lib_constant(get_rn_value(rnUpperBoundSa))
603
    rangeSo = sollya_lib_range(lowerBoundSo, upperBoundSo)
604
    currentPrecSo = sollya_lib_get_prec(None)
605
    if maxPrecSa > currentSollyaPrecSa:
606
        sollya_lib_set_prec(currentPrecSo)
607
        sollya_lib_clear_obj(currentPrecSo)
608
        sollya_lib_clear_obj(newPrecSo)
609
    sollya_lib_clear_obj(lowerBoundSo)
610
    sollya_lib_clear_obj(upperBoundSo)
611
    return(rangeSo)
612

    
613
def pobyso_range_so_sa(rangeSo, realIntervalField = None):
614
    if realIntervalField is None:
615
        precSa = pobyso_get_prec_of_range_so_sa(rangeSo)
616
        realIntervalField = RealIntervalField(precSa)
617
    intervalSa = \
618
        pobyso_get_interval_from_range_so_sa(rangeSo, realIntervalField) 
619
    return(intervalSa)
620

    
621
def pobyso_remez_canonical_sa_sa(func, \
622
                                 degree, \
623
                                 lowerBound, \
624
                                 upperBound, \
625
                                 weight = None, \
626
                                 quality = None):
627
    """
628
    All arguments are Sage/Python.
629
    The functions (func and weight) must be passed as expressions or strings.
630
    Otherwise the function fails. 
631
    The return value is a pointer is a Sage polynomial.
632
    """
633
    var('zorglub')    # Dummy variable name for type check only.
634
    polySo = pobyso_remez_canonical_sa_so(func, \
635
                                 degree, \
636
                                 lowerBound, \
637
                                 upperBound, \
638
                                 weight = None, \
639
                                 quality = None)
640
    if parent(func) == parent("string"):
641
        functionSa = eval(func)
642
    # Expression test.
643
    elif type(func) == type(zorglub):
644
        functionSa = func
645
    maxPrecision = 0
646
    if polySo is None:
647
        return(None)
648
    maxPrecision = pobyso_get_max_prec_of_exp_so_sa(polySo)
649
    RRRR = RealField(maxPrecision)
650
    polynomialRing = RRRR[functionSa.variables()[0]]
651
    expSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, RRRR)
652
    polySa = polynomial(expSa, polynomialRing)
653
    return(polySa)
654
    
655
def pobyso_remez_canonical(func, \
656
                           degree, \
657
                           lowerBound, \
658
                           upperBound, \
659
                           weight = "1", \
660
                           quality = None):
661
    """ Legacy function. See pobyso_remez_canonical_sa_so. """
662
    return(pobyso_remez_canonical_sa_so(func, \
663
                                        degree, \
664
                                        lowerBound, \
665
                                        upperBound, \
666
                                        weight, \
667
                                        quality))
668
def pobyso_remez_canonical_sa_so(func, \
669
                                 degree, \
670
                                 lowerBound, \
671
                                 upperBound, \
672
                                 weight = None, \
673
                                 quality = None):
674
    """
675
    All arguments are Sage/Python.
676
    The functions (func and weight) must be passed as expressions or strings.
677
    Otherwise the function fails. 
678
    The return value is a pointer to a Sollya function.
679
    """
680
    var('zorglub')    # Dummy variable name for type check only.
681
    currentVariableName = None
682
    # The func argument can be of different types (string, 
683
    # symbolic expression...)
684
    if parent(func) == parent("string"):
685
        functionSo = sollya_lib_parse_string(func)
686
    # Expression test.
687
    elif type(func) == type(zorglub):
688
        # Until we are able to translate Sage expressions into Sollya 
689
        # expressions : parse the string version.
690
        currentVariableName = func.variables()[0]
691
        sollya_lib_name_free_variable(str(currentVariableName))
692
        functionSo = sollya_lib_parse_string(func._assume_str())
693
    else:
694
        return(None)
695
    if weight is None:
696
        weightSo = pobyso_constant_1_sa_so()
697
    elif parent(weight) == parent("string"):
698
        weightSo = sollya_lib_parse_string(func)
699
    elif type(weight) == type(zorglub): 
700
        functionSo = sollya_lib_parse_string_sa_so(weight._assume_str())
701
    else:
702
        return(None)
703
    degreeSo = pobyso_constant_from_int(degree)
704
    rangeSo = pobyso_range_sa_so(lowerBound, upperBound)
705
    if not quality is None:
706
        qualitySo= pobyso_constant_sa_so(quality)
707
    else:
708
        qualitySo = None
709
    return(sollya_lib_remez(functionSo, \
710
                            degreeSo, \
711
                            rangeSo, \
712
                            weightSo, \
713
                            qualitySo, \
714
                            None))
715
    
716
def pobyso_remez_canonical_so_so(funcSo, \
717
                                 degreeSo, \
718
                                 rangeSo, \
719
                                 weightSo = pobyso_constant_1_sa_so(),\
720
                                 qualitySo = None):
721
    """
722
    All arguments are pointers to Sollya objects.
723
    The return value is a pointer to a Sollya function.
724
    """
725
    if not sollya_lib_obj_is_function(funcSo):
726
        return(None)
727
    return(sollya_lib_remez(funcSo, degreeSo, rangeSo, weightSo, qualitySo, None))
728
    
729
def pobyso_set_canonical_off():
730
    sollya_lib_set_canonical(sollya_lib_off())
731

    
732
def pobyso_set_canonical_on():
733
    sollya_lib_set_canonical(sollya_lib_on())
734

    
735
def pobyso_set_prec(p):
736
    """ Legacy function. See pobyso_set_prec_sa_so. """
737
    return( pobyso_set_prec_sa_so(p))
738

    
739
def pobyso_set_prec_sa_so(p):
740
    a = c_int(p)
741
    precSo = c_void_p(sollya_lib_constant_from_int(a))
742
    sollya_lib_set_prec(precSo)
743

    
744
def pobyso_supnorm_so_so(polySo, funcSo, intervalSo, errorTypeSo, accuracySo):
745
    return(sollya_lib_supnorm(polySo, funcSo, intervalSo, errorTypeSo, \
746
                              accuracySo))
747

    
748
def pobyso_taylor_expansion_no_change_var_so_so(functionSo, degreeSo, rangeSo, \
749
                                                errorTypeSo, \
750
                                                sollyaPrecSo=None):
751
    # No global change of the working precision.
752
    if not sollyaPrecSo is None:
753
        initialPrecSo = sollya_lib_get_prec(None)
754
        sollya_lib_set_prec(sollyaPrecSo)
755
    #
756
    intervalCenterSo = sollya_lib_mid(rangeSo)
757
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, \
758
                                         intervalCenterSo, \
759
                                         rangeSo, errorTypeSo, None)
760
    (taylorFormListSo, numElements, isEndElliptic) = \
761
        pobyso_get_list_elements_so_so(taylorFormSo)
762
    polySo = taylorFormListSo[0]
763
    errorRangeSo = taylorFormListSo[2]
764
    maxErrorSo = sollya_lib_sup(errorRangeSo)
765
    # If changed, reset the Sollya working precision.
766
    if not sollyaPrecSo is None:
767
        sollya_lib_set_prec(initialPrecSo)
768
        sollya_lib_clear_obj(initailPrecSo)
769
    return([polySo, intervalCenterSo, maxErrorSo])
770
# end pobyso_taylor_expansion_no_change_var_so_so
771

    
772
def pobyso_taylor(function, degree, point):
773
    """ Legacy function. See pobysoTaylor_so_so. """
774
    return(pobyso_taylor_so_so(function, degree, point))
775

    
776
def pobyso_taylor_so_so(functionSo, degreeSo, pointSo):
777
    return(sollya_lib_taylor(functionSo, degreeSo, pointSo))
778
    
779
def pobyso_taylorform(function, degree, point = None, interval = None, errorType=None):
780
    """ Legacy function. See ;"""
781
    
782
def pobyso_taylorform_sa_sa(functionSa, \
783
                            degree, \
784
                            point, \
785
                            precision, \
786
                            interval=None, \
787
                            errorType=None):
788
    """
789
    Compute the Taylor form of 'degree' for 'functionSa' at 'point' 
790
    for 'interval' with 'errorType'. 
791
    point: must be a Real or a Real interval.
792
    return the Taylor form as an array
793
    TODO: take care of the interval and of point when it is an interval;
794
          when errorType is not None;
795
          take care of the other elements of the Taylor form (coefficients errors and
796
          delta.
797
    """
798
    # Absolute as the default error.
799
    if errorType is None:
800
        errorTypeSo = sollya_lib_absolute()
801
    else:
802
        #TODO: deal with the other case.
803
        pass
804
    varSa = functionSa.variables()[0]
805
    pointBaseRingString = str(point.base_ring())
806
    if not re.search('Real', pointBaseRingString):
807
        return None
808
    # Call Sollya but first "sollyafy" the arguments.
809
    pobyso_name_free_variable_sa_so(str(varSa))
810
    #pobyso_set_prec_sa_so(300)
811
    # Sollyafy the function.
812
    functionSo = pobyso_parse_string_sa_so(functionSa._assume_str())
813
    if sollya_lib_obj_is_error(functionSo):
814
        print "pobyso_tailorform: function string can't be parsed!"
815
        return None
816
    # Sollyafy the degree
817
    degreeSo = sollya_lib_constant_from_int(int(degree))
818
    # Sollyafy the point
819
    if not re.search('Interval', pointBaseRingString):
820
        pointSo  = pobyso_constant_sa_so(point)
821
    else:
822
        # TODO: deal with the interval case.
823
        pass
824
    # Call Sollya
825
    taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, errorTypeSo,\
826
                                         None)
827
    (tfsAsList, numElements, isEndElliptic) = \
828
            pobyso_get_list_elements_so_so(taylorFormSo)
829
    polySo = tfsAsList[0]
830
    maxPrecision = pobyso_get_max_prec_of_exp_so_sa(polySo)
831
    polyRealField = RealField(maxPrecision)
832
    expSa = pobyso_get_sage_exp_from_sollya_exp_so_sa(polySo, polyRealField)
833
    sollya_lib_close()
834
    polynomialRing = polyRealField[str(varSa)]
835
    polySa = polynomial(expSa, polynomialRing)
836
    taylorFormSa = [polySa]
837
    return(taylorFormSa)
838
# End pobyso_taylor_form_sa_sa
839

    
840
def pobyso_taylorform_so_so(functionSo, degreeSo, pointSo, intervalSo=None, \
841
                            errorTypeSo=None):
842
    createdErrorType = False
843
    if errorTypeSo is None:
844
        errorTypeSo = sollya_lib_absolute()
845
        createdErrorType = True
846
    else:
847
        #TODO: deal with the other case.
848
        pass
849
    if intervalSo is None:
850
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
851
                                         errorTypeSo, None)
852
    else:
853
        resultSo = sollya_lib_taylorform(functionSo, degreeSo, pointSo, \
854
                                         intervalSo, errorTypeSo, None)
855
    if createdErrorType:
856
        sollya_lib_clear_obj(errorTypeSo)
857
    return(resultSo)
858
        
859

    
860
def pobyso_univar_polynomial_print_reverse(polySa):
861
    """ Legacy function. See pobyso_univar_polynomial_print_reverse_sa_sa. """
862
    return(pobyso_univar_polynomial_print_reverse_sa_sa(polySa))
863

    
864
def pobyso_univar_polynomial_print_reverse_sa_sa(polySa):
865
    """
866
    Return the string representation of a univariate polynomial with
867
    monomials ordered in the x^0..x^n order of the monomials.
868
    Remember: Sage
869
    """
870
    polynomialRing = polySa.base_ring()
871
    # A very expensive solution:
872
    # -create a fake multivariate polynomial field with only one variable,
873
    #   specifying a negative lexicographical order;
874
    mpolynomialRing = PolynomialRing(polynomialRing.base(), \
875
                                     polynomialRing.variable_name(), \
876
                                     1, order='neglex')
877
    # - convert the univariate argument polynomial into a multivariate
878
    #   version;
879
    p = mpolynomialRing(polySa)
880
    # - return the string representation of the converted form.
881
    # There is no simple str() method defined for p's class.
882
    return(p.__str__())
883
#
884
print pobyso_get_prec()  
885
pobyso_set_prec(165)
886
print pobyso_get_prec()  
887
a=100
888
print type(a)
889
id(a)
890
print "Max arity: ", pobyso_max_arity
891
print "Function tripleDouble (43) as a string: ", pobyso_function_type_as_string(43)
892
print "Function None (44) as a string: ", pobyso_function_type_as_string(44)
893
print "...Pobyso check done"