Algorithm for Generating Odd Primes (OEIS A065091) Without Primality Testing
rui.ferreira
Hello SeqFan members,
I would like to share a recent discovery of mine that generates the sequence of odd primes (OEIS A065091) without the need for any explicit primality tests. If this method has already been explored, I would greatly appreciate any references to prior work on this topic.
The algorithm I’ve developed combines a simple algebraic construction with a maximality filter. When executed up to a bound n, it produces all odd primes up to approximately 2n.
Algorithm Overview:
For each pair (a, b) where 0 ≤ a < b < n, the following two quantities are defined:
x = |a^2 - b^2|, y = x * |b - a|
For each value of x, we compute the maximal value of y. Remarkably, when the maximal value of y equals x, the corresponding values of x correspond to the odd primes.
This construction leads to an efficient method for generating odd primes without needing to check divisibility explicitly.
Attachments:
I’ve included two Python code implementations of this algorithm for testing. I hope you find this approach interesting and look forward to your feedback.
Best regards,
Rui Ferreira
1st code:
import numpy as np
from numba import njit
# -----------------------------------------------------------------------------
#
# This code explores an interesting mathematical conjecture involving pairs of integers.
# For two integers `a` and `b` where `0 <= a < b < n`, we define:
#
# x = b^2 - a^2 (the absolute difference of their squares)
# y = x * (b - a) (x multiplied by the absolute difference of `b` and `a`)
#
# The goal is to find the maximum value of `y` for each possible `x` by iterating over
# all integer pairs `(a, b)` such that `0 <= a < b < n`.
#
# The conjecture that underpins this code can be summarized as:
#
# **Main Conjecture (Theorem)**:
# - If `x` is an odd integer, then `max_y(x) = x` if and only if `x` is an odd prime number.
#
# This means that when we compute `max_y(x)` for any odd integer `x` and find that `max_y(x)` equals `x`,
# we can conclude that `x` is an odd prime number.
#
# **OEIS A065091**:
# - The sequence of numbers where `max_y(x) = x` matches exactly the sequence of odd prime numbers,
# as listed in the Online Encyclopedia of Integer Sequences (OEIS) as sequence A065091.
#
# By exploring this extremal principle, we can generate odd primes without directly testing for primality.
# This approach doesn't use conventional prime-checking algorithms, sieves, or divisibility tests.
#
# -----------------------------------------------------------------------------
def compute_max_y(n: int, max_x: int) -> np.ndarray:
"""
Compute the maximal y value for each x obtained from:
x = |a^2 - b^2|
y = x * |b - a|
for all integer pairs (a, b) with 0 <= a < b < n.
Args:
n (int): The upper bound for a and b.
max_x (int): The maximum value of x to consider.
Returns:
np.ndarray: An array where the index corresponds to x and the value at that index
is the maximal y for that x. If x doesn't appear in any pair (a, b),
the value at index x will be 0.
"""
# Initialize an array to store the maximal y value for each x
max_y_per_x = np.zeros(max_x + 1, dtype=np.int64)
# Iterate over all pairs of integers (a, b) where 0 <= a < b < n
for a in range(n - 1):
for b in range(a + 1, n):
# Calculate x and y as per the given formulas
x = abs(a**2 - b**2) # x = |a^2 - b^2|
y = x * abs(b - a) # y = x * |b - a|
# Update the maximal y for the current x if a larger y is found
if y > max_y_per_x[x]:
max_y_per_x[x] = y
return max_y_per_x
def run_without_prime_tests(n: int) -> None:
"""
Prints all x values where max_y(x) == x. This sequence empirically matches
the sequence of odd primes, as described in OEIS A065091, without explicitly
performing prime checks.
Args:
n (int): The upper bound for generating integer pairs (a, b) to compute the sequence.
"""
# The maximum possible value for x is 2 * (n - 1) * (n - 1)
max_possible_x = 2 * (n - 1) * (n - 1)
# Compute the maximal y values for each possible x up to max_possible_x
max_y_per_x = compute_max_y(n, max_possible_x)
results = []
idx = 0
# Print a header for the result report with column descriptions
print("\nOEIS A065091 Odd primes (https://oeis.org/A065091)")
print()
print("[Index ] = [Prime Number]")
# Iterate through all possible values of x from 2 to max_possible_x
for x in range(2, max_possible_x + 1):
y = max_y_per_x[x]
# Skip values of x where no valid y (maximal value) is found
if y == 0:
continue
# Only consider cases where the maximal y equals x (i.e., x == y)
if y != x:
continue
# Increment the index and store the result
idx += 1
results.append((idx, x, y))
# Print the prime index and corresponding prime number in the specified format
print(f"[{idx:6d}] = {x:6d}")
def main():
"""
Main function to prompt the user for input and initiate the sequence computation.
It handles user input for the value of n and calls the run_without_prime_tests
function to display the sequence.
"""
while True:
try:
# Prompt the user to input an integer n
n_input = input("Enter n (integer > 2, e.g., 1000): ").strip()
n = int(n_input)
# Ensure that n is greater than 2
if n <= 2:
raise ValueError
break
except ValueError:
print("Invalid input. Please enter an integer greater than 2.")
# Call the function to compute and print the results
run_without_prime_tests(n)
if __name__ == "__main__":
main()
import numpy as np
# Function to generate maximal values based on differences of squares
def generate_extremal_values(n):
"""
Generates a sequence of maximal y values corresponding to x values obtained
from the difference of squares formula for all integer pairs (a, b) with
0 <= a < b < n.
For each pair (a, b), x is calculated as |a^2 - b^2|, and y is computed
as y = x * |b - a|. The function returns the list of maximal y values for each unique x.
Args:
n (int): The upper bound for a and b.
Returns:
list: A list of maximal y values for each unique x.
"""
max_y = {} # Dictionary to store the maximal y value for each x
# Iterate over all possible pairs of integers a and b
for a in range(n):
for b in range(a, n + a):
x = abs(a**2 - b**2) # Calculate x = |a^2 - b^2|
y = x * abs(a - b) # Calculate y = x * |b - a|
# Update the dictionary with the maximal y for each x
if x not in max_y or y > max_y[x]:
max_y[x] = y
# Return the list of maximal y values
return list(max_y.values())
# Function to identify segments where the sequence decreases
def identify_segments(sequence):
"""
Identifies segments of the sequence where the value decreases. A segment
is defined as a continuous subsequence where each number is not smaller
than the previous number.
Args:
sequence (list): The list of values to analyze.
Returns:
list: A list of segments, where each segment is a list of consecutive values
in the sequence that do not decrease.
"""
segments = [] # List to store the identified segments
current_segment = [sequence[0]] # Start with the first value as the initial segment
# Iterate over the sequence to group consecutive values into segments
for i in range(1, len(sequence)):
if sequence[i] < sequence[i - 1]: # If the sequence decreases
segments.append(current_segment) # End the current segment and start a new one
current_segment = [sequence[i]] # Start a new segment with the current value
else:
current_segment.append(sequence[i]) # Add the current value to the current segment
segments.append(current_segment) # Append the last segment
return segments
# Example usage
n = 100 # Define the upper limit for generating integer pairs (a, b)
sequence = generate_extremal_values(n) # Generate the sequence of maximal y values
segments = identify_segments(sequence) # Identify the decreasing segments in the sequence
# Initialize index counter for reporting primes
index = 1
# Print a header for the result report with column descriptions
print("\nOEIS A065091 Odd primes (https://oeis.org/A065091)") # OEIS reference for the sequence of odd primes
print()
print("[Index ] = [Prime Number]") # Descriptive header for the columns
# Iterate through each segment to report the prime index and the first prime in each segment
for k, segment in enumerate(segments, 1):
first_prime = segment[0] # Take the first number of the segment as a candidate prime
if first_prime == 2 * k - 1: # Check if the first prime corresponds to the expected odd prime
# Print the prime index and the corresponding prime number in the specified format
print(f"[{index:6d}] = {first_prime:6d}")
index += 1 # Increment the prime index after reporting the prime
Gareth McCaughan
Hello SeqFan members,
I would like to share a recent discovery of mine that generates the sequence of odd primes (OEIS A065091) without the need for any explicit primality tests. If this method has already been explored, I would greatly appreciate any references to prior work on this topic.
The algorithm I’ve developed combines a simple algebraic construction with a maximality filter. When executed up to a bound n, it produces all odd primes up to approximately 2n.
Algorithm Overview:
For each pair (a, b) where 0 ≤ a < b < n, the following two quantities are defined:
x = |a^2 - b^2|, y = x * |b - a|
For each value of x, we compute the maximal value of y. Remarkably, when the maximal value of y equals x, the corresponding values of x correspond to the odd primes.
The comments in your code call this both a "conjecture" and a "theorem". In case you don't already have a proof, here is one.
I'll adopt a convention opposite to yours and take a>b rather
than b>a. In this case, we can write a^2-b^2 and a-b instead of
|a^2-b^2| and |a-b|.
Obviously all the values of y are >= x, so the maximal y equals x <=> the _only_ possible value of y is x <=> the only possible value of a-b is 1.
Well, first of all, if a-b=1 is possible then x must be odd since then a+b is also odd and so x = (a-b)(a+b) is odd. So we only need to consider odd values of x.
Choices of a,b such that x=(a-b)(a+b) correspond exactly to pairs of odd numbers whose product is x, with a-b being the smaller and a+b being the larger. So "the only possible value of a-b = 1" is equivalent to "if x is a product of two odd numbers then the smaller is 1". Which, for odd x, is exactly the same thing as x being prime.
(Well, except that x=1 is possible, corresponding to a=1 and b=0,
and 1 is not usually considered to be a prime number.)
This construction leads to an efficient method for generating odd primes without needing to check divisibility explicitly.
It seems to me that it doesn't so much not need to check divisibility explicitly, as _hide_ its divisibility checks.
If instead of considering pairs (a,b) you considered the pairs (a-b,a+b) then your algorithm would look like this:
For each pair (m,n) where m<n and m,n have the same parity, define the quantities x = mn and y=mx. Then x is an odd prime iff, considering all such m,n with mn=x, the largest value of mx equals x.
Or, equivalently, like this:
For each pair (m,n) where m<n and m,n have the same parity, define the quantities x = mn and y=mx. Then x is an odd prime iff, considering all such m,n with mn=x, the largest value of m is 1.
That looks like divisibility-checking to me!
--
g
Rui Ferreira
Hello SeqFan members,
For anyone interested, I’ve created a GitHub repository at https://github.com/rf-iasys/OEIS where I’ll place my number theory and OEIS-related projects. This will help avoid sending frequent code updates to the forum.
Best regards,
Rui Ferreira
Charles Greathouse
--
You received this message because you are subscribed to the Google Groups "SeqFan" group.
To unsubscribe from this group and stop receiving emails from it, send an email to seqfan+un...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/seqfan/50dced9c-b7b8-4ba2-aff2-cdb89f9541c1n%40googlegroups.com.
Robert Israel
To view this discussion visit https://groups.google.com/d/msgid/seqfan/CANXmBjznfN8t5QdO18pk01cD4Qz5Wiz8%3Dr9C4pJd%3DnSX8n7qrQ%40mail.gmail.com.