
=============
Giacpy
=============

:Name: giacpy
:Summary: A Cython frontend to the c++ library giac. (Computer Algebra System)
:Author: Frederic Han
:Author-email: han@math.jussieu.fr
:Copyright: 2012 Frederic Han
:License:  GPL v2 or above
:Home-page: http://www.math.jussieu.fr/~han/xcas/giacpy/



Access from python to the Computer Algebra System giac via libgiac

------------
Introduction
------------


This is the sage version of giacpy.

    * 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:

    * The class Pygen is the main tool to interact from python/sage with the c++ library giac via cython. 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.

    * 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.

    ::

        sage: from sage.libs.giac import libgiac,giacsettings
        // Giac share root-directory:/usr/local/sage-5.10/local/share/giac/
        // Using keyword file /usr/local/sage-5.10/local/share/giac/doc/fr/keywords
        // Giac share root-directory:/usr/local/sage-5.10/local/share/giac/
        Help file /usr/local/sage-5.10/local/share/giac/doc/fr/aide_cas not found
        Added 0 synonyms
        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)

    * To obtain more hints on giacpy consider the help of the libgiac function.

    ::

        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


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 in several languages (en, fr, el) at:

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


EXAMPLES:

    ::

        sage: from sage.libs.giac import *
        sage: x=libgiac('x');f=1/(2+sin(5*x))
        sage: f.int()
        2/5/sqrt(3)*(atan((2*tan(5*x/2)+1)/sqrt(3))+pi*floor(5*x/2/pi+1/2))
        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


Graphics 2D Output via qcas (the qt frontend to giac) is disabled in the sage version of giacpy.

This is the sage version of giacpy: an interface to be able to use from Python the Giac features. 

   * 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: 
        PARI, NTL (arithmetic), CoCoA (Groebner basis), 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.


   * giacpy is an interface to this library. It is built with cython. Graphic output is obtained with qcas by Loic Lecoq: http://git.tuxfamily.org/?p=qcas/qcas.git


-------------------------------
Short Usage of the sage version
-------------------------------

Example::

    >>> from giacpy import libgiac
    >>> x,y,z=libgiac('x,y,z')
    >>> f=(x+y+z+1)**15+1
    >>> g=(f*(f+1)).normal()
    >>> print g.nops()
    >>> print g.factor().nops()
    >>> f.diff()
    >>> f.diff(z)

Help::

    >>> help("giacpy")
    >>> from giacpy import *
    >>> normal.help() ; # to have help on some giac keyword
    >>> solve.htmlhelp('fr') ; # to have detailled help on some giac keyword
    >>> htmlhelp()  ; # to have help the global help pages.


    * Graphics 2D Output: (cf. help('giacpy') for examples)
     If your version of giacpy has qt support, you can send graphics to qcas with the .qcas() method. For experimental interactive geometry see: help(qcas)


   * For binaries of giacpy: http://www.math.jussieu.fr/~han/xcas/giacpy/



-----------------------------------
Short Tutorial on the libgiac function
-----------------------------------

    The libgiac instruction is function evaluate a python object with the giac library.

    * It creates in python a Pygen element and evaluate it with the c++ library giac:

 
     >>> from giacpy import libgiac
     >>> x,y=libgiac('x,y');type(x)
     <type 'giacpy.Pygen'>
     >>> (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.

    >>> from giacpy import *  
    >>> x,y,z=libgiac('x,y,z')
    >>> 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)


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

    >>> type(y);giac('y:=1');y
    <type 'giacpy.Pygen'>
    1
    y

   * There are some natural coercion to Pygen elements:

    >>> pi>3.14 ; pi >3.15 ; giac(3)==3
    True
    False
    True


