Lean Prover Zulip Chat Archive

Stream: general

Topic: happy pi day!


view this post on Zulip Kenny Lau (Mar 14 2019 at 11:59):

happy pi day!

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 12:41):

https://www.bbc.co.uk/news/technology-47524760

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 12:54):

But did they formally verify it? And if not, should we believe the calculations? They should be reproduced using a completely different system.

view this post on Zulip Moses Schönfinkel (Mar 14 2019 at 13:27):

They use https://en.wikipedia.org/wiki/Bellard%27s_formula to verify digits they compute with https://en.wikipedia.org/wiki/Chudnovsky_algorithm

view this post on Zulip Moses Schönfinkel (Mar 14 2019 at 13:30):

Also it would appear that they are simply using Alexander Yee's pi cruncher to simply run their computation, this is its homepage: http://www.numberworld.org/y-cruncher/

view this post on Zulip Koundinya Vajjha (Mar 14 2019 at 13:32):

But did they formally verify it? And if not, should we believe the calculations? They should be reproduced using a completely different system.

https://arxiv.org/abs/1709.01743

view this post on Zulip Moses Schönfinkel (Mar 14 2019 at 13:38):

There's also an extra verification algorithm they employ it seems, so you can trust it about as much as you can trust the (un-verified C++) implementation of Lean that runs on un-verified hardware.

view this post on Zulip Kenny Lau (Mar 14 2019 at 13:40):

Can we prove that pi starts with 3.14? i.e. 3.14 <= pi < 3.15?

view this post on Zulip Kenny Lau (Mar 14 2019 at 13:40):

I think a Chinese mathematician 1000 years ago used 3•2^n-gons to estimate pi up to 7 digits?

view this post on Zulip Chris Hughes (Mar 14 2019 at 14:04):

The easiest way using current mathlib, is to just use some lemma I proved about approximating e^x for x le 1, as the first n terms of the power series. I think tactics to do this stuff would be good.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:23):

Can we prove that pi starts with 3.14? i.e. 3.14 <= pi < 3.15?

Do we have any theorem of the form "pi is the following infinite sum" or "pi is the limit of the following sequence"?

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:26):

Maybe today is a good day to talk about pi. We know cos(pi/2)=0 so we can probably deduce any reasonable trig assertion about pi. Now the derivative of arctan is 1/(1+x^2) and...oh dear now I want to integrate term by term to deduce the power series for arctan, and then evaluate at 1 to deduce pi/4=1-1/3+1/5-1/7+..., for which one can get some reasonable approximations to pi. But we still need derivative of a power series :-/

view this post on Zulip Kenny Lau (Mar 14 2019 at 14:43):

that power series is notoriously slow to converge

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:44):

I'm well aware of that. Conversely it's an alternating series so it provides us with a rigorous proof that, for example, pi>1-1/3+1/5-1/7+1/9-1/11 and pi<1-1/3+1/5-1/7+1/9-1/11+1/13.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:45):

Which is more than we know currently in mathlib.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:47):

The 2^n-gon trick converges much more quickly and we have the trig in place, but it relies on the fact that pi = circumference over diameter, which we can't formalise yet.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:47):

Claim: pi > 3. Proof: draw a regular hexagon and chop it into 6 triangles. Then draw a circle round the hexagon.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:48):

Circumference of circle > circumference of hexagon = 6, diameter = 2.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:48):

Can we turn this idea into a mathlib proof that pi > 3? If so then the polygon trick starts to kick in.

view this post on Zulip Abhimanyu Pallavi Sudhir (Mar 14 2019 at 14:48):

You'll need pi as an integral, correct?

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:49):

I don't know. We've defined pi "in terms of trig", and we're "doing trig" in some vague way above -- but it looks to me like we're somehow doing a different trig.

view this post on Zulip Kenny Lau (Mar 14 2019 at 14:49):

just define it as the limit of the 3*2^n-gons

view this post on Zulip Abhimanyu Pallavi Sudhir (Mar 14 2019 at 14:50):

I think the easiest way to prove the circumference of that limit is pi still uses Riemann sums, though.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:50):

ha ha

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:50):

there's a problem with the polygon approach

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:51):

you end up defining 2pi as the limit, as n tends to infinity, of nsin(2π/n)n\sin(2\pi/n)

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:52):

(that's using n-gons; using 2^n-gons just change n to 2^n)

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:52):

OTOH if one could prove that nsin(2π/n)n\sin(2\pi/n) was increasing

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:53):

and that the derivative of sin was cos

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:53):

then setting n=6 you get pi>3

view this post on Zulip Reid Barton (Mar 14 2019 at 14:53):

For lower bounds, I guess you just need to prove sinx<x\sin x < x for x>0x > 0 which seems feasible

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:54):

aah so this is the point; we can compute sin(2 pi / (3*2^m)) as a solution to a quadratic

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:54):

For lower bounds, I guess you just need to prove $\sin x < x$ for $x > 0$ which seems feasible

yeah, just draw the graph :P

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:54):

that crazy quote bug again

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:55):

sin(x)<x should be easy by that alternating series thing

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:55):

Maybe Chris even proved sin(x)>x-x^3/3 which might even have been used to prove the current known bounds for pi, which are something like 2<pi<4

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:56):

good to see that google are on 33 trillion digits and we're still working on the 0th one :-)

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:56):

the thing is, when we get it, we'll have a proof

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:56):

not just some blather about cloud computing

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 14:57):

It's Xena tonight, maybe we should have a pi day celebration and try to prove pi > 3

view this post on Zulip Kenny Lau (Mar 14 2019 at 15:01):

firstly how do you define pi

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:02):

it's in mathlib

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:02):

so it is defined by official decree already

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:02):

I think they went for 22/7 but I can't remember

view this post on Zulip Kenny Lau (Mar 14 2019 at 15:03):

it's in mathlib?

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:03):

no wait, that was something else.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:03):

Isn't it in mathlib??

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:04):

yes it's there

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:04):

It's in analysis.exponential

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:04):

local notation π := 22 / 7

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:04):

no wait

view this post on Zulip Kenny Lau (Mar 14 2019 at 15:05):

oh right they proved sin(pi) = 0 right

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:05):

local notation π := real.pi

view this post on Zulip Kenny Lau (Mar 14 2019 at 15:05):

oh

view this post on Zulip Kenny Lau (Mar 14 2019 at 15:05):

what's the definition?

view this post on Zulip Abhimanyu Pallavi Sudhir (Mar 14 2019 at 15:05):

Two times the zero of cos x between 0 and 2. I don't see why they didn't just prove exp periodic and define it as the half-period. Maybe that's hard without knowing what pi is, but it feels more natural.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:06):

I think they were pressed for time and wanted to get it out for my birthday :-)

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:07):

I gave a talk yesterday to 100 people at my university interested in teaching, and told them the story of how I got pi for my birthday, and it got the biggest laugh of the talk.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:07):

You weren't around at that point Abhi, so you couldn't keep them in order :-)

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:08):

noncomputable def pi : ℝ := 2 * classical.some exists_cos_eq_zero

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:08):

holy crap, it could be twice any zero at all!

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:08):

Maybe it's not less than 4 after all!

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:08):

pesky axiom of choice

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:10):

lemma exists_cos_eq_zero :  x, 1  x  x  2  cos x = 0 :=

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:10):

Proof uses cos_one_pos and cos_two_neg and IVT

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:10):

so indeed currently we should be able to prove that pi is between 2 and 4

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:11):

but for proving it's bigger than 3 it does not suffice to prove that cos(1.5) is positive, because there might be some other zero of cos in the range (1,1.5).

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:12):

It would be good to know that cos was decreasing. That would follow from MVT perhaps, if we knew the derivative of cos was sin, but we might be a long way from that.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 15:12):

Do we know the derivative of exp is exp yet?

view this post on Zulip Abhimanyu Pallavi Sudhir (Mar 14 2019 at 15:14):

I don't. Maybe Jeremy Avigad does (if he's specialised his derivative for basic functions).

view this post on Zulip Reid Barton (Mar 14 2019 at 15:14):

I think we do know cos is decreasing, or at least I saw something similar when looking through the file for sin x < x

view this post on Zulip Reid Barton (Mar 14 2019 at 15:15):

cos_lt_cos_of_nonneg_of_le_pi_div_two

view this post on Zulip Reid Barton (Mar 14 2019 at 15:15):

or actually cos_lt_cos_of_nonneg_of_le_pi

view this post on Zulip Reid Barton (Mar 14 2019 at 15:16):

It's proved using angle addition formulas

view this post on Zulip Floris van Doorn (Mar 14 2019 at 15:49):

