The Möbius Function (and squarefree numbers)

by Ed Pegg Jr, www.wolfram.com, www.mathpuzzle.com, www.maa.org
ed@mathpuzzle.com, edp@wolfram.com
Useful sites for this notebook include the following:
http://mathworld.wolfram.com/MertensFunction.html
http://mathworld.wolfram.com/MoebiusFunction.html
http://functions.wolfram.com/NumberTheoryFunctions/MoebiusMu/

The Möbius Function

In 1831, August Ferdinand Möbius put numbers into 3 bins, as a new type of function. His procedure is represented by the Greek letter mu (µ), and is called the Möbius function in his honor.

Table[MoebiusMu[n], {n, 1, 100}]

In the "zero" bin, Möbius placed multiples of square numbers (other than 1).

ZeroBin  =  
      Flatten[Position[Table[MoebiusMu[n], {n, 1, 100}], 0]]

Five square-multiples in a row.

Flatten[
      Position[Partition[Table[MoebiusMu[n], {n, 1, 10000}], 5, 1], {0, 0, 0, 
          0, 0}]]

MatrixForm[FactorInteger[844 + Range[0, 4]]]

Differences between the actual and predicted density of square-divisible numbers.

ListPlot[
        FoldList[Plus, 0, 
            Table[If[MoebiusMu[n] == 0, 1, 0], {n, 1, 100000 - 1}]] - 
          Table[n (1 - 6/Pi^2), {n, 1, 100000}], 
        PlotJoined -> True, AspectRatio -> 1/7];

In the "negative one" bin, Möbius placed any number that factored into an odd number of distinct primes.

NegativeOneBin  =  
      Flatten[Position[Table[MoebiusMu[n], {n, 1, 100}], -1]]

The number 665499549999999999 was found by using FromDigits and IntegerDigits, and selecting the numbers that satisfied MoebiusMu[n]=-1.  When I noticed that some numbers could end in arbitrarily long strings of nines, I tried to find the best one using only 4 different digits.

Table[
      MoebiusMu[FromDigits[Take[IntegerDigits[665499549999999999], n]]], {n, 
        2, 18}]

MatrixForm[
      Table[Map[First, 
          FactorInteger[
            FromDigits[Take[IntegerDigits[665499549999999999], n]]]], {n, 2, 
          18}]]

Differences between the real and predicted density of µ(n) = -1.

ListPlot[
        FoldList[Plus, 0, 
            Table[If[MoebiusMu[n] == -1, 1, 0], {n, 1, 
                100000 - 1}]] - Table[n (3/Pi^2), {n, 1, 100000}], 
        PlotJoined -> True, AspectRatio -> 1/7];

In the "positive one" bin, Möbius put all the numbers that factor into an even number of distinct primes. For completeness, he put 1 into this bin.

PositiveOneBin  =  
      Flatten[Position[Table[MoebiusMu[n], {n, 1, 100}], 1]]

IntegerDigits[69999069090909090990909090] // Length

I found 69999069090909090990909090 with the same method as above.  This time, I restricted the search to numbers with only three digits.  I believe this is the longest possible.  No infinite

Table[
      MoebiusMu[
        FromDigits[Take[IntegerDigits[69999069090909090990909090], n]]], {n, 
        2, 26}]

MatrixForm[
      Table[Map[First, 
          FactorInteger[
            FromDigits[
              Take[IntegerDigits[69999069090909090990909090], n]]]], {n, 2, 
          26}]]

For amazing factorability demonstrations, use PowerMod.

PowerMod[10, 10^100, 316912650057057350374175801344000001]

More Riemann Zeta function zeroes can be found at Andrew Odlyzko: Tables of zeros of the Riemann zeta function.