Lists of Pygen and Giac lists:
------------------------------

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

    >>> l1=giac(range(10)); l2=[1/(i**2+1) for i in l1]
    >>> 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.

    >>> l1+l1
    [0,2,4,6,8,10,12,14,16,18]
    >>> 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.
 
    >>> V=[ [x[i]**j for i in range(9)] for j in range(9)]
    >>> giac(V).dim()
    [9,9]
    >>> det_minor(V).factor()
    (x[7]-(x[8]))*(x[6]-(x[8]))*(x[6]-(x[7]))*(x[5]-(x[8]))*(x[5]-(x[7]))*(x[5]-(x[6]))*(x[4]-(x[8]))*(x[4]-(x[7]))*(x[4]-(x[6]))*(x[4]-(x[5]))*(x[3]-(x[8]))*(x[3]-(x[7]))*(x[3]-(x[6]))*(x[3]-(x[5]))*(x[3]-(x[4]))*(x[2]-(x[8]))*(x[2]-(x[7]))*(x[2]-(x[6]))*(x[2]-(x[5]))*(x[2]-(x[4]))*(x[2]-(x[3]))*(x[1]-(x[8]))*(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[8]))*(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 %

    >>> V=ranm(5,5) % 2;
    >>> ker(V).rowdim()+V.rank()
    5
    >>> a=giac(7)%3;a;a%0;7%3
    1 % 3
    1
    1

   Do not confuse with the full python integers:
   
    >>> type(7%3);type(a)
    <type 'int'>
    <type 'giacpy.Pygen'>

Syntaxes with reserved or unknown Python symbols:
-------------------------------------------------

   * In general equations needs symbols such as = < > or that have another meaning in Python. So those objects must be quoted.  

    >>> from giacpy import *
    >>> x=giac('x')
    >>> (1+2*sin(3*x)).solve(x)
    list[-pi/3/6,7*pi/18]

    >>> solve('sin(3*x)>2*sin(x)',x)
    Traceback (most recent call last):
    ...
    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:

    >>> assume('x>-pi && x<pi')
    x
    >>> 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

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

 
   * Same problems with the ..

    >>> from giacpy import *
    >>> x=giac('x')
    >>> f=1/(5+cos(4*x));f.int(x)
    1/2/(2*sqrt(6))*(atan(2*tan(4*x/2)/sqrt(6))+pi*floor(4*x/2/pi+1/2))
    >>> fMax(f,'x=-0..pi').simplify()
    pi/4,3*pi/4
    >>> fMax.help()
    "Returns the abscissa of the maximum of the expression.
    Expr,[Var]
    fMax(-x^2+2*x+1,x)
    fMin"
    >>> sum(1/(1+x**2),'x=0..infinity').simplify()
    (pi*exp(pi)**2+pi+exp(pi)**2-1)/(2*exp(pi)**2-2)





---------
Changelog
---------


   * Version 0.2: 
      - Add a comparison function to Pygen. (with coersion) 
      - Add a basic definition for most giac functions.
      - Add some help.

   * Version 0.2.1: 
      - Add __neg__ and __pos__ support for Pygen. (Ex: -pi)
      - Change __repr__ to hide too long outputs.
      - Make ** be the default printing for powers in giac.

   * Version 0.2.2: 
      - Change Pygen() to Pygen('NULL'). (Ex: rand())
      - Add direct acces to the python double value of a Pygen: a._double
      - Add  conversion to giac modulars via the operator %
      - Add  ctrl-c support during list initialisation and iteration
      - Modification of __getitem__ to allow formal variables with indexes.
      - Add htmlhelp method for Pygen objects.
      - Improve the giac initialisation of Python long integers. (basic Horner method instead of strings)
      - Improve  help(giac) and doctests
      - Add support for the slice notation with giac lists

   * Version 0.2.3: 
      - Fix Pygen() None initialisation. Add crash test and improve speed in _wrap_gen
      - Add a small Makefile
      - Add a GiacSettings class with some frontends to the cas settings.
      - Add French keywords

   * Version 0.2.4:
      - Update giac 1.1 keywords.

   * Version 0.3:
      - Add a qt output for 2d graphics via qcas.
      - Fixes for giac 1.1

   * Version 0.4:
      - Fixes for Python 3 compatibility
      - Qt/qcas can be disabled at compilation. (cf setup.py)
      - 0.4.1:
	  + add some giac keywords.
      	  + add proba_epsilon in GiacSetting. 
	  + test if the html doc is present locally, otherwise open the web doc.