If we know that sin x < x here is a proof sketch to show that pi > 3:

The doubling formula cos(2*x) = 2*cos^2(x) - 1 has been formalized, so we can compute cos^2(pi/4) and cos^2(pi/8). Now we can also compute sin^2(pi/8). This turns out to be 1/2 - 1/sqrt(8). Now we know that pi > 8sin(pi/8) = 8sqrt(1/2-1/sqrt(8)) > 3.

view this post on Zulip Kenny Lau (Mar 14 2019 at 15:56):

nice!

view this post on Zulip Chris Hughes (Mar 14 2019 at 15:57):

I can imagine this being a long annoying proof.

view this post on Zulip Reid Barton (Mar 14 2019 at 15:57):

sin_bound should be enough to prove sin x < x as well

view this post on Zulip Floris van Doorn (Mar 14 2019 at 15:58):

I'm starting now to compute sin(pi/8).

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 16:02):

Happy pi day :-)

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 16:03):

Chris this proof doesn't look too bad to me at all. I was manipulating 1+sqrt(5)/2 without too much trouble last week when we were doing some analysis example sheet question.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 16:03):

Or is there something else you're worried about?

view this post on Zulip Floris van Doorn (Mar 14 2019 at 16:03):

How do I prove (2 : ℂ) ≠ 0?

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 16:03):

These things should go in pi.lean

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 16:03):

norm_num proves 2 ne 0

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 16:03):

or if it doesn't then take real parts

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 16:03):

because it will definitely prove (2 : real) \ne 0

view this post on Zulip Chris Hughes (Mar 14 2019 at 16:03):

two_ne_zero'

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 16:04):

rofl

view this post on Zulip Chris Hughes (Mar 14 2019 at 16:04):

It needs the prime, because I think the other version is for linear ordered things.

view this post on Zulip Floris van Doorn (Mar 14 2019 at 16:04):

yes, I found the other one. Thanks

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 16:04):

linear ordered -> char 0

view this post on Zulip Chris Hughes (Mar 14 2019 at 16:05):

The linear order version is in core I think.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 16:07):

https://leanprover.zulipchat.com/#narrow/stream/116395-maths/topic/golden.20ratio.20calculation

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 16:07):

@Floris van Doorn there's an example of basic manipulations with sqrt(5)

view this post on Zulip Floris van Doorn (Mar 14 2019 at 19:06):

import analysis.exponential

namespace real
variable (x : ℝ)

lemma sqrt_le_left {x y : ℝ} (hy : 0 ≤ y) : sqrt x ≤ y ↔ x ≤ y ^ 2 :=
begin
  rw [mul_self_le_mul_self_iff (sqrt_nonneg _) hy, pow_two],
  cases le_total 0 x with hx hx,
  { rw [mul_self_sqrt hx] },
  { have h1 : 0 ≤ y * y := mul_nonneg hy hy,
    have h2 : x ≤ y * y := le_trans hx h1,
    simp [sqrt_eq_zero_of_nonpos, hx, h1, h2] }
end

lemma le_sqrt {x y : ℝ} (hx : 0 ≤ x) (hy : 0 ≤ y) : x ≤ sqrt y ↔ x ^ 2 ≤ y :=
by rw [mul_self_le_mul_self_iff hx (sqrt_nonneg _), pow_two, mul_self_sqrt hy]

lemma div_lt_self {x y : ℝ} (hx : x > 0) (hy : y > 1) : x / y < x :=
by simpa using div_lt_div' (le_refl x) hy hx zero_lt_one

lemma cos_square : cos x ^ 2 = (1 + cos (2 * x)) / 2 :=
by simp [cos_two_mul, mul_div_cancel_left, two_ne_zero]

lemma sin_square : sin x ^ 2 = 1 - cos x ^ 2 :=
by { rw [←sin_pow_two_add_cos_pow_two x], simp }

lemma cos_square_pi_div_four : cos (pi / 4) ^ 2 = 1 / 2 :=
by { have : 2 * (pi / 4) = pi / 2, ring, simp [cos_square, this] }

lemma cos_pi_div_four : cos (pi / 4) = sqrt (1 / 2) :=
begin
  symmetry, rw [sqrt_eq_iff_sqr_eq, cos_square_pi_div_four], norm_num,
  apply le_of_lt, apply cos_pos_of_neg_pi_div_two_lt_of_lt_pi_div_two,
  { transitivity (0 : ℝ), rw neg_lt_zero, apply pi_div_two_pos,
    exact div_pos pi_pos (by norm_num) },
  apply div_lt_div' (le_refl pi) _ pi_pos _,
  norm_num, norm_num
end

/- why is division_ring.inv_inj in the file "field"?-/
lemma cos_square_pi_div_eight : cos (pi / 8) ^ 2 = 1 / 2 + sqrt (1 / 8) :=
begin
  have h1 : 2 * (pi / 8) = pi / 4, ring,
  have h2 : sqrt 4 = 2, rw [sqrt_eq_iff_sqr_eq]; norm_num,
  have h3 : 2 * sqrt 2 ≠ 0, norm_num,
  have h4 : sqrt 8 ≠ 0, norm_num,
  simp [cos_square, h1, cos_pi_div_four],
  rw [←div_add_div_same], simp [div_eq_mul_inv, (mul_inv' _ _).symm, division_ring.inv_inj, h3, h4],
  transitivity sqrt 4 * sqrt 2, rw [h2],
  rw [←sqrt_mul]; norm_num
end

lemma cos_pi_div_eight : cos (pi / 8) = sqrt (1 / 2 + sqrt (1 / 8)) :=
begin
  symmetry, rw [sqrt_eq_iff_sqr_eq, cos_square_pi_div_eight],
  apply add_nonneg, norm_num, apply sqrt_nonneg,
  apply le_of_lt, apply cos_pos_of_neg_pi_div_two_lt_of_lt_pi_div_two,
  { transitivity (0 : ℝ), rw neg_lt_zero, apply pi_div_two_pos,
    exact div_pos pi_pos (by norm_num) },
  apply div_lt_div' (le_refl pi) _ pi_pos _,
  norm_num, norm_num
end

lemma sin_square_pi_div_eight : sin (pi / 8) ^ 2 = 1 / 2 - sqrt (1 / 8) :=
by { simp [sin_square, cos_square_pi_div_eight], rw [←add_assoc], congr' 1, norm_num }

lemma sin_pi_div_eight : sin (pi / 8) = sqrt (1 / 2 - sqrt (1 / 8)) :=
begin
  symmetry, rw [sqrt_eq_iff_sqr_eq, sin_square_pi_div_eight],
  { rw [sub_nonneg, sqrt_le_left], norm_num, norm_num },
  apply le_of_lt, apply sin_pos_of_pos_of_lt_pi,
  { exact div_pos pi_pos (by norm_num) },
  exact div_lt_self pi_pos (by norm_num)
end

lemma sin_lt {x : ℝ} (h : 0 < x) : sin x < x :=
begin
  cases le_or_gt x 1 with h' h',
  { have hx : abs x = x := abs_of_nonneg (le_of_lt h),
    have : abs x ≤ 1, rwa [hx],
    have := sin_bound this, rw [abs_le] at this,
    have := this.2, rw [sub_le_iff_le_add', hx] at this,
    apply lt_of_le_of_lt this, rw [sub_add], apply lt_of_lt_of_le _ (le_of_eq (sub_zero x)),
    apply sub_lt_sub_left, rw sub_pos, apply mul_lt_mul',
    { rw [pow_succ x 3], refine le_trans _ (le_of_eq (one_mul _)),
      rw mul_le_mul_right, exact h', apply pow_pos h },
    norm_num, norm_num, apply pow_pos h },
  exact lt_of_le_of_lt (sin_le_one x) h'
end

lemma pt_gt_three : pi > 3 :=
begin
  have : pi > sin (pi / 8) * 8,
  { apply mul_lt_of_lt_div, norm_num, apply sin_lt, apply div_pos pi_pos, norm_num },
  apply lt_of_le_of_lt _ this, rw [sin_pi_div_eight],
  apply le_mul_of_div_le, norm_num,
  rw [le_sqrt, le_sub_iff_add_le, ←le_sub_iff_add_le', sqrt_le_left],
  norm_num, norm_num, norm_num,
  rw [sub_nonneg, sqrt_le_left], norm_num, norm_num
end

end real

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 19:08):

This doesn't compile for me

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 19:08):

solve1 tactic failed, focused goal has not been solved
state:
⊢ pi / 8 ≤ 1

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 19:08):

We need pi <= 8

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 19:09):

lemma pt_gt_three : pi > 3 :=
begin
  have : pi > sin (pi / 8) * 8,
  { apply mul_lt_of_lt_div, norm_num, apply sin_lt, apply div_pos pi_pos, norm_num,
  -- ⊢ pi / 8 ≤ 1
    sorry
  },
  apply lt_of_le_of_lt _ this, rw [sin_pi_div_eight],
  apply le_mul_of_div_le, norm_num,
  rw [le_sqrt, le_sub_iff_add_le, ←le_sub_iff_add_le', sqrt_le_left],
  norm_num, norm_num, norm_num,
  rw [sub_nonneg, sqrt_le_left], norm_num, norm_num
end

view this post on Zulip Floris van Doorn (Mar 14 2019 at 19:09):

oops, I missed that. One minute

view this post on Zulip Floris van Doorn (Mar 14 2019 at 19:11):

I updated the previous post, it should compile now.

view this post on Zulip Reid Barton (Mar 14 2019 at 19:11):

You can also prove sin x < x for all x > 0 by handling the case x > 1 by sin x <= 1 < x

view this post on Zulip Floris van Doorn (Mar 14 2019 at 19:14):

@Reid Barton That is indeed the nicer solution. I was planning to do that, but I was a bit too eager to post it here and forgot :)

view this post on Zulip Floris van Doorn (Mar 14 2019 at 19:15):

I updated the code once more.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 19:29):

This should all go in pi.lean, right?

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 19:30):

