GSoC 2013 Bernoulli Numbers

classic Classic list List threaded Threaded
5 messages Options
Reply | Threaded
Open this post in threaded view
|

GSoC 2013 Bernoulli Numbers

David Joy
Hi,

I'm an Open University student completing my Open degree this summer (CompSci / Pure- & Applied Maths / Classical History) and am very interested in the Bernoulli numbers project for the Google Summer of Code. As well as the OU I've done quite a bit of mathematics elsewhere (Durham and KCL). Enclosed is my (naive) response to your challenge to implement the Akiyama–Tanigawa algorithm.

I use the Boost framework for rational numbers, so as not to lose precision (though haven't yet  adapted the solution for arbitrary precision numbers). [It would be more aesthetically pleasing if I could write the main algorithm in a STL-style so that it iterates of the vector, rather than explicitly using an index.]

I have yet to get to grips with Boost.Test, but the enclosed shell script can generate a list of results, which matches the table found on the Bernoulli Number Page out to B_20.

All code criticism gratefully accepted,

David

Bernoulli.cpp
Bernoulli_test.sh
Reply | Threaded
Open this post in threaded view
|

Re: GSoC 2013 Bernoulli Numbers

Christopher Kormanyos

Hi David,


> I'm an Open University student completing my Open degree this summer
> (CompSci / Pure- & Applied Maths / Classical History) and am very interested
> in the Bernoulli numbers project for the Google Summer of Code. As well as
> the OU I've done quite a bit of mathematics elsewhere (Durham and KCL).
> Enclosed is my (naive) response to your challenge to implement the
> Akiyama–Tanigawa algorithm.

> I use the Boost framework for rational numbers, so as not to lose precision
> (though haven't yet  adapted the solution for arbitrary precision numbers).
> [It would be more aesthetically pleasing if I could write the main algorithm
> in a STL-style so that it iterates of the vector, rather than explicitly
> using an index.]

Code looks pretty good, and the use of rational
is adept and creative. It's a good solution.

> I have yet to get to grips with Boost.Test, but the enclosed shell script
> can generate a list of results, which matches the table found on the
> Bernoulli Number Page <http://www.bernoulli.org/>   out to B_20.

Or go to Wolfram's Alpha which understands the language
of Mathematica(R) to a great extent and type in something
like:

Table[N[BernoulliB[i], 51], {i, 0, 20, 1}]

> All code criticism gratefully accepted,

The code is good enough. So there's no more
work needed here.

> Bernoulli.cpp
> <http://boost.2283326.n4.nabble.com/file/n4646219/Bernoulli.cpp
> Bernoulli_test.sh
> <http://boost.2283326.n4.nabble.com/file/n4646219/Bernoulli_test.sh

You need to start working on your proposal.
Start the proposal according to this link
https://svn.boost.org/trac/boost/wiki/SoCSubmissionTemplate<https://svn.boost.org/trac/boost/wiki/SoCSubmissionTemplate>.
This link may also help https://svn.boost.org/trac/boost/wiki/SoCHints

We have more applications for the Boost.Math project and
fewer for the Multiprecision project. Are you interested
in the Multiprecision project? Do you have any experience
with those kinds of algorithms.

Sincerely, Chris.

_______________________________________________
Unsubscribe & other changes: http://lists.boost.org/mailman/listinfo.cgi/boost
Reply | Threaded
Open this post in threaded view
|

[GSoC 2013] Draft Proposal for Bernoulli Numbers project

David Joy
Hi Chris,

Thanks for your reply. I hadn't looked too closely at the Multiprecision project, and it was only after checking out the material on high-precision integer arithmetic at MIT Open Courseware that I realised I had seen something similar once before: a fascinating lecture by the (Fields Medallist) Timothy Gowers given at Gresham College.

My explorations of this are must wait, though; I feel much more at home in the world of infinite sums, the Riemann Zeta function, the Gamma function and the Chinese Remainder Theorem so have decided only to submit a proposal for the Bernoulli numbers project. Here follows my first draft:


PROJECT PROPOSAL

Objectives
* Add support for Bernoulli numbers along with thread-safe caching.
* Add support for polygamma of positive integer order (requires Bernoulli numbers).
* Improve tgamma/lgamma for multiprecision types based on Stirling's approx.
* Optional: Add support for the Hurwitz zeta function.


Preliminaries
* Main reference: (Richard P. Brent, David Harvey) Fast computation of Bernoulli, Tangent and Secant numbers, plus formulae gleaned from Wikipedia pages (to be cross-checked for accuracy)

* Bernoulli numbers can be found most useful in finite- or infinite sums, rather than on their own, e.g.

  - polygamma:
  \psi^{(m)}(z) = (-1)^{m+1}\sum_{k=0}^{\infty}\frac{(k+m-1)!}{k!}\frac{B_k}{z^{k+m}} \qquad m \ge 1
  \psi^{(0)}(z) = \ln(z) - \sum_{k=1}^\infty \frac{B_k}{k z^k}

  - Stirling's approximation to (the natural logarithm of) the gamma function:
  \ln (\Gamma (z)) \sim \left(z-\tfrac{1}{2}\right)\ln(z) -z + \tfrac{1}{2}\ln(2 \pi) + \sum_{n=1}^\infty \frac{B_{2n}} {2n(2n-1) z^{2n-1}}

  - Hurwitz zeta function (for negative n):
  \zeta(-n,x) = - \frac{B_{n+1}(x)}{n+1}
  where B_{n+1}(x) is the Bernoulli polynomial of degree n+1: B_n(x) = \sum_{k=0}^n {n \choose k} B_{n-k} x^k
  (for positive n it is related to the polygamma function, above)

therefore the principal aim should be to return a suitable data structure containing the first n Bernoulli numbers.

* The above paper shows that the first n Bernoulli numbers can be computed with O( n^2 (log n)^{2 + o(1)} ) bit-operations and provides an algorithm, FastTangentNumbers, by which the first n Tangent numbers can be computed with that complexity, whence the first n Bernoulli numbers are obtained via the relation:

T_n = (-1)^{n-1} 2^{2n} \frac{B_{2n}}{2n}

* The paper also provides another algorithm, TangentNumbers, similar to the Akiyama-Tanigawa algorithm for Bernoulli numbers but only involving integer arithmetic, and the suggestion that it may be more efficient than FastTangentNumbers for n < 1000.


Intentions
* Implement the algorithms TangentNumbers and FastTangentNumbers from the reference, using the Boost.Multiprecision library.
* (Trivially) extend these to calculate the first n Bernoulli numbers (and the first n Bernoulli polynomials?).
* Using m Bernoulli numbers (each of precision k - for some k and m to be determined) implement the gamma function, polygamma function(s) and Hurwitz zeta function, all to a given precision n, using the relations given above.


Miscelleous Thoughts
* It would make sense for the Bernoulli numbers to be returned as rational numbers, using the Boost.Rational library, so that precision is not needlessly lost by dividing.
* By extension, should the infinite sums needed in the implementation of the other functions be returned as polynomials which can be evaluated similarly to the Boost.Polynomial library? (So that we return a single object representing the function, which the user can then evaluate at multiple points.)


That's as far as I've got at present, but I'll be filling in the remaining gaps before Friday, of course.

Regards,

David
Reply | Threaded
Open this post in threaded view
|

Re: [GSoC 2013] Draft Proposal for Bernoulli Numbers project

Christopher Kormanyos


> I feel much more at home in the world of infinite sums,
> the Riemann Zeta function, the Gamma function

> and the Chinese Remainder Theorem so have decided
> only to submit a proposal for the Bernoulli numbers
> project. Here follows my first draft:

Good work on the proposal so far David.

Try to add a brief section about potential design
choices and how the proposed solution will fit
into Boost.Math.

Please note that the proposal must be submitted
to Googles melange site. Some candidates have
been sending the proposals directly to us, and
this is incomplete.

You can send the proposal to us for comments,
but it's really close to the deadline. So please
try to get your proposal up on Google melange.

Some students put the equations in another file
linked in the proposal because melange is
predominantly text-based.

Keep up the good work.

Sincerely, Chris.

_______________________________________________
Unsubscribe & other changes: http://lists.boost.org/mailman/listinfo.cgi/boost
Reply | Threaded
Open this post in threaded view
|

Re: [GSoC 2013] Draft Proposal for Bernoulli Numbers project

David Joy