# Python code for the OEIS Sequence A089580 # Coder: Karl-Heinz Hofmann, Germany, mailto@kallimusic.de # Date: Sep 17 2023 # Special feature: Pure counting (with 10^x jumps if m and m+10^x have the same "A060298 signature") # and without using a formula. import numpy from math import isqrt import time #██████ Main parameter █████████████████████████████████████████████████ # ██ nupto = 28 # 28: 1 second, 60: 3 seconds, 100: 19 seconds, ██ # 200: 285 seconds, 300: 1376 seconds, 400: 4541 seconds ██ # ██ #███████████████████████████████████████████████████████████████████████ t1 = time.time() A060298 = numpy.zeros(nupto,dtype=object) m = 2 m_limit = isqrt(10**nupto) while m <= m_limit: signature = numpy.zeros(nupto, dtype=object) k = 2 while m**k < 10**nupto: signature[len(str(m**k))-1] += 1 k += 1 pot = nupto//2 while pot > 0: signature_pot = numpy.zeros(nupto, dtype=object) k = 2 while (m + 10**pot)**k < 10**nupto: signature_pot[len(str((m + 10**pot)**k))-1] += 1 k += 1 if all(signature == signature_pot): A060298 += signature*(10**pot) m += 10**pot #print("Now jumping: ",10**pot, " m´s with same "'"A060298 signature"'"") # Slow!! if switched on pot = 1 pot -= 1 if pot == 0: A060298 += signature m += 1 A089580 = list(numpy.cumsum(A060298)) print(A089580) t2 = time.time() print("Time needed: < ",int(t2-t1)+1, " seconds")