giacpy_sage.pyx 186 KB
Newer Older
1
2
3
4
5
6
7
8
#*****************************************************************************
#    AUTHOR:  Han Frederic <frederic.han@imj-prg.fr>
#  2012
#  Distributed under the terms of the GNU General Public License (GPL)
#  version 2 or higher.
#                  http://www.gnu.org/licenses/
#*****************************************************************************
r"""
9
:Name: giacpy_sage
10
11
12
13
14
15
16
:Summary: A Cython frontend to the c++ library giac. (Computer Algebra System)
:License:  GPL v2 or above
:Home-page: http://www.math.jussieu.fr/~han/xcas/giacpy/
:Author: Frederic Han
:Author-email: frederic.han@imj-prg.fr


17
This is the sage version of giacpy. Since 0.6 it is named giacpy_sage to avoid confusions
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37

- Giacpy is an interface to the c++ giac library. This interface is built with cython, and the resulting speed is similar to giac's one.

- Giac is a general purpose Computer algebra system by Bernard Parisse released under GPLv3.

      * http://www-fourier.ujf-grenoble.fr/~parisse/giac.html
      * It is build on C and C++ libraries:
        NTL (arithmetic), GSL (numerics), GMP (big integers), MPFR (bigfloats)
      * It  provides fast  algorithms  for multivariate polynomial operations (product, GCD, factorisation) and
      * symbolic  computations: solver, simplifications, limits/series, integration, summation...
      * Linear Algebra with numerical or symbolic coefficients.


AUTHORS:

- Frederic Han (2013-09-23): initial version


EXAMPLES:

Frederic HAN's avatar
Frederic HAN committed
38
- The class Pygen is the main tool to interact from python/sage with the c++ library giac via cython.
39
40
41
42
43
44
45
46

  The initialisation of a Pygen just create an object in giac, but the mathematical computation  is not done. This class is mainly for cython users.


  Here A is a Pygen element, and it is ready for any giac function.

  ::

47
        sage: from giacpy_sage import *    # random
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
        //...
        sage: A=Pygen('2+2');A
        2+2
        sage: A.eval()
        4

  In general, you may prefer to directly create a Pygen and execute the evaluation in giac. This is exactly  the meaning of the libgiac() function.

  ::

        sage: a=libgiac('2+2');a;isinstance(a,Pygen)
        4
        True

- Most common usage of this package in sage will be with the libgiac() function. This function is just the composition of the Pygen initialisation and the evaluation of this object in giac.

  ::

66
        sage: from giacpy_sage import libgiac,giacsettings
67
68
69
70
71
72
73
        sage: x,y,z=libgiac('x,y,z');  # add some giac objects
        sage: f=(x+3*y)/(x+z+1)^2 -(x+z+1)^2/(x+3*y)
        sage: f.factor()
        (3*y-x^2-2*x*z-x-z^2-2*z-1)*(3*y+x^2+2*x*z+3*x+z^2+2*z+1)/((x+z+1)^2*(3*y+x))
        sage: f.normal()
        (-x^4-4*x^3*z-4*x^3-6*x^2*z^2-12*x^2*z-5*x^2+6*x*y-4*x*z^3-12*x*z^2-12*x*z-4*x+9*y^2-z^4-4*z^3-6*z^2-4*z-1)/(x^3+3*x^2*y+2*x^2*z+2*x^2+6*x*y*z+6*x*y+x*z^2+2*x*z+x+3*y*z^2+6*y*z+3*y)

74
- To obtain more hints on giacpy_sage consider the help of the :func:`libgiac<_giac>` function.
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94

  ::

        sage: libgiac?             # doctest: +SKIP

- Some settings of giac are available via the ``giacsettings`` element. (Ex: maximal number of threads in computations, allowing probabilistic algorithms or not...

  ::

        sage: R=PolynomialRing(QQ,8,'x')
        sage: I=sage.rings.ideal.Katsura(R,8)
        sage: giacsettings.proba_epsilon=1e-15;
        sage: Igiac=libgiac(I.gens());
        sage: time Bgiac=Igiac.gbasis([R.gens()],'revlex')  # doctest: +SKIP
        Running a probabilistic check for the reconstructed Groebner basis. If successfull, error probability is less than 1e-15 and is estimated to be less than 10^-109. Use proba_epsilon:=0 to certify (this takes more time).
        Time: CPU 0.46 s, Wall: 0.50 s
        sage: giacsettings.proba_epsilon=0;
        sage: Igiac=libgiac(I.gens());
        sage: time Bgiac=Igiac.gbasis([R.gens()],'revlex') # doctest: +SKIP
        Time: CPU 2.74 s, Wall: 2.75 s
Frederic HAN's avatar
Frederic HAN committed
95
        sage: giacsettings.proba_epsilon=1e-15;
96
97
98

  ::

99
        sage: from giacpy_sage import *
100
        sage: x=libgiac('x');f=1/(2+sin(5*x))
101
102
103
104
        sage: oldrep=2/5/sqrt(3)*(atan((2*tan(5*x/2)+1)/sqrt(3))+pi*floor(5*x/2/pi+1/2))
        sage: newrep=2/5/sqrt(3)*(atan((-sqrt(3)*sin(5*x)+cos(5*x)+2*sin(5*x)+1)/(sqrt(3)*cos(5*x)+sqrt(3)-2*cos(5*x)+sin(5*x)+2))+5*x/2)
        sage: ((f.int()-newrep)*(f.int()-oldrep())).normal()
        0
105
106
107
108
109
110
111
112
        sage: f.series(x,0,3)
        1/2-5/4*x+25/8*x^2-125/48*x^3+x^4*order_size(x)
        sage: libgiac(sqrt(5)+pi).approx(100)
        5.377660631089582934871817052010779119637787758986631545245841837718337331924013898042449233720899343


SEEALSO:

Frederic HAN's avatar
Frederic HAN committed
113
    ``libgiac``, ``giacsettings``, ``Pygen``,``loadgiacgen``
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142



GETTING HELP:

- To obtain some help on a giac keyword use the help() method. In sage the htmlhelp() method for Pygen element is disabled. Just use the ? or .help() method.

  ::

        sage: libgiac.gcd?             # doctest: +SKIP
        "Returns the greatest common divisor of 2 polynomials of several variables or of 2 integers or of 2 rationals.
        (Intg or Poly),(Intg or Poly)
        gcd(45,75);gcd(15/7,50/9);gcd(x^2-2*x+1,x^3-1);gcd(t^2-2*t+1,t^2+t-2);gcd((x^2-1)*(y^2-1)*z^2,x^3*y^3*z+(-(y^3))*z+x^3*z-z)
        lcm,euler,modgcd,ezgcd,psrgcd,heugcd,Gcd"

- You can find full html documentation about the **giac** functions  at:

      * http://www-fourier.ujf-grenoble.fr/~parisse/giac/doc/en/cascmd_en/

      * http://www-fourier.ujf-grenoble.fr/~parisse/giac/doc/fr/cascmd_fr/

      * http://www-fourier.ujf-grenoble.fr/~parisse/giac/doc/el/cascmd_el/

      * or in :doc:`$SAGE_LOCAL/share/giac/doc/en/cascmd_en/index.html`



REMARK:

143
- Graphics 2D Output via qcas (the qt frontend to giac) is removed in the sage version of giacpy.
144
145
146
147
148

"""
########################################################