So now we can prove that floor(pi) = 3

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 19:30):

So the next question is whether we can prove digit 10 pi 0 = 1

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 19:31):

ie floor(pi*10)%10=1

view this post on Zulip Johan Commelin (Mar 14 2019 at 19:31):

This should all go in pi.lean, right?

Isn't that about pi types :rolling_on_the_floor_laughing: :cake:

view this post on Zulip Floris van Doorn (Mar 14 2019 at 19:37):

Computing sin(pi/16) we can conclude that pi > 3.1, computing sin(pi/64) we can conclude that pi>3.14. Using the other side of sin_bound we can probably also get upper bounds.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 19:45):

Can you get to 32 trillion by the end of pi day?

view this post on Zulip Floris van Doorn (Mar 14 2019 at 19:49):

probably not

view this post on Zulip Floris van Doorn (Mar 14 2019 at 23:35):

Is there a better way than norm_num to prove the following?

(2 : ℝ) ≤ ((((2 - (157 / 50 / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2) ^ 2 - 2) ^ 2 - 2) ^ 2

(the RHS is approx 2.00033). norm_num is struggling on this one, and gets a deterministic timeout.

view this post on Zulip Floris van Doorn (Mar 14 2019 at 23:45):

If I can prove inequalities like that, I can approximate pi for a few decimal places, at least from below (I need to still work on approximation from above). The approximation I have of pi is of the form

π ≈ 2 ^ 7 * sqrt (2 - sqrt (2 + sqrt (2 + sqrt (2 + sqrt (2 + sqrt (2 + sqrt 2))))))

This seems to converge at an okay rate: the number of correct digits seems to be linear in the number of square roots.

In particular, I formalized

lemma pt_gt_314 : (2 : ℝ) ≤ ((((2 - (157 / 50 / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2) ^ 2 - 2) ^ 2 - 2) ^ 2 → pi > 3.14

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 23:57):

I can reduce your time-out goal to ⊢ 0 < 3402823669209384634633746074317682114560000000000000000000000000000000000000000000000000000000000000000

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 23:57):

It's not even close

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 23:58):

but I'm struggling with the last bit.

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 23:58):

oh -- these are reals

view this post on Zulip Kevin Buzzard (Mar 14 2019 at 23:59):

import data.real.basic

