difference_of_two_squares_solutions.pl

#!/usr/bin/perl

# Daniel "Trizen" Șuteu
# Date: 09 August 2017
# Edit: 26 October 2017
# https://github.com/trizen

# A simple and efficient algorithm for finding all the non-negative integer solutions to the equation:
#
#   x^2 - y^2 = n
#
# where `n` is known (along with its prime factorization).

# Blog post:
#   https://trizenx.blogspot.com/2017/10/representing-integers-as-difference-of.html

use 5.010;
use strict;
use warnings;

use ntheory qw(divisors);

sub difference_of_two_squares_solutions {
    my ($n) = @_;

    my @solutions;
    foreach my $divisor (divisors($n)) {

        last if $divisor > sqrt($n);

        my $p = $divisor;
        my $q = $n / $divisor;

        ($p + $q) % 2 == 0 or next;

        my $x = ($q + $p) >> 1;
        my $y = ($q - $p) >> 1;

        unshift @solutions, [$x, $y];
    }

    return @solutions;
}

foreach my $n (1 .. 1e2) {
    (my @solutions = difference_of_two_squares_solutions($n)) || next;

    say "$n = ", join(' = ', map { "$_->[0]^2 - $_->[1]^2" } @solutions);

    # Verify solutions
    foreach my $solution(@solutions) {
        if ($solution->[0]**2 - $solution->[1]**2 != $n) {
            die "Error for $n: (@$solution)";
        }
    }
}