149
cimport giacpy_sage
150
151
152
153
154
155
156
157
158
159
160
161

# sage includes
from sage.libs.gmp.mpz cimport mpz_t, mpz_init_set
from sage.rings.integer cimport Integer
from sage.rings.all import ZZ, QQ, IntegerModRing
from sage.symbolic.ring import SR
from sage.rings.rational cimport Rational
from sage.matrix.matrix import is_Matrix

from sage.plot.line import line
from sage.plot.scatter_plot import scatter_plot

162
163
164
165
166
from sage.libs.pynac.pynac import symbol_table
from sage.calculus.calculus import symbolic_expression_from_string
from sage.symbolic.expression import Expression
from sage.symbolic.expression_conversions import InterfaceInit
from sage.interfaces.giac import giac
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200

from sys import maxsize as Pymaxint
import os
import math


#### Python3 compatibility  #############################
from sys import version_info as Pyversioninfo
if Pyversioninfo[0]>2 :
   PythonVersion3 = True
else:
   PythonVersion3 = False

def decstring23(s):
  if PythonVersion3 :
      return s.decode()
  else:
      return s


def encstring23(s):
  if PythonVersion3 :
      return bytes(s,'UTF-8')
  else:
      return s


if PythonVersion3:
    listrange=list,range
else:
    listrange=list

####End of Python3 compatibility########################

201
202
#NB: include 'sage/ext/interrupt.pxi' is deprecated in sage 7.2.beta1
#    use cysignals instead.
203
204
205
206
207
208
#
# stdsage.pxi and cydignal.pxi are deprecated, replace by cimport
#include 'cysignals/signals.pxi'
from cysignals.signals cimport *
#include 'sage/ext/stdsage.pxi'
from sage.ext.stdsage cimport PY_NEW
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227

########################################################


########################################################
# A global context pointer. One by giac session.
########################################################
cdef context * context_ptr=new context()
#
# Some global variables for optimisation
GIACNULL=Pygen('NULL')
#



# Create a giac setting instance
giacsettings=GiacSetting()
Pygen('I:=sqrt(-1)').eval()

228
229
230
231
232
233
# A function to convert SR Expression with defined giac conversion to a string
# for giac/libgiac.
# NB: We want to do this without starting an external giac program and
# self._giac_() does.
SRexpressiontoGiac=InterfaceInit(giac)
#
234
235
236
237
238
239
240
241
242
243
244
245
246

