A Revision

In my previous post, I made a modification to the following algorithm:

Suppose we take any positive integer \(n \gt 1\) 

• if prime, double it and add 1: \(n \rightarrow 2n+1\)
• if composite, determine its number of divisors \(d\)
• if \( n \pmod d \equiv 0\) then \(n \rightarrow \dfrac{n}{d} \)
• if \( n \pmod d \not\equiv 0 \) then \(n \rightarrow n \times d\)

Keep repeating this process until a loop is reached or call a stop after a fixed number of iterations. 

The modification I made prevented the numbers generated from becoming too large too quickly. Instead of doubling a prime and adding 1, I decided to do this only if the number was a \(4k+1\) prime. If it was a \(4k+3\), I subtracted 1 and divided by 2. The new algorithm looks like this:

• if a \(4k+1\) prime, double it and add 1: \(n \rightarrow 2n+1\)
• if a \(4k+3\) prime, subtract 1 and divide by 2: \(n \rightarrow (n-1)/2\)
• if composite, determine its number of divisors \(d\)
• if \( n \pmod d \equiv 0\) then \(n \rightarrow \dfrac{n}{d} \)
• if \( n \pmod d \not\equiv 0 \) then \(n \rightarrow n \times d\) 

Here is an example using 28069 (permalink):

--- Loop detected at value 149708 ---

Divisors to Sequence:

28069, 56139, 224556, 18713, 37427, 149708, 1796496, 71859840, 561405, 8982480, 112281, 898248, 28743936, 2069563392, 10778976, 149708

------------------------------

Sequence Length: 16

Highest Value:   2069563392

This method of dealing with primes is also better suited to the sequences I mentioned earlier in my two posts: 

• Number's Factors to Sequence Algorithm 1
• Number's Factors to Sequence Algorithm 2

Here are the two new algorithms. 

Suppose we take any positive integer \(n \gt 1\) and apply the following rules to it:

• if a \(4k+1\) prime, double it and add 1: \(n \rightarrow 2n+1\)
• if a \(4k+3\) prime, subtract 1 and divide by 2: \(n \rightarrow (n-1)/2\)
• if composite, determine its number of factors \(f\) counted \( \textbf{with multiplicity}\)
• if \( n \pmod f \equiv 0\) then \(n \rightarrow \dfrac{n}{f} \)
• if \( n \pmod f \not\equiv 0 \) then \(n \rightarrow n \times f\)

Keep repeating this process until a loop is reached or call a stop after a fixed number of iterations. 

Here is an example using 28069 (permalink):

--- Loop detected at value 37427 ---

Number of factors to sequence with multiplicity:

28069, 56139, 112278, 37426, 18713, 37427, 74854, 224562, 898248, 149708, 37427

--------------------

Sequence Length: 11 

Highest Value:   898248 

Suppose we take any positive integer \(n \gt 1\) and apply the following rules to it:

• if a \(4k+1\) prime, double it and add 1: \(n \rightarrow 2n+1\)
• if a \(4k+3\) prime, subtract 1 and divide by 2: \(n \rightarrow (n-1)/2\)
• if composite, determine its number of factors \(f\) counted \( \textbf{without multiplicity}\)
• if \( n \pmod f \equiv 0\) then \(n \rightarrow \dfrac{n}{f} \)
• if \( n \pmod f \not\equiv 0 \) then \(n \rightarrow n \times f\)

Keep repeating this process until a loop is reached or call a stop after a fixed number of iterations. 

Here is an example using 28069 (permalink):

--- Loop detected at value 224562 ---

28069, 56139, 112278, 37426, 18713, 37427, 74854, 224562, 898248, 224562

--------------------

Sequence Length: 10

Highest Value:   898248

A Correction

After creating my post Number's Divisors to Sequence Algorithm, I was feeling satisfied. I'd gotten Gemini to create a neat little table for me. It showed the record lengths reached by numbers under the algorithm in the range up to one million but suffered from the fact that it was simply wrong. Here is the original table:

