Introducing operators

[GoogleCodeExporter]

Defining functions and operations is well settled down
in sympycore by now. Now we need a standard way to define operators.

Definitions:
- domain is a set of scalar values or a direct product
of such sets.
- functions (operations Add, Mul, Union, etc are also functions
in sympycore) are mappings that map domain elements to domain elements.
- operators are mappings that map functions to functions.

Applying a function to arguments either produces scalar value
or unevaluated function value. Similarly, applying an operator
to a function, should either produce a function 
(e.g. Lambda instance) or unevaluated operator.

Unapplied functions can be used in arithmetic operations. Similarly,
unapplied operators should be usable in arithmetic operations.

Different from functions, most interesting operators are additive
or/and homogenous. Should linear operators expand sums by default or not?

Applying a function always tries to evaluate. In sympycore only few
functions (Pow, for instance) have a `evaluated=True` kw argument
that prevent evaluation. The `evaluated` argument is used internally,
though.
In the case of operators, the unevaluated forms are more often used
(eg because for lots of problems algorithms are too limited to be
able to evaluate operator results). 
We need a unique way to support unevaluated forms of operators
as well as functions.

Features needed:
- differential operators, both unapplied and unevaluated applied forms.
Examples:
  x = Symbol('x')
  f = FunctionType('f') - function instance
  D1 = FD(1)             - derivative with respect to the first argument
  Dx = D(x)              - derivative with respect to x

  D1(D1) -> FD((1,2))    - 2nd derivative with respect to the first argument
  Dx(Dx) -> D((x,2))     - 2nd derivative with respect to x

  D1(f)                  - derivative function of f with respect to 1st arg
  Dx(f(x))               - derivative of f(x) with respect to x
  D1(f)(x) == Dx(f(x))   - must hold
  D1 + D1 -> 2*D1
  D1 * D1 -> D1**2 such that (D1**2)(f)(x) == Dx(f(x))**2

To represent Dx(f(x)), calculus uses Derivative(f(x), x).

- integral operators, both unapplied and unevaluated applied forms.

- sums, products, limits - not sure if they should be operators
or functions or separate operations.

Original issue reported on code.google.com by pearu.peterson on 6 Dec 2007 at 11:06

[GoogleCodeExporter]

Maybe get rid of Derivative and only use D?

Original comment by ondrej.c...@gmail.com on 6 Dec 2007 at 1:29

[GoogleCodeExporter]

> `evaluated=True`

Ondrej disagreed with me the last time here, but I'd like to restate my 
position: the
default behavior of calling Integral, Sum, etc should be to attempt to fully
evaluate, and you should have to specify evaluate=False to bypass that. In other
words, they should behave like Add, Mul, etc (for which we might also want to 
add
support for the evaluate flag).

> Should linear operators expand sums by default or not?

Probably no. One reason is that you might get infinite terms instead of an 
intended
cancellation, for example if you try to integrate 1/x - 1/(x+1) from 1 to 
infinity
and it gets split up prematurely. I would also say that, in general, f(x)+f(y) 
is a
more complex expression than f(x+y).

> - sums, products, limits - not sure if they should be operators
> or functions or separate operations.

They could be: you can define "primitive functions" for sums and products just 
as for
integrals. However, the notion is less common. The most practical reason to 
derive
Sum and Product from the same class as Integral is to share the code for 
handling the
indices, but that functionality could also be derived from a more abstract 
class or
mixin.

> Maybe get rid of Derivative and only use D?

At least in the default case, it would be very nice if you could just type 
D(sin(x), x).

Original comment by fredrik....@gmail.com on 6 Dec 2007 at 2:10

[GoogleCodeExporter]

> Ondrej disagreed with me the last time here, but I'd like to restate my 
position: the
> default behavior of calling Integral, Sum, etc should be to attempt to fully
> evaluate, and you should have to specify evaluate=False to bypass that. In 
other
> words, they should behave like Add, Mul, etc (for which we might also want to 
add
> support for the evaluate flag).

My position was based on a prior discussion and consensus in SymPy (search in 
issues
for that).

Generally I think that a construction of expression, using classes like Mul, 
Integral
etc. needs to be fast and only a quick evaluation sould occur. Integral.doit() 
is not
quick, so it imho shouldn't happen by default. If you want to integrate, use 
integrate().

Original comment by ondrej.c...@gmail.com on 6 Dec 2007 at 3:12

[GoogleCodeExporter]

> Generally I think that a construction of expression, using classes like Mul,
Integral etc. needs to be fast and only a quick evaluation sould occur.

Either you want an evaluated integral or you don't. I would bet that most 
people most
of the time want their integrals evaluated if possible (e.g. when in interactive
mode), and this should be the default behavior. I don't really see how speed 
enters
the equation; those users who need to process a lot of hard or unpredictable
integrals are probably advanced enough to figure out that they should ask for
suppressed evaluation and this should not affect the default behavior.

A possible alternative would be to specify an option that selects algorithm. The
default could be to try only polynomial integration which currently is fast, 
e.g.
0.005 seconds to integrate a cubic polynomial in sympycore, combined with a 
quick
table lookup for common functions like cos(a*x+b). Calling Integral with 
risch=True
could use the Risch algorithm when the fast method fails.

(However, I suspect that if we do this, users will eventually start complaining 
about
risch=True not being the default :-)

> If you want to integrate, use integrate().

Having both Integrate() and Integral() is redundant and confusing. Also, 
methods like
.doit() are ugly IMO.

Original comment by fredrik....@gmail.com on 6 Dec 2007 at 4:41

[GoogleCodeExporter]