(Subscript[ζ, 0]) =1/2  +  
            I {{14.134725142}, {21.022039639}, {25.010857580}, {30.424876126}, {32.935061588}, {37.586178159}, {40.918719012}, {43.327073281}, {48.005150881}, {49.773832478}, {52.970321478}, {56.446247697}, {59.347044003}, {60.831778525}, {65.112544048}, {67.079810529}, {69.546401711}, {72.067157674}, {75.704690699}, {77.144840069}, {79.337375020}, {82.910380854}, {84.735492981}, {87.425274613}, {88.809111208}, {92.491899271}, {94.651344041}, {95.870634228}, {98.831194218}, {101.317851006}, {103.725538040}, {105.446623052}, {107.168611184}, {111.029535543}, {111.874659177}, {114.320220915}, {116.226680321}, {118.790782866}, {121.370125002}, {122.946829294}, {124.256818554}, {127.516683880}, {129.578704200}, {131.087688531}, {133.497737203}, {134.756509753}, {138.116042055}, {139.736208952}, {141.123707404}, {143.111845808}, {146.000982487}, {147.422765343}, {150.053520421}, {150.925257612}, {153.024693811}, {156.112909294}, {157.597591818}, {158.849988171}, {161.188964138}, {163.030709687}, {165.537069188}, {167.184439978}, {169.094515416}, {169.911976479}, {173.411536520}, {174.754191523}, {176.441434298}, {178.377407776}, {179.916484020}, {182.207078484}, {184.874467848}, {185.598783678}, {187.228922584}, {189.416158656}, {192.026656361}, {193.079726604}, {195.265396680}, {196.876481841}, {198.015309676}, {201.264751944}, {202.493594514}, {204.189671803}, {205.394697202}, {207.906258888}, {209.576509717}, {211.690862595}, {213.347919360}, {214.547044783}, {216.169538508}, {219.067596349}, {220.714918839}, {221.430705555}, {224.007000255}, {224.983324670}, {227.421444280}, {229.337413306}, {231.250188700}, {231.987235253}, {233.693404179}, {236.524229666}};

The (Subscript[ζ, 0]) powers correlate nicely with MoebiusMu. I do not know if this is a new technique. To me, it was an obvious consequence of (∑_(n = 1)^∞)((μ(n))/((n^s)))= (1/(ζ(s))). Note that it only has a tendancy to be correct.

