Essays/Mertens Function

Basics

The Mertens function on positive integer n is the sum of the Möbius function on the integers from 1 to n .

mobius  =: */ @: - @: ~: @: q:

mertens0=: +/ @: (mobius"0) @: >: @: i.

   mertens0"0 >:i.5 10
 1  0 _1 _1 _2 _1 _2 _2 _2 _1
_2 _2 _3 _2 _1 _1 _2 _2 _3 _3
_2 _1 _2 _2 _2 _1 _1 _1 _2 _3
_4 _4 _3 _2 _1 _1 _2 _1  0  0
_1 _2 _3 _3 _3 _2 _3 _3 _3 _3

Faster and leaner computations of the Mertens function are possible.

Avoiding 0s

mertens1=: 3 : 0
 p=. i.@(1&>.)&.(p:^:_1) 1+<.%:y
 b=. -. > ([ +. #@[ $ {.&1@])~&.>/(<"0 *:p),<(1+y)$0
 +/ _1 ^ #@q: I. b
)

   (mertens1"0 -: mertens1"0) >: i.5 10
1

   mertens1"0 ]10^i.8
1 _1 1 2 _23 _48 212 1037

mertens1 avoids applying the Möbius function to exactly those numbers having a 0 value; that is, all numbers containing a square; that is, all multiples of the squares of 2 3 5 7 11 13 ... up to the largest prime less than or equal to the square root of n . The Möbius function on the remaining numbers (I. b in mertens1) is just the parity of the length of the prime factorization, 1 if even and _1 if odd.

The second line in mertens1 is equivalent to b=. -. +./ (1+y) ($ {.&1)"0 *:p but is much leaner by not creating a large intermediate matrix.

Avoiding Factoring

mq=: 0&$: : (4 : 0)
 p=. i.&.(p:^:_1) 10>.1+<.y%2
 q=. ,. p {.~ p (<:i.0:) <.%:y
 z=. (2<.y) $ (1,p:^:_1)`((i.1 0) ; [: ,. i.&.(p:^:_1))@.x >:y
 for. i. 0 >. <: (*/\p:i.20) I. >:y do.
  p=. p {.~ p (<:i.0:) <.y%*/{.q
  q=. q #~ ({:"1 q)<y%*/"1 q
  j=. 1 + p i. {:"1 q
  k=. j -~ p I. <.1+y%*/"1 q
  q=. (k#q) ,. (j,:"0 k) ;@:(<;.0) p
  z=. z, #`<@.x q
 end.
)

mertens2=: -/ @ mq

mertens2 avoids factoring altogether, but instead enumerates lists of distinct primes of each length. The sub-function mq provides counts if the left argument is 0  and the actual tables of primes if it is 1 ; the monad mq is equivalent to 0&mq .

For example:

   mq 1e6
1 78498 209867 206964 92966 18387 1235 8
   -/ mq 1e6
212
   mertens2 1e6
212

   1 mq 45
┌┬──┬────┬─────┐
││ 2│2  3│2 3 5│
││ 3│2  5│2 3 7│
││ 5│2  7│     │
││ 7│2 11│     │
││11│2 13│     │
││13│2 17│     │
││17│2 19│     │
││19│3  5│     │
││23│3  7│     │
││29│3 11│     │
││31│3 13│     │
││37│5  7│     │
││41│    │     │
││43│    │     │
└┴──┴────┴─────┘
   #&> 1 mq 45
1 14 12 2
   mq 45
1 14 12 2

Given integer n and q , the c-column table of prime factorizations of square-free numbers less than or equal to n , the dyad n nextq q generates the (1+c)-column table of such prime factorizations.

nextq=: 4 : 0
 p=. i.&.(p:^:_1) 1+<.x%*/{.y
 y=. y #~ ({:"1 y)<x%*/"1 y
 j=. (**/$y) * 1 + p i. {:"1 y
 k=. j -~ p I. <.1+x%*/"1 y
 (k#y) ,. (j,:"0 k) ;@:(<;.0) p
)

   45 nextq i.1 0
 2
 3
 5
 7
11
13
17
19
23
29
31
37
41
43
   45 nextq 45 nextq i.1 0
2  3
2  5
2  7
2 11
2 13
2 17
2 19
3  5
3  7
3 11
3 13
5  7
   45 nextq 45 nextq 45 nextq i.1 0
2 3 5
2 3 7
   45 nextq&.>^:(<4) <i.1 0
┌┬──┬────┬─────┐
││ 2│2  3│2 3 5│
││ 3│2  5│2 3 7│
││ 5│2  7│     │
││ 7│2 11│     │
││11│2 13│     │
││13│2 17│     │
││17│2 19│     │
││19│3  5│     │
││23│3  7│     │
││29│3 11│     │
││31│3 13│     │
││37│5  7│     │
││41│    │     │
││43│    │     │
└┴──┴────┴─────┘
   (1 mq 45) -: 45 nextq&.>^:(<4) <i.1 0
1

   y=: ?1e6
   p: i.20
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71
   */\ p: i.20
2 6 30 210 2310 30030 510510 9.69969e6 2.23093e8 6.46969e9 2.0056e11 7.42074e12 3.0425e14 ...

   1 + (*/\p:i.20) I. >:y
8
   (1 mq y) -: y nextq&.>^:(<8) <i.1 0
1

The last equivalence suggests an alternative computation of 1&mq :

qn =: 1 + (*/\p:i.20) I. >:
mq1=: 3 : 'y nextq&.>^:(<qn y) <i.1 0'

   (1&mq -: mq1) ?1e6
1

mqcheck is like mq but incorporates checks.

mqcheck=: 0&$: : (4 : 0)
 z=. x mq y
 assert. ($z) = ,1 + (*/\p:i.20) I. >:y
 if. 0=x do.
  assert. (0 < z) *. z = <.z
  assert. (mertens0 y) = -/ z
  assert. 1={.z
  assert. (((1<y)$1){z) = p:^:_1 >:y
 else.
  assert. (32=3!:0 z) *. 2=#@$&>z
  assert. (i.#z) -: {:@$&> z
  assert. * #&> z
  assert. 1 p: ; ,&.> z
  assert. z -: ~.  &.> z
  assert. z -: ~."1&.> z
  assert. z -: /:~  &.> z
  assert. z -: /:~"1&.> z
  assert. y >: ; */"1&.> z
  assert. (0 mq y) -: #&> z
  assert. (mertens0 y) = -/ #&> z
  assert. ({.z) = <i.1 0
  assert. (((1<y)$1){z) = <,.i.&.(p:^:_1) >:y
  assert. ({.&.>}.z) -: <\p:i.<:#z
 end.
)

Avoiding Factoring (2)

A leaner variation of mertens2 derives by observing that for present purposes, the prime factorization can be replaced by the product and the largest prime in a factorization. Thus:

mr=: 0&$: : (4 : 0)
 p=. i.&.(p:^:_1) 1+<.y%2
 r=. ,. i.&.(p:^:_1) 1+<.%:y
 z=. (2<.y) $ (1,p:^:_1)`((,1) ; i.&.(p:^:_1))@.x >:y
 for. i. 0 >. <: (*/\p:i.20) I. >:y do.
  p=. p {.~ p (<:i.0:) <.y%{.{.r
  r=. r #~ ({:"1 r)<y%{."1 r
  j=. 1 + p i. {:"1 r
  k=. j -~ p I. <.1+y%{."1 r
  r=. (k#{."1 r) (* ,. ]) (j,:"0 k) ;@:(<;.0) p
  z=. z, #`(<@:({."1))@.x r
 end.
)

mertens3=: -/ @ mr

1 mr y are the square-free integers less than or equal to y grouped by the lengths of the prime factorizations. 0 mr y (also mr y) are the number of integers in each group (#&> 1 mr y). For example:

   mr 45
1 14 12 2
   -/ mr 45
_3
   mertens3 45
_3

   1 mr 45
┌─┬─────────────────────────────────────┬──────────────────────────────────┬─────┐
│1│2 3 5 7 11 13 17 19 23 29 31 37 41 43│6 10 14 22 26 34 38 15 21 33 39 35│30 42│
└─┴─────────────────────────────────────┴──────────────────────────────────┴─────┘
   (1 mr 45) -: */"1&.> 1 mq 45
1

   (1&mr -: */"1&.>@(1&mq))"0 >: 2 10 ?@$ 1e6
1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1

mv y exploits 1 mr y to calculate the Mertens function from 1 to y . It is much more efficient than the equivalent mertens3"0 >:i.y . Thus:

mv=: 3 : '+/\ ((#&>t)#($t)$1 _1) (<:;t=. 1 mr y)} y$0'

   mv 30
1 0 _1 _1 _2 _1 _2 _2 _2 _1 _2 _2 _3 _2 _1 _1 _2 _2 _3 _3 _2 _1 _2 _2 _2 _1 _1 _1 _2 _3
   mertens3"0 >: i.30
1 0 _1 _1 _2 _1 _2 _2 _2 _1 _2 _2 _3 _2 _1 _1 _2 _2 _3 _3 _2 _1 _2 _2 _2 _1 _1 _1 _2 _3

   (mv -: mertens3"0@:>:@:i.) >:?1e4
1

   load 'plot'
   plot (>:i.1e6);mv 1e6

Benchmark

Collected Definitions

mobius  =: */ @: - @: ~: @: q:

mertens0=: +/ @: (mobius"0) @: >: @: i.

mertens1=: 3 : 0
 p=. i.@(1&>.)&.(p:^:_1) 1+<.%:y
 b=. -. > ([ +. #@[ $ {.&1@])~&.>/(<"0 *:p),<(1+y)$0
 +/ _1 ^ #@q: I. b
)

mq=: 0&$: : (4 : 0)
 p=. i.&.(p:^:_1) 10>.1+<.y%2
 q=. ,. p {.~ p (<:i.0:) <.%:y
 z=. (2<.y) $ (1,p:^:_1)`((i.1 0) ; [: ,. i.&.(p:^:_1))@.x >:y
 for. i. 0 >. <: (*/\p:i.20) I. >:y do.
  p=. p {.~ p (<:i.0:) <.y%*/{.q
  q=. q #~ ({:"1 q)<y%*/"1 q
  j=. 1 + p i. {:"1 q
  k=. j -~ p I. <.1+y%*/"1 q
  q=. (k#q) ,. (j,:"0 k) ;@:(<;.0) p
  z=. z, #`<@.x q
 end.
)

mertens2=: -/ @ mq

nextq=: 4 : 0
 p=. i.&.(p:^:_1) 1+<.x%*/{.y
 y=. y #~ ({:"1 y)<x%*/"1 y
 j=. (**/$y) * 1 + p i. {:"1 y
 k=. j -~ p I. <.1+x%*/"1 y
 (k#y) ,. (j,:"0 k) ;@:(<;.0) p
)

qn =: 1 + (*/\p:i.20) I. >:
mq1=: 3 : 'y nextq&.>^:(<qn y) <i.1 0'

mqcheck=: 0&$: : (4 : 0)
 z=. x mq y
 assert. ($z) = ,1 + (*/\p:i.20) I. >:y
 if. 0=x do.
  assert. (0 < z) *. z = <.z
  assert. (mertens0 y) = -/ z
  assert. 1={.z
  assert. (((1<y)$1){z) = p:^:_1 >:y
 else.
  assert. (32=3!:0 z) *. 2=#@$&>z
  assert. (i.#z) -: {:@$&> z
  assert. * #&> z
  assert. 1 p: ; ,&.> z
  assert. z -: ~.  &.> z
  assert. z -: ~."1&.> z
  assert. z -: /:~  &.> z
  assert. z -: /:~"1&.> z
  assert. y >: ; */"1&.> z
  assert. (0 mq y) -: #&> z
  assert. (mertens0 y) = -/ #&> z
  assert. ({.z) = <i.1 0
  assert. (((1<y)$1){z) = <,.i.&.(p:^:_1) >:y
  assert. ({.&.>}.z) -: <\p:i.<:#z
 end.
)

mr=: 0&$: : (4 : 0)
 p=. i.&.(p:^:_1) 1+<.y%2
 r=. ,. i.&.(p:^:_1) 1+<.%:y
 z=. (2<.y) $ (1,p:^:_1)`((,1) ; i.&.(p:^:_1))@.x >:y
 for. i. 0 >. <: (*/\p:i.20) I. >:y do.
  p=. p {.~ p (<:i.0:) <.y%{.{.r
  r=. r #~ ({:"1 r)<y%{."1 r
  j=. 1 + p i. {:"1 r
  k=. j -~ p I. <.1+y%{."1 r
  r=. (k#{."1 r) (* ,. ]) (j,:"0 k) ;@:(<;.0) p
  z=. z, #`(<@:({."1))@.x r
 end.
)

mertens3=: -/ @ mr

mv=: 3 : '+/\ ((#&>t)#($t)$1 _1) (<:;t=. 1 mr y)} y$0'

See also

• MathWorld

• On-line Encyclopedia of Integer Sequences A002321

• Wikipedia

• Möbius Function

Contributed by RogerHui.