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.

Previous Showing 21-30 of 81 results. Next

A018806 Sum of gcd(x, y) for 1 <= x, y <= n.

Original entry on oeis.org

1, 5, 12, 24, 37, 61, 80, 112, 145, 189, 220, 288, 325, 389, 464, 544, 593, 701, 756, 880, 989, 1093, 1160, 1336, 1441, 1565, 1700, 1880, 1965, 2205, 2296, 2488, 2665, 2829, 3028, 3328, 3437, 3621, 3832, 4152, 4273, 4621, 4748, 5040, 5373, 5597, 5736, 6168
Offset: 1

Views

Author

David W. Wilson

Keywords

Comments

a(n) is also the entrywise 1-norm of the n X n GCD matrix.

Links

Crossrefs

Cf. A000010, A018805, A064951, A344479.

Programs

  • Maple
    N:= 1000 # to get a(1) to a(N)
g:= add(numtheory:-phi(k)*x^k*(1+x^k)/((1-x^k)^2*(1-x)),k=1..N):
S:= series(g, x, N+1):
seq(coeff(S,x,j), j=1..N); # Robert Israel, Jan 14 2015
  • Mathematica
    Table[nn = n;Total[Level[Table[Table[GCD[i, j], {i, 1, nn}], {j, 1, nn}], {2}]], {n, 1, 48}] (* Geoffrey Critzer, Jan 14 2015 *)
  • PARI
    a(n)=2*sum(i=1,n,sum(j=1,i-1,gcd(i,j)))+n*(n+1)/2 \\ Charles R Greathouse IV, Jun 21 2013
  • PARI
    a(n)=sum(k=1,n,eulerphi(k)*(n\k)^2) \\ Charles R Greathouse IV, Jun 21 2013
  • Python
    from sympy import totient
