JCSBimp's Hexagonal Orthogonal

Nothing earthshaking or historic is here. This is just my latest and greatest attempt to create a Perl script that improves on Fermat factorization, with the simplest and shortest code possible.

As always, if anyone wants to try this algorithm out in a language, or an extension of Perl, that can use arbitrary-precision integers, it would be interesting to see how it behaves with products of very large primes.

The code follows:

#!/usr/bin/perl
use strict;
use warnings;
use POSIX qw/ceil/;

# Name: quickfact.pl
# Author: J. Calvin Smith
#
# Notice: Copyright 2023 J. Calvin Smith <jcsbimp@me.com>
#         <https://jcsbimp.com> Ceiling-Square-based
#         integer factorization script and algorithm.
#
#         Licensed for use under terms of GNU-GPL v3,
#         publicly available at the web link:
#         https://opensource.org/license/gpl-3-0/, as
#         follows:
#
# This program is free software: you can redistribute it
# and/or modify it under the terms of the GNU General Public
# License as published by the Free Software Foundation, either
# version 3 of the License, or (at your option) any later version.
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
# General Public License (linked above) for more details.
#
# --------
#
# Revision History:
#   28 September 2023 - Initial Release; put on jcsbimp.com blog
#   26 October 2023 - Improvement using intersection points
# 
# --------
#
# Usage: quickfact.pl <pos_int>
#
#    pos_int - the integer the user wishes to factor
#
# This Perl script finds a positive integer factor of pos_int,
# finding a non-trivial factor, if possible.
# It uses a method derived from an analytic geometry construction
# based on boundary lines of the integer grid points of f(x,y) =
# x(x+1)+y, using straight lines on the grid thus numbered between
# consecutive perfect squares to create "zone boundaries." It bases
# calculations on the first two elements of the series of remainders
# obtained by finding the Ceiling Square c^2 of pos_int, the least
# perfect square greater than pos_int that is even or odd as
# ( pos_int + 1 ) / 2 is even or odd, and computing c^2 - pos_int
# and (c+1)^2 - pos_int. Calculations of a "remainder line" that
# intersects the zone boundaries lead the algorithm to a solution
# when the points of intersection have integer x values.
#
# NOTES FOR USE: Integers here are of the standard size in Perl.
# No use of arbitrary-precision variables or arithmetic appears,
# but a developer wishing to make a future such version to test the
# algorithm on arbitrarily large integers has both the author's permis-
# sion and his blessing, this being free and proof-of-concept  software.

# Test for argument vector; give usage line if empty.

my $num_args = $#ARGV + 1;

if ( $num_args == 0 ) {
  print "Usage: quickfact.pl <pos_int>\n",
        "  pos_int = positive integer to factor\n";
  exit( 0 );
}
my $pos_int = $ARGV[0];
if ( $pos_int < 1 ) {
  print "ERROR: $pos_int is less than 1. Quitting.\n";
  exit( 0 );
}
my $pos_ceil = ceil( $pos_int );
if ( $pos_int != $pos_ceil ) {
  print "ERROR: $pos_int is not an integer. Quitting.\n";
  exit( 0 );
}

# Declare variables.

open my $out_fh, '>', "outfile_$pos_int.csv" or die $!;

my $c; # the (Adjusted) Ceiling Root of m
my $r0; # c^2-m, the first remainder
my $csr0; # Ceiling of sqrt of the remainder

my ( $meo, $ceo ); # Even/odd tests to set Ceiling Root
my ( $numer, $numinc, $denom, $the_mod ); # For the test fraction
my $the_gcd; # for the secondary test
my $num_iters = 0; # The number of iterations of the main loop
my $not_done_yet = 1; # Loop control variable
my ( $s, $t ); # the greater and lesser squares
my $succstr = 'SUCCESS';
my $primstr = 'PRIME';

# Determine if pos_int is even or odd.

if ( ( $pos_int % 2 ) == 0 ) {
  print "SUCCESS: $pos_int = 2 * ", $pos_int/2, ".\n";
  exit( 0 );
}

# pos_int is an odd positive integer.

# Calculate Ceiling Root of pos_int.

$c = ceil( sqrt( $pos_int ) );
$ceo = $c % 2;
$meo = ( ( $pos_int + 1 ) / 2 ) % 2;
if ( $ceo != $meo ) { # if unequal, make equal
  $c++;               #   by adjusting Ceiling Root
}

# Calculate the first term of the Remainder Series of pos_int.

$r0 = $c * $c - $pos_int;

# Check this remainder for easy Fermat factorization.

$csr0 = ceil( sqrt( $r0 ) );
if ( ( $csr0 * $csr0 ) == $r0 ) { # Easy answer: r0 is a perfect square.
  print "SUCCESS: $pos_int = ", ( $c - $csr0 ), " * ", ( $c + $csr0 );
  print ",\n  found by Fermat factorization on Remainder 0.\n";
  exit( 0 );
}

# Now ready to execute the main loop. Set up initial values.

$numer = 0 - $r0;
$denom = 2 * $c;
if ( ( $numer % 2 ) == 1 ) {
  $numer++;
  $denom -= 2;
  $numinc = 8;
}
else {
  $numinc = 4;
}

# The Main Loop

while ( $not_done_yet ) {
  $the_mod = $numer % $denom;
  if ( $the_mod == 0 ) { # Success!
    $not_done_yet = 0;
    $s = $c + $numer / $denom;
    $t = sqrt( $s * $s - $pos_int );
    if ( ( $s - $t ) == 1 ) {
      print $primstr;
    }
    else {
      print $succstr;
    }
    print ": $pos_int = ";
    print ( $s - $t );
    print " * ";
    print ( $s + $t );
    print ".\nFound on iteration $num_iters.\n";
  }
  else {
    # Now check for a gcd match.
    $the_gcd = gcd( $numer, $denom );
    if ( $the_gcd > 2 ) {
      $the_gcd = gcd( $the_gcd, $pos_int );
      if ( $the_gcd > 2 ) {
        print "SUCCESS: $pos_int = $the_gcd * ";
        print ( $pos_int / $the_gcd );
        print ".\nFound on iteration $num_iters.\n";
        $not_done_yet = 0;
      }
    }
  }
  $numer += $numinc;
  $numinc += 8;
  $denom -= 4;
  if ( $denom <= 0 ) {
    if ( $not_done_yet ) {
      print "No factorization found. Sorry!\n";
      print "  Especially since this isn't supposed to happen!\n";
    }
  }  
  $num_iters++;
}

sub gcd {
  my ($u, $v) = @_;
  if ($v) {
    return gcd($v, $u % $v);
  } else {
    return abs($u);
  }
}