Number     | Length     | Status
-----------------------------------
2          | 5          | New Record!     
3          | 8          | New Record!     
6          | 10         | New Record!     
11         | 21         | New Record!     
22         | 23         | New Record!     
44         | 25         | New Record!     
50         | 28         | New Record!     
99         | 32         | New Record!     
125        | 33         | New Record!     
206        | 34         | New Record!     
350        | 37         | New Record!     
463        | 44         | New Record!     
487        | 46         | New Record!     
974        | 48         | New Record!     
1375       | 51         | New Record!     
1573       | 52         | New Record!     
1625       | 60         | New Record!     
5200       | 62         | New Record!     
14157      | 63         | New Record!     
16879      | 64         | New Record!     
18747      | 67         | New Record!     
39325      | 70         | New Record!     
89237      | 71         | New Record!     
151911     | 75         | New Record!     
563553     | 77         | New Record!     
803133     | 82         | New Record!  

Fortunately I spotted the problem when looking at the trajectory of the number associated with my diurnal age: 28068. It showed a length of 201, far in excess of the 82 associated with 803133. Furthermore, the algorithm only counted up to 201 and so the sequence was longer than that.

Naturally I confronted Gemini about the problem and here was its response:

The Reason for the DiscrepancyThe sequence for 28068 exhibits a behavior known as a Cunningham Chain.Your rule if is_prime: number = number * 2 + 1 matches the definition of a Cunningham Chain (where $p \to 2p+1$).If the sequence hits a "vein" of consecutive primes (e.g., a chain of 40 primes), the number grows by a factor of $2^{40}$ (approx. 1 trillion). My previous Python code had a safety brake: if number > \(10^12\): break. 28068 likely skyrockets past  one trillion ($10^{12}$) very quickly. My code hit this safety limit and stopped counting, reporting a truncated (incorrect) length.

Gemini proposed modifications to the code but these overtaxed the Jupyter notebook and it got stuck at a point shown in Figure 1.

Figure 1

Eventually I decided to get Gemini to write the code for just a single number as input and not a range. Even this proved too much for the Jupyter notebook. I don't know whether the sequence for 28068 goes on forever or not but I do know that the original algorithm that Gemini provided me with was flawed. It was only by chance that I found this out so caution is advised in accepting anything Gemini offers up. I've commented in the past on its coding limitations.

The failure of the algorithm to reach completion is attributed to Cunningham chains so I'd thought I'd better look back at that topic. It turns out that the longest known Cunningham chain of the first kind (\(2p+1\) is 17 primes long so why does Gemini mention chains of 40 primes when none are known? Nonetheless, it's likely that the \(2p+1\) rule that I'm applying to prime numbers does lead quickly to very large number. My choice of this rule was quite arbitrary. Given that \(4k+1\) and \(4k-1\) primes are pretty much in equal abundance, I could apply an alternative rule such as the following to a prime \(p\):

• if \(p \pmod 4 \equiv 1\) then \(p \rightarrow 2 \times p +1 \)
• if \(p \pmod 4 \equiv 3\) then \(p \rightarrow (p -1 )/2 \)

This should keep the progressive numbers from growing too large. I should do this with the factors as well (see earlier posts Number's Factors to Sequence Algorithm 1 and Number's Factors to Sequence Algorithm 2). Figure 2 shows the record breaking numbers under these new rules (up to one million).

Figure 2: permalink

In summary, the record breaking numbers are 1, 2, 4, 13, 17, 34, 50, 98, 294, 650, 722, 2166, 4751, 5313, 9502, 11979, 19773, 46137, 125229, 257049, 385573, 714025.

I also got Gemini to write the code for the input of a single number and the trajectory of the number as output. Let's use 28068 as an example. Figure 2 shows the output.

Figure 3: permalink

ADDENDUM on 16th of February 2026:

It appears that this latest program to determine the record lengths under the new algorithm is faulty. If 28077 is entered the program crashes because the numbers become too large. I confronted Gemini with this discovery and it pointed out that number size is capped at $10^{25}$) and some sequences explode so quickly that they exceed the number cap before they reach a record length. Gemini modified the product to detect numbers leading to these exploding sequences. There are quite a few. Figure 4 shows the 27 of them in the range from 28004 to 28449. This gives some idea of their frequency (27 out 445 or a little over 6%).