#######################################################
# The wrapper to eval with giac
#######################################################
# in sage we don't export the giac function. We replace it by libgiac
def _giac(s):
   """
   This function evaluate a python/sage object with the giac library. It creates in python/sage a Pygen element and evaluate it with giac:

   - **First Example**:

     ::

247
       sage: from giacpy_sage import libgiac
248
249
250
251
252
253
254
255
256
257
258
259
       sage: x,y=libgiac('x,y')
       sage: (x+2*y).cos().texpand()
       cos(x)*(2*cos(y)^2-1)-sin(x)*2*cos(y)*sin(y)



   - **Coercion, Pygen and internal giac variables**:

     The most usefull objects will be the Python object of type Pygen.

     ::

260
       sage: from giacpy_sage import *
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
       sage: x,y,z=libgiac('x,y,z')
       sage: f=sum([x[i] for i in range(5)])^15/(y+z);f.coeff(x[0],12)
       (455*(x[1])^3+1365*(x[1])^2*x[2]+1365*(x[1])^2*x[3]+1365*(x[1])^2*x[4]+1365*x[1]*(x[2])^2+2730*x[1]*x[2]*x[3]+2730*x[1]*x[2]*x[4]+1365*x[1]*(x[3])^2+2730*x[1]*x[3]*x[4]+1365*x[1]*(x[4])^2+455*(x[2])^3+1365*(x[2])^2*x[3]+1365*(x[2])^2*x[4]+1365*x[2]*(x[3])^2+2730*x[2]*x[3]*x[4]+1365*x[2]*(x[4])^2+455*(x[3])^3+1365*(x[3])^2*x[4]+1365*x[3]*(x[4])^2+455*(x[4])^3)/(y+z)

     Warning: The complex number sqrt(-1) is exported in python as I. (But it may appears as i)

     ::

       sage: libgiac((1+I*sqrt(3))^3).normal(); libgiac(1+I)
       -8
       1+i

     Python integers and reals can be directly converted to giac.

     ::

277
       sage: from giacpy_sage import *
278
279
280
281
282
283
284
285
       sage: a=libgiac(2^1024);a.nextprime();(libgiac(1.234567)).erf().approx(10)
       179769313486231590772930519078902473361797697894230657273430081157732675805500963132708477322407536021120113879871393357658789768814416622492847430639474124377767893424865485276302219601246094119453082952085005768838150682342462881473913110540827237163350510684586298239947245938479716304835356329624224137859
       0.9191788641

     The Python object y defined above is of type Pygen. It is not an internal giac variable. (Most of the time you won't need to use internal giac variables).

     ::

286
       sage: from giacpy_sage import *
287
288
289
290
291
292
293
294
295
296
297
298
       sage: libgiac('y:=1');y
       1
       y
       sage: libgiac.purge('y')
       1
       sage: libgiac('y')
       y

     There are some natural coercion to Pygen elements:

     ::

299
       sage: from giacpy_sage import *
300
301
302
303
304
       sage: libgiac(pi)>3.14 ; libgiac(pi) >3.15 ; libgiac(3)==3
       True
       False
       True

Frederic HAN's avatar
Frederic HAN committed
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
   - **Linear Algebra**:

    In Giac/Xcas vectors are just lists and matrices are lists of list:

    ::

       sage: from giacpy_sage import *
       sage: x,y=libgiac('x,y')
       sage: A=libgiac([[1,2],[3,4]])  # we create a giac matrix from it lines
       sage: v=libgiac([x,y]); v   # a giac vector
       [x,y]
       sage: A*v # matrix product with a vector outputs a vector
       [x+2*y,3*x+4*y]
       sage: v*v  # dot product
       x*x+y*y

    Remark that w=giac([[x],[y]]) is a matrix of 1 column and 2 rows. It is not a vector\
    so w*w doesn't make sense:

    ::

       sage: w=libgiac([[x],[y]])
       sage: w.transpose()*w   # this matrix product makes sense and output a 1x1 matrix.
       matrix[[x*x+y*y]]

    In sage affectation doesn't create a new matrix. (cf. pointers) see also \
 the doc of  'Pygen.__setitem__':

    ::

       sage: B1=A;
       sage: B1[0,0]=43; B1 # in place affectation changes both B1 and A
       [[43,2],[3,4]]
       sage: A
       [[43,2],[3,4]]
       sage: A[0][0]=A[0][0]+1; A  # similar as A[0,0]=A[0,0]+1
       [[44,2],[3,4]]
       sage: A.pcar(x)  # compute the characteristic polynomial of A
       x^2-48*x+170
       sage: B2=A.copy() # use copy to create another object
       sage: B2[0,0]=55; B2  # here A is not modified
       [[55,2],[3,4]]
       sage: A
       [[44,2],[3,4]]

    Sparse Matrices are avaible via the table function:

    ::

       sage: from giacpy_sage import *
       sage: A=libgiac.table(()); A  # create an empty giac table
       table(
       )
       sage: A[2,3]=33; A[0,2]='2/7' # set non zero entries of the sparse matrix
       sage: A*A  # basic matrix operation are supported with sparse matrices
       table(
       (0,3) = 66/7
       )
       sage: D=libgiac.diag([22,3,'1/7']); D  # some diagonal matrix
       [[22,0,0],[0,3,0],[0,0,1/7]]
       sage: libgiac.table(D)    # to create a sparse matrix from an ordinary one
       table(
       (0,0) = 22,
       (1,1) = 3,
       (2,2) = 1/7
       )

    But many matrix functions apply only with ordinary matrices so need conversions:

    ::

       sage: B1=A.matrix(); B1 # convert the sparse matrix to a matrix, but the size is minimal
       [[0,0,2/7,0],[0,0,0,0],[0,0,0,33]]
       sage: B2=B1.redim(4,4) # so we may need to resize B1
       sage: B2.pmin(x)  
       x^3

382
383
384
385
386
387
   - **Lists of Pygen and Giac lists**:

     Here l1 is a giac list and l2 is a python list of Pygen type objects.

     ::

388
       sage: from giacpy_sage import *
389
390
391
392
393
394
395
396
       sage: l1=libgiac(range(10)); l2=[1/(i^2+1) for i in l1]
       sage: sum(l2)
       33054527/16762850

     So l1+l1 is done in giac and means a vector addition. But l2+l2 is done in Python so it is the list concatenation.

     ::

397
       sage: from giacpy_sage import *
398
399
400
401
402
403
404
405
406
       sage: l1+l1
       [0,2,4,6,8,10,12,14,16,18]
       sage: l2+l2
       [1, 1/2, 1/5, 1/10, 1/17, 1/26, 1/37, 1/50, 1/65, 1/82, 1, 1/2, 1/5, 1/10, 1/17, 1/26, 1/37, 1/50, 1/65, 1/82]

     Here V is not a Pygen element. We need to push it to giac to use a giac method like dim, or we need to use an imported function.

     ::

407
       sage: from giacpy_sage import *
408
409
410
411
412
413
414
415
416
417
418
       sage: V=[ [x[i]^j for i in range(8)] for j in range(8)]
       sage: libgiac(V).dim()
       [8,8]
       sage: libgiac.det_minor(V).factor()
       (x[6]-(x[7]))*(x[5]-(x[7]))*(x[5]-(x[6]))*(x[4]-(x[7]))*(x[4]-(x[6]))*(x[4]-(x[5]))*(x[3]-(x[7]))*(x[3]-(x[6]))*(x[3]-(x[5]))*(x[3]-(x[4]))*(x[2]-(x[7]))*(x[2]-(x[6]))*(x[2]-(x[5]))*(x[2]-(x[4]))*(x[2]-(x[3]))*(x[1]-(x[7]))*(x[1]-(x[6]))*(x[1]-(x[5]))*(x[1]-(x[4]))*(x[1]-(x[3]))*(x[1]-(x[2]))*(x[0]-(x[7]))*(x[0]-(x[6]))*(x[0]-(x[5]))*(x[0]-(x[4]))*(x[0]-(x[3]))*(x[0]-(x[2]))*(x[0]-(x[1]))


   - **Modular objects with %**

     ::

419
       sage: from giacpy_sage import *
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
       sage: V=libgiac.ranm(5,6) % 2;
       sage: V.ker().rowdim()+V.rank()
       6
       sage: a=libgiac(7)%3;a;a%0;7%3
       1 % 3
       1
       1

     Do not confuse with the  python integers:

     ::

       sage: type(7%3)==type(a);type(a)==type(7%3)
       False
       False

   - **Syntaxes with reserved or unknown Python/sage symbols**:

Frederic HAN's avatar
Frederic HAN committed
438
     In general equations needs symbols such as ``=`` ``<`` ``>`` that have another meaning in Python or Sage.
439
440
441
442
     So those objects must be quoted.

     ::

443
       sage: from giacpy_sage import *
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
       sage: x=libgiac('x')
       sage: (1+2*sin(3*x)).solve(x).simplify()
       Warning, argument is not an equation, solving 1+2*sin(3*x)=0
       list[-pi/18,7*pi/18]

       sage: libgiac.solve('sin(3*x)>2*sin(x)',x)        # doctest: +SKIP
       ...
       RuntimeError: Unable to find numeric values solving equation. For trigonometric equations this may be solved using assumptions, e.g. assume(x>-pi && x<pi) Error: Bad Argument Value


     You can also add some hypothesis to a giac symbol:

     ::

       sage: libgiac.assume('x>-pi && x<pi')
       x
       sage: libgiac.solve('sin(3*x)>2*sin(x)',x)
       list[((x>(-5*pi/6)) and (x<(-pi/6))),((x>0) and (x<(pi/6))),((x>(5*pi/6)) and (x<pi))]

     To remove those hypothesis use the giac function: purge

     ::

       sage: libgiac.purge('x')
       assume[[],[line[-pi,pi]],[-pi,pi]]
       sage: libgiac.solve('x>0')
       list[x>0]

     Same problems with the ..

     ::

476
       sage: from giacpy_sage import *
477
       sage: x=libgiac('x')
478
479
480
481
482
       sage: f=1/(5+cos(4*x))
       sage: oldrep=1/2/(2*sqrt(6))*(atan(2*tan(4*x/2)/sqrt(6))+pi*floor(4*x/2/pi+1/2))
       sage: newrep=1/2/(2*sqrt(6))*(atan((-sqrt(6)*sin(4*x)+2*sin(4*x))/(sqrt(6)*cos(4*x)+sqrt(6)-2*cos(4*x)+2))+4*x/2)
       sage: ((f.int(x)-newrep)*(f.int(x)-oldrep)).normal()
       0
483
484
485
486
487
488
489
490
491
492
493
494
       sage: libgiac.fMax(f,'x=-0..pi').simplify()
       pi/4,3*pi/4
       sage: libgiac.fMax.help()  # doctest: +SKIP
       "Returns the abscissa of the maximum of the expression.
       Expr,[Var]
       fMax(-x^2+2*x+1,x)
       fMin"
       sage: libgiac.sum(1/(1+x^2),'x=0..infinity').simplify()
       (pi*exp(pi)^2+pi+exp(pi)^2-1)/(2*exp(pi)^2-2)

   - **From giac to sage**:

Frederic HAN's avatar
Frederic HAN committed
495
     One can convert a Pygen element to sage with the sage method.
496
497
498
499
500
501
502
503
     Get more details with:

     ::

       sage: Pygen.sage?          # doctest: +SKIP

     ::

504
       sage: from giacpy_sage import *
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
       sage: L=libgiac('[1,sqrt(5),[1.3,x]]')
       sage: L.sage()       # All entries are converted recursively
       [1, sqrt(5), [1.3, x]]

     To obtain matrices and vectors, use the :meth:`matrix<Pygen._matrix_>` and :meth:`vector<Pygen._vector_>` commands. Get more details with:

     ::

       sage: Pygen._matrix_?          # doctest: +SKIP
       sage: Pygen._vector_?          # doctest: +SKIP

     ::

       sage: n=var('n');A=matrix([[1,2],[-1,1]])
       sage: B=libgiac(A).matpow(n)    # We compute the symbolic power on A via libgiac
Frederic HAN's avatar
Frederic HAN committed
520
       sage: C=matrix(SR,B); C         # We convert B to sage
521
522
523
524
525
526
527
528
529
530
       [                     1/2*(I*sqrt(2) + 1)^n + 1/2*(-I*sqrt(2) + 1)^n -1/2*I*sqrt(2)*(I*sqrt(2) + 1)^n + 1/2*I*sqrt(2)*(-I*sqrt(2) + 1)^n]
       [ 1/4*I*sqrt(2)*(I*sqrt(2) + 1)^n - 1/4*I*sqrt(2)*(-I*sqrt(2) + 1)^n                      1/2*(I*sqrt(2) + 1)^n + 1/2*(-I*sqrt(2) + 1)^n]
       sage: (C.subs(n=3)-A^3).expand()
       [0 0]
       [0 0]


   **MEMENTO of usual GIAC functions**:

   - *Expand with simplification*
Frederic HAN's avatar
Frederic HAN committed
531

532
533
534
535
536
537
538
539
540
541
542
543
544
545
         * ``ratnormal``, ``normal``, ``simplify``   (from the fastest to the most sophisticated)

         *  NB: ``expand`` function doesn't regroup nor cancel terms, so it could be slow. (pedagogical purpose only?)

   - *Factor/Regroup*

         * ``factor``, ``factors``, ``regroup``, ``cfactor``, ``ifactor``

   - *Misc*

         * ``unapply``, ``op``, ``subst``

   - *Polynomials/Fractions*

Frederic HAN's avatar
Frederic HAN committed
546
         * ``coeff``,  ``gbasis``, ``greduce``, ``lcoeff``, ``pcoeff``, ``canonical_form``,
547
548
549

         * ``proot``,  ``poly2symb``,  ``symb2poly``, ``posubLMQ``, ``poslbdLMQ``, ``VAS``, ``tcoeff``,  ``valuation``

Frederic HAN's avatar
Frederic HAN committed
550
         * ``gcd``, ``egcd``, ``lcm``, ``quo``, ``rem``, ``quorem``, ``abcuv``, ``chinrem``,
551
552
553
554
555
556
557
558
559
560

         * ``peval``, ``horner``, ``lagrange``, ``ptayl``, ``spline``,  ``sturm``,  ``sturmab``

         * ``partfrac``, ``cpartfrac``

   - *Memory/Variables*

         * ``assume``, ``about``, ``purge``, ``ans``

   - *Calculus/Exact*
Frederic HAN's avatar
Frederic HAN committed
561

562
563
564
565
         * ``linsolve``,  ``solve``,  ``csolve``,  ``desolve``,  ``seqsolve``, ``reverse_rsolve``, ``matpow``

         * ``limit``, ``series``, ``sum``, ``diff``, ``fMax``, ``fMin``,

Frederic HAN's avatar
Frederic HAN committed
566
         * ``integrate``, ``subst``, ``ibpdv``, ``ibpu``, ``preval``
567
568
569
570
571
572
573

   - *Calculus/Exp, Log, powers*

         * ``exp2pow``, ``hyp2exp``, ``expexpand``, ``lin``, ``lncollect``, ``lnexpand``, ``powexpand``, ``pow2exp``

   - *Trigo*

Frederic HAN's avatar
Frederic HAN committed
574
         * ``trigexpand``, ``tsimplify``, ``tlin``, ``tcollect``,
575
576
577

         * ``halftan``, ``cos2sintan``, ``sin2costan``, ``tan2sincos``, ``tan2cossin2``, ``tan2sincos2``, ``trigcos``, ``trigsin``, ``trigtan``, ``shift_phase``

Frederic HAN's avatar
Frederic HAN committed
578
         * ``exp2trig``, ``trig2exp``
579
580
581
582
583

         * ``atrig2ln``, ``acos2asin``, ``acos2atan``, ``asin2acos``, ``asin2atan``, ``atan2acos``, ``atan2asin``

   - *Linear Algebra*

Frederic HAN's avatar
Frederic HAN committed
584
         * ``identity``, ``matrix``, ``makemat``, ``syst2mat``, ``matpow``, ``table``, ``redim``
Frederic HAN's avatar
Frederic HAN committed
585

586
587
588
589
590
591
         * ``det``,  ``det_minor``, ``rank``, ``ker``, ``image``, ``rref``, ``simplex_reduce``,

         * ``egv``, ``egvl``,  ``eigenvalues``, ``pcar``, ``pcar_hessenberg``, ``pmin``,

         * ``jordan``, ``adjoint_matrix``, ``companion``, ``hessenberg``, ``transpose``,

Frederic HAN's avatar
Frederic HAN committed
592
         * ``cholesky``, ``lll``,  ``lu``, ``qr``, ``svd``, ``a2q``, ``gauss``, ``gramschmidt``,
593
594
595
596
597
598
599
600
601
602
           ``q2a``, ``isom``, ``mkisom``


   - *Finite Fieds*

         * ``%``, ``% 0``, ``mod``, ``GF``, ``powmod``


   - *Integers*

Frederic HAN's avatar
Frederic HAN committed
603
         * ``gcd``, ``iabcuv``, ``ichinrem``, ``idivis``, ``iegcd``,
604

Frederic HAN's avatar
Frederic HAN committed
605
         * ``ifactor``, ``ifactors``, ``iquo``, ``iquorem``, ``irem``,
606
607
608
609
610

         * ``is_prime, is_pseudoprime``, ``lcm``, ``mod``, ``nextprime``, ``pa2b2``, ``prevprime``,
           ``smod``, ``euler``, ``fracmod``

   - *List*
Frederic HAN's avatar
Frederic HAN committed
611

612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
         * ``append``, ``accumulate_head_tail``, ``concat``, ``head``, ``makelist``, ``member``, ``mid``, ``revlist``, ``rotate``, ``shift``, ``size``, ``sizes``, ``sort``, ``suppress``, ``tail``

   - *Set*

         * ``intersect``, ``minus``, ``union``, ``is_element``, ``is_included``

   """
   return Pygen(s).eval()


