Opened 12 years ago

Closed 11 years ago

# add a compose function to sage

Reported by: Owned by: was AlexGhitza major sage-4.7 misc sage-4.7.alpha2 Christopher Olah, Felix Lawrence Paul Zimmermann N/A

```def compose(f, n, a):
"""
Return f(f(...f(a)...)), where the composition occurs n times.

INPUT:
- `f` -- anything that is callable
- `n` -- a nonnegative integer
- `a` -- any input for `f`

OUTPUT:
result of composing `f` with itself `n` times and applying to `a`.

EXAMPLES::

sage: def f(x): return x^2 + 1
sage: x = var('x')
sage: compose(f, 3, x)
((x^2 + 1)^2 + 1)^2 + 1
"""
n = Integer(n)
if n <= 0: return a
a = f(a)
for i in range(n-1): a = f(a)
return a

```

To the release manager: apply only the last patch (trac_7742-add-compose-etc_v2.1.patch).

### comment:1 in reply to: ↑ description ; follow-up: ↓ 2 Changed 12 years ago by hivert

n = Integer(n)

Converting to int should be sufficient... Or the implementation below should'nt be very useful :-).

if n <= 0: return a

I'd rather sage to raise a `ValueError` if `n < 0` than returning silently a wrong value...

Florent

### comment:2 in reply to: ↑ 1 Changed 12 years ago by colah

n = Integer(n)

Converting to int should be sufficient... Or the implementation below should'nt be very useful :-).

if n <= 0: return a

I'd rather sage to raise a `ValueError` if `n < 0` than returning silently a wrong value...

Florent

We probably want to return an inverse function, if possible, when n<0... I'm also wondering if there is a generalisation of this to all real numbers.

### comment:3 Changed 12 years ago by was

```Hi,

There is a function in Mathematica that does exactly what you want and it's called Nest. I've reimplemented it in python as general effort to port some functional tools that I like to python. Please check this link if you are interested on it:

http://bitbucket.org/juliusc/functional/src/tip/functional/functional_mathematica.py

Hope this helps,
Carlos
```

### comment:4 Changed 12 years ago by was

Hi,

There is a function in Mathematica that does exactly what you want and it's called Nest. I've reimplemented it in python as general effort to port some functional tools that I like to python. Please check this link if you are interested on it:

Hope this helps, Carlos

### comment:5 Changed 12 years ago by was

We probably want to return an inverse function, if possible, when n<0..

Note that usually the inverse isn't available.

### comment:6 follow-up: ↓ 7 Changed 12 years ago by ncohen

What about doing it in log(n) compositions, like when you try to compute x12 ? :-)

Nathann

### comment:7 in reply to: ↑ 6 Changed 12 years ago by was

What about doing it in log(n) compositions, like when you try to compute x12 ? :-)

That's an interesting idea. For a Python function that will make absolutely no difference. For a symbolic callable expression it would make a noticeable difference though. So we should definitely think about it.

I think the most important thing though is to implement something with a clear interface and definition now and get it into sage. Then we can make it faster whenever we want.

### comment:8 Changed 12 years ago by ncohen

I completely agree ! That's what we are doing in graph theory : first implement, then try to think about how to optimize it :-)

I this case it should not be too different from these few lines ( if I make none of my usual mistakes )

```if n==0:
return "identity function"( no idea how it is called in Sage )

if n%2 == 1:
ff = compose(f, (n-1)/2, a)
return f(ff(ff(a)))
else:
ff = compose(f, n/2, a)
return ff(ff(a))
```

But perhaps the function "compose" should be a method of these symbolic callable expressions you mentionned, which would be called by this "compose" function : it may be interesting to compute symbolically the nth composed of a function without evaluating it :-)

### comment:9 Changed 12 years ago by colah

I ran into some problems producing the patch. On mhansen's suggestion, I added {{from sage.rings.integer import Integer}} at the beginning of the function instead of the beginning of misc.py and it worked... But I don't know why and am concerned that it adds unnecessary overhead...

### comment:10 Changed 12 years ago by colah

• Status changed from new to needs_review

### comment:11 Changed 12 years ago by zimmerma

• Status changed from needs_review to needs_work

I have two comments: (1) the patch contains twice the new function, and (2) it would be nice to be able to write `compose(f,x,n)` for a symbolic variable n (or any expression). Only when n is an integer would the "nesting" be evaluated. Comment (1) needs more work, comment (2) could be a separate ticket, but if easy it would be nice to do it.

