Re: [math-fun] Compute a Function From its Fourier Series?
GFee>See the following paper.
Authors: Gottlieb and Shu Title: On the Gibb's phenomenon and its resolution Journal: SIAM Review Issue: Vol.39, No. 4, pages 644-668, December 1997
Also see http://www.tweedledum.com/rwg/gibbs.htm --rwg --------------------------------- Tonight's top picks. What will you watch tonight? Preview the hottest shows on Yahoo! TV.
Aren't you just reducing the "slew" ? I seem to recall that there are theorems relating the height of the "ears" and the steepness of the ramp. --- Didn't Gibbs sing for the Bee Gees ? At 12:30 AM 9/30/2007, Bill Gosper wrote:
GFee>See the following paper.
Authors: Gottlieb and Shu Title: On the Gibb's phenomenon and its resolution Journal: SIAM Review Issue: Vol.39, No. 4, pages 644-668, December 1997
Also see http://www.tweedledum.com/rwg/gibbs.htm --rwg
No time right now to compute, but I think my ears are *way* smaller than from a ramp of as wide as the interpolated step. It's a resonance effect. Widening the step will make it *worse*. --rwg Henry Baker <hbaker1@pipeline.com> wrote: Aren't you just reducing the "slew" ? I seem to recall that there are theorems relating the height of the "ears" and the steepness of the ramp. --- Didn't Gibbs sing for the Bee Gees ? At 12:30 AM 9/30/2007, Bill Gosper wrote:
GFee>See the following paper.
Authors: Gottlieb and Shu Title: On the Gibb's phenomenon and its resolution Journal: SIAM Review Issue: Vol.39, No. 4, pages 644-668, December 1997
Also see http://www.tweedledum.com/rwg/gibbs.htm --rwg
--------------------------------- Yahoo! oneSearch: Finally, mobile search that gives answers, not web links.
Joshua Zucker plugged http://sanjosemathcircle.org . Would they be willing guinea pigs for my latest disk puzzle? (http://gosper.org/Arnold.rtf .) I'm anxious to know if it has multiple solutions. Regrettably, Octave (http://gosper.org/octav2.rtf) and Twubblesome Twelve (home page) do. http://gosper.org/7^3leg.rtf has 7.9 million solutions (thanks for counting, Mike Beeler) but is nevertheless difficult (and more instructive than the disk puzzles, which are just hard). This invitation extends to all funsters handy to the South Bay. --rwg PS, I assume the Math Circle leaders have checked out www.tweedledum.com/rwg/ for useful material. --------------------------------- Shape Yahoo! in your own image. Join our Network Research Panel today!
Here is part of my private response to Alan:
It's not a strictly kosher solution, however, because there are definitely two non-triangular interstices.
You have correctly guessed that I intended only tricuspid voids, and assumed that no other solution would fit. [...] This puzzle is probably more evil than you think. The "oval" has no axis of symmetry, and the tolerances are very close. In the intended solution, the pieces slide in easily, but it is difficult to insert a strip of paper between a piece and the edge of the cavity. [The puzzle is nearly 1' in diameter.] So I had to try all four reflections of your solution, and one very nearly fit, but not quite. Sadly, however, I woke up today realizing how I had wasted space, and found several bogus solutions inside five minutes. Back to the drawing board. So far, my only disk puzzle where the only solution is tricuspid is the circular http://gosper.org/Huddle.RTF . This was only Arnold Laidanegger. The real Dozenegger is still in the gym. --rwg --------------------------------- Looking for a deal? Find great prices on flights and hotels with Yahoo! FareChase.
Last night Alan Schoen sent me a solution so much more efficient than my tricuspid that it worked in all four orientations! And here I thought line thickness and the weirdness of the cavity would prevent even copying, let alone solution, given just that lo-res drawing . I'm tempted to ask Alan to sit on his arrangement while I try to "tricuspidize" it. --rwg --------------------------------- Yahoo! oneSearch: Finally, mobile search that gives answers, not web links.
Bill: Be my guest! I'm more than happy to leave the tricuspidizing to you. I'll bet you succeed. Alan Bill Gosper wrote:
Last night Alan Schoen sent me a solution so much more efficient than my tricuspid that it worked in all four orientations! And here I thought line thickness and the weirdness of the cavity would prevent even copying, let alone solution, given just that lo-res drawing . I'm tempted to ask Alan to sit on his arrangement while I try to "tricuspidize" it. --rwg
--------------------------------- Yahoo! oneSearch: Finally, mobile search that gives answers, not web links. _______________________________________________ math-fun mailing list math-fun@mailman.xmission.com http://mailman.xmission.com/cgi-bin/mailman/listinfo/math-fun
Bill Gosper <rwmgosper@yahoo.com> wrote: Last night Alan Schoen sent me a solution so much more efficient than my tricuspid that it worked in all four orientations! And here I thought line thickness and the weirdness of the cavity would prevent even copying, let alone solution, given just that lo-res drawing . I'm tempted to ask Alan to sit on his arrangement while I try to "tricuspidize" it. --rwg Partial success: A huge effort produced a tricuspid packing, but with undesirable radii. The good news is that there are three degrees of freedom in the equations. Small matter of additional effort. The main difficulty is with arnoldp(x,y,r), a 3367 term polynomial which vanishes if a circle at (x,y) of radius r is tangent to the ovoid. Unfortunately, not iff. The only way I could eliminate extraneous parameters and reduce arnoldp to three variables was to numerically approximate the coefficients. Then the elimination process (aided by Mathematica, since my copy of Macsyma was compiled for a 286 in NT) introduced huge spurious factors. The false solutions are raising hell with Newton's method, several of whose 36 equations are arnoldps. It's taking forever, and screwing up. The 3367 floating point coefficients preclude symbolically extracting the bogus factors. But DUH! With a fish dinner and a shot of chai, I just realized an obvious way for Arnold to visit Jenny Craig. If any funsters chirp up with the answer, I'll kick myself for not requesting help earlier. If not, stay tuned for more arnage. --rwg PS, Cool mental factorization of 3367 (one third, two thirds): 10101/3 = 111111/3/11 = 111/3 1001/11 = 37 7 13. --------------------------------- Looking for last minute shopping deals? Find them fast with Yahoo! Search.
I said Partial success: A huge effort produced a tricuspid packing, but with undesirable radii. The good news is that there are three degrees of freedom in the equations. Small matter of additional effort. The main difficulty is with arnoldp(x,y,r), a 3367 term polynomial which vanishes if a circle at (x,y) of radius r is tangent to the ovoid. Unfortunately, not iff. The only way I could eliminate extraneous parameters and reduce arnoldp to three variables was to numerically approximate the coefficients. Then the elimination process (aided by Mathematica, since my copy of Macsyma was compiled for a 286 in NT) introduced huge spurious factors. The false solutions are raising hell with Newton's method, many of whose 36 equations are arnoldps. It's taking forever, and screwing up. The 3367 floating point coefficients preclude symbolically extracting the bogus factors. But DUH! With a fish dinner and a shot of chai, I just realized an obvious way for Arnold to visit Jenny Craig. If any funsters chirp up with the answer, I'll kick myself for not requesting help earlier. If not, stay tuned for more arnage. No chirpers. It so happened that when given x and y, the desired r was always the least positive real root of arnoldp(x,y,r). This is of no help to the Newton problem where x, y, and r are all unknown, but it lets you to construct a couple of hundred valid (x,y,r) triples, to which you fit successively larger polynomials until they start successfully predic(a)ting new cases. It turns out that 161 terms, of degree 10 in x and y, and 5 in r^2, do the trick. (Although I am disconcerted by the large remainder I get from dividing the old arnoldp by the new one.) The resulting hundredfold speedup in the search process led to a promising design. http://gosper.org/Dozenegger.rtf (2MB)-: But instead of running off to laser-cut fancy new Dozeneggers, I'm just making a set of pieces from thin scrap to test in the existing Arnolds, having been cooked once too often by unintended solutions. Alan Schoen will no doubt find a solution, which I pray is the unique, intended one. "This time for sure." --rwg --------------------------------- Looking for last minute shopping deals? Find them fast with Yahoo! Search.
Simply sumswapping the constant term in the Weierstrass P double series, 1 inf y + - ==== y \ y 1 1 ----- - > (------------- + -------------) = --. 6 / 2 2 m pi pi ==== sinh (m pi y) sinh (----) y m = 1 y E.g., for y=1, inf ==== 1 \ 1 1 - - 2 > ----------- = --. 3 / 2 pi ==== sinh (m pi) m = 1 Combining with y=1+i, inf ==== \ 1 1 2 > ----------------- = -- / 2 1 pi ==== cosh ((m - -) pi) m = 1 2 With y=sqrt 2 inf ==== 1 \ cosh(sqrt(2) m pi) + 2 sqrt(2) - - 2 > ---------------------- = -------. 2 / 2 pi ==== sinh (sqrt(2) m pi) m = 1 With y=1 and 2, inf ==== \ 2 cosh(m pi) 1 1 1 > (------------ - -----------) = - - ----. / 2 2 3 2 pi ==== sinh (m pi) cosh (m pi) m = 1 With y=cis x, inf ==== i x cos(x) \ e 1 ------ - 2 > Re ---------------- = --. 3 / 2 i x pi ==== sinh (pi m e ) m = 1 With x= pi/6, inf ==== \ 1 1 1 > ------------- = - - ----------. / 2 6 sqrt(3) pi ==== sinh (m pi y) m = 1 (There are doubtless better ways to take log(e^pi).) --rwg --------------------------------- Never miss a thing. Make Yahoo your homepage.
On 17 Oct 93 on this list I mentioned getting B_n(x) from inverting a matrix made from Pascal's triangle. Here it is with some wrinkles. The Pascal matrix: (c1039) block([gcd:false], genmatrix(lambda([i,j],if j>i then 0 else binomial(i,j-1)/i/x^(i-1)),8)) [ 1 0 0 0 0 0 0 0 ] [ ] [ 1 2 ] [ --- --- 0 0 0 0 0 0 ] [ 2 x 2 x ] [ ] [ 1 3 3 ] [ ---- ---- ---- 0 0 0 0 0 ] [ 2 2 2 ] [ 3 x 3 x 3 x ] [ ] [ 1 4 6 4 ] [ ---- ---- ---- ---- 0 0 0 0 ] [ 3 3 3 3 ] [ 4 x 4 x 4 x 4 x ] [ ] [ 1 5 10 10 5 ] (d1039) [ ---- ---- ---- ---- ---- 0 0 0 ] [ 4 4 4 4 4 ] [ 5 x 5 x 5 x 5 x 5 x ] [ ] [ 1 6 15 20 15 6 ] [ ---- ---- ---- ---- ---- ---- 0 0 ] [ 5 5 5 5 5 5 ] [ 6 x 6 x 6 x 6 x 6 x 6 x ] [ ] [ 1 7 21 35 35 21 7 ] [ ---- ---- ---- ---- ---- ---- ---- 0 ] [ 6 6 6 6 6 6 6 ] [ 7 x 7 x 7 x 7 x 7 x 7 x 7 x ] [ ] [ 1 8 28 56 70 56 28 8 ] [ ---- ---- ---- ---- ---- ---- ---- ---- ] [ 7 7 7 7 7 7 7 7 ] [ 8 x 8 x 8 x 8 x 8 x 8 x 8 x 8 x ] The first approximation to its inverse: (c1040) (t^(1..8),coefmatrix(%%*x^(0..7),%%)) [ 1 0 0 0 0 0 0 0 ] [ ] [ 0 x 0 0 0 0 0 0 ] [ ] [ 2 ] [ 0 0 x 0 0 0 0 0 ] [ ] [ 3 ] [ 0 0 0 x 0 0 0 0 ] [ ] (d1040) [ 4 ] [ 0 0 0 0 x 0 0 0 ] [ ] [ 5 ] [ 0 0 0 0 0 x 0 0 ] [ ] [ 6 ] [ 0 0 0 0 0 0 x 0 ] [ ] [ 7 ] [ 0 0 0 0 0 0 0 x ] Newton step #1 for inverse: (c1041) 2*%-% . d1039 . % [ 1 0 0 0 0 0 0 0 ] [ ] [ 1 ] [ - - x 0 0 0 0 0 0 ] [ 2 ] [ ] [ 1 2 ] [ - - - x x 0 0 0 0 0 ] [ 3 ] [ ] [ 2 ] [ 1 3 x 3 ] [ - - - x - ---- x 0 0 0 0 ] [ 4 2 ] [ ] (d1041) [ 1 2 3 4 ] [ - - - x - 2 x - 2 x x 0 0 0 ] [ 5 ] [ ] [ 2 3 4 ] [ 1 5 x 10 x 5 x 5 ] [ - - - x - ---- - ----- - ---- x 0 0 ] [ 6 2 3 2 ] [ ] [ 1 2 3 4 5 6 ] [ - - - x - 3 x - 5 x - 5 x - 3 x x 0 ] [ 7 ] [ ] [ 2 4 6 ] [ 1 7 x 3 35 x 5 7 x 7 ] [ - - - x - ---- - 7 x - ----- - 7 x - ---- x ] [ 8 2 4 2 ] #2: (c1042) 2*%-% . D1039 . % (d1042) [ 1 0 0 0 0 0 0 0 ] [ ] [ 1 ] [ - - x 0 0 0 0 0 0 ] [ 2 ] [ ] [ 1 2 ] [ - - x x 0 0 0 0 0 ] [ 6 ] [ ] [ 2 ] [ x 3 x 3 ] [ 0 - - ---- x 0 0 0 0 ] [ 2 2 ] [ ] [ 23 2 3 4 ] [ - -- 0 x - 2 x x 0 0 0 ] [ 15 ] [ ] [ 3 4 ] [ 25 23 x 5 x 5 x 5 ] [ - -- - ---- 0 ---- - ---- x 0 0 ] [ 4 3 3 2 ] [ ] [ 4 ] [ 1573 75 x 2 5 x 5 6 ] [ - ---- - ---- - 23 x 0 ---- - 3 x x 0 ] [ 84 2 2 ] [ ] [ 2 3 5 6 ] [ 301 1573 x 525 x 161 x 7 x 7 x 7 ] [ - --- - ------ - ------ - ------ 0 ---- - ---- x ] [ 6 12 4 3 2 2 ] Finally, #3: (c1043) 2*%-% . d1039 . % [ 1 0 0 0 0 0 0 0 ] [ ] [ 1 ] [ - - x 0 0 0 0 0 0 ] [ 2 ] [ ] [ 1 2 ] [ - - x x 0 0 0 0 0 ] [ 6 ] [ ] [ 2 ] [ x 3 x 3 ] [ 0 - - ---- x 0 0 0 0 ] [ 2 2 ] [ ] [ 1 2 3 4 ] (d1043) [ - -- 0 x - 2 x x 0 0 0 ] [ 30 ] [ ] [ 3 4 ] [ x 5 x 5 x 5 ] [ 0 - - 0 ---- - ---- x 0 0 ] [ 6 3 2 ] [ ] [ 2 4 ] [ 1 x 5 x 5 6 ] [ -- 0 - -- 0 ---- - 3 x x 0 ] [ 42 2 2 ] [ ] [ 3 5 6 ] [ x 7 x 7 x 7 x 7 ] [ 0 - 0 - ---- 0 ---- - ---- x ] [ 6 6 2 2 ] The row sums are the first eight Bernoulli polynomials. This isn't a bad way to compute them, especially if you use the iteration x: x + (x - x.a.x), where 2^k subdiagonals of the matrix in parens vanish on the kth iteration. And there is no need to compute the unconverged diagonals: (c1112) genmatrix(lambda([i,j],if i-j>3 then [this,crud,doesn\'t,matter,one,bit,believe,it,"or","not"][binomial(i-4,2)+j] else d1042[i,j]),8) (d1112) [ 1 0 0 0 0 0 0 0 ] [ ] [ 1 ] [ - - x 0 0 0 0 0 0 ] [ 2 ] [ ] [ 1 2 ] [ - - x x 0 0 0 0 0 ] [ 6 ] [ ] [ 2 ] [ x 3 x 3 ] [ 0 - - ---- x 0 0 0 0 ] [ 2 2 ] [ ] [ 2 3 4 ] [ this 0 x - 2 x x 0 0 0 ] [ ] [ 3 4 ] [ 5 x 5 x 5 ] [ crud doesn't 0 ---- - ---- x 0 0 ] [ 3 2 ] [ ] [ 4 ] [ 5 x 5 6 ] [ matter one bit 0 ---- - 3 x x 0 ] [ 2 ] [ ] [ 5 6 ] [ 7 x 7 x 7 ] [ believe it or not 0 ---- - ---- x ] [ 2 2 ] (c1113) expand(2*%-% . d1039 . %) [ 1 0 0 0 0 0 0 0 ] [ ] [ 1 ] [ - - x 0 0 0 0 0 0 ] [ 2 ] [ ] [ 1 2 ] [ - - x x 0 0 0 0 0 ] [ 6 ] [ ] [ 2 ] [ x 3 x 3 ] [ 0 - - ---- x 0 0 0 0 ] [ 2 2 ] [ ] [ 1 2 3 4 ] (d1113) [ - -- 0 x - 2 x x 0 0 0 ] [ 30 ] [ ] [ 3 4 ] [ x 5 x 5 x 5 ] [ 0 - - 0 ---- - ---- x 0 0 ] [ 6 3 2 ] [ ] [ 2 4 ] [ 1 x 5 x 5 6 ] [ -- 0 - -- 0 ---- - 3 x x 0 ] [ 42 2 2 ] [ ] [ 3 5 6 ] [ x 7 x 7 x 7 x 7 ] [ 0 - 0 - ---- 0 ---- - ---- x ] [ 6 6 2 2 ] Further hack: only compute the next 2^k terms of the last row, and partially fill in the preceding rows by (trivial) differentiation. Finally, one can double the matrix size each iteration by filling in n diagonals via trivial integration, then running one Newton step. Two days earlier I gave [ n n ] m [ k - m (k - m) (x + k - 1) (x + k - 1) ] /===\ [ ----- -------------------- ------------ ] | | [ k + 1 k + 1 k ] (d241) ( | | [ ]) m = | | [ k - m 1 ] k = 1 [ 0 ----- - ] [ k + 1 k ] [ ] [ 0 0 1 ] [ 0 0 bernpoly(x, n) ] [ ] [ 0 0 harmonic(m) ] [ ] [ 0 0 m ] where m is any integer > n. Then, From conway@math.Princeton.EDU Sat Oct 16 18:55:31 1993 Message-Id: <9310162255.AA14049@broccoli.princeton.edu> To: math-fun@cs.arizona.edu, rwg@NEWTON.macsyma.com Subject: Re: gamma function Those coefficients are just the partial sums of the series (nC0)/n - (nC1)/(n-1) + (nC2)/(n-2) - ... with alternating signs. So this is a FANTASTIC explicit formula [equivalent to k - 1 ==== k \ n (- 1) binomial(n + 1, k) > (x + j) n + 1 / ==== ==== \ j = 0 bernpoly(x, n) = - > ---------------------------------------- / k ==== k = 1 --rwg] for the Bernoulli polynomials. May I have partial credit, or did you get it first? JHC To which I sniffed It's all yours. Only, why is a single product of [matrices of] rational functions less fantastic than a double sum of binomials? it simplifies further! jhc Yes...? ----- But he didn't say. The answer to my sniff is clearly that there are more uses for formulas than mere recipes for calculation. However, it turns out "It's all yours" was too generous. His double sum is term-for-term the m=n+1 case of my matrix product. The product of upper triangular 3x3s n /===\ [ a(k) b(k) c(k) ] | | [ ] | | [ 0 d(k) e(k) ] | | [ ] k = 1 [ 0 0 1 ] computes in its upper right element the triangular double sum j - 1 /===\ | | a(i) n k - 1 k - 1 ( | | ----) b(j) k - 1 ==== /===\ ==== | | d(i) /===\ \ | | a(r) \ i = 1 | | (*) > (c(k) | | ---- + > ----------------- e(k)) | | d(r). / | | d(r) / d(j) | | ==== r = 1 ==== r = 1 k = 1 j = 1 (Notice the pun on "triangular". They're unconnected, since the product of an upper triangular 4x4 computes a tetrahedral triple sum.) When fed the assignments of a(k),...,e(k) in (d241) above, the leftmost product in (*) merges with the inner sum as the j=0 term. There was also one of those ironic email misunderstandings: JHC> The definition I prefer is (in the exponential notation, rwg>Umbrage! JHC> which it's worth using) B^n = (B-1)^n for n > 1. Putting n = 2, and expanding formally (which is what the notation means): B^2 = B^2 - 2B^1 + 1, so B^1 = 1/2, [...], My interjection was just playful recognition of the umbral calculus, but Conway took it to be sincere umbrage to umbrage, and, unfortunately, went somewhat on the defensive. --rwg --------------------------------- Never miss a thing. Make Yahoo your homepage.
participants (3)
-
Alan Schoen -
Bill Gosper -
Henry Baker