#######################################
# A class to adjust giac configuration
#######################################
cdef class GiacSetting(Pygen):
     """
     A class to customise the Computer Algebra  System settings

     * property threads:

631
632
633
634
635
636
637
638
639
640
641
642
         Maximal number of allowed theads in giac:

         ::

            sage: from giacpy_sage import giacsettings
            sage: import os
            sage: try:
            ....:     ncpu = int(os.environ['SAGE_NUM_THREADS'])
            ....: except KeyError:
            ....:     ncpu =1
            sage: giacsettings.threads == ncpu
            True
643
644
645
646
647
648
649

     * property digits:

         Default digits number used for approximations:

         ::

650
             sage: from giacpy_sage import giacsettings,libgiac
651
652
653
654
655
656
657
658
659
660
661
662
             sage: giacsettings.digits=20;giacsettings.digits
             20
             sage: libgiac.approx('1/7')
             0.14285714285714285714
             sage: giacsettings.digits=12;

     * property sqrtflag:

         Flag to allow sqrt extractions during solve and factorizations:

         ::

663
             sage: from giacpy_sage import libgiac,giacsettings
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
             sage: giacsettings.sqrtflag=False;giacsettings.sqrtflag
             False
             sage: libgiac('x**2-2').factor()
             x^2-2
             sage: giacsettings.sqrtflag=True;
             sage: libgiac('x**2-2').factor()
             (x-sqrt(2))*(x+sqrt(2))


     property complexflag:

         Flag to allow complex number in solving equations or factorizations:

         ::

679
            sage: from giacpy_sage import libgiac,giacsettings
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
            sage: giacsettings.complexflag=False;giacsettings.complexflag
            False
            sage: libgiac('x**2+4').factor()
            x^2+4
            sage: giacsettings.complexflag=True;
            sage: libgiac('x**2+4').factor()
            (x+2*i)*(x-2*i)


     * property eval_level:

         Recursive level of substitution of variables during an evaluation:

         ::

695
             sage: from giacpy_sage import giacsettings,libgiac
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
             sage: giacsettings.eval_level=1
             sage: libgiac("purge(a):;b:=a;a:=1;b")
             "Done",a,1,a
             sage: giacsettings.eval_level=25; giacsettings.eval_level
             25
             sage: libgiac("purge(a):;b:=a;a:=1;b")
             "Done",a,1,1

     * property proba_epsilon:

         Maximum probability of a wrong answer with a probabilist algorithm.
         Set this number to 0 to disable probabilist algorithms (slower):

         ::

711
             sage: from giacpy_sage import giacsettings,libgiac
712
713
714
715
716
717
718
719
720
721
722
723
             sage: giacsettings.proba_epsilon=0;libgiac('proba_epsilon')
             0.0
             sage: giacsettings.proba_epsilon=10^(-13)
             sage: libgiac('proba_epsilon')<10^(-14)
             False

     * property epsilon:

         Value used by the epsilon2zero function:

         ::

724
             sage: from giacpy_sage import giacsettings,libgiac
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
             sage: giacsettings.epsilon=1e-10
             sage: P=libgiac('1e-11+x+5')
             sage: P==x+5
             False
             sage: (P.epsilon2zero()).simplify()
             x+5

     """

     def __repr__(self):
        return "Giac Settings"

     def _sage_doc_(self):
        return GiacSetting.__doc__


     property digits:
         r"""
         Default digits number used for approximations:

         ::

747
            sage: from giacpy_sage import giacsettings,libgiac
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
            sage: giacsettings.digits=20;giacsettings.digits
            20
            sage: libgiac.approx('1/7')
            0.14285714285714285714
            sage: giacsettings.digits=12;

         """
         def __get__(self):
           return (self.cas_setup()[6])._val


         def __set__(self,value):
           l=Pygen('cas_setup()').eval()
           pl=[ i for i in l ]
           pl[6]=value
           Pygen('cas_setup(%s)'%(pl)).eval()


     property sqrtflag:
         r"""
         Flag to allow square roots in solving equations or factorizations:
         """
         def __get__(self):
           return (self.cas_setup()[9])._val == 1


         def __set__(self,value):
           l=Pygen('cas_setup()').eval()
           pl=[ i for i in l ]
           if value:
              pl[9]=1
           else:
              pl[9]=0
           Pygen('cas_setup(%s)'%(pl)).eval()


     property complexflag:
         r"""
         Flag to allow complex number in solving equations or factorizations:

         ::

790
            sage: from giacpy_sage import libgiac,giacsettings
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
            sage: giacsettings.complexflag=False;giacsettings.complexflag
            False
            sage: libgiac('x**2+4').factor()
            x^2+4
            sage: giacsettings.complexflag=True;
            sage: libgiac('x**2+4').factor()
            (x+2*i)*(x-2*i)
         """
         def __get__(self):
           return (self.cas_setup()[2])._val == 1


         def __set__(self,value):
           l=Pygen('cas_setup()').eval()
           pl=[ i for i in l ]
           if value:
              pl[2]=1
           else:
              pl[2]=0
           Pygen('cas_setup(%s)'%(pl)).eval()


     property eval_level:
         r"""
         Recursive level of substitution of variables during an evaluation:

         ::

819
            sage: from giacpy_sage import giacsettings,libgiac
820
821
822
823
824
825
826
            sage: giacsettings.eval_level=1
            sage: libgiac("purge(a):;b:=a;a:=1;b")
            "Done",a,1,a
            sage: giacsettings.eval_level=25; giacsettings.eval_level
            25
            sage: libgiac("purge(a):;b:=a;a:=1;b")
            "Done",a,1,1
827
828
            sage: libgiac.purge('a,b')
            1,a
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849

         """
         def __get__(self):
           return (self.cas_setup()[7][3])._val


         def __set__(self,value):
           l=Pygen('cas_setup()').eval()
           pl=[ i for i in l ]
           pl[7]=[l[7][0],l[7][1],l[7][2],value]
           Pygen('cas_setup(%s)'%(pl)).eval()



     property proba_epsilon:
         r"""
         Maximum probability of a wrong answer with a probabilist algorithm.
         Set this number to 0 to disable probabilist algorithms (slower):

         ::

850
            sage: from giacpy_sage import giacsettings,libgiac
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
            sage: giacsettings.proba_epsilon=0;libgiac('proba_epsilon')
            0.0
            sage: giacsettings.proba_epsilon=10^(-13)
            sage: libgiac('proba_epsilon')<10^(-14)
            False

         """
         def __get__(self):
           return (self.cas_setup()[5][1])._double


         def __set__(self,value):
           l=Pygen('cas_setup()').eval()
           pl=[ i for i in l ]
           pl[5]=[l[5][0],value]
           Pygen('cas_setup(%s)'%(pl)).eval()


     property epsilon:
         r"""
         Value used by the epsilon2zero function:

         ::

875
            sage: from giacpy_sage import giacsettings,libgiac
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
            sage: giacsettings.epsilon=1e-10
            sage: P=libgiac('1e-11+x+5')
            sage: P==x+5
            False
            sage: (P.epsilon2zero()).simplify()
            x+5

         """
         def __get__(self):
           return (self.cas_setup()[5][0])._double


         def __set__(self,value):
           l=Pygen('cas_setup()').eval()
           pl=[ i for i in l ]
           pl[5]=[value,l[5][1]]
           Pygen('cas_setup(%s)'%(pl)).eval()

     property threads:
         """
         Maximal number of allowed theads in giac
         """
         def __get__(self):
           return (self.cas_setup()[7][0])._val

         def __set__(self,value):
           Pygen('threads:=%s'%(str(value))).eval()