### comment:12 follow-up: ↓ 16 Changed 12 years ago by flawrence

• Status changed from needs_work to needs_review

The new patch addresses many of the above concerns (other than symbolic n), and splits the functionality into three functions:

• compose(f,g) which composes functions such that compose(f,g)(x) = f(g(x))
• self_compose(f, n) which composes a function with itself such that self_compose(f,n)(x) = f(...(f(x))...)
• nest(f, n, a) which effectively self_composes and evaluates (but does it slightly faster than a self_compose)

Questions:

• compose(f,g)(x) = f(g(x)) - is this the most common convention among mathematicians? or should the RHS be g(f(x))?
• Is (f, n, a) or (f, a, n) the best order of parameters for nest?
• Is there somewhere more suitable than misc/misc.py for these functions?

### comment:13 Changed 12 years ago by jrp

• Status changed from needs_review to needs_work

What version of Sage are you targeting? I tried applying both of these (individually) to 4.6.0 and to 4.5.3 and couldn't apply them.

### Changed 12 years ago by flawrence

adds compose, self_compose and nest functions. actually replaces previous patch this time

### comment:14 Changed 12 years ago by flawrence

• Status changed from needs_work to needs_review

Sorry, the second patch accidentally relied on the first (which needed rebasing). I've fixed it and the second patch should now apply and work for both 4.6.0 and 4.6.1alpha0.

### comment:15 Changed 12 years ago by zimmerma

• Authors set to Christopher Olah, Felix Lawrence
• Reviewers set to Paul Zimmermann
• Status changed from needs_review to needs_work

This needs more work, since the function `self_compose` fails for n=0:

```sage: function('f');
sage: h=self_compose(f,0)
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)

/localdisk/tmp/sage-4.6/<ipython console> in <module>()

/localdisk/tmp/sage-4.6/local/lib/python2.6/site-packages/sage/misc/misc.pyc in self_compose(f, n)
890     # all necessary powers of two (f, f^2, f^4), storing them in a list
891     fs = [f]
--> 892     for i in xrange(int(floor(log(n, 2)))):
893         fs.append(compose(fs[i], fs[i]))
894

ValueError: math domain error
```

Paul

### comment:16 in reply to: ↑ 12 Changed 12 years ago by zimmerma

• compose(f,g)(x) = f(g(x)) - is this the most common convention among mathematicians? or should the RHS be g(f(x))?

compose(f,g) seems good to me.

• Is (f, n, a) or (f, a, n) the best order of parameters for nest?

I prefer (f,n,a) since this is coherent with `self_compose(f,n)`.

• Is there somewhere more suitable than misc/misc.py for these functions?

I don't know.

Also if a symbolic n is not allowed (yet) in `self_compose`, some type-checking should be done:

```sage: function('f')
sage: var('m')
sage: self_compose(f,m)(x)
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)

/localdisk/tmp/sage-4.6/<ipython console> in <module>()

/localdisk/tmp/sage-4.6/local/lib/python2.6/site-packages/sage/misc/misc.pyc in self_compose(f, n)
890     # all necessary powers of two (f, f^2, f^4), storing them in a list
891     fs = [f]
--> 892     for i in xrange(int(floor(log(n, 2)))):
893         fs.append(compose(fs[i], fs[i]))
894

/localdisk/tmp/sage-4.6/local/lib/python2.6/site-packages/sage/symbolic/expression.so in sage.symbolic.expression.Expression.__float__ (sage/symbolic/expression.cpp:5398)()

TypeError: unable to simplify to float approximation
```

Paul

### comment:17 follow-up: ↓ 18 Changed 12 years ago by jrp

At the moment the powers-of-two speed up isn't helping even in the symbolic case, because the implementation is to return lambda x: f(g(x)). If both inputs to compose are symbolic, then I think we want to just return f(g(x)) instead.

### comment:18 in reply to: ↑ 17 ; follow-ups: ↓ 19 ↓ 21 Changed 12 years ago by jrp

At the moment the powers-of-two speed up isn't helping even in the symbolic case, because the implementation is to return lambda x: f(g(x)). If both inputs to compose are symbolic, then I think we want to just return f(g(x)) instead.

Given that symbolic functions have named variables, is this the correct behavior?