example : (2 : )  ((((2 - (157 / 50 / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2) ^ 2 - 2) ^ 2 - 2) ^ 2 :=
begin
  rw (show 4 + 1 = 5, from rfl),
  rw (show (2 : ) ^ 5 = 32, by norm_num),
  rw (show (157 : ) / 50 / 32 = 157 / 1600, by norm_num),
  rw (show ((157 : ) / 1600) ^ 2 = 24649 / 2560000, by norm_num),
  rw (show (2 : ) - 24649 / 2560000 = 5095351/2560000, by norm_num),
  rw (show ((5095351 : ) /2560000) ^ 2 = 25962601813201/6553600000000, by norm_num),
  -- let's try skipping steps
  rw (show ((((25962601813201 : ) / 6553600000000 - 2) ^ 2 - 2) ^ 2 - 2) ^ 2 =
    6806775565993596917798714352237438749189222725565123913061124308255143227446400872965401873904861225601/3402823669209384634633746074317682114560000000000000000000000000000000000000000000000000000000000000000,
    by norm_num),
  -- worked!
  --  2 
  --  6806775565993596917798714352237438749189222725565123913061124308255143227446400872965401873904861225601 /
  --    3402823669209384634633746074317682114560000000000000000000000000000000000000000000000000000000000000000
  rw le_div_iff,
    norm_num,
  --  0 < 3402823669209384634633746074317682114560000000000000000000000000000000000000000000000000000000000000000
  sorry
end

view this post on Zulip Chris Hughes (Mar 15 2019 at 00:05):

norm_num just works for that last goal for me.

example : (2 : )  ((((2 - (157 / 50 / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2) ^ 2 - 2) ^ 2 - 2) ^ 2 :=
begin
  rw (show 4 + 1 = 5, from rfl),
  rw (show (2 : ) ^ 5 = 32, by norm_num),
  rw (show (157 : ) / 50 / 32 = 157 / 1600, by norm_num),
  rw (show ((157 : ) / 1600) ^ 2 = 24649 / 2560000, by norm_num),
  rw (show (2 : ) - 24649 / 2560000 = 5095351/2560000, by norm_num),
  rw (show ((5095351 : ) /2560000) ^ 2 = 25962601813201/6553600000000, by norm_num),
  -- let's try skipping steps
  rw (show ((((25962601813201 : ) / 6553600000000 - 2) ^ 2 - 2) ^ 2 - 2) ^ 2 =
    6806775565993596917798714352237438749189222725565123913061124308255143227446400872965401873904861225601/3402823669209384634633746074317682114560000000000000000000000000000000000000000000000000000000000000000,
    by norm_num),
  -- worked!
  -- ⊢ 2 ≤
  --  6806775565993596917798714352237438749189222725565123913061124308255143227446400872965401873904861225601 /
  --    3402823669209384634633746074317682114560000000000000000000000000000000000000000000000000000000000000000
  rw le_div_iff,
  { norm_num },
  { norm_num }
  -- ⊢ 0 < 3402823669209384634633746074317682114560000000000000000000000000000000000000000000000000000000000000000

end

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 00:08):

import data.real.basic

--set_option pp.all true
example : (2 : )  ((((2 - (157 / 50 / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2) ^ 2 - 2) ^ 2 - 2) ^ 2 :=
begin
  rw (show 4 + 1 = 5, from rfl),
  rw (show (2 : ) ^ 5 = 32, by norm_num),
  rw (show (157 : ) / 50 / 32 = 157 / 1600, by norm_num),
  rw (show ((157 : ) / 1600) ^ 2 = 24649 / 2560000, by norm_num),
  rw (show (2 : ) - 24649 / 2560000 = 5095351/2560000, by norm_num),
  rw (show ((5095351 : ) /2560000) ^ 2 = 25962601813201/6553600000000, by norm_num),
  -- let's try skipping steps
  rw (show ((((25962601813201 : ) / 6553600000000 - 2) ^ 2 - 2) ^ 2 - 2) ^ 2 =
    6806775565993596917798714352237438749189222725565123913061124308255143227446400872965401873904861225601/3402823669209384634633746074317682114560000000000000000000000000000000000000000000000000000000000000000,
    by norm_num),
  -- worked!
  -- ⊢ 2 ≤
  --  6806775565993596917798714352237438749189222725565123913061124308255143227446400872965401873904861225601 /
  --    3402823669209384634633746074317682114560000000000000000000000000000000000000000000000000000000000000000
  rw le_div_iff,
    norm_num,
  -- ⊢ 0 < 3402823669209384634633746074317682114560000000000000000000000000000000000000000000000000000000000000000
  rw (show (3402823669209384634633746074317682114560000000000000000000000000000000000000000000000000000000000000000 : )
           = 40 ^ 64, by norm_num),
  apply pow_pos,
  norm_num,
end

This worked for me.

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 00:08):

Oh, I think that I upset Lean with some overflow and then failed to restart, and had problems.

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 00:20):

Only 31999999999998 digits to go and we have the record!

view this post on Zulip Kenny Lau (Mar 15 2019 at 00:31):

oh boy

view this post on Zulip Floris van Doorn (Mar 15 2019 at 01:49):

I have proven

lemma pi_gt_314 : pi > 3.14
lemma pi_lt_315 : pi < 3.15

assuming the following 12 inequalities:

⊢ 2 ≤ ((((2 - (157 / 50 / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2) ^ 2 - 2) ^ 2 - 2) ^ 2
⊢ 0 ≤ (((2 - (157 / 50 / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2) ^ 2 - 2) ^ 2 - 2
⊢ 0 ≤ ((2 - (157 / 50 / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2) ^ 2 - 2
⊢ 0 ≤ (2 - (157 / 50 / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2
⊢ 0 ≤ 2 - (157 / 50 / 2 ^ (4 + 1)) ^ 2
⊢ 0 < 157 / 50 / 2 ^ (4 + 1)
⊢ ((((2 - ((63 / 20 - 1 / 4 ^ 4) / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2) ^ 2 - 2) ^ 2 - 2) ^ 2 ≤ 2
⊢ 0 < (((2 - ((63 / 20 - 1 / 4 ^ 4) / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2) ^ 2 - 2) ^ 2 - 2
⊢ 0 < ((2 - ((63 / 20 - 1 / 4 ^ 4) / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2) ^ 2 - 2
⊢ 0 < (2 - ((63 / 20 - 1 / 4 ^ 4) / 2 ^ (4 + 1)) ^ 2) ^ 2 - 2
⊢ 0 < 2 - ((63 / 20 - 1 / 4 ^ 4) / 2 ^ (4 + 1)) ^ 2
⊢ 0 ≤ (63 / 20 - 1 / 4 ^ 4) / 2 ^ (4 + 1)

view this post on Zulip Kenny Lau (Mar 15 2019 at 01:54):

I've verified all of your inequalities on WolframAlpha

view this post on Zulip Mario Carneiro (Mar 15 2019 at 01:56):

Oh boy I've missed a lot

view this post on Zulip Mario Carneiro (Mar 15 2019 at 01:57):

For these kinds of numeric calculations you can simplify a lot by decreasing the precision in the middle of the calculation

view this post on Zulip Mario Carneiro (Mar 15 2019 at 01:58):

Are all the inequalities implied by one of them?

view this post on Zulip Kenny Lau (Mar 15 2019 at 01:58):

the first one is exceptionally tight

view this post on Zulip Kenny Lau (Mar 15 2019 at 01:58):

the 2nd to 5th looks like a more general theorem?

view this post on Zulip Kenny Lau (Mar 15 2019 at 01:58):

wait are those approximations to sqrt(2)?

view this post on Zulip Kenny Lau (Mar 15 2019 at 01:59):

@Floris van Doorn I think it would be better if you tell us more about what's actually going on

view this post on Zulip Kenny Lau (Mar 15 2019 at 01:59):

if we want approximations to sqrt(2) maybe we should turn to Pell theory

view this post on Zulip Mario Carneiro (Mar 15 2019 at 01:59):

I think they are approximations to 2

view this post on Zulip Mario Carneiro (Mar 15 2019 at 01:59):

fixed points of x |-> x^2-2

view this post on Zulip Kenny Lau (Mar 15 2019 at 02:00):

interesting...?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:01):

The first one isn't that tight: it is saying that

   1128227574827648531222203602074520069222725565123913061124308255143227446400872965401873904861225601 /
3402823669209384634633746074317682114560000000000000000000000000000000000000000000000000000000000000000
>=0

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:01):

Notice that there are like 60 digits of precision there and only 4 of them matter

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:02):

so you can vastly simplify the calculation by precision reduction in the middle

view this post on Zulip Kenny Lau (Mar 15 2019 at 02:02):

aha

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:05):

π ≈ 2 ^ 7 * sqrt (2 - sqrt (2 + sqrt (2 + sqrt (2 + sqrt (2 + sqrt (2 + sqrt 2))))))

Where did this come from?

view this post on Zulip Floris van Doorn (Mar 15 2019 at 02:06):

The file I'm working on is here: https://github.com/fpvandoorn/mathlib/blob/pi/src/data/real/pi.lean

view this post on Zulip Kenny Lau (Mar 15 2019 at 02:08):

/- the series sqrt(2+sqrt(2+ ... ))-/
noncomputable def sqrt_two_add_series :   
| 0     := 0
| (n+1) := sqrt (2 + sqrt_two_add_series n)

lemma pi_gt_sqrt_two_add_series (n : ) : pi > 2 ^ (n+1) * sqrt (2 - sqrt_two_add_series n) :=

lemma pi_lt_sqrt_two_add_series (n : ) : pi < 2 ^ (n+1) * sqrt (2 - sqrt_two_add_series n) + 1 / 4 ^ n :=

lemma pi_gt_314 : pi > 3.14 :=
begin
  refine lt_of_le_of_lt _ (pi_gt_sqrt_two_add_series 4),

lemma pi_lt_315 : pi < 3.15 :=
begin
  refine lt_of_lt_of_le (pi_lt_sqrt_two_add_series 4) _,

view this post on Zulip Kenny Lau (Mar 15 2019 at 02:08):

extract ^

view this post on Zulip Kenny Lau (Mar 15 2019 at 02:09):

anyway @Floris van Doorn remarkable work!

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:11):

The initial starting point sqrt_two_add_series 0 is arbitrary? What constraints does it require?

view this post on Zulip Floris van Doorn (Mar 15 2019 at 02:12):

The first and 7th inequalities are relatively tight (it compares 2 with 2.0003 and 1.994). The other inequalities are not tight at all, the numbers on the right are all between 1 and 2.
What's going on is that I proved that

2 ^ 5 * sqrt (2 - sqrt (2 + sqrt (2 + sqrt (2 + sqrt (2)))) < π < 2 ^ 5 * sqrt (2 - sqrt (2 + sqrt (2 + sqrt (2 + sqrt (2)))) + 1 / 4 ^ 4

Now we can check the bounds by replacing π by 3.14 and 3.15 and then moving all the square roots of the inequality to the other side. To do that, we need to check that stuff is positive (that gives all the non-tight inequalities), and the final inequality are the two tight ones.

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:12):

If we let sqrt_two_add_series a n be the series starting at a, then my hope is to get an inequality chain based on sqrt_two_add_series a (n+1) <= sqrt_two_add_series b n when a and b are suitably related rational numbers

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:14):

is sqrt_two_add_series a n monotonic in a? I think it is

view this post on Zulip Floris van Doorn (Mar 15 2019 at 02:14):

The defining equation of sqrt_two_add_series is

cos (pi / 2 ^ (n+1)) = sqrt_two_add_series n / 2

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:14):

aha, you inverted sin^2+cos^2 = 1

view this post on Zulip Floris van Doorn (Mar 15 2019 at 02:15):

yes, the main equality is lemma cos_square : cos x ^ 2 = 1 / 2 + cos (2 * x) / 2

view this post on Zulip Kenny Lau (Mar 15 2019 at 02:16):

how do we formalize "precision reduction"? @Mario Carneiro

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:17):

blasty lean users aren't going to like the answer

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:17):

you give intermediate small rational numbers as certificates

view this post on Zulip Floris van Doorn (Mar 15 2019 at 02:18):

yes, it might be nice to let sqrt_two_add_series start at an arbitrary number, and yes, sqrt_two_add_series a n is monotone in both a and n.

view this post on Zulip Kenny Lau (Mar 15 2019 at 02:19):

fortunately I'm not a blasty lean user

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:19):

Is it true that sqrt_two_add_series a (n+1) = sqrt_two_add_series (sqrt a + 2) n or something like that?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:19):

sorry playing catchup

view this post on Zulip Floris van Doorn (Mar 15 2019 at 02:21):

almost: sqrt_two_add_series a (n+1) = sqrt_two_add_series (sqrt (2 + a)) n

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:21):

Okay, so with monotonicity we get the desired result: if 2+a <= b^2 and b >= 0 then sqrt_two_add_series a (n+1) <= sqrt_two_add_series b n

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:22):

so now you apply this lemma 5 times, and supply 5 small rational numbers that satisfy the inequalities (check in mathematica or something)

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:23):

and those numbers are the certificate

view this post on Zulip Floris van Doorn (Mar 15 2019 at 02:23):

that sounds good

view this post on Zulip Floris van Doorn (Mar 15 2019 at 02:23):

(of course I hoped for a lazier solution :) )

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:24):

Oh the lazy solution is also possible, you can just tell lean "round to n digits at each stage"

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:24):

and pick n so it works

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:25):

but the really really important part is don't work out all the algebra with the original number, that causes exponential blowup

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:26):

