## Stream: new members

### Topic: Can't define Chebyshev polynomials

#### Abhimanyu Pallavi Sudhir (Jan 12 2019 at 17:34):

Why doesn't this work?

def chebyshev : ℕ → polynomial ℝ
| 0 := polynomial.C 1
| 1 := polynomial.X
| (n + 2) := 2 * polynomial.X * chebyshev (n + 1) - chebyshev n


Lean tells me the definition relies on classical.choice. Indeed putting noncomputable at the front fixes things, but why is it noncomputable?

#### Kevin Buzzard (Jan 12 2019 at 17:52):

The definition of polynomial X is "functions from nat to X which vanish outside a finite set"; "vanish" means "equals zero", and equality is undecidable in the reals. This might be something to do with it. Try polynomial int and see if it fixes things!

#### Kevin Buzzard (Jan 12 2019 at 17:54):

NB import analysis.polynomial to make it work, but indeed int fixes things, and is probably the "correct" thing to do.

#### Kevin Buzzard (Jan 12 2019 at 18:00):

import analysis.polynomial

def chebyshev : ℕ → polynomial ℤ
| 0 := polynomial.C 1
| 1 := polynomial.X
| (n + 2) := 2 * polynomial.X * chebyshev (n + 1) - chebyshev n

#eval chebyshev 17 -- just about works on my desktop
/-
C (65536) * X ^ 17 + C (-278528) * X ^ 15 + C (487424) * X ^ 13 + C (-452608) * X ^ 11 + C (239360) * X ^ 9 + C (-71808) * X ^ 7 + C (11424) * X ^ 5 + C (-816) * X ^ 3 + C (17) * X
-/

#eval chebyshev 19 -- (deterministic) timeout


#### Mario Carneiro (Jan 12 2019 at 18:06):

you can probably get farther with a list based representation of the polynomials

#### Kevin Buzzard (Jan 12 2019 at 18:07):

I guess the reason lists aren't used in general is that the leading term might be 0, which causes problems in general; however it will not cause problems here.

#### Mario Carneiro (Jan 12 2019 at 18:07):

well, it would make it harder to print if the leading term was 0

#### Mario Carneiro (Jan 12 2019 at 18:07):

and get the degree and other things

#### Kevin Buzzard (Jan 12 2019 at 18:07):

In general I guess one could define some more computationally efficient but slightly broken "non-zero polynomials" with addition only defined if you can prove that the degrees are unequal, and multiplication OK over an integral semi-domain :-)

#### Kevin Buzzard (Jan 12 2019 at 18:08):

But even lists are inefficient, right, because they're linked lists in reality?

#### Mario Carneiro (Jan 12 2019 at 18:09):

If you want to forgo the equality checks, you can define a polynomial as a quotient of representations with some zeros at the end

#### Mario Carneiro (Jan 12 2019 at 18:09):

then addition and multiplication are easy and degree needs equality

#### Mario Carneiro (Jan 12 2019 at 18:09):

yes lists are linked lists which aren't so great

#### Mario Carneiro (Jan 12 2019 at 18:10):

you could use buffer for really good VM performance

#### Mario Carneiro (Jan 12 2019 at 18:10):

but they are better than functions which is what is currently used

#### Mario Carneiro (Jan 12 2019 at 18:10):

I think the polynomials become big if else chains

#### Mario Carneiro (Jan 12 2019 at 18:12):

in fact, they probably aren't even reduced, you get these big thunks for calculating the coefficients

#### Mario Carneiro (Jan 12 2019 at 18:12):

the advantage of lists is they are a strict data structure, you calculate all the values before recursing

#### Mario Carneiro (Jan 12 2019 at 18:13):

also -- should have mentioned this before -- chebyshev there has an exponential time implementation

#### Mario Carneiro (Jan 12 2019 at 18:13):

just like the naive fib algorithm

#### Daniel Keys (Feb 13 2020 at 18:02):

It seems some things may have changed since this discussion took place. Trying to get started working with polynomials; analysis.polynomial is not a valid import any longer. In the following I try to reproduce Kevin's piece of code listed above in this same topic, but it doesn't work on my computer because of the VM not having code for classical.choice. Am I missing something? Also, the definition of p1 is non-computable; do we have a way to compute (like #eval) a value from it (or any other polynomial on the reals, or even nats) at a given point?

import data.polynomial
import data.real.basic

-- this is noncomputable
def p1 : polynomial ℝ := polynomial.C (2:ℝ)
-- these check but I can't produce a value
#check polynomial.eval ( 5 : ℝ ) ( polynomial.C (2:ℝ) )
#check polynomial.C 1

def chebyshev : ℕ → polynomial ℤ -- fails,  VM: no classical.choice
| 0 := polynomial.C 1
| 1 := polynomial.X
| (n + 2) := 2 * polynomial.X * chebyshev (n + 1) - chebyshev n

#eval chebyshev 3 --fails, VM: no classical.choice


#### Johan Commelin (Feb 13 2020 at 18:09):

choice will certainly be uncomputable. That's a fact of life that we'll need to deal with.

#### Johan Commelin (Feb 13 2020 at 18:09):

So you will want

noncomputable theory


at the top of your file, just after the imports

#### Johan Commelin (Feb 13 2020 at 18:11):

Next, #eval won't work, so you can't simply compute the value. But you could write some simp-lemmas that can do some work for you.

#### Johan Commelin (Feb 13 2020 at 18:11):

Which is almost as good.

#### Johan Commelin (Feb 13 2020 at 18:13):

It is unfortunate that it is currently not easy to run the simplifier against an expression like chebyshev 3. You need to write example : chebyshev 3 = the_answer := by simp [chebyshev], which defeats the purpose because you will have to precompute the answer. And then Lean will tell you that it agrees :expressionless:

#### Rob Lewis (Feb 13 2020 at 18:13):

Maybe we should define a command #simp that behaves like #eval except using the simplifier.

#### Daniel Keys (Feb 13 2020 at 18:16):

Just wondering how Kevin could #eval chebyshev 17, as he shows above (that is, the older messages on the same topic).

#### Rob Lewis (Feb 13 2020 at 18:18):

#1391 removed this. It worked on toy examples, but isn't very helpful in practice and causes problems elsewhere.

#### Kevin Buzzard (Feb 13 2020 at 19:58):

Chris is sitting next to me talking about perhaps having some attribute for computing and then you just unfold a lot until you can use norm_num or something...

#### Rob Lewis (Feb 13 2020 at 20:27):

I think that's tricky. Do you unfold has_add.add? Sometimes yes, because you want to unfold addition on whatever type you're working on. Sometimes no, because you want norm_num to deal with it.

Last updated: May 09 2021 at 19:11 UTC