I'm sure this isn't what you are looking for but I had some time today, so I
calculated Nivens constant using the zeta function and bernoullis. My guess is
that it is computationally less expensive than using q:
Step 1: Generate Bernoulli Numbers. We get a list of (even) Bernoulli numbers
up to the given number.
NB. Generates all even bernoulli numbers up to y
generatebernoulli=: 3 : 0
arr=. (y + 2) $ 0
arr=. (1) 1} arr
bernoulliList=. 1
flag=. 1 NB. true
a=. 1
b=. 1
c=. _2
d=. 1
size=. (2 * y) - 2
for_j. i. size do.
if. flag = 1 do.
a=. a + 1
b=. 4 * b
c=. _1 * c
d=. c * (b - 1)
k=. a
while. k > 0 do.
tmp0=. k { arr
tmp1=. (k+1) { arr
arr=. (tmp0 + tmp1) k} arr
k=. k - 1
end.
else.
k=. 1
while. k < a do.
tmp0=. k { arr
tmp1=. (k-1) { arr
arr=. (tmp0 + tmp1) k} arr
k=. >: k
end.
bernoulliList=. bernoulliList, (x: ((a-1){arr) % d)
end.
flag=. 2 | >: flag
end.
bernoulliList
)
STEP 2: Use the relation between bernoulli numbers and zeta to write a verb
that computes the zeta (2 *y), for argument y. left argument, x, should be the
previously
calculated bernoulli number list.
NB. Calculates Even Zeta, i.e. Zeta (2 *y).
NB. Uses the relation between even zeta values and bernoulli numbers.
calcevenzeta =: 4 : 0
bernoullis =. x
m =. 2 * y
(_1^(y+1)) * ((y{bernoullis) * (2 * o. 1)^m) % (2 * ! m)
)
STEP 3: For odd numbers we have to do raw computation. So compute zeta(y) for
200 terms.
NB. Calculates zeta (y), but should be used for odd y
calcoddzeta =: 3 : 0
arr =.%>: i. 200 NB. 200 is enough
+/arr ^ y
)
STEP 4: calculate Nivens constant
bernoullilist =: generatebernoulli 100 NB. we don't even need this many
zetaapproximations =:(bl&calcevenzeta@:-:)`calcoddzeta@.(1&=(2&|)) NB. two
cases, even or odd integers.
NIVENS_CONSTANT =: 1+ +/ 1- % zetaapproximations"(0) 2 + i. 100 NB. we will
approximate with the first 100 integers, but I doubt this many is necessary.
this gives 1.7052
I wrapped all the above in a single verb and ran timespacex
timespacex 'allcalc 100' NB. calculate first 100 terms only.
0.044685 34240
Doing raw calcs with q: is way more expensive:
timespacex 'mean >./@(0&,)@(_&q:)"0 }.i.1e5'
1.5124 9.04083e6
--------------------------------------------
On Wed, 7/6/16, jonghough via Programming <[email protected]> wrote:
Subject: Re: [Jprogramming] Niven's constant
To: [email protected]
Date: Wednesday, July 6, 2016, 12:40 AM
Hi,
I don't have my PC at hand now, but what is
the problem with your calculation?
I had never heard of Nivens constant, but
looking at Wikipedia there is a formula for it containing
the zeta function.
Why not
use this formula?
Zeta 2 3
4 are well known constants, and you can calculate or
hardcode the values for other even arguments,
https://en.m.wikipedia.org/wiki/Particular_values_of_Riemann_zeta_function
You might need to
calculate bernoulli numbers.
For odd values I think you
need to calculate from scratch, but since zeta 14 is almost
1 you don't need to do too many calculations, just from
5 to 15, say.
This is just thinking out loud. But my main
point is you only probably need the first 15 or so zeta
values before convergence is 'good enough'.
Regards,
Jon
On Wed, Jul
6, 2016 at 12:12 AM +0900, "mikel paternain"
<[email protected]>
wrote:
________________________________
De: mikel paternain
Enviado:
domingo, 3 de julio de 2016 12:07:37
Para:
[email protected]
Asunto: Niven's constant
Hi,
I need improve the code for calculating
Niven's constant. Could you help me?
A first test is (see "Computo ergo
sum" section in JoJ):
exp=:+/"1@=@q:
max=:>./
a=.max@exp(2+1.1000)
C=.(+/%#)a
Thanks in advance,
Mikel Paternain
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