I really like certificate proofs, they show the added power of ITP

view this post on Zulip Kenny Lau (Mar 15 2019 at 02:26):

didn't you write a floating point library or something?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:26):

you can solve NP problems in polynomial time, people! That's big news

view this post on Zulip Kenny Lau (Mar 15 2019 at 02:26):

would that be useful?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:27):

this is why metamath is fast and every tactic based prover is 5 orders of magnitude slower, because certificates matter

view this post on Zulip Kenny Lau (Mar 15 2019 at 02:28):

so we still use norm_num in each stage or something like that?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:28):

yes

view this post on Zulip Kenny Lau (Mar 15 2019 at 02:28):

it's definitely the wrong time for me to try, so maybe later

view this post on Zulip Kenny Lau (Mar 15 2019 at 02:28):

it's 2:28 AM

view this post on Zulip Mario Carneiro (Mar 15 2019 at 02:29):

I will come up with appropriate certificates in a bit, once I understand the situation right

view this post on Zulip Floris van Doorn (Mar 15 2019 at 02:40):

I'm currently adding a second argument and showing the monotonicity

view this post on Zulip Floris van Doorn (Mar 15 2019 at 02:51):

(updated the file on Github)

view this post on Zulip Floris van Doorn (Mar 15 2019 at 03:17):

are you working on the certificates @Mario Carneiro ? Otherwise I'll start on it now.

view this post on Zulip Mario Carneiro (Mar 15 2019 at 03:18):

I have this so far, for the lower bound:

{32 a5 >= 157/50, 2 - a5^2 >= a4, a4^2 >= 2 + a3, a3^2 >= 2 + a2,  a2^2 >= 2 + a1, a1^2 >= 2}

view this post on Zulip Mario Carneiro (Mar 15 2019 at 03:26):

Oh, I didn't expect this would work:

{a1 -> 2, a2 -> 2, a3 -> 2, a4 -> -(2687/1350), a5 -> 21/214}

view this post on Zulip Mario Carneiro (Mar 15 2019 at 03:27):

the last few (first few? my numbering is weird) are exactly 2, which works out nicely with the algebra. You probably want to special case that since it makes everything easy for going to higher values

view this post on Zulip Floris van Doorn (Mar 15 2019 at 03:27):

I'm still a bit confused what your variables a1-a5 correspond to.

view this post on Zulip Mario Carneiro (Mar 15 2019 at 03:28):

Here's the abstract argument:

pi >= 2^5 Sqrt[2 - f[0, 4]] >= 2^5 a5 >= 157/50
2 - a5^2 >= a4 >= f[a3, 1] >= f[a2, 2] >= f[a1, 3] >= f[0, 4]

view this post on Zulip Mario Carneiro (Mar 15 2019 at 03:28):

where f is your series function

view this post on Zulip Mario Carneiro (Mar 15 2019 at 03:31):

Note that there are no side conditions on this proof, because a4^2 >= 2 + a3 implies a4 >= Sqrt[2 + a3] regardless of whether a4 is positive

view this post on Zulip Mario Carneiro (Mar 15 2019 at 03:31):

for the upper bound I think there will be side conditions

view this post on Zulip Floris van Doorn (Mar 15 2019 at 03:37):

a4^2 >= 2 + a3 implies a4 >= Sqrt[2 + a3]

uhm... no: (-10) ^ 2 >= 2 + 0 but not -10 >= sqrt(2 + 0).

view this post on Zulip Floris van Doorn (Mar 15 2019 at 03:37):

For the upper bound there won't be side conditions, I indeed realized that a couple of minutes ago.

view this post on Zulip Mario Carneiro (Mar 15 2019 at 03:45):

Oh boy, I keep getting myself confused here. You are right, each of the a's must be positive, and I messed up a sign somewhere so my cert doesn't work

view this post on Zulip Floris van Doorn (Mar 15 2019 at 03:46):

The signs are pretty confusing indeed. But if you give me the certificate for the inequalities

f[a4, 0] >= f[a3, 1] >= f[a2, 2] >= f[a1, 3] >= f[0, 4]

then I'm ready to plug them in my Lean code now.

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:02):

Okay, here are the values for the lower bound:

{a1 -> 577/408, a2 -> 619/335, a3 -> 2144/1093, a4 -> 2687/1350, a5 -> 89/907}

view this post on Zulip Floris van Doorn (Mar 15 2019 at 04:04):

it works!

view this post on Zulip Floris van Doorn (Mar 15 2019 at 04:07):

lemma pi_gt_314 : pi > 3.14 :=
begin
  refine lt_of_le_of_lt _ (pi_gt_sqrt_two_add_series 4),
  rw [mul_comm],
  apply le_mul_of_div_le, norm_num, apply le_sqrt_of_sqr_le,
  rw [le_sub],
  refine le_trans (sqrt_two_add_series_step_up (577/408) (by norm_num) (by norm_num)) _,
  refine le_trans (sqrt_two_add_series_step_up (619/335) (by norm_num) (by norm_num)) _,
  refine le_trans (sqrt_two_add_series_step_up (2144/1093) (by norm_num) (by norm_num)) _,
  refine le_trans (sqrt_two_add_series_step_up (2687/1350) (by norm_num) (by norm_num)) _,
  rw [sqrt_two_add_series_zero],
  norm_num
end

view this post on Zulip Floris van Doorn (Mar 15 2019 at 04:19):

lemma pi_lt_315 : pi < 3.15 :=
begin
  refine lt_of_lt_of_le (pi_lt_sqrt_two_add_series 4) _,
  apply add_le_of_le_sub_right, rw [mul_comm], apply mul_le_of_le_div, apply pow_pos, norm_num,
  rw [sqrt_le_left, sub_le], swap, norm_num,
  refine le_trans _ (sqrt_two_add_series_step_down (41/29) (by norm_num)),
  refine le_trans _ (sqrt_two_add_series_step_down (194/105) (by norm_num)),
  refine le_trans _ (sqrt_two_add_series_step_down (13513/6889) (by norm_num)),
  refine le_trans _ (sqrt_two_add_series_step_down (2271/1141) (by norm_num)),
  norm_num
end

view this post on Zulip Floris van Doorn (Mar 15 2019 at 04:22):

The compiling proofs can now be found here: https://github.com/fpvandoorn/mathlib/blob/pi/src/data/real/pi.lean

They do each take ~5 seconds to compile. Maybe we can find smaller certificates (I found mine in a very naive way).

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:23):

I do want to figure out a general mechanism for optimizing all the certificates simultaneously; right now I'm optimizing them one by one and tweaking the results

view this post on Zulip Johan Commelin (Mar 15 2019 at 04:23):

Can't a tactic do that for you?

view this post on Zulip Johan Commelin (Mar 15 2019 at 04:23):

This must be a well-known idea in ITP, right?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:23):

yes and no; this is stuff that should not be re-run every time

view this post on Zulip Floris van Doorn (Mar 15 2019 at 04:23):

Yes, I'm just using Mathematica's Rationalize myself, until Rationalize finds a rational number which is on the correct side of the number I'm approximating

view this post on Zulip Johan Commelin (Mar 15 2019 at 04:24):

The tactic could trace its result

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:24):

me too (you can ask for a number in a range with Rationalize[(a+b)/2,(a-b)/2])

view this post on Zulip Floris van Doorn (Mar 15 2019 at 04:25):

oh yes, good point

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:26):

there isn't much point doing this in lean, unless lean gets as good at numbers as a CAS

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:27):

this is really the sort of thing CASs are designed for: unverified computation / symbolic manipulation

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:27):

If I was to do this sort of thing a lot I would want to automate some parts of it

view this post on Zulip Johan Commelin (Mar 15 2019 at 04:34):

So we need a Sage bridge...

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:38):

sure, that would do it

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:39):

or code extraction

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:39):

in this case, it's just pure computation, although it's nice having mathematica with builtin functions like Rationalize so we don't have to write it on the spot

view this post on Zulip Johan Commelin (Mar 15 2019 at 04:52):

@Mario Carneiro Do you think meta code could at some point be just as fast a proper CAS? I don't see why we shouldn't strive for world domination there as well...

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:52):

sure, we just need more efficiency

view this post on Zulip Johan Commelin (Mar 15 2019 at 04:52):

Is there hope that :four_leaf_clover: will bring this efficiency?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:52):

I've used lean plenty for this kind of thing, and it's a big bottleneck

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:53):

Yes, there is hope... I'm not putting all my eggs in that basket though

view this post on Zulip Johan Commelin (Mar 15 2019 at 04:53):

You haven't even put all your ITP eggs in Leans basket...

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:54):

I don't believe in world dominion

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:55):