########################################################
#                                                      #
#    The python class that points to a cpp giac gen    #
#                                                      #
########################################################
cdef class Pygen:

     cdef gen * gptr   #pointer to the corresponding C++ element of type giac::gen

     def __cinit__(self, s=None):

         #NB: the  != here gives problems with  the __richcmp__ function
         #if (s!=None):
         # so it's better to use isinstance
         if (isinstance(s,None.__class__)):
             # Do NOT replace with: self=GIACNULL  (cf the doctest in __repr__
             sig_on()
             self.gptr = new gen ((<Pygen>GIACNULL).gptr[0])
             sig_off()
             return

         if isinstance(s,int):       # This looks 100 faster than the str initialisation
           if PythonVersion3:
             #in python3 int and long are int
             if s.bit_length()< Pymaxint.bit_length():
                sig_on()
                self.gptr = new gen(<long long>s)
                sig_off()
             else:
                sig_on()
                self.gptr = new gen(pylongtogen(s))
                sig_off()
           else:
                sig_on()
                self.gptr = new gen(<long long>s)
                sig_off()

         #
Frederic HAN's avatar
Frederic HAN committed
943
         elif isinstance(s,Integer):       # for sage int (gmp)
944
945
946
947
948
949
950
951
             sig_on()
             if (abs(s)>Pymaxint):
                 self.gptr = new gen((<Integer>s).value)
             else:
                self.gptr = new gen(<long long>s)   #important for pow to have a int
             sig_off()