Not Integrate() but a function integrate(). Those users who want to calculate 
the
integral immediately will use integrate. That's the default in SymPy and
Maple/Mathematica too imho. 

It's like Derivative and diff. Or Limit and limit. The upper case should 
represent
the operation (without doing it), the lower case (function) will do it.

This was the consensus in SymPy so far. But I am open to discuss it again and 
make a
different consensus. But I, speaking for myself, prefer the current behavior.

Original comment by ondrej.c...@gmail.com on 6 Dec 2007 at 4:59

[GoogleCodeExporter]

In Mathematica there is only one function: Integrate, D, Limit.

Original comment by fredrik....@gmail.com on 6 Dec 2007 at 5:15

[GoogleCodeExporter]

In sympycore all functions are CamelCase (this is the difference with sympy)
and they should return expressions in canonical form so that any algorithms
applied to expressions later one, can assume that expressions are in certain 
form.

Since differentation algorithm is simple and in many cases doing
full differentation produces simpler expressions (representing
canonical form) and therefore there is almost no need for
unevaluated derivative operator. The unevaluated derivative operator
is only needed for functions that do not define fdiff? method(s).

This situation is completely different for integrals and limits.
The canonical form for integral expressions, for instance, would be
in most cases unevaluated form (may be only apply homogenous property).
Only integrals of polynomials should be evaluated for obtaining
the corresponding canonical form.

Similar to the problem of simplifying expressions with methods like
expand, factor, collect, etc (and the corresponding functions Expand, Factor,
..), evaluating difficult integrals should be probably carried out via
special methods (and the corresponding functions) using Risch algorithm
or partial integration or some other trick that usage cannot be decided
algorithmically.

I think we need a visitor concept to walk through expressions that
applies procedures that are attached to given visitor. This would
allow us to apply different algorithms to the same expression. Note
that using keyword arguments here would not be sufficent as we cannot
attach all possible algorithm keywords/methods to all classes.
Eg if one has
  expr = 1 + Integrate(<integrand>)
then forcing to Integrate to use, say, risch algorithm, one cannot
do expr.try_risch(..) as Add does not define try_risch method (but Integrate 
does).
With visitor concept, one could do

  expr.visit(apply_risch)

where

def apply_risch(e):
  if e.is_Integrate:
    return e.try_risch()
  return

and visit() methods walks recursively through expressions and
replaces any part of the expression with the result of apply_risch
when it returns not None.

Original comment by pearu.peterson on 6 Dec 2007 at 6:50

[GoogleCodeExporter]

I think having unevaluated derivative is needed, because the expression tends to
grow, i.e.:

In [5]: (x**10*cos(x)).diff(x, 6)
Out[5]: 
   10                  5                 6              9                8    
- x  *cos(x) - 181440*x *sin(x) - 75600*x *cos(x) - 60*x *sin(x) + 1350*x *cos

             7                  4       
(x) + 14400*x *sin(x) + 151200*x *cos(x)

So being able to manipulate the expression In [5] is essential imho. The same 
about
integrals and limits. I think we just need both evaluated and unevaluated.

Visitor concept is nice, if it simplifies code.

Original comment by ondrej.c...@gmail.com on 6 Dec 2007 at 9:20

[GoogleCodeExporter]

sympycore now implements a base class `Operator` for operators
in sympy/arithmetic/operator.py. It differs from `Function`
class by the fact that Operator does not have methods
related to differentation and integration; and Operator
constructor expects one or more of the following arguments:

 - an expression that will be mapped to tuple. For example,
   FD(1) -> FD((1,1)) - derivative with respect to 1st argument
   D(x) -> D((x,1)) - derivative with respect to x variable
   AD(x) -> AD((x,1)) - antiderivative with respect to x variable

 - a tuple of an expressions, see also the examples above.
   For example,
     AD((x,a,b)) - defined integral operator with respect to
                   x variable over the range (a,b).

Currently the only defined operators are FD, D, AD. However,
Operator is a suitable base class also to Sum and Product
or any operator containing index information. For example:

  S((i,1,n)) - a sum operator with index i running from 1 to n
  S((i,1,n),(j,i,n)) - a double sum operator where the second
                sum depends on the index of the first sum
  P((i,1,n)) - a product operator

When applying operator to an expression, it either computes
the result or returns an unevaluated result as an instance
of the corresponding class in the following table:

  Operator  Cls of uneval.result    Examples
  FD        FDerivative           FD(1)(f) -> FDerivative(f,(1,1))
  D         Derivative            D(x)(f(x)) -> Derivative(f(x),(x,1))
  AD        Integral              AD(x)(f(x)) -> Integral(f(x),(x,1))
                                  AD((x,a,b))(f(x)) -> Integral(f(x),(x,a,b))
  S         Sum                   S((i,1,n))(a(i)) -> Sum(a(i),(i,1,n))
  P         Product               P((i,1,n))(a(i)) -> Product(a(i),(i,1,n))
  L         Limit                 L((x,oo))(f(x)) -> Limit(f(x),(x,oo))

The canonize methods of FDerivative, Derivative, etc functions
construct relevant operators and apply them to given expressions.
To skip this (eg when operator cannot evaluate the expression),
that is, to obtain an unevaluated instance,
an `is_canonical=False` kw argument is introduced to all 
sympycore functions.

Original comment by pearu.peterson on 9 Dec 2007 at 3:19

[GoogleCodeExporter]

Btw, I am open to changing `is_canonical=False` kw to something more appropiate.
Eg `evaluated=False`, `skip_canonize=False`. The last one might be closest
to what action is actually taken when True.

Original comment by pearu.peterson on 9 Dec 2007 at 3:23