Map[{First[#], Length[#]} &, 
      Split[Sort[
          Table[{MoebiusMu[
                n], Sign[Re[n^(Subscript[ζ, 0])[[1]]]][[1]]}, {n, 1, 
              1000}]]]]

Table[
      If[MoebiusMu[n] == 1, 
        Sum[Re[n^(Subscript[ζ, 0])[[k]]], {k, 1, 100}], 0], {n, 1, 
        100}]

Table[
      If[MoebiusMu[n] == -1, 
        Sum[Re[n^(Subscript[ζ, 0])[[k]]], {k, 1, 100}], 0], {n, 1, 
        100}]

ListPlot[
        FoldList[Plus, 0, 
          Flatten[Table[
              Sign[Re[1/n^(Subscript[ζ, 0])[[1]]]] MoebiusMu[n], {n, 1, 
                100000}]]], PlotJoined -> True, 
        AspectRatio -> 1/7];

ListPlot[
        FoldList[Plus, 0, 
          Flatten[Table[
              Sign[Re[1/n^(Subscript[ζ, 0])[[45]]]] MoebiusMu[n], {n, 1, 
                100000}]]], PlotJoined -> True, 
        AspectRatio -> 1/7];

ListPlot[
        FoldList[Plus, 0, 
          Flatten[Table[
              Sign[Re[1/n^(Subscript[ζ, 0])[[92]]]] MoebiusMu[n], {n, 1, 
                100000}]]], PlotJoined -> True, 
        AspectRatio -> 1/7];

For n = 1 to 20, µ(n) = {1, -1, -1, 0, -1, 1, -1, 0, 0, 1, -1, 0, -1, 1, 1, 0, -1, 0, -1, 0}. The cumulative sum is thus {1, 0, -1, -1, -2, -1, -2, -2, -2, -1, -2, -2, -3, -2, -1, -1, -2, -2, -3, -3}. This is known as Mertens Function, or M(x). In 1885, Stieltjes claimed to have a proof that |M(x)/(Sqrt[x])| < 1 for all x. In 1897, Mertens said the same thing, and now gets all the credit. The first time |M(x)/(Sqrt[x])| exceeds ½ is at M(7725030629) = 43947.

ListPlot[
        Drop[FoldList[Plus, 0, Table[MoebiusMu[n], {n, 1, 100000}]], 1]/
          Sqrt[Range[100000]], PlotJoined -> True, 
        AspectRatio -> 1/7, 
        PlotRange -> {- .5,  .5}];

ComplexCircle  =  
        Take[Sort[
            Flatten[Table[a + b I, {a, -100, 100}, {b, -100, 100}]], 
            Abs[#1] < Abs[#2] &], 10000];

The density of the squarefree Gaussian integers is (6/((π^2) C)), where C is the Catalan constant.

N[(6/((π^2) C)),50]

Here's a plot showing a convergence to that value.

ListPlot[
        FoldList[Plus, 0, 
            Table[If[
                Max[Map[Last, 
                      FactorInteger[ComplexCircle[[n]]]]]  ==  1, 
                1, 0], {n, 1, 10000}]]/Range[10001], PlotJoined -> True, 
        AspectRatio ->  1/7];

There is a pattern to how Re[1/n^(14.1347251417347 I+1/2)] changes sign.

Plot[Re[1/n^(14.1347251417347 I + 1/2)], {n, 1, 2000}, 
        AspectRatio ->  1/7,  
        PlotRange ->  {- .1,  .1}];

First, we find the sign changes.

signchanges  = 
      Flatten[Position[
          Table[Sign[
              Re[(1/(n^(14.1347251417347 I + (1/2))))] 
                Re[(1/(((n + 1))^(14.1347251417347 I + (1/2))))]], {n, 1, 2000}], -1]]

Next, we note that there is a ratio between successive terms.

N[Drop[signchanges, 1]/Drop[signchanges, -1]]

Once the ratio is determined, a constant can be sought.

signchanges/Table[1.2489^k, {k, 1, 30}]

Table[Round[2.17607  1.2489^k], {k, 1, 30}]  -  
      signchanges

Good enough.  Now, we can use logarithms to tackle (10^((10^100)))+1.  So far as I know, no-one has ever predicted the number of factors for a number this large before.  Note that k represents the number of sign changes.

Solve[2.17607  1.2489^k  ==  10^x, k]

Round[4.49917 Log[0.459544]  +  4.49917 Log[10] 10^100]

We are interested in the number of sign changes.  Out near a Googleplex, the sign has changed an even number of times, so there is a slight suggestion that googolplex+1 has an even number of factors.  Much more accuracy would be desirable.

RasterGraphics code by Oyvind Tafjord of Wolfram Research.

RasterGraphics[state_, colors_:  2, size_:  1] := 
      With[{dim = Reverse[Dimensions[state]]}, 
        Graphics[Raster[Reverse[1 - state/(colors - 1)]], 
          AspectRatio -> Automatic, 
          PlotRange -> {{0, dim[[1]]}, {0, dim[[2]]}}, 
          ImageSize -> size*dim + 1]]

GaussianMu[a_]  :=  
      With[{jj = FactorInteger[a, GaussianIntegers ->  True]}, 
        With[{kk = 
              If[jj == {}, {}, 
                If[Abs[jj[[1, 1]]] == 1, Drop[jj, 1], jj]]}, 
          If[Max[Map[Last, kk]] > 1, 0, 
            2  Mod[Length[kk] + 1, 2] - 1]]]

Nice pictures of the behavior of the Moebius function on Gaussian integers can be made.

Show[
        RasterGraphics[
          Table[If[GaussianMu[a + b I] == -1, 1, 0], {a, -200, 
              200}, {b, -200, 200}]]];

Integrate[ArcTan[x]/x, {x, 0, 1}]

(Zeta[2, 1/4]  -  Zeta[2, 3/4])/16

16/(Zeta[2] (Zeta[2, 1/4]  -  Zeta[2, 3/4])) // 
      FullSimplify

FactorInteger[11 + 3  I]

Gauss observed that for numbers one less than a prime, MoebiusMu[p-1] is equal to the sum of the primitive roots of p, mod p.

\(PrimitiveRootQ[a_Integer, p_Integer] := 
      Block[{fac, res}, 
        fac = FactorInteger[p - 1]; 
res = 
          Table[PowerMod[a, \((p - 1)\)/fac[\([i, 1]\)], p], {i, 
              Length[fac]}]; \(! MemberQ[res, 
            1]\)]
\)
\(PrimitiveRoots[p_Integer] := 
      Select[Range[p - 1], 
        PrimitiveRootQ[#, p] &]
\)
Table[{Prime[n] - 1, MoebiusMu[Prime[n] - 1], 
        Mod[Total[PrimitiveRoots[Prime[n]]], Prime[n], -1]}, {n, 2, 
        100}]