it's not a healthy mindset. Better to support open communication and interaction between systems, because there will always be people that prefer another language

view this post on Zulip Johan Commelin (Mar 15 2019 at 04:56):

It also means you have to support n^2 communications...

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:57):

the common answer to that is "not if it's a star graph", i.e. world domination

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:57):

but other sparse graphs work too

view this post on Zulip Mario Carneiro (Mar 15 2019 at 04:57):

right now it's just an almost completely disconnected graph

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 08:39):

I swear that was the happiest pi day of my life.

view this post on Zulip Mario Carneiro (Mar 15 2019 at 08:42):

this was a pretty good pi day activity

view this post on Zulip Johan Commelin (Mar 15 2019 at 08:57):

It would be cool to see some end result go into mathlib. Does that make sense?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 09:02):

yes, I think it should certainly go in mathlib

view this post on Zulip Scott Morrison (Mar 15 2019 at 09:02):

Seems pretty reasonable. We can chuck out old estimates that have been superseded (both in terms of precision of the bounds and ideas), but the best available bounds sound like a useful thing.

view this post on Zulip Mario Carneiro (Mar 15 2019 at 09:03):

we might want to revisit the methods later, but something is better than nothing (and this method is pretty good already)

view this post on Zulip Mario Carneiro (Mar 15 2019 at 09:03):

I don't necessarily even think that we should only provide the best bounds. Tighter bounds are harder to calculate and harder to use

view this post on Zulip Mario Carneiro (Mar 15 2019 at 09:04):

In fact, as I argued in this thread, the best way to prove these kinds of inequalities is by controlled relaxation of the bounds to something simpler

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 09:31):

https://xenaproject.wordpress.com/2019/03/15/happy-pi-day/

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 09:31):

Comments welcome as ever.

view this post on Zulip Johan Commelin (Mar 15 2019 at 09:51):

@Kevin Buzzard Here is a link to the archive: https://leanprover-community.github.io/archive/113488general/89755happypiday.html
It doesn't require login

view this post on Zulip Kenny Lau (Mar 15 2019 at 09:53):

it's hard to generate good certs...

view this post on Zulip Kenny Lau (Mar 15 2019 at 09:53):

import math

# pi > 2 ^ (n+1) * sqrt (2 - sqrt_two_add_series 0 n)
# pi < 2 ^ (n+1) * sqrt (2 - sqrt_two_add_series 0 n) + 1 / 4 ^ n
# sqrt_two_add_series(x,n+1) = sqrt(2+sqrt_two_add_series(x,n))
def sqrt_two_add_series(x,n):
    while n: x=math.sqrt(2+x); n-=1
    return x

for n in range(10):
    print("n=%d:\t%.10f < pi < %.10f"%(n, 2**(n+1)*math.sqrt(2-sqrt_two_add_series(0,n)),
        2**(n+1)*math.sqrt(2-sqrt_two_add_series(0,n))+1/4**n))

"""OUTPUT:
n=0:    2.8284271247 < pi < 3.8284271247
n=1:    3.0614674589 < pi < 3.3114674589
n=2:    3.1214451523 < pi < 3.1839451523
n=3:    3.1365484905 < pi < 3.1521734905
n=4:    3.1403311570 < pi < 3.1442374070
n=5:    3.1412772509 < pi < 3.1422538134
n=6:    3.1415138011 < pi < 3.1417579418
n=7:    3.1415729404 < pi < 3.1416339755
n=8:    3.1415877253 < pi < 3.1416029841
n=9:    3.1415914215 < pi < 3.1415952362
"""

# to prove: 3.14 <=  2**5 * math.sqrt(2 - sqrt_two_add_series(0,4))
# to prove: 157/1600 <=  math.sqrt(2 - sqrt_two_add_series(0,4))
# to prove: 24649/2560000 <=  2 - sqrt_two_add_series(0,4)
# to prove: sqrt_two_add_series(0,4) <= 5095351/2560000
# rule: 2+x<=y^2 |- sqrt_two_add_series(x,n+1) <= sqrt_two_add_series(y,n)

print(sqrt_two_add_series(0,4))
print(sqrt_two_add_series(577/408,3))
print(sqrt_two_add_series(2537/1373,2))
print(sqrt_two_add_series(15163/7730,1))
print(sqrt_two_add_series(1447/727,0))
print(5095351/2560000)

"""OUTPUT:
1.9903694533443939
1.9903694901455364
1.9903707033912401
1.9903713890865138
1.9903713892709767
1.990371484375
"""

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 09:58):

Use continued fractions

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 09:58):

That is a cheap yet powerful method to approximate reals by rationals

view this post on Zulip Kenny Lau (Mar 15 2019 at 10:01):

that's what I used

view this post on Zulip Kenny Lau (Mar 15 2019 at 10:03):

here are some smaller certs:

print(sqrt_two_add_series(0,4))
print(sqrt_two_add_series(577/408,3))
print(sqrt_two_add_series(619/335,2))
print(sqrt_two_add_series(2297/1171,1))
print(sqrt_two_add_series(1447/727,0))
print(5095351/2560000)

"""OUTPUT:
1.9903694533443939
1.9903694901455364
1.9903695896706308
1.9903696406887783
1.9903713892709767
1.990371484375
"""

view this post on Zulip Kenny Lau (Mar 15 2019 at 10:08):

here are even smaller certs:

print(sqrt_two_add_series(0,4))
print(sqrt_two_add_series(577/408,3))
print(sqrt_two_add_series(619/335,2))
print(sqrt_two_add_series(1940/989,1))
print(sqrt_two_add_series(1447/727,0))
print(5095351/2560000)

"""OUTPUT:
1.9903694533443939
1.9903694901455364
1.9903695896706308
1.9903711590704518
1.9903713892709767
1.990371484375
"""

view this post on Zulip Kenny Lau (Mar 15 2019 at 10:18):

@Floris van Doorn I've made a PR to your file

view this post on Zulip Kenny Lau (Mar 15 2019 at 10:48):

even smaller:

print(sqrt_two_add_series(0,4))
print(sqrt_two_add_series(338/239,3))
print(sqrt_two_add_series(704/381,2))
print(sqrt_two_add_series(1940/989,1))
print(sqrt_two_add_series(1447/727,0))
print(2-(3.14/2**5)**2)

"""OUTPUT:
1.9903694533443939
1.990369667837741
1.9903700913568452
1.9903711590704518
1.9903713892709767
1.990371484375
"""

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:05):

semi-brute-forced the (locally) smallest certs:

print(sqrt_two_add_series(0,4))
print(sqrt_two_add_series(99/70,3))
print(sqrt_two_add_series(874/473,2))
print(sqrt_two_add_series(1940/989,1))
print(sqrt_two_add_series(1447/727,0))
print(2-(3.14/2**5)**2)

"""OUTPUT:
1.9903694533443939
1.9903707035225453
1.9903708019895043
1.9903711590704518
1.9903713892709767
1.990371484375
"""

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:05):

fitness = sum of all the numerators and denominators

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:09):

def sqrt_two_add_series(x,n):
    while n: x=(2+x)**0.5; n-=1
    return x

def gcd(a,b):
    while a: a,b=b%a,a
    return b

def approx(lo,hi,N):
    res = []
    for denom in range(N):
        num = int(hi*denom)
        if gcd(num,denom) == 1 and lo <= num/denom:
            res.append((num,denom))
    return res

N=5000

res = [[L] for L in approx(sqrt_two_add_series(0,4), 2-(3.14/2**5)**2, N)]

for i in [3,2,1]:
    lo = sqrt_two_add_series(0,i)
    new = []
    for L in res:
        size = sum(num+denom for num,denom in L)
        hi = (L[0][0]/L[0][1])**2-2
        for num,denom in approx(lo,hi,N):
            if num+denom+size < 7000:
                new.append([(num,denom)] + L)
    res = new

min_size = 7000
for L in res:
    size = sum(num+denom for num,denom in L)
    if size < min_size:
        print(size, L)
        min_size = size

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:12):

however the bottleneck is actually proving:

sqrt_two_add_series (1447 / 727) 0 ≤ 2 - (157 / 50 / 2 ^ (4 + 1)) ^ 2

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:12):

we didn't certify the right hand side

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 11:13):

I think a Chinese mathematician 1000 years ago used 3•2^n-gons to estimate pi up to 7 digits?

@Kenny Lau do you think it's feasible to get up to 7 digits in Lean?

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:13):

yes

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:22):

@Mario Carneiro how would you prove ⊢ 1447 / 727 ≤ 2 - (157 / 50 / 2 ^ (4 + 1)) ^ 2?

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 11:24):

Are these rationals or reals?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:24):

