cp's OEIS Frontend

This is a front-end for the Online Encyclopedia of Integer Sequences, made by Christian Perfect. The idea is to provide OEIS entries in non-ancient HTML, and then to think about how they're presented visually. The source code is on GitHub.

A100613 Number of elements in the set {(x,y): 1 <= x,y <= n, gcd(x,y) > 1}.

Original entry on oeis.org

0, 1, 2, 5, 6, 13, 14, 21, 26, 37, 38, 53, 54, 69, 82, 97, 98, 121, 122, 145, 162, 185, 186, 217, 226, 253, 270, 301, 302, 345, 346, 377, 402, 437, 458, 505, 506, 545, 574, 621, 622, 681, 682, 729, 770, 817, 818, 881, 894, 953, 990, 1045, 1046, 1117, 1146, 1209
Offset: 1

Views

Author

Douglas Stones (dssto1(AT)student.monash.edu.au), Dec 02 2004

Keywords

Crossrefs

Programs

  • Haskell
    a100613 n = length [()| x <- [1..n], y <- [1..n], gcd x y > 1]
    -- Reinhard Zumkeller, Jan 21 2013
    
  • Mathematica
    f[n_] := Table[ #^2 &[m], {m, 1, n + 1}] - FoldList[Plus, 1, 2 Array[EulerPhi, n, 2]] (* Gregg K. Whisler, Jun 25 2008 *)
  • PARI
    a(n) = sum(i=1, n, sum(j=1, n, gcd(i,j)>1)); \\ Michel Marcus, Jan 30 2017
    
  • Python
    from functools import lru_cache
    @lru_cache(maxsize=None)
    def A100613(n): # based on second formula in A018805
        if n == 0:
            return 0
        c, j = 1, 2
        k1 = n//j
        while k1 > 1:
            j2 = n//k1 + 1
            c += (j2-j)*(k1**2-A100613(k1))
            j, k1 = j2, n//j2
        return n+c-j # Chai Wah Wu, Mar 24 2021

Formula

a(n) = A000290(n) - A018805(n) = A185670(n) + A063985(n). - Reinhard Zumkeller, Jan 21 2013
a(n) = Sum_{k = 2..n} A242114(n,k). - Reinhard Zumkeller, May 04 2014
a(n) ~ kn^2, where k = 1 - 6/Pi^2 = 0.39207... (A229099). - Charles R Greathouse IV, Mar 29 2024