On changing Bernoulli(1) to +½

[Jeremy Tan]

My name is Jeremy Tan, or Parcly Taxel in the furry/MLP art scene. As of this post I am a recent graduate from the National University of Singapore with two degrees in maths and computer science.

Over the past month I had a good read of Peter Luschny's Bernoulli Manifesto (http://luschny.de/math/zeta/The-Bernoulli-Manifesto.html) and was thoroughly convinced that B_1 (the first Bernoulli number) has to be +½, not -½. (Much of Luschny's argument centres on being able to (1) interpolate the Bernoulli numbers when B_1 = +½ with an entire function intimately related to the zeta function, and (2) extend the range of validity of or simplify several important equations like the Euler–Maclaurin formula. Have a read yourself though – it is close to divine truth.)

So I went to SymPy – one of SageMath's dependencies, and where a discussion on this topic was open (https://github.com/sympy/sympy/issues/23866) – and successfully merged several PRs there (https://github.com/sympy/sympy/pull/23926) implementing both that change and some functions in Luschny's "An introduction to the Bernoulli function" (https://arxiv.org/abs/2009.06743).

I thought I was also done with changing B_1 = +½ for SageMath, but then someone pointed out that the latter currently uses other libraries that all have B_1 = -½. I have already opened a PR for one such library, FLINT, to change B_1 = +½ there (https://github.com/wbhart/flint2/pull/1179). However Fredrik Johansson has advised me that I take the discussion right here, to sage-devel, because (in his words)

> if FLINT and Arb change their definitions but the Sage developers decide that they don't like it, they will just treat the new behavior as a bug and add a special case in the wrapper to return B_1 = -½.

So my proposal is to special-case it the other way – before the backend selection in Sage's Bernoulli code (https://github.com/sagemath/sage/blob/08202bc1ba7caea46327908db8e3715d1adf6f9a/src/sage/arith/misc.py#L349), add a check for argument 1 and immediately return +½ if that is the case. This also has the advantage of bypassing libraries that haven't or don't want to change.

What do you think?

Jeremy Tan / Parcly Taxel

[William Stein]

It could be done via the "1 year deprecation policy". I.e., return the
current value by default with a warning message
(and note about an option to change it) for the next year, then when
there is a release in late 2023 (?), the default would change. This
would give people time to update their code.

I have no comment on the pros and cons of this personally, though I'm
curious if the change breaks any code anywhere else in Sage (e.g.,
maybe for computing q-expansions of modular forms?)...

>
> Jeremy Tan / Parcly Taxel
>

> --
> You received this message because you are subscribed to the Google Groups "sage-devel" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.
> To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/CAGYgO94gF%3DBKo7gRnUj8c3H0bJyuLp_Apr%3D8Y9NC%2BFM%2BSZHNOg%40mail.gmail.com.



--
William (http://wstein.org)

[davida...@gmail.com]

>  I'm curious if the change breaks any code anywhere else in Sage (e.g., maybe for computing q-expansions of modular forms?)...

You guessed right. I did a quick local change to the bernoulli function and it indeed breaks some tests in sage/modular/modform:

    ~/sage$ ./sage -t src/sage/modular/modform

    sage -t --random-seed=279226112023210448433794639443228726052 src/sage/modular/modform/ambient.py  # 1 doctest failed
    sage -t --random-seed=279226112023210448433794639443228726052 src/sage/modular/modform/element.py  # 11 doctests failed
    sage -t --random-seed=279226112023210448433794639443228726052 src/sage/modular/modform/ambient_g1.py  # 1 doctest failed
    sage -t --random-seed=279226112023210448433794639443228726052 src/sage/modular/modform/eisenstein_submodule.py  # 4 doctests failed
    sage -t --random-seed=279226112023210448433794639443228726052 src/sage/modular/modform/ring.py  # 3 doctests failed
    sage -t --random-seed=279226112023210448433794639443228726052 src/sage/modular/modform/constructor.py  # 1 doctest failed

However, I would be in favor for this change. I would also be glad to lend a hand for fixing those doctests.

[William Stein]

On Sat, Sep 10, 2022 at 10:04 AM davida...@gmail.com
<davida...@gmail.com> wrote:
>
> > I'm curious if the change breaks any code anywhere else in Sage (e.g., maybe for computing q-expansions of modular forms?)...
>
> You guessed right. I did a quick local change to the bernoulli function and it indeed breaks some tests in sage/modular/modform:

I created all of the files listed below. My guess is that code for
computing q-expansions of Eisenstein series assume B(1) is what it is,
and one would just need to change that code by changing a sign
somewhere.

William

> To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/da3c85ec-1442-4c60-a714-fba2f14f908bn%40googlegroups.com.



--
William (http://wstein.org)

[Oscar Benjamin]

On Sat, 10 Sept 2022 at 18:49, William Stein <wst...@gmail.com> wrote:
>
> On Sat, Sep 10, 2022 at 10:04 AM davida...@gmail.com
> <davida...@gmail.com> wrote:
> >
> > > I'm curious if the change breaks any code anywhere else in Sage (e.g., maybe for computing q-expansions of modular forms?)...
> >
> > You guessed right. I did a quick local change to the bernoulli function and it indeed breaks some tests in sage/modular/modform:
>
> I created all of the files listed below. My guess is that code for
> computing q-expansions of Eisenstein series assume B(1) is what it is,
> and one would just need to change that code by changing a sign
> somewhere.

The two different versions of the Bernoulli sequence are related very
simply (B-(n) = (-1)^n B+(n)) so it should be easy to update any code
that depends on the current behaviour. This is why SymPy isn't
providing any direct replacement for its previous B- behaviour.


Oscar

[Jeremy Tan]

Considering the subsequent replies to my proposal, I'm perfectly happy with implementing B_1 = +½ in Sage with the 1-year deprecation policy. During the deprecation period I would expect an extra boolean keyword argument to bernoulli() enforcing B_1 = +½ if true and relying on the backend libraries if false; at the end of said period the extra argument will be removed and B_1 = +½ will be permanently enforced.

If you want me to, I can fork Sage and make a trac ticket for this, but then I'd need your help in identifying other files/tests to change for the new behaviour.

[davida...@gmail.com]

When you open a trac ticket, post the link here, we will continue the discussion there.

[Jeremy Tan]

https://trac.sagemath.org/ticket/34521

Jeremy Tan / Parcly Taxel

[Jeremy Tan]

Perhaps I need your help, William. The associated trac ticket https://trac.sagemath.org/ticket/34521 was marked invalid offhand by the release manager Frédéric Chapoton despite me managing to fix all the doctests simply and the ticket itself not conflicting with any Sage policies.

He claims that "Mr Luschny has not published any mathematical article" – but his manifesto is indistinguishable from an article in letter from (it is a reply to Donald Knuth) and his introduction to the Bernoulli function is even more of an article. I believe I have responded in kind to his claim that B_1 = -½ is "standard".

There must be a way to get Chapoton out of the way here so we can effect the deprecation.

[Dima Pasechnik]

On Mon, Sep 12, 2022 at 2:21 PM Jeremy Tan <redde...@gmail.com> wrote:
>
> Perhaps I need your help, William. The associated trac ticket https://trac.sagemath.org/ticket/34521 was marked invalid offhand by the release manager Frédéric Chapoton despite me managing to fix all the doctests simply and the ticket itself not conflicting with any Sage policies.
>
> He claims that "Mr Luschny has not published any mathematical article" – but his manifesto is indistinguishable from an article in letter from (it is a reply to Donald Knuth) and his introduction to the Bernoulli function is even more of an article. I believe I have responded in kind to his claim that B_1 = -½ is "standard".
>
> There must be a way to get Chapoton out of the way here so we can effect the deprecation.

Frederic has not been helpful here, unfortunately - sorry, Jeremy.
I've made a comment to this effect on the ticket. People communicated
here on the topic, it was deemed to have a merit, and closing tickets
like this is just putting people off.

Dima

> To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/3fe77fe4-266e-4650-ae15-cc0f436b1b28n%40googlegroups.com.

[David Joyner]

On Mon, Sep 12, 2022 at 9:49 AM Dima Pasechnik <dim...@gmail.com> wrote:
>
> On Mon, Sep 12, 2022 at 2:21 PM Jeremy Tan <redde...@gmail.com> wrote:
> >
> > Perhaps I need your help, William. The associated trac ticket https://trac.sagemath.org/ticket/34521 was marked invalid offhand by the release manager Frédéric Chapoton despite me managing to fix all the doctests simply and the ticket itself not conflicting with any Sage policies.
> >
> > He claims that "Mr Luschny has not published any mathematical article" – but his manifesto is indistinguishable from an article in letter from (it is a reply to Donald Knuth) and his introduction to the Bernoulli function is even more of an article. I believe I have responded in kind to his claim that B_1 = -½ is "standard".
> >
> > There must be a way to get Chapoton out of the way here so we can effect the deprecation.
>

I'm probably too late on this, but I'd like to suggest a sort of compromise:

According to wikipedia, there is no standard. There's a B+ (favored by
NIST and "modern authors") and a B- (favored by Donald Knuth some others).
Why don't we create a B+ and a B-? This will make the OP happy.
To make William Stein happy (or happiest), we (temporarily?) set B=B+, so
he doesn't have to change any code.

> To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/CAAWYfq2CbfNa8d_sw9BaWPH_-vFEERq7Mqy5vSqoFnpryBv_-A%40mail.gmail.com.

[davida...@gmail.com]

>  Why don't we create a B+ and a B-?

This was one of the idea of the ticket https://trac.sagemath.org/ticket/34521.  A new option to the bernoulli function was added ("plus=False"), giving the option to the user to choose their prefered value. However, I think that the ticket was set to invalid and closed prematurely without any possibility for a meaningful discussion. I think that the debate here would be: what is the default value of this new option. There is two choices:

1. "plus=False": bernoulli(1) = -1/2, nothing changes.
   
2. "plus=True": bernoulli(1) = +1/2 (new value). Deprecation needs to be implemented.

The choice two was chosen in ticket #34521. If two is chosen, I believe that, as mentioned by William, there should be a 1-year deprecation period before the change actually takes place. However, if the community prefer (1), I believe this can be done also (discussion is useful). In any case, I think that something should be done about this, because FLINT has been changed to return +1/2 (and FLINT is the default for small Bernoulli numbers)

David A.

[David Joyner]

On Mon, Sep 12, 2022 at 11:24 AM davida...@gmail.com
<davida...@gmail.com> wrote:
>
> > Why don't we create a B+ and a B-?
>
> This was one of the idea of the ticket https://trac.sagemath.org/ticket/34521. A new option to the bernoulli function was added ("plus=False"), giving the option to the user to choose their prefered value. However, I think that the ticket was set to invalid and closed prematurely without any possibility for a meaningful discussion. I think that the debate here would be: what is the default value of this new option. There is two choices:
>
> "plus=False": bernoulli(1) = -1/2, nothing changes.
> "plus=True": bernoulli(1) = +1/2 (new value). Deprecation needs to be implemented.
>
> The choice two was chosen in ticket #34521. If two is chosen, I believe that, as mentioned by William, there should be a 1-year deprecation period before the change actually takes place. However, if the community prefer (1), I believe this can be done also (discussion is useful). In any case, I think that something should be done about this, because FLINT has been changed to return +1/2 (and FLINT is the default for small Bernoulli numbers)
>

My vote is for (1) for 2 reasons. (a) It makes it clear the SageMath
community acknowledges there is no clear standard. (b) The ticket is
relatively easy to implement since it won't break anything.

Once implemented, debate on the default behavior can be moved to a new
ticket, if desired.

> To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/3b526211-160d-4135-9c53-0c2bc33c7990n%40googlegroups.com.

[Oscar Benjamin]

On Mon, 12 Sept 2022 at 16:24, davida...@gmail.com

<davida...@gmail.com> wrote:
>
> > Why don't we create a B+ and a B-?
>
> This was one of the idea of the ticket https://trac.sagemath.org/ticket/34521. A new option to the bernoulli function was added ("plus=False"), giving the option to the user to choose their prefered value. However, I think that the ticket was set to invalid and closed prematurely without any possibility for a meaningful discussion. I think that the debate here would be: what is the default value of this new option. There is two choices:
>
> "plus=False": bernoulli(1) = -1/2, nothing changes.
> "plus=True": bernoulli(1) = +1/2 (new value). Deprecation needs to be implemented.

If there is a need to provide both behaviours then I would suggest
using two new functions bernoulli_plus and bernoulli_minus rather than
a boolean flag. The bernoulli function can continue to return
bernoulli_minus but can refer to the other functions in its
documentation. That way the conventions are clearly distinguished and
in future it can be possible to make it so that bernoulli emits a
warning (prompting the caller to use bernoulli_plus or bernoulli_minus
explicitly instead) or perhaps is eventually changed to call
bernoulli_plus after some time of giving a warning.

The reason I didn't want to go that way for SymPy is that
bernoulli_plus admits a natural generalisation to real and complex
numbers but bernoulli_minus does not. It would be awkward to implement
both bernoulli_plus and bernoulli_minus as symbolic functions where
the rules and domains, numerical evaluation etc would all be
different: one would be a complex analytic function while the other an
integer sequence.

> In any case, I think that something should be done about this, because FLINT has been changed to return +1/2 (and FLINT is the default for small Bernoulli numbers)

Has FLINT been changed? I see a PR here but it is still open:
https://github.com/wbhart/flint2/pull/1179

--
Oscar

[davida...@gmail.com]

>  Has FLINT been changed? I see a PR here but it is still open

You're totally right. I thought that it was already closed there! Thanks for pointing this out. 

Note that I don't have a very strong opinion on this matter, I'm only worrying about preventing potential bugs.

David A.

[Fredrik Johansson]

I'm pretty neutral about this change, but I've received PRs for FLINT and mpmath (presumably there will be one for Arb as well) so I suppose I will need to make a decision about merging or closing them sooner or later :-)

The sign convention for B_1 is fairly arbitrary, and the downside of changing it is that this introduces ambiguity and inconsistency where there was none before. On the other hand, you can make the case that "patching" deficient conventions is better in the long run... at least if the new convention eventually sees near-universal adoption (which is not at all guaranteed here).

An advantage of B+ is that it would agree with the definition used in Sage for generalized Bernoulli numbers when restricted to the trivial character; see https://doc.sagemath.org/html/en/reference/modmisc/sage/modular/dirichlet.html#sage.modular.dirichlet.DirichletCharacter.bernoulli

The claim "bernoulli_plus admits a natural generalisation to real and complex numbers but bernoulli_minus does not" (made elsewhere in this thread) seems a bit hyperbolic. For B+ this natural generalization is -n*zeta(1-n); for B- one can just use -n*zeta(1-n)*cos(pi*n). OK, one is a bit simpler than the other, but both are perfectly fine entire functions.

For Sage and mpmath, I suppose the least intrusive option is to introduce a keyword argument to select the plus/minus convention, with B- as default until there is community consensus to change it. There is not really a strong need to change low-level libraries like FLINT at this point since it's trivial to handle both conventions in a wrapper.

Fredrik

[edgardi...@gmail.com]

The choice of the sign is arbitrary. So why make this change? What is the benefit?
Most recent papers, algebra systems (Maple/Mathematica/Magma/Matlab/Oscar), and libraries (Pari/Flint/Mpmath/ARB) seemed to have picked B_1 = -1/2.
Thus why put work into changing the default value and go against the current norm?

By doing this you will generate a lot of unnecessary work to go from one arbitrary choice (that most people use) to another that few seem to use.

[William Stein]

On Mon, Sep 12, 2022 at 7:00 PM edgardi...@gmail.com
<edgardi...@gmail.com> wrote:
>
> The choice of the sign is arbitrary. So why make this change? What is the benefit?

The answer to that question is what the website

http://luschny.de/math/zeta/The-Bernoulli-Manifesto.html

is all about. I think you have to read that website. Here's a random
quote from that website that might peak your interest:

"By now, hundreds of books that use "the minus-one-half" convention
have unfortunately been written. Even the major software systems for
symbolic mathematics have that 20th-century aberration deeply
embedded. Yet Luschny convinced me that we have all been wrong, and
that it's high time to change back to the correct definition before
the situation gets even worse. Therefore I changed the definition of
$B_1$ in all printings of The Art of Computer Programming..."

- Donald Knuth




--
William (http://wstein.org)

[Jeremy Tan]

On Tuesday, 13 September 2022 at 10:00:20 UTC+8 edgardi...@gmail.com wrote:

The choice of the sign is arbitrary. So why make this change? What is the benefit?
Most recent papers, algebra systems (Maple/Mathematica/Magma/Matlab/Oscar), and libraries (Pari/Flint/Mpmath/ARB) seemed to have picked B_1 = -1/2.
Thus why put work into changing the default value and go against the current norm?

By doing this you will generate a lot of unnecessary work to go from one arbitrary choice (that most people use) to another that few seem to use.

The people and CAS's who use B_1 = -½ are good, but wrong. The central thesis of Luschny's manifesto is that B_1 = +½ shows the central place of the Bernoulli numbers in all mathematics much better than B_1 = -½. As he puts it in the closing paragraphs:

"I believe in the unity of mathematics as described by Barry Mazur… In any case mathematics is not a bunch of disparate conventions which must be maintained in its present form just to protect some mathematicians from a possible confusion."

Every general-purpose CAS is an embodiment of said unity of mathematics, and should strive to capture – make real – that unity as much as possible. Yes, the Four Heavenly Kings (Magma, Maple, Mathematica, Matlab) may provide B_1 = -½, but this is SageMath, an open-source alternative to all of them. We don't have to follow them or any convention.

As well as that, SageMath originally stood for System for Algebra and Geometry Experimentation. Given the almost excessive number of relations between the Bernoulli numbers with B_1 = +½ and other mathematical entities, I am confident that a secondary school/high school student/undergraduate first introduced to the Bernoulli numbers with B_1 = +½, and given SageMath with B_1 = +½, would have a lot more fun than if B_1 = -½ was initially given instead.

Finally I would like to drop a quip from the manifesto: "My impression is that the number of modern authors who use B_1 = +½ is steadily increasing."

[Vincent Delecroix]

PARI/GP actually has a better convention : only even Bernoulli numbers exist

? bernvec(5)
%1 = [1, 1/6, -1/30, 1/42, -1/30, 5/66]

And the two conventions can be recovered as evaluations of Bernoulli
polynomials at 0 and 1 respectively

? [subst(bernpol(n), x, 0) | n <- [1..6]]
%2 = [-1/2, 1/6, 0, -1/30, 0, 1/42]
? [subst(bernpol(n), x, 1) | n <- [1..6]]
%3 = [1/2, 1/6, 0, -1/30, 0, 1/42]

I second Edgar that I don't see much of the point on doing this
change. If we had to do any change I would suggest to raise an error
when bernoulli(n) is called with n=1 so that nobody could complain
about sage convention choices.

Vincent

On Tue, 13 Sept 2022 at 04:00, edgardi...@gmail.com

> --
> You received this message because you are subscribed to the Google Groups "sage-devel" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.

> To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/9c61b1ee-d2c4-4c6c-b539-f89ad6d45941n%40googlegroups.com.

[Jeremy Tan]

A simpleton's way of getting out of the problem indeed. PARI/GP's documentation says:

? ?bernvec
bernvec(n): returns a vector containing, as rational numbers, the Bernoulli numbers B_0, B_2, ..., B_{2n}.
? bernfrac(3)
%2 = 0
? bernfrac(5)
%3 = 0

So not only B_1 but also B_3, B_5, etc. have values.

Jeremy Tan / Parcly Taxel

[David Joyner]

On Tue, Sep 13, 2022 at 3:43 AM Jeremy Tan <redde...@gmail.com> wrote:
>
> A simpleton's way of getting out of the problem indeed. PARI/GP's documentation says:
>

Let's play nice here, okay?

> --
> You received this message because you are subscribed to the Google Groups "sage-devel" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.

> To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/57e2bd6c-98eb-4617-928b-feea43fce978n%40googlegroups.com.

[Oscar Benjamin]

On Mon, 12 Sept 2022 at 22:09, Fredrik Johansson

<fredrik....@gmail.com> wrote:
>
> The claim "bernoulli_plus admits a natural generalisation to real and complex numbers but bernoulli_minus does not" (made elsewhere in this thread) seems a bit hyperbolic. For B+ this natural generalization is -n*zeta(1-n); for B- one can just use -n*zeta(1-n)*cos(pi*n). OK, one is a bit simpler than the other, but both are perfectly fine entire functions.

I chose the word "natural" carefully because to me it conveys the fact
that it is a subjective statement rather than a formal one. I assume
it is formally possible to make infinitely many complex analytic
functions to interpolate any function on the nonnegative integers (is
that a named theorem?).

The question is whether the factor of cos(pi*z) that you've included
is something that is meaningful for noninteger values of z or
something included only artifically. Did you have any reason to choose
it besides fixing bernoulli(1)=-1/2 (which can also be done in
infinitely many other ways)?

Most mathematical functions will never gain the status of being a
named function in any context so would defining a complex analytic
bernoulli_minus(z) = -z*zeta(1-z)*cos(pi*z) have any particular value?
Any formula involving bernoulli_plus would apply equally to
bernoulli_minus but with some factors of cos(pi*z) lying around. Would
any of the formulas have some other factor of cos(pi*z) for those to
cancel with in order to produce something more "natural"? Would there
be any value in SymPy, mpmath, Arb, Sage etc implementing such a
bernoulli_minus?

--
Oscar

[Dima Pasechnik]

On Tue, Sep 13, 2022 at 10:35 AM Oscar Benjamin
<oscar.j....@gmail.com> wrote:
>
> On Mon, 12 Sept 2022 at 22:09, Fredrik Johansson
> <fredrik....@gmail.com> wrote:
> >
> > The claim "bernoulli_plus admits a natural generalisation to real and complex numbers but bernoulli_minus does not" (made elsewhere in this thread) seems a bit hyperbolic. For B+ this natural generalization is -n*zeta(1-n); for B- one can just use -n*zeta(1-n)*cos(pi*n). OK, one is a bit simpler than the other, but both are perfectly fine entire functions.
>
> I chose the word "natural" carefully because to me it conveys the fact
> that it is a subjective statement rather than a formal one. I assume
> it is formally possible to make infinitely many complex analytic
> functions to interpolate any function on the nonnegative integers (is
> that a named theorem?).

you can always add sin(x)/pi - which vanishes on 0,1,2,..., - to your
interpolant.
So that's an excercise :-)

>
> The question is whether the factor of cos(pi*z) that you've included
> is something that is meaningful for noninteger values of z or
> something included only artifically. Did you have any reason to choose
> it besides fixing bernoulli(1)=-1/2 (which can also be done in
> infinitely many other ways)?
>
> Most mathematical functions will never gain the status of being a
> named function in any context so would defining a complex analytic
> bernoulli_minus(z) = -z*zeta(1-z)*cos(pi*z) have any particular value?
> Any formula involving bernoulli_plus would apply equally to
> bernoulli_minus but with some factors of cos(pi*z) lying around. Would
> any of the formulas have some other factor of cos(pi*z) for those to
> cancel with in order to produce something more "natural"? Would there
> be any value in SymPy, mpmath, Arb, Sage etc implementing such a
> bernoulli_minus?
>
> --
> Oscar
>

> --
> You received this message because you are subscribed to the Google Groups "sage-devel" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.

> To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/CAHVvXxSXzmOXcNnO3scCZOm7pZG36tGB9G-GEavqQi2hKRq-Xw%40mail.gmail.com.

[Jeremy Tan]

On Tue, 13 Sept 2022, 17:00 David Joyner, <wdjo...@gmail.com> wrote:

Let's play nice here, okay?

Let me explain what I mean in a nicer way. Not defining B_1 looks good on the surface given the current discussion, but is really a strictly worse option than defining B_1 = +½ or -½ because then the n = 1 case has to be special-cased out of formulae where B_1 is normally used.

That B_n for odd n at least 3 is zero may provide computational advantages, but is a non-trivial fact. It is more parsimonious to put the even, one and greater-odd instances of the Bernoulli numbers into one bag.

[John H Palmieri]

I have no opinions about what B_1 should be, but I am concerned about potential confusion: some users will expect one value for B_1, others will expect a different value, and so one group or other will end up being confused when answers don't come out the way they expect. The safest course might be to drop "bernoulli" altogether and have two functions for the two different conventions. Or keep "bernoulli" but have it just print a message saying that there are two conventions, so instead use either "bernoulli_plus" or "bernoulli_minus" (or whatever they might be called). I don't know if the safest course is the right course, though.

(This is how I feel about the "range" of a function in mathematics: some people use this word to mean the image, some people mean it to use the codomain. Since there is serious possibility of confusion, I personally just never use the term at all.)

--

John

[Peter Luschny]

tl;dr, 80 lines, sorry.

I think there is a better alternative than changing the current
implementation of the Bernoulli numbers.

Fredrik: "The sign convention for B_1 is fairly arbitrary, ..."

Calling the question a 'convention' sets a wrong framing from
the start. Conventions are treated differently than bugs.
And this is a bug. Setting B- leads to inconsistency in several
identities, which I have described on my known web-page.

Fredrik, honestly, do you really think Knuth rewrites TAOCP
and CM just to change an arbitrary sign convention? Knuth
expressly considers it a mistake.

As long as the discussion is viewed as one about conventions,
it is uninteresting and does not improve our understanding of
the situation. People who believe this can equally well flip
a coin when in doubt.

Best leave the whole discussion behind, the question of + or -,
the question of convention or bug, the question of decision
or compromise, or what people on Twitter think (Fredrik asked :))
All answers to such questions must be profoundly unsatisfactory.

The matter should be solved mathematically in a way Euler
would have liked, without resorting to this absurd fixation
on a single value.

And by that, I mean to define a suitable complex function,
call it the 'Bernoulli function', and then say: The 'Bernoulli 

numbers' are the values of this function on the non-negative integers.

Everyone knows a possible way to do this: B(z) = -z*zeta(1-z).
Of course, you could add a sin or cos factor to it (Fredrik
mentioned that, Knuth too, by the way), but why should you do
that? Occam's razor speaks against it. Maximum simplicity has
always belonged to the beauty criteria in mathematics.

Similar approaches have already been done; Helmut Hasse, for
example, has constructed the corresponding function without
recourse to the zeta function. And there are also other
representations without reference to the zeta function.

Another pleasing possibility is via the Jensen-Johansson-Blagouchine
formulas, first given by Jensen and first proved by Johansson
and Blagouchine ("Computing Stieltjes constants using complex
integration").

These formulas are also the starting point in my arXiv paper.
There I try to show that this function is also essential in
other contexts in the theory of special numbers and offers
some technical advantages.

* My suggestion: Offer this function in Sage. This may not mean
much more than a three-line wrapper around the zeta function.
Then the current code does not have to be changed, no deprecation
warnings, no keywords for alternatives are required.

The decision is then solely up to the user whether he wants
to continue the current usage and use bernoulli(n) or find
pleasure in bernoulli_function(n).

Fredrik's fear that new "ambiguity and inconsistency" could
creep in is then unfounded: The names and the definitions are
too different; here polynomials, there a zeta-like function,
both with different applications.

If you want to describe this approach in a highfalutin way,
it is Grothendieck-ish. It does not answer the question which
of the two values is the 'correct' one; it shows how this problem
disappears when put into the proper conceptual framework.

This approach also makes sense in all the other cases that
Fredrik has to decide. With it, the burden of making a sensible
decision is turned into the freedom for the user to explore
a fascinating function.

 --
 Peter
 

[Jeremy Tan]

Thank you for coming, Luschny.

I not only wholeheartedly believe B_1 = +½ and that there is no convention about it, but also that I believe in the reality and usefulness of the Bernoulli, Euler and other functions you defined in your 2020 paper.

When I implemented those functions in SymPy there was already a base – two-argument bernoulli() and euler() with polynomial support were there, as well as genocchi(). So I could easily generalise those functions into yours.

On the other hand, until writing this reply I did not know that SageMath exposed a hurwitz_zeta() function by default from which I can define the generalised Bernoulli function. Hold on; I'll have a new ticket ready.

Jeremy Tan / Parcly Taxel

--
You received this message because you are subscribed to a topic in the Google Groups "sage-devel" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/sage-devel/G4sbKoibXK4/unsubscribe.
To unsubscribe from this group and all its topics, send an email to sage-devel+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/19d2d0c8-556b-4bbf-a9d1-e87e7e8ed4c9n%40googlegroups.com.

[Jeremy Tan]

https://trac.sagemath.org/ticket/34536

There, Luschny!

On why I initially only changed B_1 = +½, I have a rebuttable demand of my code: it must be as beautiful as the maths it implements. Having bernoulli_gen() (your generalised Bernoulli function as implemented in the SageMath ticket) alongside an integer-only bernoulli() is not pretty to me, and especially not pretty when they return inconsistent values at 1.

Hence my first idea was to overload bernoulli() like I did with SymPy, but seeing as that was not very feasible, my second idea was to simply change B_1 = +½. Now that you prodded me into finding hurwitz_zeta() at top-level, I have the new ticket only implementing bernoulli_gen().

Jeremy Tan / Parcly Taxel

On Thu, 15 Sept 2022, 02:01 Peter Luschny, <peter....@gmail.com> wrote:

[Jeremy Tan]

There is now an Arb equivalent of the Sage ticket that implements the generalised Bernoulli function but nothing else:

https://github.com/fredrik-johansson/arb/pull/438

Any further equivalents (e.g. on mpmath) will be edited into the description of the Sage ticket rather than being mailed here.

Parcly Taxel

[William Stein]

On Mon, Sep 12, 2022 at 7:58 PM Jeremy Tan <redde...@gmail.com> wrote:
[...]

> As well as that, SageMath originally stood for System for Algebra and Geometry Experimentation.

I just happened to stumble again on the original scrap of paper just
now where I made up the name
Sage back in 2004. The name was "System for **Arithmetic** Geometry
Experimentation" for the first year or two.
It was supposed to be better than "MANIN", which was the original
name. I also considered "pymath",
which looks like it never got used for anything seriously [1].

[1] https://pypi.org/project/pymath/



--
William (http://wstein.org)

[Kwankyu Lee]

On Friday, September 16, 2022 at 10:46:09 AM UTC+9 wst...@gmail.com wrote:

I just happened to stumble again on the original scrap of paper just
now where I made up the name
Sage back in 2004. 

 

A good photo of it deserves to be permanently placed in the history section of the sagemath website.  But I could not find the history section.

 

[William Stein]

Wow a history section would be a great idea. It could also highlight how much work people put into the sage notebook back in the day, and the extent to which sage builds on pari, singular, gap, etc…

 

--
You received this message because you are subscribed to the Google Groups "sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/155dc41c-6d30-4792-923f-2ff85d831de7n%40googlegroups.com.

--

-- William Stein