There was a fifth cert in my version that floris omitted

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 11:24):

I mean, I guess they're reals. Would it make life easier if they were rationals?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:25):

you want to insert an intermediate: 1447 / 727 ≤ 2 - a^2 and a <= 157 / 50 / 2 ^ (4 + 1)

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:26):

good work on the cert optimization tho

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:26):

I mean, I guess they're reals. Would it make life easier if they were rationals?

no, you almost never want to use the decidable instance

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:26):

you want to insert an intermediate: 1447 / 727 ≤ 2 - a^2 and a <= 157 / 50 / 2 ^ (4 + 1)

ok lemme work on it

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 11:27):

a=0 works great

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:29):

wait what? Something must be wrong then

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:30):

oh, a >= 157 / 50 / 2 ^ (4 + 1)

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:30):

like I mentioned before, sign headaches are real here

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 11:31):

246/2507

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 11:32):

@Kenny Lau

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:33):

(157/50/2^5 = 157/1600 which works better, but actually this is not the intermediate we want)

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:33):

lemma pi_gt_314 : pi > 3.14 :=
begin
  refine lt_of_le_of_lt _ (pi_gt_sqrt_two_add_series 4), rw [mul_comm],
  apply le_mul_of_div_le, norm_num, apply le_sqrt_of_sqr_le, rw [le_sub],
  refine le_trans (sqrt_two_add_series_step_up (99/70) (by norm_num) (by norm_num)) _,
  refine le_trans (sqrt_two_add_series_step_up (874/473) (by norm_num) (by norm_num)) _,
  refine le_trans (sqrt_two_add_series_step_up (1940/989) (by norm_num) (by norm_num)) _,
  refine le_trans (sqrt_two_add_series_step_up (1447/727) (by norm_num) (by norm_num)) _,
  refine le_trans (by norm_num : (1447 / 727 : )  2 - 141 / 14644) _,
  rw sub_le_sub_iff_left, norm_num
end

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:34):

this is 2 seconds faster on my computer (6s vs 8s)

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:35):

@Mario Carneiro do you see how to make this faster?

⊢ (157 / 50 / 2 ^ (4 + 1)) ^ 2 ≤ 141 / 14644

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:35):

same trick

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:35):

prove a^2 <= bla and bla <= a

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:37):

By the way, the actual optimization criterion is the number of theorem applications in the norm_num proofs, although this is crazy hard to predict

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:39):

@Mario Carneiro but like bla would just be 141/14644

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:39):

I already optimized this

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:39):

no, you are optimizing a

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:40):

157 / 50 / 2 ^ (4 + 1) <= magic and magic^2 <= 141/14644

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:40):

the magic would just be 157/50/2^5 which is 157/1600

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:40):

I already reached the limit I think

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:41):

confirmed

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 11:41):

(deleted)

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:41):

yeah okay then it's just a lot of working out

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:42):

If this is really expensive compared to the other steps, you could try making 141/14644 a bit larger and see if you can compensate with the others

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 11:44):

Is clearing denominators expensive? There's a theorem which says you can do that

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:44):

Also I'm not sure but it might help to prepare the division cancellation steps since all of the magic numbers are fractions

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:44):

yeah that

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:44):

what does that mean?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:45):

meaning your intermediate theorems, the ones that you stick the magic numbers into, have side goals that talk only about multiplication and addition of natural numbers

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:46):

what goal in particular?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:47):

uh, I'm not sure what the exact form of the theorems being applied are

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:47):

but you can have the magic numbers be \u a / \u b instead of a real fraction

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:49):

bottlenecks:

⊢ 2 + 0 ≤ (99 / 70) ^ 2
⊢ 2 + 99 / 70 ≤ (874 / 473) ^ 2
⊢ 2 + 874 / 473 ≤ (1940 / 989) ^ 2
⊢ 2 + 1940 / 989 ≤ (1447 / 727) ^ 2
⊢ 1447 / 727 ≤ 2 - 141 / 14644
⊢ (157 / 50 / 2 ^ (4 + 1)) ^ 2 ≤ 141 / 14644

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:49):

and then the side goals are things like a^2 * d <= (c + 2 d) * b^2 instead of (a/b)^2 <= c/d + 2

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:49):

why would this be faster?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:50):

Depending on how norm_num unfolds the things, it might work out larger intermediates than necessary

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:50):

hmm

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:51):

In theory, norm_num should be doing the right thing, and so this won't actually gain much

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:51):

You might consider looking at the proof to profile it

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:52):

(split off the norm_num subproofs to lemmas)

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:54):

the bottlenecks are essentially uncertifiable right?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:54):

Well, at some point you actually have to prove a tight bound

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 11:54):

Kenny how do you time how long a command takes?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:55):

Generally speaking, the difficulty of the proof (size of the numbers) is inversely proportional to the size of the gap

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:55):

and we have placed 5 intermediate points in the original gap

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:55):

so if we shift the points around a bit we can equalize the load

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:55):

for ⊢ 2 + 99 / 70 ≤ (874 / 473) ^ 2, the largest number norm_num produces in the proof is 361313348

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:56):

that sounds about right, for those numbers

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:56):

(norm_num.subst_into_div (763876 / 473) 473 (874 * (874 / 473)) 473 (763876 / 223729)
                       (norm_num.div_eq_div_helper (763876 / 473) 473 763876 223729 361313348
                          (norm_num.subst_into_prod (763876 / 473) 223729 (763876 / 473) 223729 361313348
                             (norm_num.subst_into_div 763876 473 763876 473 (763876 / 473)
                                (norm_num.div_eq_div_helper 763876 473 763876 473 361313348
                                   (norm_num.subst_into_prod 763876 473 763876 473 361313348 (eq.refl 763876)
                                      (eq.refl 473)
                                      (norm_num.mul_bit1_helper 763876 236 180274736 361313348
                                         (norm_num.mul_bit0_helper 763876 118 90137368

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:56):

it should be 874^2 * 70

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:56):

well it isn't

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 11:56):

It's something times 473

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 11:57):

it's 874^2 * 473

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:57):

wat

view this post on Zulip Kenny Lau (Mar 15 2019 at 11:57):

how did that number get there

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:58):

I should note this is the one part of norm_num I didn't write

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:58):

I offload + - * / to leo's impl

view this post on Zulip Mario Carneiro (Mar 15 2019 at 11:59):

but it's not like there are that many ways to write a (correct) arithmetic prover

view this post on Zulip Kenny Lau (Mar 15 2019 at 12:00):

anyway I'm going now, this was fun

view this post on Zulip Mario Carneiro (Mar 15 2019 at 12:00):

good work

view this post on Zulip Mario Carneiro (Mar 15 2019 at 12:04):

aha:

variables {α : Type*} [field α]
#check @norm_num.div_eq_div_helper α _ (763876 / 473) 473 763876 223729 361313348
-- norm_num.div_eq_div_helper (763876 / 473) 473 763876 223729 361313348 :
--   763876 / 473 * 223729 = 361313348 →
--   763876 * 473 = 361313348 → 473 ≠ 0 → 223729 ≠ 0 → 763876 / 473 / 473 = 763876 / 223729

It is cross-multiplying to prove that (874 * (874 / 473)) / 473 = (763876 / 473) / 473 = 763876 / 223729

view this post on Zulip Mario Carneiro (Mar 15 2019 at 12:04):

that's... not a good idea for simplifying divisions of divisions

view this post on Zulip Kevin Buzzard (Mar 15 2019 at 12:06):

PR to core Kenny!

view this post on Zulip Floris van Doorn (Mar 15 2019 at 13:38):

At the bottom of the file https://github.com/fpvandoorn/mathlib/blob/pi/src/data/real/pi.lean I also computed pi up to 5 decimal digits. It does take 20-25 seconds to compile for me though.

view this post on Zulip Floris van Doorn (Mar 15 2019 at 13:39):

It's definitely possible to get better approximations, but the compilation time will just be bigger for that. So I don't think we want the tightest bounds in mathlib. But soon I'll PR the first 3 digits of pi to mathlib.

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:16):

what is the profiler result? Is it norm_num or the kernel taking up the time?

view this post on Zulip Floris van Doorn (Mar 15 2019 at 14:22):

Not sure how to read the output, here it is for pi_gt_31415:

elaboration of pi_gt_31415 took 11.9s
elaboration: tactic execution took 8.56s
num. allocated objects:  8023
num. allocated closures: 9813
num. allocated big nums: 96
 8560ms   100.0%   tactic.step
 8560ms   100.0%   scope_trace
 8560ms   100.0%   tactic.istep
 8560ms   100.0%   _interaction._lambda_2
 8560ms   100.0%   tactic.istep._lambda_1
 5972ms    69.8%   tactic.refine
 5945ms    69.5%   tactic.to_expr
 2583ms    30.2%   interaction_monad.orelse'
 2580ms    30.1%   tactic.seq
 1641ms    19.2%   interaction_monad_orelse
 1621ms    18.9%   tactic.interactive.norm_num1
 1620ms    18.9%   tactic.replace_at._lambda_4
 1620ms    18.9%   tactic.alternative._lambda_3
 1620ms    18.9%   tactic.replace_at
 1615ms    18.9%   interaction_monad.monad._lambda_9
 1614ms    18.9%   tactic.ext_simplify_core
 1614ms    18.9%   norm_num.derive
 1607ms    18.8%   tactic.repeat
 1607ms    18.8%   _private.3096317357.repeat_aux._main._lambda_1
 1607ms    18.8%   tactic.all_goals
 1607ms    18.8%   _private.3321893233.all_goals_core._main._lambda_2
 1607ms    18.8%   tactic.try_core
 1607ms    18.8%   all_goals_core
 1607ms    18.8%   repeat_aux
 1548ms    18.1%   norm_num.derive._main._lambda_3
 1364ms    15.9%   tactic.norm_num
  959ms    11.2%   tactic.interactive.simp_core
  937ms    10.9%   tactic.join_user_simp_lemmas
  937ms    10.9%   tactic.mk_simp_set_core
  937ms    10.9%   simp_lemmas.mk_default
  937ms    10.9%   tactic.mk_simp_set
[...]

view this post on Zulip Floris van Doorn (Mar 15 2019 at 14:24):

So most of the tactic time is in refine, so building the proof term(?). But there is also 3-4 seconds of time which is not used to execute tactics. Is that just to type-check the final proof?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:30):

I think you missed out on some part... there should be a line for kernel typechecking

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:30):

it's next to the elaboration of foo took n s line

view this post on Zulip Floris van Doorn (Mar 15 2019 at 14:30):

Yeah, I thought so too. I cannot find it anywhere in my Lean messages

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:31):

maybe that just means it was below the threshold

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:31):

but damn that's a lot of to_expr

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:32):

and istep too... the tactic overhead is significant

view this post on Zulip Floris van Doorn (Mar 15 2019 at 14:32):

Oh no, I know what is going on: all the tactic.refine are the nested norm_nums

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:33):

The way you read this is to look for nodes that come just before a percentage drop

view this post on Zulip Floris van Doorn (Mar 15 2019 at 14:33):

When we execute

refine sqrt_two_add_series_step_up (11482/8119)  (by norm_num) (by norm_num) _,

the inner norm_num reports that the inner norm_num takes up almost a second:

713ms    98.5%   norm_num.derive

but for the outer tactic block, this is all part of the refine tactic.

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:34):

refine is quick, it is just calling another tactic that is slow

view this post on Zulip Floris van Doorn (Mar 15 2019 at 14:34):

ok, sure. But the to_expr just calls the nested norm_num tactics, right?

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:35):

yes, but that should show up in the total profile

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:35):

1614ms    18.9%   norm_num.derive
1607ms    18.8%   tactic.repeat
...

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:36):

compared to 8560ms for all the tactics, that's not as large a fraction as I would have hoped

view this post on Zulip Floris van Doorn (Mar 15 2019 at 14:37):

Mario Carneiro: yes, but that should show up in the total profile

I don't think it does. The sum of the times spend in norm_num.derive in the nested by norm_num tactics is much larger than the reported time on norm_num.derive when I hover over the end of the tactic block

view this post on Zulip Floris van Doorn (Mar 15 2019 at 14:38):

If I add everything up, about 5.3s is spend in norm_num.derive

view this post on Zulip Floris van Doorn (Mar 15 2019 at 14:39):

A bit more actually, but between 5s and 6s.

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:44):

It's surprising how inconclusive this test is, but it looks like the outer block does in fact reflect slow downs caused by inner blocks

set_option profiler true

example : true :=
begin
  tactic.sleep 10000,
  refine begin
    tactic.sleep 10000,
    trivial
  end
end

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:44):

either the profiler is way off, or sleep is not very accurate

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:46):

elaboration: tactic compilation took 17.8ms
test.lean:10:0: information tactic profile data
elaboration: tactic execution took 14.4s
num. allocated objects:  15
num. allocated closures: 46
14404ms   100.0%   tactic.istep
14404ms   100.0%   _interaction._lambda_2
14404ms   100.0%   scope_trace
14404ms   100.0%   tactic.istep._lambda_1
14404ms   100.0%   tactic.step
 7470ms    51.9%   tactic.to_expr
 7470ms    51.9%   tactic.refine
 6934ms    48.1%   tactic.sleep

view this post on Zulip Floris van Doorn (Mar 15 2019 at 14:46):

That seems pretty conclusive to me... The bottom end reports 1.7s spend at sleep and 5.2s at refine, the top end reports 5.2s spend at sleep. But yeah, those numbers are not very close to 2000 and 4000...

view this post on Zulip Floris van Doorn (Mar 15 2019 at 14:47):

(for future reference: the sleep numbers in Mario's post were originally 2000 and 4000)

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:47):

even when I crank the numbers up they are still not so close

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:47):

I think the units must be off

view this post on Zulip Floris van Doorn (Mar 15 2019 at 14:48):

yeah, there must be something off with sleep. But it's clear that the inner sleep time only counts towards refine, not towards sleep.

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:49):

yeah. To make it even more obvious:

example :=
begin
  have := begin
    tactic.sleep 2000,
    trivial
  end,
  refine begin
    tactic.sleep 2000,
    trivial
  end
end
elaboration: tactic execution took 2.8s
num. allocated objects:  21
num. allocated closures: 52
 2798ms   100.0%   tactic.to_expr
 2798ms   100.0%   tactic.istep
 2798ms   100.0%   _interaction._lambda_2
 2798ms   100.0%   scope_trace
 2798ms   100.0%   tactic.istep._lambda_1
 2798ms   100.0%   tactic.step
 1429ms    51.1%   tactic.interactive.have._lambda_1
 1369ms    48.9%   tactic.refine

view this post on Zulip Mario Carneiro (Mar 15 2019 at 14:52):

There is an easy solution to this for your situation though: don't use nested tactics, use a refine with more underscores and kill the subgoals with norm_num in the main tactic block

view this post on Zulip Floris van Doorn (Mar 15 2019 at 15:00):

If I do that, we see that indeed most of the time is spend at norm_num:

elaboration: tactic execution took 7.48s
num. allocated objects:  29602
num. allocated closures: 30223
num. allocated big nums: 322
 7482ms   100.0%   tactic.istep._lambda_1
 7482ms   100.0%   _interaction._lambda_2
 7482ms   100.0%   tactic.istep
 7482ms   100.0%   tactic.step
 7482ms   100.0%   scope_trace
 7442ms    99.5%   interaction_monad.orelse'
 7439ms    99.4%   tactic.seq
 7420ms    99.2%   interaction_monad.monad._lambda_9
 7388ms    98.7%   tactic.all_goals
 7388ms    98.7%   _private.3321893233.all_goals_core._main._lambda_2
 7388ms    98.7%   all_goals_core
 6640ms    88.7%   interaction_monad_orelse
 6619ms    88.5%   tactic.replace_at._lambda_4
 6619ms    88.5%   tactic.replace_at
 6619ms    88.5%   tactic.alternative._lambda_3
 6619ms    88.5%   tactic.interactive.norm_num1
 6563ms    87.7%   norm_num.derive
 6563ms    87.7%   tactic.ext_simplify_core
 6499ms    86.9%   norm_num.derive._main._lambda_3
 5696ms    76.1%   tactic.norm_num
 2613ms    34.9%   norm_num.eval_ineq
 2608ms    34.9%   norm_num.prove_lt
 1357ms    18.1%   _private.3096317357.repeat_aux._main._lambda_1
[...]

view this post on Zulip Floris van Doorn (Mar 15 2019 at 15:02):

Interestingly, there is a drop after all_goals_core (I'm now accumulating all the inequalities as goals, and killing them all with all_goals {norm_num} at the end)

view this post on Zulip Floris van Doorn (Mar 15 2019 at 15:06):

oh, in the other proof (pi_lt_31416) the drop-offs from 100% to 90% are at very different tactics. So maybe that is just noise/semi-randomness.

view this post on Zulip Floris van Doorn (Mar 15 2019 at 15:28):

Since people requested it, I also computed the first 7 digits of pi at the bottom of https://github.com/fpvandoorn/mathlib/blob/pi/src/data/real/pi.lean