Figure 4: permalink

Why Is 313131 An Interesting Number?

Recently a friend of mine who was staying at a hotel revealed that the six digit security code for the room was 313131. This looked like an easy code to crack and I was reminded of a post that I'd made recently titled Passcodes and Repeated Digits in November of 2025. Figure 1 shows the probabilities for the number of distinct digits chosen:

Figure 1

So the code 313131 is indeed easy to crack, requiring only a maximum of 62 attempts or 31 attempts on average. However, what was of most interest to me regarding 313131 was its prime factorisation:$$ 313131=3 \times 7 \times 13 \times 31 \times 37$$If we rearrange the order of multiplication we get the following:$$ 313131=7 \times 3 \times 13 \times 31 \times 37$$Concatenating these digits we get the number \(73133137\) which is palindromic.

Naturally, I wondered how many other numbers share this property and so I set Gemini to investigate this with the following prompt:

Write a program in Python, tailored for insertion into SageMathCell or a Jupyter notebook, that identifies all composite numbers from 4 up to a user chosen upper limit with the property that the prime factors, when arranged in a suitable order and then concatenated, produce a palindromic number. The default upper limit can be set to 40000. An example of such a number would be 313131 = 3 x 7 x 13 x 31 x 37 that can written as 7 x 3 x 13 x 31 x 37 to produce the palindromic number 73133137 when concatenated. The output of the program should be a table showing number and factorisation followed by a comma separated list of the qualifying numbers.

Naturally this a processor intensive process and my Jupyter notebook struggled mightily to generate the list of suitable numbers. After some time, I decided to restrict the range to between 28000 and 29000 and after further modifications by Gemini, I was able to achieve the following list of numbers:

28072, 28125, 28194, 28224, 28242, 28273, 28308, 28322, 28332, 28416, 28431, 28448, 28585, 28589, 28593, 28601, 28602, 28609, 28620, 28672, 28685, 28692, 28750, 28800, 28812, 28847, 28951

Figure 2 shows the details:

Figure 2

I've incorporated this algorithm into my daily search program that has become somewhat impressive if I do say so myself. See Figure 3 for the results of 28067, my diurnal age today.

Figure 3

The output is generated by the same program but I've just restricted the range to my daily number, rather lazy but it does the job. Remember the prime factorisation takes multiplicity into account. For completeness I'll include the Python code:

import sys

def get_prime_factors(n):
    """Returns a list of prime factors of n."""
    factors = []
    d = 2
    temp = n
    while d * d <= temp:
        while temp % d == 0:
            factors.append(d)
            temp //= d
        d += 1
    if temp > 1:
        factors.append(temp)
    return factors

def is_palindrome(s):
    """Checks if a string is a palindrome."""
    return s == s[::-1]

def distinct_permutations(iterable):
    """
    Yields unique permutations of items in iterable.
    This handles repeated elements efficiently (e.g., [2, 2, 3])
    without generating all n! redundant combinations.
    """
    # Sort the list to start with the lexicographically first permutation
    s = sorted(iterable)
    yield tuple(s)
    
    n = len(s)
    while True:
        # 1. Find the largest index i such that s[i] < s[i+1]
        i = n - 2
        while i >= 0 and s[i] >= s[i+1]:
            i -= 1
        
        if i == -1:
            return # All permutations generated
        
        # 2. Find the largest index j such that s[i] < s[j]
        j = n - 1
        while s[j] <= s[i]:
            j -= 1
        
        # 3. Swap s[i] and s[j]
        s[i], s[j] = s[j], s[i]
        
        # 4. Reverse the sequence from s[i+1] up to the end
        s[i+1:] = s[i+1:][::-1]
        
        yield tuple(s)

