experimental-projects/pseudoprimes/generators/erdos_generate_carmichael_with_prefix.pl

#!/usr/bin/perl

# Erdos construction method for Carmichael numbers:
#   1. Choose an even integer L with many prime factors.
#   2. Let P be the set of primes d+1, where d|L and d+1 does not divide L.
#   3. Find a subset S of P such that prod(S) == 1 (mod L). Then prod(S) is a Carmichael number.

# Alternatively:
#   3. Find a subset S of P such that prod(S) == prod(P) (mod L). Then prod(P) / prod(S) is a Carmichael number.

use 5.020;
use warnings;
use ntheory qw(:all);
use experimental qw(signatures);

# Modular product of a list of integers
sub vecprodmod ($arr, $mod) {
    my $prod = 1;
    foreach my $k (@$arr) {
        $prod = mulmod($prod, $k, $mod);
    }
    $prod;
}

# Primes p such that p-1 divides L and p does not divide L
sub lambda_primes ($L) {
    grep { $L % $_ != 0 } grep { $_ > 2 and is_prime($_) } map { $_ + 1 } divisors($L);
}

sub method_1 ($L) {     # smallest numbers first

    my @P = lambda_primes($L);
    #my $prefix = vecprod(3, 5, 17, 23);
    my $prefix = vecprod(3, 5, 17, 23, 29, 53, 83, 89);

    if (gcd(vecprod(@P), $prefix) != $prefix) {
        say "Skipping L = $L...";
        return;
    }

    @P = grep { gcd($_, $prefix) == 1 } @P;

    my $n = scalar(@P);

    foreach my $k (1 .. @P>>1) {
        next if (binomial($n, $k) > 1e6);

        forcomb {
            if (mulmod($prefix, vecprodmod([@P[@_]], $L), $L) == 1) {
                say vecprod(@P[@_]) * $prefix;
            }
        } scalar(@P), $k;
    }

    #~ my $B = vecprodmod(\@P, $L);
    #~ my $T = vecprod(@P);

    #~ foreach my $k (1 .. @P>>1) {

        #~ next if (binomial($n, $k) > 1e6);

        #~ forcomb {
            #~ if (vecprodmod([@P[@_]], $L) == $B) {
                #~ my $S = vecprod(@P[@_]);
                #~ say ($T / $S) if ($T != $S);
            #~ }
        #~ } $n, $k;
    #~ }
}

foreach my $L(

#24, 36, 48, 60, 120, 180, 240, 360, 720, 840, 1260, 1680, 2520, 5040, 7560, 10080, 15120, 20160, 25200, 27720, 45360, 50400, 55440, 83160, 110880, 166320, 221760, 277200, 332640, 498960, 554400, 665280, 720720, 1081080, 1441440, 2162160

#10080, 15120, 20160, 25200, 27720, 45360, 50400, 55440, 83160, 110880, 166320, 221760, 277200, 332640, 498960, 554400, 665280, 720720, 1081080, 1441440, 2162160, 2882880, 3603600, 4324320, 6486480, 7207200, 8648640, 10810800, 14414400, 17297280, 21621600, 32432400, 36756720, 43243200, 61261200, 73513440, 110270160, 122522400, 147026880, 183783600, 245044800, 294053760, 367567200, 551350800, 698377680, 735134400, 1102701600, 1396755360, 2095133040, 2205403200, 2327925600, 2793510720, 3491888400, 4655851200, 5587021440, 6983776800, 10475665200, 13967553600, 20951330400, 27935107200, 41902660800, 48886437600, 64250746560, 73329656400, 80313433200, 97772875200, 128501493120, 146659312800, 160626866400, 240940299600, 293318625600, 321253732800, 481880599200, 642507465600, 963761198400, 1124388064800, 1606268664000, 1686582097200, 1927522396800, 2248776129600, 3212537328000, 3373164194400, 4497552259200, 6746328388800, 8995104518400, 9316358251200

#34082048, 64352288, 73545472, 128704576, 147090944, 257409152, 294181888, 374902528, 514818304, 1176727552, 1231886656, 2059273216

#11531520, 80720640, 126846720, 149909760, 329801472, 887927040, 1049368320, 1233872640, 1649007360, 2410087680, 11543051520, 16040344320, 31331139840, 40971490560, 48420852480, 106525875456, 219317978880, 286800433920, 473034481920, 532629377280, 745681128192, 3728405640960, 3979054759680, 5203379301120, 6149448264960, 67643930914560, 473507516401920

#144144, 720720, 1585584, 2306304, 15135120, 25369344, 31279248, 157549392, 1049368320, 1233872640, 1533692160, 1815061248, 2048142096, 6146300160

#5040, 7560, 10080, 21420, 57960, 59616, 109200, 110160, 163800, 178020, 308448, 540540, 811440, 3908520, 7592400, 304045280
#3780, 4032, 5040, 10080, 22176, 37800, 75600, 648000

#73545472, 128704576, 147090944, 202250048, 257409152, 294181888, 374902528, 514818304, 1029636608, 1231886656

#2130128, 2434432, 2818816, 3587584, 4260256, 4868864, 8520512, 10506496, 18386368, 19731712, 34082048, 36644608, 36772736, 53557504, 57209152, 63295232, 73545472, 87030944, 128704576, 147090944, 174061888, 202250048, 294181888, 348123776, 374902528, 618345728

#73545472, 128704576, 147090944, 202250048, 257409152, 294181888, 374902528, 514818304, 1029636608, 1231886656

#5040, 75600, 5040, 3780, 37800, 4032, 10080, 22176, 648000
#1240624, 2626624, 7771456, 10506496, 15542912, 139678448

#map{ $_ * 16 } 1188, 2024, 285824, 5789168, 7201568, 15542912, 4352299952

#map { $_ * 7 } (2130128, 2434432, 2818816, 3587584, 4260256, 4868864, 8520512, 10506496, 18386368, 19731712, 34082048, 36644608, 36772736, 53557504, 57209152, 63295232, 73545472, 87030944, 128704576, 147090944, 174061888, 202250048, 257409152, 294181888, 348123776, 374902528, 514818304, 618345728, 1029636608, 1231886656)

#divisors(2237587312658553856)

#map {$_} (285824, 655424, 1350272, 2618528, 2700544, 5789168, 7201568, 15542912, 34105456, 127359232, 251075968, 1426704224, 4352299952, 83510730688, 1189651277584, 38988702524464)
(1240624, 4868864, 8520512, 10506496, 18386368, 27419392, 34082048, 128704576, 139678448, 147090944, 202250048, 294181888, 374902528, 34818399616)

) {
    method_1($L);
    #method_2($L);
}