Very interesting question, thanks!
Spent whole day programming C++ code from scratch that implements very fast Elliptic Curve ECM Factorization Method and Amicable Numbers search.
It is known that ECM method is very very fast. This ECM method, if optimized well, can easily within one second find factors even of 100-bit number (consisting of two 50-bit primes). It is slower only than Quadratic Sieve and GNFS.
As you said you're searching for 25-digit Amicable Numbers, this is 80-bit numbers.
As I measured in my below C++ program, on my old 2-core laptop, it can factor 150 random numbers 80-bit in size (25 digits) per second out of easy factorable set (I set a timer to cutoff difficult to factor numbers). And it can factor 20_000 numbers if their size is 27-30 bits (9 digits).
I said above easy factorable set. It appears that 60-70% of all random numbers 80-bit in size (25 digit) are easy factorable, meaning that it takes just 7ms to factor single such number which results in 150 numbers per second.
I'll search Amicable Pairs only among such easy-to-factor numbers. There are as many amicable numbers among easy-to-factor set as there are in all the rest. Amicable Numbers don't suffer from such easy-to-factor filtering.
Below is console output on how my program searches for 27-bit amicable numbers using ECM factorization:
.46501 _43919 *41268 , 20512 nums/sec, tdiv 99.34% 7mcs, ecm 25.3% 39mcs, total 37mcs, done 7%
found (97041735, 97945785)
.321639 _303557 *285333 , 20389 nums/sec, tdiv 99.35% 7mcs, ecm 25.4% 40mcs, total 38mcs, done 48%
found (34765731, 36939357)
.29933 _28247 *26517 , 20265 nums/sec, tdiv 99.35% 7mcs, ecm 25.4% 40mcs, total 38mcs, done 4%
found (9773505, 11791935)
.529659 _499736 *470070 , 19312 nums/sec, tdiv 99.35% 7mcs, ecm 25.5% 41mcs, total 40mcs, done 79%
found (5357625, 5684679)
.25937 _24411 *22905 , 19340 nums/sec, tdiv 99.35% 7mcs, ecm 25.4% 41mcs, total 40mcs, done 4%
found (998104, 1043096)
.26367 _24902 *23462 , 19360 nums/sec, tdiv 99.35% 7mcs, ecm 25.5% 41mcs, total 40mcs, done 4%
found (1184, 1210)
.130374 _122896 *115508 , 19648 nums/sec, tdiv 99.35% 7mcs, ecm 25.5% 40mcs, total 39mcs, done 20%
found (35390008, 39259592)
found (133178325, 133471275)
.32032 _30226 *28368 , 19753 nums/sec, tdiv 99.35% 7mcs, ecm 25.5% 40mcs, total 39mcs, done 5%
found (898216, 980984)
You can see above that program outputs found and amicable pair as soon as it finds something. Also it show some detailed statistics about speed and timings. It takes several seconds to find one amicable number of 27-bit size if you search 20_000 numbers per second as I do above.
If you think above 80-bit number (25 digits) then although I factor 150 numbers per second, still there SO FEW amicable numbers that it will take trillions of seconds to find at least one Amicable Number of such 80-bit size. Below is statistics of this 80-bit factoring. You can see that within several seconds it has done only 0.0000000000064% percentage of finding just single amicable pair, very tiny amount of whole work.
.4218 _2261 *1699 , 123 nums/sec, tdiv 43.81% 133mcs,
ecm 2.8% 1886mcs, total 8144mcs, done 0.0000000000064%
You may have heard about project Amicable Numbers for BOINC - it is great distributed search of such numbers. Millions of computer search for them. Of course if there are millions of computers then they can find some numbers much faster. In fact they have found already 1_227_817_736 such numbers, really great amount, see all numbers here.
In my C++ program I also present my ECM factorizer, as an example at program start I do a very fast factoring of 100-bit number, random one. It shows console output similar to following:
Initial N to factor with ECM: 1044877821881837432173434582227 (100-bit)
Number of threads 2.
Factors TrialDivA: [2953]
Factoring 35383603856479425437736059
Curves [ 0, 8), bound 2^ 9.000, 0.110 sec
Curves [ 8, 16), bound 2^ 10.000, 0.222 sec
Curves [ 16, 24), bound 2^ 10.585, 0.418 sec
Curves [ 24, 32), bound 2^ 11.000, 0.416 sec
Curves [ 32, 40), bound 2^ 11.322, 0.401 sec
Factors ECM: [21788197859]
Fully factored. N 1044877821881837432173434582227: [2953, 21788197859, 16239802890289801]
Time 1.692 sec.
As you see it takes just one second to factor quite complex 100-bit number. See above that first I do Trial Division stage, this is simplest factorization algorithm working in O(Sqrt(N)) time. After Trial Division I create many curves with growing Bound, each new curve has a bigger chance to factor a number, but takes more time.
My ECM algorithm according to ECM in Wiki is as following:
Try to find small prime factors, by doing Trial Division method. Spend around 5-10% of total time in this step. Remove found factors from a number.
Check if number is probably prime with high condifence, for that I use Fermat Probability Test with 32 trials. To overcome Carmichael numbers you may also use Miller Rabin Test instead. If number is primes return it as only factor.
Generate curve parameters A, X, Y randomly and derive B from curve equation Y^2 = X^3 + AX + B (mod N). Check if curve is OK, value 4 * A ** 3 - 27 * B ** 2 should be non-zero.
Generate small primes through Sieve of Eratosthenes, primes below our Bound. Each prime raise to some small power, this raised prime would be called K. I do raising into power while it is smaller than some Bound2, which is Sqrt(Bound) in my case.
Compute elliptic point multiplication P = k * P, where K taken from previous step and P is (X, Y). Compute according to Elliptic Curve Multiplication Wiki.
Point multiplication uses Modular Inverse, which computes GCD(SomeValue, N) (greatest common divisor) according to Extended Euclidean Algorithm Wiki. If this GCD is not 1, then it gives non-1 factor of N, hence we found a factor and can advance ECM algorithm to search for remaining factors.
If all primes till Bound were multiplied and gave no factor then re-run ECM factorization algorithm (3.-6. above) again with another random curve and bigger Bound. In my code I take new bound by adding 512 to old bound.
If you want to tweak my code for computing some other bit sizes, then look into Test() function. Modify bits_factor = 100 to some other values, this sets to factor 100 bit random number as an example. This bits_factor variable controls only example of factorization. Modify also bits_amicable_max = 27 if you want to set what Amicable Number bit size to search for, 27 stands for 27-bit numbers, you may set it to 80 to find 80-bit amicable numbers.
In my code I support 3 different sources of Big Integer implementation, they allow you to have any big bit size, even thousands of bits. First source is CLang's/GCC's built-in __int128, this gives 128-bit integer arithmetics for free. If you have other compiler, you may disable macro SUPPORT_STD_U128. Secondly I used Boost.Multiprecision library, it allows you to do big integer arithmetics efficiently till 1024 bits, disable macro SUPPORT_BOOST if you don't want Boost. Afterwards I used GMP Library, this allows you to have even bitsize of Millions bits, disable SUPPORT_GMP if you don't want GMP.
Under Linux it is easy to install Boost and GMP, sudo apt install libgmp-dev libboost-all-dev. In Windows you have to install VCPKG package manager first, afterwards do vcpkg install boost gmp inside its installation folder.
Click on Try it online! link below if you want to see my code running online on GodBolt servers. But NOTE that GodBolt version of code is different from the one below on Github Gist, this GodBolt version has reduced sizes of parameters specifically to run faster, because their server has very had limit for program to run, probably 5-10 seconds, not more.
SOURCE CODE HERE. Below Try it online! link there is Github Gist link. I placed my code there because it is 36 KB in size, and StackOverflow has a post size limit of 30_000 symbols in total, so together post text and source code will be a way to much.
Try it online on GodBolt!
Github Gist Full Source Code
i <= 1 + boundof theforloop? – Minuet