sage: f(x) = x+1 sage: g(y) = y+1 sage: compose(f,g) y + 2 sage: compose(g,f) x + 2

### comment:19 in reply to: ↑ 18 ; follow-ups: ↓ 20 ↓ 27 Changed 12 years ago by flawrence

At the moment the powers-of-two speed up isn't helping even in the symbolic case, because the implementation is to return lambda x: f(g(x)). If both inputs to compose are symbolic, then I think we want to just return f(g(x)) instead.

Given that symbolic functions have named variables, is this the correct behavior?

sage: f(x) = x+1 sage: g(y) = y+1 sage: compose(f,g) y + 2 sage: compose(g,f) x + 2

I don't actually know anything about symbolic functions in Sage. As such, all doctests and examples in the patch are for regular python functions not symbolic functions. The compose functions in the patch use lambda functions, which may or may not break symbolic functions.

With normal python functions (such as those in the doctest), the powers-of-two provides a major speedup for large n over the other techniques I tried, both for the function composition and for the evaluation of the resulting function (benchmarked using "timeit"). In fact, without a powers-of-two-type algorithm, the self_compose function regularly hit the recursion limit and crashed for large n (above 500 or 1000). But I don't know if anyone is interested in large n cases.

### comment:20 in reply to: ↑ 19 Changed 12 years ago by jrp

With normal python functions (such as those in the doctest), the powers-of-two provides a major speedup for large n over the other techniques I tried, both for the function composition and for the evaluation of the resulting function (benchmarked using "timeit"). In fact, without a powers-of-two-type algorithm, the self_compose function regularly hit the recursion limit and crashed for large n (above 500 or 1000).

I ran the following:

```sage: def f(x):
....:     return x + 1
....:
sage: g = self_compose(f,10000)
sage: g
<function <lambda> at 0xb0aa80c>
sage: timeit('g(0)')
25 loops, best of 3: 9.41 ms per loop
sage: timeit('nest(f,10000,0)')
125 loops, best of 3: 6.35 ms per loop
```

Since nest just loops through xrange and doesn't use powers of two, maybe the slowdown before was due to the recursion. But in any case, it definitely helps with the symbolic ones.

The compose functions in the patch use lambda functions, which may or may not break symbolic functions.

I'm attaching a patch to check for symbolic functions and (hopefully) do the right thing. I also make self_compose(f,0) return the identity.

### comment:21 in reply to: ↑ 18 Changed 12 years ago by zimmerma

Given that symbolic functions have named variables, is this the correct behavior?

```sage: f(x) = x+1
sage: g(y) = y+1
sage: compose(f,g)
y + 2
sage: compose(g,f)
x + 2
```

[please use { { { and } } } to quote Sage examples]

I can't reproduce this. I get with the patch from Felix:

```sage: f(x) = x+1
sage: g(y) = y+1
sage: compose(f,g)
<function <lambda> at 0x1787c80>
sage: compose(f,g)(x)
x + 2
```

which seems ok to me.

Paul

### Changed 12 years ago by jrp

Special treatment of symbolic functions.

### comment:22 Changed 12 years ago by jrp

• Status changed from needs_work to needs_review

I realized that what really matters is the type of the first function in the composition - i.e., the type of g in compose(f,g).

### comment:23 Changed 12 years ago by zimmerma

jrp, you didn't answer comment 21. What was wrong with Felix patch?

Paul

### comment:24 Changed 12 years ago by jrp

Paul - I wanted it to automatically do the right thing when g is a symbolic function. I agree that this isn't very important since we can already just call compose(f,g)(x). More important is that with symbolics, now the powers-of-two optimization will work. (I don't think it did before, but will think more on this).

### comment:25 follow-up: ↓ 26 Changed 12 years ago by zimmerma

• Status changed from needs_review to needs_work

I am not very happy with the fact that when f is an expression, the code assumes it has exactly one symbolic variable:

```sage: f = 2
sage: self_compose(f,0)
<function <lambda> at 0x41031b8>
sage: self_compose(f,1)
2
sage: self_compose(f,2)
<function <lambda> at 0x7fa8ac2b0d70>
```

or

```sage: f = x+y
sage: self_compose(f,0)
x
sage: self_compose(f,1)
x + y
sage: self_compose(f,2)
x + 2*y
```

Consider also:

```sage: f = pi+1
sage: self_compose(f,0)
---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)

/usr/local/sage-4.6/sage/local/lib/python2.6/site-packages/sage/misc/misc.pyc in self_compose(f, n)
911     if n == 0:
912         if type(f) == sage.symbolic.expression.Expression:
--> 913             return f.parent()(f.variables()[0])
914         else:
915             return lambda x: x

IndexError: tuple index out of range
```

I would prefer that symbolic expressions such as `f=x+1` are not allowed, and we should write instead `f(x)=x+1` or `f = lambda x : x+1`.

Paul

### comment:26 in reply to: ↑ 25 Changed 11 years ago by jrp

I would prefer that symbolic expressions such as `f=x+1` are not allowed, and we should write instead `f(x)=x+1` or `f = lambda x : x+1`.

I now think it's best to leave compose,self_compose, and nest agnostic to the type of f; the user can supply anything callable and we won't attempt to turn it back into a symbolic function.

Then as a separate issue, symbolic functions can grow an f.self_compose(n) function that plays nice with types.

### comment:27 in reply to: ↑ 19 ; follow-up: ↓ 28 Changed 11 years ago by jrp

With normal python functions (such as those in the doctest), the powers-of-two provides a major speedup for large n over the other techniques I tried, both for the function composition and for the evaluation of the resulting function (benchmarked using "timeit"). In fact, without a powers-of-two-type algorithm, the self_compose function regularly hit the recursion limit and crashed for large n (above 500 or 1000). But I don't know if anyone is interested in large n cases.

With normal python functions, I don't think there is any gain from powers of two.

```sage: f = lambda x: x + 1
sage: g = x + 1
sage: f1 = self_compose(f,10000)
sage: g1 = self_compose(g,10000)
sage: timeit('self_compose(f,10000)')
625 loops, best of 3: 66.6 µs per loop
sage: timeit('self_compose(g,10000)')
625 loops, best of 3: 66.9 µs per loop
sage: timeit('f1(0)')
25 loops, best of 3: 9.32 ms per loop
sage: timeit('g1(0)')
5 loops, best of 3: 375 ms per loop
sage: timeit('nest(f,10000,0)')
125 loops, best of 3: 6.46 ms per loop
sage: timeit('nest(g,10000,0)')
5 loops, best of 3: 379 ms per loop
sage: def new_self_compose(f,n):
....:     return lambda a: nest(f,n,a)
....:
sage: timeit('new_self_compose(f,10000)')
625 loops, best of 3: 1.17 µs per loop
sage: timeit('new_self_compose(g,10000)')
625 loops, best of 3: 2.68 µs per loop
```

Since self_compose is implemented via `lambda x: f(g(x))`, we have to eventually pass the argument through `f` the correct number of times.

Let's move the powers of two, and the symbolics handling, into a method of symbolic functions. Then the compose tools for python functions can have fast, dirt-simple implementations.

### Changed 11 years ago by flawrence

adds compose, self_compose and nest functions. self_compose is a wrapper for nest. a self-contained patch.

### comment:28 in reply to: ↑ 27 Changed 11 years ago by flawrence

• Component changed from basic arithmetic to misc
• Status changed from needs_work to needs_review
• Type changed from defect to enhancement

With normal python functions, I don't think there is any gain from powers of two.

(snip)

Since self_compose is implemented via `lambda x: f(g(x))`, we have to eventually pass the argument through `f` the correct number of times.

Yes, `nest` is always faster than `self_compose` if you actually want to evaluate the function. However the functions do different things: `nest` isn't really a compose function, since it returns fn(x) not fn. `self_compose` is more useful than `nest` if you want to forget about f and only want to apply `f^n`. In hindsight, `self_compose` should actually look like this, which I agree would be faster to construct and to evaluate than the powers of two method:

```def self_compose(f,n):
return lambda x: nest(f,n,x)
```

Let's move the powers of two, and the symbolics handling, into a method of symbolic functions. Then the compose tools for python functions can have fast, dirt-simple implementations.

I agree that we should scrap the powers of two, and not go out of our way to support symbolics in this ticket. I'm uploading a patch that makes the change detailed above. I'm not going to have any more time to work on this for at least a month, so someone else will have to finish the ticket if further changes are required.

The new patch works for n = 0, but still does not check whether n is symbolic. I don't think this is the worst thing in the world: the documentation does state that n should be a nonnegative integer, and if it is symbolic then the error raised is appropriate IMO:

```sage: _ = var('x y')
sage: nest(sin,x,y)
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)

/Applications/sage/devel/sage-dev/<ipython console> in <module>()

/Applications/sage/local/lib/python2.6/site-packages/sage/misc/misc.pyc in nest(f, n, x)
939
940     """
--> 941     for i in xrange(n):
942         x = f(x)
943     return x

/Applications/sage/local/lib/python2.6/site-packages/sage/symbolic/expression.so in sage.symbolic.expression.Expression.__int__ (sage/symbolic/expression.cpp:4375)()

/Applications/sage/local/lib/python2.6/site-packages/sage/symbolic/expression.so in sage.symbolic.expression.Expression.n (sage/symbolic/expression.cpp:17185)()

TypeError: cannot evaluate symbolic expression numerically
```

### comment:29 Changed 11 years ago by zimmerma

• Status changed from needs_review to positive_review

### comment:30 Changed 11 years ago by jdemeyer

• Milestone changed from sage-4.6.1 to sage-4.6.2

### comment:31 Changed 11 years ago by jdemeyer

• Status changed from positive_review to needs_info

Please state clearly in the ticket description which patches have to be applied.

### comment:32 Changed 11 years ago by zimmerma

• Description modified (diff)
• Status changed from needs_info to needs_review

### comment:33 Changed 11 years ago by zimmerma

• Status changed from needs_review to positive_review

### comment:34 Changed 11 years ago by jdemeyer

• Status changed from positive_review to needs_work
• Work issues set to rebase

This needs to be rebased to sage-4.6.2.alpha2. You might also want to check for conflicts with #8456 (so it might be safer to base your patch to sage-4.6.2.alpha2 + #8456).

### comment:35 follow-up: ↓ 36 Changed 11 years ago by zimmerma

does anybody know why this patch with positive review was not included in 4.6.1?

Paul

### comment:36 in reply to: ↑ 35 Changed 11 years ago by jdemeyer

does anybody know why this patch with positive review was not included in 4.6.1?

I know because I was the one making that decision.

Most likely, the 4.6.1 release cycle was in bug-fix-only status at the time when this patch got positive_review. Unfortunately, 4.6.1 stayed in that status for quite a while (much longer than usual Sage releases) because there were a lot of changes in the Sage build scripts giving rise to many system-specific issues which had to be tracked down.

### comment:37 Changed 11 years ago by jdemeyer

• Milestone changed from sage-4.6.2 to sage-4.7

Anybody wants to look at this?...

### Changed 11 years ago by flawrence

Rebased for 4.6.2.rc0

### comment:38 Changed 11 years ago by flawrence

• Description modified (diff)
• Status changed from needs_work to needs_review
• Work issues rebase deleted

Only change since the positive review was updating the change to sage/misc/all.py

### comment:39 follow-up: ↓ 40 Changed 11 years ago by zimmerma

• Status changed from needs_review to needs_work
• Work issues set to check argument type

a minor point: it would be good to check that n is indeed a nonnegative integer in `self_compose` and `nest`, otherwise we can get strange results:

```sage: self_compose(sin, -3)(x)
x
sage: self_compose(sin, x)(x)
...
TypeError: cannot evaluate symbolic expression numerically
```

Paul

### Changed 11 years ago by flawrence

checks that n is a nonnegative integer

### comment:40 in reply to: ↑ 39 Changed 11 years ago by flawrence

• Description modified (diff)
• Status changed from needs_work to needs_review

a minor point: it would be good to check that n is indeed a nonnegative integer in `self_compose` and `nest`, otherwise we can get strange results

Done:

```sage: self_compose(sin, -3)
...
ValueError: n must be a nonnegative integer, not -3.
sage: self_compose(sin, x)
...
TypeError: n (=x) must be of type (<type 'int'>, <type 'long'>, <type 'sage.rings.integer.Integer'>).
```

### comment:41 Changed 11 years ago by zimmerma

• Status changed from needs_review to positive_review
• Work issues check argument type deleted

positive review, thank you Christoph and Felix for that nice work!

Paul

### comment:42 Changed 11 years ago by zimmerma

and thank you for the extra rebase work, I hope this time the patch will be included soon...

Paul

### comment:43 Changed 11 years ago by jdemeyer

• Merged in set to sage-4.7.alpha2
• Resolution set to fixed
• Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.