__END__
99937 = 721^2 - 648^2 = 1369^2 - 1332^2 = 49969^2 - 49968^2
99939 = 2390^2 - 2369^2 = 7142^2 - 7135^2 = 16658^2 - 16655^2 = 49970^2 - 49969^2
99940 = 358^2 - 168^2 = 1334^2 - 1296^2 = 5002^2 - 4992^2 = 24986^2 - 24984^2
99941 = 429^2 - 290^2 = 49971^2 - 49970^2
99943 = 2948^2 - 2931^2 = 49972^2 - 49971^2
99944 = 465^2 - 341^2 = 837^2 - 775^2 = 987^2 - 935^2 = 1935^2 - 1909^2 = 12495^2 - 12491^2 = 24987^2 - 24985^2
99945 = 1133^2 - 1088^2 = 3339^2 - 3324^2 = 5557^2 - 5548^2 = 9997^2 - 9992^2 = 16659^2 - 16656^2 = 49973^2 - 49972^2
99947 = 606^2 - 517^2 = 49974^2 - 49973^2
99948 = 8332^2 - 8326^2 = 24988^2 - 24986^2
99949 = 457^2 - 330^2 = 49975^2 - 49974^2
99951 = 16660^2 - 16657^2 = 49976^2 - 49975^2
99952 = 6251^2 - 6243^2 = 12496^2 - 12492^2 = 24989^2 - 24987^2
99953 = 447^2 - 316^2 = 513^2 - 404^2 = 7143^2 - 7136^2 = 49977^2 - 49976^2
99955 = 9998^2 - 9993^2 = 49978^2 - 49977^2
99956 = 24990^2 - 24988^2
99957 = 331^2 - 98^2 = 421^2 - 278^2 = 1301^2 - 1262^2 = 1531^2 - 1498^2 = 3851^2 - 3838^2 = 4549^2 - 4538^2 = 16661^2 - 16658^2 = 49979^2 - 49978^2
99959 = 2640^2 - 2621^2 = 49980^2 - 49979^2
99960 = 317^2 - 23^2 = 329^2 - 91^2 = 343^2 - 133^2 = 347^2 - 143^2 = 353^2 - 157^2 = 379^2 - 209^2 = 427^2 - 287^2 = 541^2 - 439^2 = 559^2 - 461^2 = 637^2 - 553^2 = 749^2 - 679^2 = 769^2 - 701^2 = 863^2 - 803^2 = 1211^2 - 1169^2 = 1487^2 - 1453^2 = 1681^2 - 1651^2 = 1799^2 - 1771^2 = 2509^2 - 2489^2 = 3577^2 - 3563^2 = 4171^2 - 4159^2 = 5003^2 - 4993^2 = 8333^2 - 8327^2 = 12497^2 - 12493^2 = 24991^2 - 24989^2
99961 = 49981^2 - 49980^2
99963 = 322^2 - 61^2 = 618^2 - 531^2 = 1738^2 - 1709^2 = 5558^2 - 5549^2 = 16662^2 - 16659^2 = 49982^2 - 49981^2
99964 = 440^2 - 306^2 = 24992^2 - 24990^2
99965 = 9999^2 - 9994^2 = 49983^2 - 49982^2
99967 = 7144^2 - 7137^2 = 49984^2 - 49983^2
99968 = 318^2 - 34^2 = 372^2 - 196^2 = 423^2 - 281^2 = 612^2 - 524^2 = 813^2 - 749^2 = 1158^2 - 1114^2 = 1578^2 - 1546^2 = 2283^2 - 2261^2 = 3132^2 - 3116^2 = 6252^2 - 6244^2 = 12498^2 - 12494^2 = 24993^2 - 24991^2
99969 = 425^2 - 284^2 = 1087^2 - 1040^2 = 16663^2 - 16660^2 = 49985^2 - 49984^2
99971 = 49986^2 - 49985^2
99972 = 2786^2 - 2768^2 = 8334^2 - 8328^2 = 24994^2 - 24992^2
99973 = 323^2 - 66^2 = 49987^2 - 49986^2
99975 = 340^2 - 125^2 = 400^2 - 245^2 = 452^2 - 323^2 = 584^2 - 491^2 = 704^2 - 629^2 = 1184^2 - 1141^2 = 1628^2 - 1597^2 = 2012^2 - 1987^2 = 3340^2 - 3325^2 = 10000^2 - 9995^2 = 16664^2 - 16661^2 = 49988^2 - 49987^2
99976 = 12499^2 - 12495^2 = 24995^2 - 24993^2
99977 = 2949^2 - 2932^2 = 49989^2 - 49988^2
99979 = 410^2 - 261^2 = 850^2 - 789^2 = 4550^2 - 4539^2 = 49990^2 - 49989^2
99980 = 5004^2 - 4994^2 = 24996^2 - 24994^2
99981 = 345^2 - 138^2 = 359^2 - 170^2 = 391^2 - 230^2 = 759^2 - 690^2 = 825^2 - 762^2 = 1865^2 - 1838^2 = 2185^2 - 2162^2 = 2391^2 - 2370^2 = 5559^2 - 5550^2 = 7145^2 - 7138^2 = 16665^2 - 16662^2 = 49991^2 - 49990^2
99983 = 3852^2 - 3839^2 = 49992^2 - 49991^2
99984 = 2095^2 - 2071^2 = 4172^2 - 4160^2 = 6253^2 - 6245^2 = 8335^2 - 8329^2 = 12500^2 - 12496^2 = 24997^2 - 24995^2
99985 = 10001^2 - 9996^2 = 49993^2 - 49992^2
99987 = 16666^2 - 16663^2 = 49994^2 - 49993^2
99988 = 3578^2 - 3564^2 = 24998^2 - 24996^2
99989 = 49995^2 - 49994^2
99991 = 49996^2 - 49995^2
99992 = 489^2 - 373^2 = 891^2 - 833^2 = 12501^2 - 12497^2 = 24999^2 - 24997^2
99993 = 16667^2 - 16664^2 = 49997^2 - 49996^2
99995 = 1446^2 - 1411^2 = 7146^2 - 7139^2 = 10002^2 - 9997^2 = 49998^2 - 49997^2
99996 = 680^2 - 602^2 = 1936^2 - 1910^2 = 8336^2 - 8330^2 = 25000^2 - 24998^2
99997 = 319^2 - 42^2 = 2641^2 - 2622^2 = 49999^2 - 49998^2
99999 = 320^2 - 49^2 = 468^2 - 345^2 = 1240^2 - 1199^2 = 5560^2 - 5551^2 = 16668^2 - 16665^2 = 50000^2 - 49999^2
100000 = 325^2 - 75^2 = 350^2 - 150^2 = 550^2 - 450^2 = 665^2 - 585^2 = 1025^2 - 975^2 = 1270^2 - 1230^2 = 2510^2 - 2490^2 = 3133^2 - 3117^2 = 5005^2 - 4995^2 = 6254^2 - 6246^2 = 12502^2 - 12498^2 = 25001^2 - 24999^2