def find_palindromic_composites(start_n, end_n):
    results = []
    actual_start = max(4, start_n)

    for n in range(actual_start, end_n + 1):
        factors = get_prime_factors(n)
        
        # Skip primes
        if len(factors) < 2:
            continue
            
        found_property = False
        winning_perm = None
        winning_concat = None
        
        # Use the optimized generator instead of itertools.permutations
        for p in distinct_permutations(factors):
            concat_str = "".join(map(str, p))
            
            if is_palindrome(concat_str):
                found_property = True
                winning_perm = p
                winning_concat = concat_str
                break 
        
        if found_property:
            fact_str = " x ".join(map(str, winning_perm))
            results.append((n, fact_str, winning_concat))

    return results

# --- Configuration ---
LOWER_BOUND = 28000
UPPER_BOUND = 29000

# --- Main Execution ---
print(f"Searching for composite numbers between {LOWER_BOUND} and {UPPER_BOUND}...")
print("-" * 80)
print(f"{'Number':<8} | {'Valid Factor Ordering':<35} | {'Resulting Palindrome'}")
print("-" * 80)

qualifying_numbers = find_palindromic_composites(LOWER_BOUND, UPPER_BOUND)

if not qualifying_numbers:
    print(f"No qualifying numbers found in the range {LOWER_BOUND} to {UPPER_BOUND}.")
else:
    for num, factor_str, pal_str in qualifying_numbers:
        print(f"{num:<8} | {factor_str:<35} | {pal_str}")

    print("-" * 80)
    print("\nComma separated list of qualifying numbers:")
    print(", ".join(str(r[0]) for r in qualifying_numbers))

New Project

As I've written in earlier posts, the site \( \textbf{Numbers Aplenty}\) has been unreliable for quite some time now, often going down for days. For this reason, I developed a SageMath program that produced a similar but even more extensive output. It does the job but my coding skills are somewhat primitive and I made the mistake of getting Gemini to revise my code which it struggled with. It seems that Gemini is more comfortable coding in Python rather than SageMath that is built on top of Python. Sometime the SageMath code malfuntioned. I think it's better to get Gemini to stick with Python when coding.

I asked Gemini the following question:

Numbers Aplenty is a website that I regularly make use of but of late it has become increasingly unreliable because of down time. Is it possible to create a Python program that replicates the output of Numbers Aplenty and also allows me to add additional output of my own. I’m just asking about the feasibility of such a program, you don’t have to create any code. The URL is numbersaplenty.com.

 Yes, it is absolutely feasible to create a Python program that replicates the functionality of Numbers Aplenty. In fact, this is an excellent project for Python because the language excels at the exact type of tasks the website performs: handling large integers and running logical checks.

Here is a breakdown of why this is feasible and how you would likely approach it:

1. Python is Built for This