def A018806(n): return sum(totient(k)*(n//k)**2 for k in range(1,n+1)) # Chai Wah Wu, Aug 05 2024

Formula

Sum_{k=1..n} phi(k)*(floor(n/k))^2. - Vladeta Jovovic, Nov 10 2002
a(n) ~ kn^2 log n, with k = 6/Pi^2. - Charles R Greathouse IV, Jun 21 2013
G.f.: Sum_{k >= 1} phi(k)*x^k*(1+x^k)/((1-x^k)^2*(1-x)). - Robert Israel, Jan 14 2015

A276187 Number of subsets of {1,..,n} of cardinality >= 2 such that the elements of each counted subset are pairwise coprime.

Original entry on oeis.org

0, 1, 4, 7, 18, 21, 48, 63, 94, 105, 220, 235, 482, 529, 600, 711, 1438, 1501, 3020, 3211, 3594, 3849, 7720, 7975, 11142, 11877, 14628, 15459, 30946, 31201, 62432, 69855, 76126, 80221, 89820, 91611, 183258, 192601, 208600, 214231, 428502, 431573, 863188, 900563
Offset: 1

Views

Author

Robert C. Lyons, Aug 23 2016

Keywords

Comments

n is prime if and only if a(n) = 2*a(n-1)+n-1. - Robert Israel, Aug 24 2016

Examples

			From _Gus Wiseman_, May 08 2021: (Start)
The a(2) = 1 through a(6) = 21 sets:
  {1,2}   {1,2}    {1,2}     {1,2}      {1,2}
          {1,3}    {1,3}     {1,3}      {1,3}
          {2,3}    {1,4}     {1,4}      {1,4}
         {1,2,3}   {2,3}     {1,5}      {1,5}
                   {3,4}     {2,3}      {1,6}
                  {1,2,3}    {2,5}      {2,3}
                  {1,3,4}    {3,4}      {2,5}
                             {3,5}      {3,4}
                             {4,5}      {3,5}
                            {1,2,3}     {4,5}
                            {1,2,5}     {5,6}
                            {1,3,4}    {1,2,3}
                            {1,3,5}    {1,2,5}
                            {1,4,5}    {1,3,4}
                            {2,3,5}    {1,3,5}
                            {3,4,5}    {1,4,5}
                           {1,2,3,5}   {1,5,6}
                           {1,3,4,5}   {2,3,5}
                                       {3,4,5}
                                      {1,2,3,5}
                                      {1,3,4,5}
(End)

Links

Crossrefs

Cf. A000010, A002088, A018805.
The case of pairs is A015614.
The indivisible instead of coprime version is A051026(n) - n.
Allowing empty sets and singletons gives A084422.
The relatively prime instead of pairwise coprime version is A085945(n) - 1.
Allowing all singletons gives A187106.
Allowing only the singleton {1} gives A320426.
Row sums of A320436, each minus one.
The maximal case is counted by A343659.
The version for sets of divisors is A343655(n) - 1.
A000005 counts divisors.
A186972 counts pairwise coprime k-sets containing n.
A186974 counts pairwise coprime k-sets.
A326675 ranks pairwise coprime non-singleton sets.
Cf. A007360, A063647, A186971, A302696, A305713, A320423, A320430, A326077, A327516, A337485.

Programs

  • Maple
    f:= proc(S) option remember;
    local s, Sp;
    if S = {} then return 1 fi;
    s:= S[-1];
    Sp:= S[1..-2];
    procname(Sp) + procname(select(t -> igcd(t,s)=1, Sp))
end proc:
seq(f({$1..n}) - n - 1, n=1..50); # Robert Israel, Aug 24 2016
  • Mathematica
    f[S_] := f[S] = Module[{s, Sp}, If[S == {}, Return[1]]; s = S[[-1]]; Sp = S[[1;;-2]]; f[Sp] + f[Select[Sp, GCD[#, s] == 1&]]];
Table[f[Range[n]] - n - 1, {n, 1, 50}] (* Jean-François Alcover, Sep 15 2022, after Robert Israel *)
  • PARI
    f(n,k=1)=if(n==1, return(2)); if(gcd(k,n)==1, f(n-1,n*k)) + f(n-1,k)
a(n)=f(n)-n-1 \\ Charles R Greathouse IV, Aug 24 2016
  • Sage
    from sage.combinat.subsets_pairwise import PairwiseCompatibleSubsets
def is_coprime(x, y): return gcd(x, y) == 1
max_n = 40
seq = []
for n in range(1, max_n+1):
    P = PairwiseCompatibleSubsets(range(1,n+1), is_coprime)
    a_n = len([1 for s in P.list() if len(s) > 1])
    seq.append(a_n)
print(seq)

Formula

a(n) = A320426(n) - 1. - Gus Wiseman, May 08 2021

Extensions

Name and example edited by Robert Israel, Aug 24 2016

A343978 Number of ordered 6-tuples (a,b,c,d,e,f) with gcd(a,b,c,d,e,f)=1 (1<= {a,b,c,d,e,f} <= n).

Original entry on oeis.org

1, 63, 727, 4031, 15559, 45863, 116855, 257983, 526615, 983583, 1755143, 2935231, 4776055, 7407727, 11256623, 16498719, 23859071, 33434063, 46467719, 62949975, 84644439, 111486599, 146142583, 187854119, 240880239, 303814503, 382049919, 473813703, 586746719
Offset: 1

Views

Author

Karl-Heinz Hofmann, May 06 2021

Keywords

References

  • Joachim von zur Gathen and Jürgen Gerhard, Modern Computer Algebra, Cambridge University Press, Second Edition 2003, pp. 53-54.

Links

Crossrefs

Column k=6 of A344527.
Cf. A343359, A013664, A018805, A071778, A082540, A082544.
Cf. A342632, A342586, A342935, A342841, A343527, A343193.

Programs

  • PARI
    a(n)={sum(k=1, n+1, moebius(k)*(n\k)^6)} \\ Andrew Howroyd, May 08 2021
  • Python
    from labmath import mobius
def A343978(n): return sum(mobius(k)*(n//k)**6 for k in range(1, n+1))
  • Python
    from functools import lru_cache
@lru_cache(maxsize=None)
def A343978(n):
    if n == 0:
        return 0
    c, j, k1 = 1, 2, n//2
    while k1 > 1:
        j2 = n//k1 + 1
        c += (j2-j)*A343978(k1)
        j, k1 = j2, n//j2
    return n*(n**5-1)-c+j # Chai Wah Wu, May 17 2021

Formula

a(n) = Sum_{k=1..n} mu(k)*floor(n/k)^6.
Lim_{n->infinity} a(n)/n^6 = 1/zeta(6) = A343359 = 945/Pi^6.
a(n) = n^6 - Sum_{k=2..n} a(floor(n/k)). - Seiichi Manyama, Sep 13 2024

Extensions

Edited by N. J. A. Sloane, Jun 13 2021

A063985 Partial sums of cototient sequence A051953.

Original entry on oeis.org

0, 1, 2, 4, 5, 9, 10, 14, 17, 23, 24, 32, 33, 41, 48, 56, 57, 69, 70, 82, 91, 103, 104, 120, 125, 139, 148, 164, 165, 187, 188, 204, 217, 235, 246, 270, 271, 291, 306, 330, 331, 361, 362, 386, 407, 431, 432, 464, 471, 501, 520, 548, 549, 585, 600, 632, 653, 683
Offset: 1

Views

Author

Labos Elemer, Sep 06 2001

Keywords

Comments

Number of elements in the set {(x,y): 1 <= x <= y <= n, 1 = gcd(x,y)}; a(n) = A000217(n) - A002088(n) = A100613(n) - A185670(n). - Reinhard Zumkeller, Jan 21 2013
8*a(n) is the number of dots not in direct reach via a straight line from the center of a 2*n+1 X 2*n+1 array of dots. - Kiran Ananthpur Bacche, May 25 2022

Links

Crossrefs

Cf. A000010, A000217, A002088, A048290, A051953.

Programs

  • Haskell
    a063985 n = length [()| x <- [1..n], y <- [x..n], gcd x y > 1]
-- Reinhard Zumkeller, Jan 21 2013
  • Java
    // Save the file as A063985.java to compile and run
import java.util.stream.IntStream;
import java.util.*;
public class A063985 {
  public static int getInvisiblePoints(int n) {
    Set slopes = new HashSet();
    IntStream.rangeClosed(1, n).forEach(i ->
      {IntStream.rangeClosed(1, n).forEach(j ->
        slopes.add(Float.valueOf((float)i/(float)j))); });
    return (n * n - slopes.size() + n - 1) / 2;
  }
  public static void main(String args[]) throws Exception {
    IntStream.rangeClosed(1, 30).forEach(i ->
      System.out.println(getInvisiblePoints(i)));
  }
} // Kiran Ananthpur Bacche, May 25 2022
  • Mathematica
    f[n_] := n(n + 1)/2 - Sum[ EulerPhi@i, {i, n}]; Array[f, 58] (* Robert G. Wilson v *)
Accumulate[Table[n-EulerPhi[n],{n,1,60}]] (* Harvey P. Dale, Aug 19 2015 *)
  • PARI
    { a=0; for (n=1, 1000, write("b063985.txt", n, " ", a+=n - eulerphi(n)) ) } \\ Harry J. Smith, Sep 04 2009
  • Python
    from sympy.ntheory import totient
def a(n): return sum(x - totient(x) for x in range(1,n + 1))
[a(n) for n in range(1, 51)] # Indranil Ghosh, Mar 18 2017
  • Python
    from functools import lru_cache
@lru_cache(maxsize=None)
def A063985(n): # based on second formula in A018805
    if n == 0:
        return 0
    c, j = 0, 2
    k1 = n//j
    while k1 > 1:
        j2 = n//k1 + 1
        c += (j2-j)*(k1*(k1+1)-2*A063985(k1)-1)
        j, k1 = j2, n//j2
    return (2*n+c-j)//2 # Chai Wah Wu, Mar 24 2021

Formula

a(n) = Sum_{x=1..n} (x - phi(x)) = Sum(x) - Sum(phi(x)) = A000217(n) - A002088(n), phi(n) = A000010(n), cototient(n) = A051953(n).
a(n) = n^2 - A091369(n). - Enrique Pérez Herrero, Feb 25 2012
G.f.: x/(1 - x)^3 - (1/(1 - x))*Sum_{k>=1} mu(k)*x^k/(1 - x^k)^2. - Ilya Gutkovskiy, Mar 18 2017
a(n) = (1/2 - 3/Pi^2)*n^2 + O(n*log(n)). - Amiram Eldar, Jul 26 2022

Extensions

Corrected by Robert G. Wilson v, Dec 13 2006

A082544 Number of ordered quintuples (a,b,c,d,e) with gcd(a,b,c,d,e)=1 (1<= {a,b,c,d,e} <= n).

Original entry on oeis.org

1, 31, 241, 991, 3091, 7501, 16531, 31711, 57781, 96601, 157651, 240031, 362491, 519961, 739201, 1012441, 1383721, 1822711, 2409241, 3091441, 3966301, 4974751, 6257461, 7680781, 9481681, 11474941, 13916191, 16610371, 19911151, 23435191
Offset: 1

Views

Author

Benoit Cloitre, May 11 2003

Keywords

Links

Crossrefs

Column k=5 of A344527.
Cf. A018805 (pairs), A071778 (triples), A082540 (quadruples), A343978.
Cf. A015650.

Programs

  • PARI
    a(n)=sum(k=1,n,moebius(k)*floor(n/k)^5)
  • Python
    from functools import lru_cache
@lru_cache(maxsize=None)
def A082544(n):
    if n == 0:
        return 0
    c, j = 1, 2
    k1 = n//j
    while k1 > 1:
        j2 = n//k1 + 1
        c += (j2-j)*A082544(k1)
        j, k1 = j2, n//j2
    return n*(n**4-1)-c+j # Chai Wah Wu, Mar 29 2021

Formula

a(n) = Sum_{k=1..n} mu(k)*floor(n/k)^5; a(n) is asymptotic to c*n^5 with c=0.9643....
Lim_{n->infinity} a(n)/n^5 = 1/zeta(5) = A343308. - Karl-Heinz Hofmann, Apr 11 2021
Lim_{n->infinity} n^5/a(n) = zeta(5) = A013663. - Karl-Heinz Hofmann, Apr 11 2021
a(n) = n^5 - Sum_{k=2..n} a(floor(n/k)). - Seiichi Manyama, Sep 13 2024

A213207 Number of distinct products i*j*k over all triples (i,j,k) with |i| + |j| + |k| <= n.

Original entry on oeis.org

1, 1, 1, 3, 5, 9, 13, 19, 25, 35, 43, 55, 65, 79, 91, 111, 127, 149, 167, 193, 217, 249, 273, 311, 339, 383, 419, 463, 501, 551, 591, 643, 693, 751, 799, 869, 925, 995, 1057, 1133, 1199, 1281, 1347, 1439, 1515, 1615, 1697, 1801, 1883, 2001, 2101, 2219, 2313
Offset: 0

Views

Author

Robert Price, Mar 01 2013

Keywords

Comments

This sequence is in reply to an extension request made in A100450.

Links

Crossrefs

Cf. A018805, A027430, A100448, A100449, A100450, A213208, A213212, A213213.

Programs

  • Maple
    h:= proc() true end:
b:= proc(n) local c, i, j, p;
      c:=0;
      for i to iquo(n, 3) do
        for j from i to iquo(n-i, 2) do p:= i*j*(n-i-j);
          if h(p) then h(p):= false; c:=c+1 fi
        od
      od; c
    end:
a:= proc(n) a(n):= `if`(n=0, 1, a(n-1) +2*b(n)) end:
seq(a(n), n=0..60);  # Alois P. Heinz, Mar 02 2013
  • Mathematica
    f[n_] := Length[ Union[ Flatten[ Table[ If[ Abs[i] + Abs[j] + Abs[k] <= n, {i*j*k}, {0}], {i, -n, n}, {j, -n, n}, {k, -n, n}], 2]]]; Table[ f[n], {n, 0, 100}]

A213208 Number of distinct products i*j*k over all triples (i,j,k) with |i| + |j| + |k| <= n and gcd(i,j,k) <= 1.

Original entry on oeis.org

1, 1, 1, 3, 5, 9, 11, 19, 23, 33, 39, 51, 57, 75, 87, 103, 117, 143, 155, 187, 207, 235, 259, 297, 319, 363, 395, 441, 473, 525, 555, 615, 659, 721, 765, 831, 875, 959, 1017, 1091, 1147, 1239, 1291, 1397, 1467, 1553, 1631, 1743, 1813, 1937, 2023, 2141, 2233, 2379, 2465
Offset: 0

Views

Author

Robert Price, Mar 01 2013

Keywords

Comments

This sequence is in reply to an extension request made in A100450.
Note that gcd(0,m) = m for any m.

Links

Crossrefs

Cf. A018805, A027430, A100448, A100449, A100450, A213207, A213212, A213213.

Programs

  • Maple
    h:= proc() true end:
b:= proc(n) local c, i, j, p;
      c:=0;
      for i to iquo(n, 3) do
        for j from i to iquo(n-i, 2) do
          if igcd(i, j, n-i-j)=1 then p:= i*j*(n-i-j);
            if h(p) then h(p):= false; c:=c+1 fi
          fi
        od
      od; c
    end:
a:= proc(n) a(n):= `if`(n=0, 1, a(n-1) +2*b(n)) end:
seq(a(n), n=0..60);  # Alois P. Heinz, Mar 01 2013
  • Mathematica
    f[n_] := Length[ Union[ Flatten[ Table[ If[ Abs[i] + Abs[j] + Abs[k] <= n&& GCD[i, j, k] <= 1, i*j*k, 0], {i, -n, n}, {j, -n, n}, {k, -n, n}], 2]]]; Table[ f[n], {n, 0, 100}]

A342586 a(n) is the number of pairs (x,y) with 1 <= x, y <= 10^n and gcd(x,y)=1.

Original entry on oeis.org

1, 63, 6087, 608383, 60794971, 6079301507, 607927104783, 60792712854483, 6079271032731815, 607927102346016827, 60792710185772432731, 6079271018566772422279, 607927101854119608051819, 60792710185405797839054887, 6079271018540289787820715707, 607927101854027018957417670303
Offset: 0

Views

Author

Karl-Heinz Hofmann, Mar 16 2021

Keywords

References

  • Joachim von zur Gathen and Jürgen Gerhard, Modern Computer Algebra, Cambridge University Press, Second Edition 2003, pp. 53-54. (See link below.)

Links

Crossrefs

a(n) = 2*A064018(n) - 1. - Hugo Pfoertner, Mar 16 2021
a(n) = A018805(10^n). - Michel Marcus, Mar 16 2021
Cf. A000010, A059956 (6/Pi^2), A342632, A342841, A343193, A343282.
Related counts of k-tuples:
pairs: A018805, A342632, A342586;
triples: A071778, A342935, A342841;
quadruples: A082540, A343527, A343193;
5-tuples: A343282;
6-tuples: A343978, A344038. - N. J. A. Sloane, Jun 13 2021

Programs

  • PARI
    a342586(n)=my(s, m=10^n); forfactored(k=1,m,s+=eulerphi(k)); s*2-1 \\ Bruce Garner, Mar 29 2021
  • PARI
    a342586(n)=my(s, m=10^n); forsquarefree(k=1,m,s+=moebius(k)*(m\k[1])^2); s \\ Bruce Garner, Mar 29 2021
  • Python
    import math
for n in range (0,10):
     counter = 0
     for x in range (1, pow(10,n)+1):
        for y in range(1, pow(10,n)+1):
            if math.gcd(y,x) ==  1:
                counter += 1
     print(n, counter)
  • Python
    from functools import lru_cache
@lru_cache(maxsize=None)
def A018805(n):
  if n == 1: return 1
  return n*n - sum(A018805(n//j) for j in range(2, n//2+1)) - (n+1)//2
print([A018805(10**n) for n in range(8)]) # Michael S. Branicky, Mar 18 2021

Formula

Lim_{n->infinity} a(n)/10^(2*n) = 6/Pi^2 = 1/zeta(2).

Extensions

a(10) from Michael S. Branicky, Mar 18 2021
More terms using A064018 from Hugo Pfoertner, Mar 18 2021
Edited by N. J. A. Sloane, Jun 13 2021

A342632 Number of ordered pairs (x, y) with gcd(x, y) = 1 and 1 <= {x, y} <= 2^n.

Original entry on oeis.org

1, 3, 11, 43, 159, 647, 2519, 10043, 39895, 159703, 637927, 2551171, 10200039, 40803219, 163198675, 652774767, 2611029851, 10444211447, 41776529287, 167106121619, 668423198491, 2673693100831, 10694768891659, 42779072149475, 171116268699455, 684465093334979, 2737860308070095
Offset: 0

Views

Author

Karl-Heinz Hofmann, Mar 17 2021

Keywords

Examples

			Only fractions with gcd(numerator, denominator) = 1 are counted. E.g.,
  1/2 counts, but 2/4, 3/6, 4/8 ... do not, because they reduce to 1/2;
  1/1 counts, but 2/2, 3/3, 4/4 ... do not, because they reduce to 1/1.
.
For n=0, the size of the grid is 1 X 1:
.
    | 1
  --+--
  1 | o      Sum:  1
.
For n=1, the size of the grid is 2 X 2:
.
    | 1 2
  --+----
  1 | o o          2
  2 | o .          1
                  --
             Sum:  3
.
For n=2, the size of the grid is 4 X 4:
.
    | 1 2 3 4
  --+--------
  1 | o o o o      4
  2 | o . o .      2
  3 | o o . o      3
  4 | o . o .      2
                  --
             Sum: 11
.
For n=3, the size of the grid is 8 X 8:
.
    | 1 2 3 4 5 6 7 8
  --+----------------
  1 | o o o o o o o o     8
  2 | o . o . o . o .     4
  3 | o o . o o . o o     6
  4 | o . o . o . o .     4
  5 | o o o o . o o o     7
  6 | o . . . o . o .     3
  7 | o o o o o o . o     7
  8 | o . o . o . o .     4
                         --
                    Sum: 43

Links

Crossrefs

a(n) = A018805(2^n).
Cf. A000010, A002088, A059956 (6/Pi^2), A064018, A342586.

Programs

  • PARI
    for(n=0,24,my(j=2^n);print1(2*sum(k=1,j,eulerphi(k))-1,", ")) \\ Hugo Pfoertner, Mar 17 2021
  • Python
    import math
for n in range (0, 21):
     counter = 0
     for x in range (1, pow(2, n)+1):
        for y in range(1, pow(2, n)+1):
            if math.gcd(y, x) ==  1:
                counter += 1
     print(n, counter)
  • Python
    from sympy import sieve
def A342632(n): return 2*sum(t for t in sieve.totientrange(1,2**n+1)) - 1 # Chai Wah Wu, Mar 23 2021
  • Python
    from functools import lru_cache
@lru_cache(maxsize=None)
def A018805(n):
  if n == 1: return 1
  return n*n - sum(A018805(n//j) for j in range(2, n//2+1)) - (n+1)//2
print([A018805(2**n) for n in range(25)]) # Michael S. Branicky, Mar 23 2021

Formula

Lim_{n->infinity} a(n)/2^(2*n) = 6/Pi^2 = 1/zeta(2).

Extensions

Edited by N. J. A. Sloane, Jun 13 2021

A171503 Number of 2 X 2 integer matrices with entries from {0,1,...,n} having determinant 1.

Original entry on oeis.org

0, 3, 7, 15, 23, 39, 47, 71, 87, 111, 127, 167, 183, 231, 255, 287, 319, 383, 407, 479, 511, 559, 599, 687, 719, 799, 847, 919, 967, 1079, 1111, 1231, 1295, 1375, 1439, 1535, 1583, 1727, 1799, 1895, 1959, 2119, 2167, 2335, 2415, 2511, 2599
Offset: 0

Views

Author

Jacob A. Siehler, Dec 10 2009

Keywords

Comments

Number of distinct solutions to k*x+h=0, where |h|<=n and k=1,2,...,n. - Giovanni Resta, Jan 08 2013.
Number of reduced rational numbers r/s with |r|<=n and 0Juan M. Marquez, Apr 13 2015

Links

Crossrefs

Cf. A062801, A000010, A018805. Differences are A002246.
See A326354 for an essentially identical sequence.

Programs

  • Maple
    with(numtheory):
a:= proc(n) option remember;
       `if`(n<2, [0, 3][n+1], a(n-1) + 4*phi(n))
    end:
seq(a(n), n=0..60);
  • Mathematica
    a[n_]:=Count[Det/@(Partition[ #,2]&/@Tuples[Range[0,n],4]),1]
(* Second program: *)
a[0] = 0; a[1] = 3; a[n_] := a[n] = a[n-1] + 4*EulerPhi[n];
Table[a[n], {n, 0, 60}] (* Jean-François Alcover, Jun 16 2018 *)
  • PARI
    a(n)=(n>0)+2*sum(k=1, n, moebius(k)*(n\k)^2) \\ Charles R Greathouse IV, Apr 20 2015
  • Python
    from functools import lru_cache
@lru_cache(maxsize=None)
def A171503(n): # based on second formula in A018805
    if n == 0:
        return 0
    c, j = 0, 2
    k1 = n//j
    while k1 > 1:
        j2 = n//k1 + 1
        c += (j2-j)*(A171503(k1)-1)//2
        j, k1 = j2, n//j2
    return 2*(n*(n-1)-c+j) - 1 # Chai Wah Wu, Mar 25 2021

Formula

Recursion: a(n) = a(n - 1) + 4*phi(n) for n > 1, with phi being Euler's totient function. - Juan M. Marquez, Jan 19 2010
a(n) = 4 * A002088(n) - 1 for n >= 1. - Robert Israel, Jun 01 2014

Extensions

Edited by Alois P. Heinz, Jan 19 2011
Previous Showing 21-30 of 81 results. Next

Started in 1964 by Neil J. A. Sloane | Maintained by The OEIS Foundation Inc.

Content is available under The OEIS End-User License Agreement.