Frederic HAN's avatar
Frederic HAN committed
952
         elif isinstance(s,Rational):       # for sage rat (gmp)
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
             #self.gptr = new gen((<Pygen>(Pygen(s.numerator())/Pygen(s.denominator()))).gptr[0])
             # FIXME: it's slow
             sig_on()
             self.gptr = new gen(GIAC_rdiv(gen((<Integer>(s.numerator())).value),gen((<Integer>(s.denominator())).value)))
             sig_off()


         elif is_Matrix(s):
             s = Pygen(s.list()).list2mat(s.ncols())
             sig_on()
             self.gptr = new gen((<Pygen>s).gptr[0])
             sig_off()


         elif isinstance(s,float):
             sig_on()
             self.gptr = new gen(<double>s)
             sig_off()


         elif isinstance(s,long):
             #s=str(s)
             #self.gptr = new gen(<string>s,context_ptr)
             sig_on()
             self.gptr = new gen(pylongtogen(s))
             sig_off()


         elif isinstance(s,Pygen):
             #in the test: x,y=Pygen('x,y');((x+2*y).sin()).texpand()
             # the y are lost without this case.
             sig_on()
             self.gptr = new gen((<Pygen>s).gptr[0])
             sig_off()

         elif isinstance(s,listrange):
             sig_on()
             self.gptr = new gen(_wrap_pylist(<list>s),<short int>0)
             sig_off()

         elif isinstance(s,tuple):
             sig_on()
             self.gptr = new gen(_wrap_pylist(<tuple>s),<short int>1)
             sig_off()

         # Other types are converted with strings.
         else:
           sig_on()
For faster browsing, not all history is shown. View entire blame