• Arbitrary Precision Integers: Unlike many other languages, Python handles massive integers automatically. If you want to check a number with 50 digits (far beyond Numbers Aplenty's 15-digit limit), Python can handle it natively without needing special code.

• No Downtime: Since the code runs locally on your machine, you never have to worry about the website being offline.

2. Libraries Do the Heavy Lifting

You would not need to write the complex math from scratch. The Python ecosystem has powerful libraries that already "know" number theory:

• SymPy: This is the most critical library for your goal. It is a symbolic mathematics library that includes built-in functions for:

  • is_prime(n) (Checks if a number is prime)

  • factorint(n) (Finds prime factorization)

  • divisors(n) (Lists all divisors)

  • is_square(n), is_perfect(n), etc.

• NumPy: Useful if you decide to generate generating arrays or grids of numbers for comparison.

3. Replicating the "Families" (The "Aplenty" Part)

The website categorizes numbers into families (e.g., "Happy Numbers," "Fibonacci Numbers," "Palindromic Numbers").

• Standard Families: Checking if a number is a Fibonacci number or a Palindrome is trivial in Python (often 1–3 lines of code).

• Exotic Families: For the more obscure types (like "Vampire numbers" or "Zigzag numbers"), you can look up their definitions (often found on the OEIS - On-Line Encyclopedia of Integer Sequences) and write a small function for each.

  • Example logic: "If the sum of the digits squared equals 1 (eventually), mark as Happy Number."

4. Customization (Your Requirement)

You mentioned wanting to add your own output. This is the biggest advantage of a custom script. You could design your program to:

• Output data in a specific format (JSON, CSV, or a readable text report).

• Run specific checks that Numbers Aplenty doesn't have.

• Batch process hundreds of numbers at once (which would take forever on a web interface).

Summary of Feasibility

Feature	Feasibility in Python	Difficulty
Basic Properties (Odd, Even, Digits)	Very High	Trivial
Number Theory (Primes, Factors)	Very High	Easy (using SymPy)
Replicating 170+ Families	High	Moderate (requires time to write checks for each)
Speed	High	Fast for <20 digits; slower for massive numbers
Adding Custom Output	Very High	Easy

Conceptual Workflow

If you were to commission or build this, the program structure would look like this:

1. Input: User enters a number (e.g., 1729).

2. The "Checker" Engine: The number is passed through a list of functions:

   • check_prime(1729) -> False

   • check_taxicab(1729) -> True

   • check_harshad(1729) -> True

   • ...plus your custom checks.

3. Output: The program prints a report combining the standard results with your custom metrics.

Verdict: Highly feasible. The "hard" part is just the manual labor of coding the definitions for the hundreds of specific number types you might want to track.

Getting Gemini to code in Python from scratch is much preferable to its revising code that I've written. I can gradually build up a list of functions that I want the program to execute in a Jupyter notebook. It will produce output similar but more extensive than Numbers Aplenty, it will not require Internet access, it will be efficiently coded and the code can be easily modified whenever necessary.

Revisiting the Odd (+) and Even (-) Algorithm

I've written extensively about this algorithm in previous posts and even uploaded a PDF to Academia (link). For some reason I decided to reread this PDF and this motivated me to get Gemini Pro 3.0 to try its hand at writing some Python code to implement this algorithm across a chosen range of numbers. Of course, I'd already done this previously using SageMath but my algorithm timed out on SageMathCell above 100,000 and I thought that any code that Gemini created would be far more efficient than any code that I could write. I tested it out on SageMathCell for a range up to one million but it still timed out. However, on a range up to 100,000, it only took a few seconds. I tried to run the code in my Jupyter notebook using the range up to one million but it spat the dummy. No problem, I'm mainly interested in the range up to 40,000 given my focus on my diurnal age. Here is a permalink to the algorithm on SageMathCell,

Here was the prompt that I gave Gemini:

Implement the following program in Python. Here are the details:

\( \textbf{Odd Even Algorithm} \)

\( \textbf{The Basic Algorithm:}\)

Let’s describe the basic algorithm first. It takes 0 or any positive integer as input, computes the sum of the number’s odd digits and the sum of the number’s even digits. The sum of the number’s odd digits is added to the number while the sum of the number’s even digits is subtracted. This process is repeated until a stable number is reached, meaning the sums of odd and even digits are equal, OR a loop is entered. 

\( \textbf{Nomenclature:}\)

I’m choosing to call stable numbers (sums of odd and even digits are equal) ATTRACTORS because under the algorithm the trajectory of many numbers will lead to such an attractor. 0 is the first such attractor and 112 is the next.

Numbers that lead back to themselves I’m calling VORTICALS. An example is 11 because its trajectory is 11, 13, 17, 25, 28, 18, 11. Similarly 13 is a vortical because its trajectory is 13, 17, 25, 28, 18, 11, 13. All these numbers lead back to themselves and collectively I call this collection of verticals a VORTEX. It can be represented as [11, 13, 17, 25, 28, 18] but any vortical in it can be placed first and only the cyclic order needs to be preserved.

Numbers whose trajectories lead to a vortex or an attractor are called CAPTIVES. 9 is a captive of a vortex because its trajectory is 9, 18, 11, 13, 17, 25, 28, 18 leads it to the vortex [18, 11, 13, 17, 25, 28]. 

\( \textbf{Applying the algorithm to a range of numbers:}\)

I’m interested in selecting a range of numbers (let’s say from 0 to 100,000) and applying the algorithm to each number in this range. What I want to keep track of are:

• Display list of attractors and their total number

• Display list of the captives of each attractor and the number of these captives for each

• Display a ranking of the attractors in order of number of captives (highest to lowest)

• Display list of vortices (plural of vortex) together with the vorticals that comprise them

• Display the captives of each each vortex and how many captives each has

• Display of ranking of vortices in order of number of captives (highest to lowest)Display overall statistics: number of attractors, number of vorticals, number of captives of attractors, number of captives of vortices.

Gemini carried these instructions out perfectly as the implemented code revealed. Here is a link to its response. I ran the program using a restricted range up to 40,000 and copied the output to a Google document (link). I can now search this document to find details concerning an attractor or vortex. For example, consider these forthcoming attractors: 28019, 28037, 28055, 28073 and 28091. Here are the results concerning their number of captives:

• 28019 has three captives (27951, 27971, 27993)

• 28037 has no captives

• 28055 has no captives

• 28073 has no captives

• 28091 has 138 captives (28100, 28102, 28104, 28106, 28108, 28110, 28111, 28112, 28113, 28114, 28115, 28116, 28117, 28118, 28121, 28123, 28127, 28129, 28139, 28140 ... and 118 more)

This of course is a very useful tool and I'll be making use of it for my diurnal number investigations. This was the summary generated for numbers in the range up to 40000:

• Total Attractors Found: 1527

• Total Vortices Found: 428

• Total Unique Vorticals: 1451

• Total Captives of Attractors: 20417

• Total Captives of Vortices: 16610

Claude to the Rescue

The number 27983, associated with my diurnal age today, has a property than confers it membership in OEIS A335752: the lexicographically earliest sequence of primes, each two of which sum to a triprime. The sequence begins: 3, 5, 47, 239, 3023, 27983, 53231, ... so let's confirm that these numbers satisfy using 3 as a starting point:

• \(3 + 5 = 8 =2^3\)
• \(3 + 47 = 50 = 2 \times 5^2\)
• \(3 + 239 = 242 = 2 \times 11^2 \)
• \(3 + 3023 = 3026 = 2 \times 17 \times 89 \)
• \(3+27983 = 27986 = 2 \times 7 \times 1999\)
  
  

The sums above are all triprimes and we will find that any number in the list paired with any other number will produce a triprime. The code for this was included in the OEIS comments but was in Mathematica. I got Gemini to convert it to Python code and then Grok but in both cases the primeomega function was misunderstood with both claiming that primeomega(12) =2 and not 3. 

I didn't fully understand the Python code and was looking to Gemini and then Grok to explain it to me. It was only Claude that succeeded to my satisfaction. Here is the code followed by Claude's explanation. I think Claude will be my go to resource for anything to do with coding.

from sympy import nextprime, primeomega

ps = [3]

for _ in range(10):

    p = nextprime(ps[-1]) 

    while not all(primeomega(x + p) == 3 for x in ps):

        p = nextprime(p)

    ps.append(p)

print(ps)

I'd be happy to explain this code line by line! This code generates a sequence of prime numbers with a special property.

Code Breakdown

from sympy import nextprime, primeomega

Imports two functions from SymPy:

• nextprime(\(n\)): Returns the next prime number greater than \(n\)
• primeomega(\(n\)): Returns the total number of prime factors of \(n\) (counting multiplicity). For example, primeomega(12) = 3 because 12 = 2² × 3 (two 2's and one 3)

ps = [3]

Initializes the sequence with the first prime number 3. This list will store our special sequence of primes.

for _ in range(10):

Loops 10 times to find 10 more primes (so we'll have 11 total including the initial 3). The underscore _ is used when we don't need the loop variable.

p = nextprime(ps[-1])

Gets the next prime after the last element in our list. 

• ps[-1] accesses the last element of the list
• Initially, this gives nextprime(3) = 5
• This is our candidate prime to add to the sequence

while not all(primeomega(x + p) == 3 for x in ps):

    p = nextprime(p)

This is the key part! It checks if the candidate prime \(p\) satisfies a special condition:

• For every number \(x\) already in our list ps
• The sum \(x + p\) must have exactly 3 prime factors (counting multiplicity)
• If this condition fails for any \(x\), we move to the next prime candidate
• This continues until we find a prime where all sums have exactly 3 prime factors

ps.append(p)

Adds the qualified prime to our sequence once it passes the test.

print(ps)

Prints the final sequence** of primes.

Example Walkthrough

Let's trace the first iteration:

1. Start with ps = [3]

2. Try p = 5: Check if 3 + 5 = 8 = 2³ has 3 prime factors ✓ (yes, three 2's)

3. Add 5 to list: ps = [3, 5]

4. Try p = 7: Check 3 + 7 = 10 = 2 × 5 (only 2 factors) ✗

5. Try p = 11: Check 3 + 11 = 14 = 2 × 7 (only 2 factors) ✗

6. Continue until finding a prime that works with both 3 and 5...

This creates a sequence where each new prime maintains the property that when added to any previous prime in the sequence, the result always has exactly 3 prime factors!

*************************

Now if we change the condition to primeomega = 2, we get the following sequence: 3, 7, 19, 139, 859, 8179, 173059, ... where it can be seen there is a big jump from 8179 to 173059. Testing this out using 3 again we see that:

• \(3+7=10=2 \times 5\)
• \(3+19=22=2 \times 11\)
• \(3+139=142=2 \times 71\)
• \(3+859=862=2 \times 431\)
• \(3+8179=8182=2 \times 4091\)
• \(3+173059=173062=2 \times 86531\)

The sequence is included in the OEIS as A113875: slowest growing sequence of primes having the prime-pairwise-average property: if i<j, (a(i)+a(j))/2 is prime. This makes sense in that 2 seems to be the common factor is all of the sums. 

The sequence with additional (very large) terms looks like this:

3, 7, 19, 139, 859, 8179, 173059, 1026199, 1827139, 15828679, 13187242759, 18732483199, 912492556939, 9130567625119

When we set primeomega = 4 we get the following sequence: 3, 13, 137, 1193, 4373, 26017, 33889, 495433, ... which is NOT listed in the OEIS. Again let's use 3 as our test number:

• \(3+13=16=2^4\)
• \(3+137=140=2^2 \times 5 \times 7\)
• \(3+1193=1196=2^2 \times 13 \times 23\)
• \(3+4373=4376=2^3 \times 547\)
• \(3+26017=26020=2^2 \times 5 \times 1301\)
• \(3+33889=33892=2^2 \times 37 \times 229\)
• \(3+495433=495436=2^2 \times 29 \times 4271\)

That's probably a good place to stop.

Coding and Integers

The thought struck me to write a program that looked at all base 10 numbers in the range from 10 to 40000 that, when converted to a number base between 2 and 9, returned a number with the same digits as the base 10 number but in a different order. It's a simple enough program to write in SageMath and this is what I created originally:

T=[("number","permutation","base")]

T.append(("","",""))

for n in [10..40000]:

    for b in [2..9]:

        number=int(n.str(base=b))

        if sorted(n.digits())==sorted(number.digits()):

            L.append(n)

            T.append((n,n.str(base=b),b))

print(L)

print(len(L))

print()

print(table(T))

The program will not run because there is a problem with "sorted(number.digits())". I get the error message "AttributeError: 'int' object has no attribute 'digits' ". I know from experience how to fix the problem. I simply replace "int" with "Integer" but I've never fully understood why this change is necessary. I queried Gemini and it came up with the following explanation:

Of course. Your code requires SageMath's Integer type because it uses special mathematical methods that Python's standard int type doesn't have.

The specific line of code that causes the issue is if sorted(n.digits())==sorted(number.digits()):.

---

SageMath's Integer vs. Python's int

Think of Python's built-in int as a basic calculator and SageMath's Integer as a powerful scientific calculator. While both can store whole numbers, the SageMath version comes with a lot of extra functions built specifically for number theory and advanced mathematics.

Python int: A standard data type for storing whole numbers. It's great for general programming, loops, and basic arithmetic. However, it does not have a built-in method called .digits().
SageMath Integer: A more specialized object. It's designed for mathematical computation and has many useful methods, including:
.digits(): Returns a list of the number's digits.
.str(base=b): Returns the string representation of the number in a given base.
.is_prime(): Checks if the number is prime.
.divisors(): Returns a list of the number's divisors.

---

Why Your Code Fails with int

In your loop, the variable n is a SageMath Integer because it's created by SageMath's range notation ([10..40000]).
When you write number = int(n.str(base=b)), you successfully convert the base-b string representation of n into a standard Python int and store it in the number variable.
The program then fails at the next line: sorted(number.digits()). Since number is a plain Python int, it doesn't have the .digits() method, and SageMath throws an AttributeError.

Why Your Code Works with Integer

When you use number = Integer(n.str(base=b)), you are using SageMath's Integer constructor.
This creates a SageMath Integer object and stores it in the number variable.
Now, when the program reaches sorted(number.digits()), the call is valid because SageMath Integer objects do have the .digits() method. 👍

In short, you must use Integer to ensure that both n and number are SageMath objects that understand the special mathematical methods your algorithm relies on.

So now I know. So what is the output from my program when it is run? It turns out that there are 76 numbers that satisfy and they are (permalink):

13, 23, 46, 158, 227, 265, 316, 445, 1030, 1045, 1135, 1234, 1236, 1273, 1366, 1380, 1431, 1454, 1653, 2027, 2060, 2116, 2154, 2315, 2534, 3160, 3161, 3162, 3163, 3164, 3165, 3166, 3167, 3226, 5270, 5567, 5637, 5783, 10144, 10235, 10342, 10453, 10542, 11425, 11750, 12415, 12450, 12564, 12651, 13045, 13245, 13260, 13402, 13534, 13620, 14610, 15226, 15643, 16255, 16273, 16546, 16633, 21322, 21753, 21763, 21835, 23568, 26804, 30576, 31457, 32348, 34582, 35001, 35081, 35228, 37465

These numbers are the initial members of OEIS A090144 (I've discounted the trivial numbers from 1 to 8). Here are the details:

  number   permutation   base

  13       31            4

  23       32            7

  46       64            7

  158      185           9

  227      272           9

  265      526           7

  316      631           7

  445      544           9

  1030     3001          7

  1045     4501          6

  1135     5131          6

  1234     3412          7

  1236     1623          9

  1273     2371          8

  1366     3661          7

  1380     1803          9

  1431     4113          7

  1454     4145          7

  1653     3165          8

  2027     2702          9

  2060     6002          7

  2116     6112          7

  2154     4152          8

  2315     3152          9

  2534     3425          9

  3160     6130          8

  3161     6131          8

  3162     6132          8

  3163     6133          8

  3164     6134          8

  3165     6135          8

  3166     6136          8

  3167     6137          8

  3226     6232          8

  5270     7205          9

  5567     7565          9

  5637     7653          9

  5783     7835          9

  10144    41401         7

  10235    15032         9

  10342    42103         7

  10453    15304         9

  10542    42510         7

  11425    45211         7

  11750    17105         9

  12415    51124         7

  12450    51204         7

  12564    51426         7

  12651    51612         7

  13045    53014         7

  13245    53421         7

  13260    20163         9

  13402    20341         9

  13534    54313         7

  13620    20613         9

  14610    60411         7

  15226    62251         7

  15643    63415         7

  16255    65251         7

  16273    37621         8

  16546    66145         7

  16633    66331         7

  21322    32221         9

  21753    52371         8

  21763    32761         9

  21835    32851         9

  23568    35286         9

  26804    40682         9

  30576    73560         8

  31457    75341         8

  32348    48332         9

  34582    52384         9

  35001    53010         9

  35081    53108         9

  35228    53282         9

  37465    56347         9