Zulip Chat Archive

Stream: new members

Topic: Can't define Chebyshev polynomials


view this post on Zulip 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?

view this post on Zulip 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!

view this post on Zulip 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.

view this post on Zulip 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

view this post on Zulip Mario Carneiro (Jan 12 2019 at 18:06):

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

view this post on Zulip 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.

view this post on Zulip Mario Carneiro (Jan 12 2019 at 18:07):

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

view this post on Zulip Mario Carneiro (Jan 12 2019 at 18:07):

and get the degree and other things

view this post on Zulip 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 :-)

view this post on Zulip Kevin Buzzard (Jan 12 2019 at 18:08):

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

view this post on Zulip 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

view this post on Zulip Mario Carneiro (Jan 12 2019 at 18:09):

then addition and multiplication are easy and degree needs equality

view this post on Zulip Mario Carneiro (Jan 12 2019 at 18:09):

yes lists are linked lists which aren't so great

view this post on Zulip Mario Carneiro (Jan 12 2019 at 18:10):

you could use buffer for really good VM performance

view this post on Zulip Mario Carneiro (Jan 12 2019 at 18:10):

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

view this post on Zulip Mario Carneiro (Jan 12 2019 at 18:10):

I think the polynomials become big if else chains

view this post on Zulip 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

view this post on Zulip 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

view this post on Zulip Mario Carneiro (Jan 12 2019 at 18:13):

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

view this post on Zulip Mario Carneiro (Jan 12 2019 at 18:13):

just like the naive fib algorithm

view this post on Zulip 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

view this post on Zulip 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.

view this post on Zulip Johan Commelin (Feb 13 2020 at 18:09):

So you will want

noncomputable theory

at the top of your file, just after the imports

view this post on Zulip 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.

view this post on Zulip Johan Commelin (Feb 13 2020 at 18:11):

Which is almost as good.

view this post on Zulip 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:

view this post on Zulip Rob Lewis (Feb 13 2020 at 18:13):

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

view this post on Zulip 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).

view this post on Zulip 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.

view this post on Zulip 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...

view this post on Zulip 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