A381873 a(1) = 1, a(2) = 2; for n > 2, a(n) is the smallest unused positive number that shares a factor with a(n-1) while containing at most two distinct prime factors.

1, 2, 4, 6, 3, 9, 12, 8, 10, 5, 15, 18, 14, 7, 21, 24, 16, 20, 22, 11, 33, 27, 36, 26, 13, 39, 45, 25, 35, 28, 32, 34, 17, 51, 48, 38, 19, 57, 54, 40, 44, 46, 23, 69, 63, 49, 56, 50, 52, 58, 29, 87, 72, 62, 31, 93, 75, 55, 65, 80, 64, 68, 74, 37, 111, 81

Offset: 1

Author

Scott R. Shannon, Mar 09 2025

Keywords

• nonn
• look

Comments

Like the EKG sequence A064413 for the terms studied the primes appear in their natural order, although unlike A064413 some primes p are preceded by 3*p and followed by 2*p. The first time this occurs is a(682) = 1401, a(683) = 467, a(684) = 934, although as n increases this becomes more common.

The primes are all contained in the lowest line of values which has up upward curvature - see the attached image. This leads to it crossing the line a(n) = n and creating the fixed point a(15527). The only other fixed points are 1, 2, 8 and 40, and it is likely no more exist.

No further fixed points through a(8*10^5). - Michael S. Branicky, Mar 10 2025

Examples

			a(23) = 36 = 2^2*3^3 as a(22) = 27 and 36 is unused and shares a factor with 27 while containing two distinct prime factors. Note that 30 = 2*3*5 cannot be chosen as it contains three distinct prime factors; this is the first term to differ from A064413.
		

Links

• Scott R. Shannon, Table of n, a(n) for n = 1..10000
• Scott R. Shannon, Image of the first 30000 terms. The green line is a(n) = n.

Crossrefs

Cf. A070915, A027748, A000977, A064413.

Programs

Python

from math import gcd
from sympy import factorint
from functools import cache
from itertools import count, islice
@cache
def omega(n): return len(factorint(n))
def agen(): # generator of terms
    yield 1
    aset, an, m = {1}, 2, 3
    while True:
        yield an
        aset.add(an)
        an = next(k for k in count(m) if k not in aset and gcd(an, k) > 1 and omega(k) <= 2)
        while m in aset or omega(m) > 2: m += 1
print(list(islice(agen(), 66))) # Michael S. Branicky, Mar 09 2025