Carmichael lambda function

Carmichael lambda function is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

Background

The Carmichael function, or Carmichael lambda function, is a function in number theory. The Carmichael lambda λ(n) of a positive integer n is the smallest positive integer m such that

${\displaystyle a^{m}\equiv 1{\pmod {n}}}$

holds for every integer coprime to n.

The Carmichael lambda function can be iterated, that is, called repeatedly on its result. If this iteration is performed k times, this is considered the k-iterated Carmichael lambda function. Thus, λ(λ(n)) would be the 2-iterated lambda function. With repeated, sufficiently large but finite k-iterations, the iterated Carmichael function will eventually compute as 1 for all positive integers n.

Task

• Write a function to obtain the Carmichael lambda of a positive integer. If the function is supplied by a core library of the language, this also may be used.

• Write a function to count the number of iterations k of the k-iterated lambda function needed to get to a value of 1. Show the value of λ(n) and the number of k-iterations for integers from 1 to 25.

• Find the lowest integer for which the value of iterated k is i, for i from 1 to 15.

Stretch task (an extension of the third task above)

• Find, additionally, for i from 16 to 25, the lowest integer n for which the number of iterations k of the Carmichael lambda function from n to get to 1 is i.

See also

• Wikipedia entry for Carmichael function
•   OEIS lambda(lambda(n)) function sequence

FreeBASIC

Translation of: Phix

Type PrimePower
    prime As Uinteger
    power As Uinteger
End Type

Function GCD(n As Uinteger, d As Uinteger) As Uinteger
    Return Iif(d = 0, n, GCD(d, n Mod d))
End Function

Function lcm(n As Uinteger, d As Uinteger) As Uinteger
    Return (n * d) \ GCD(n, d)
End Function

Function getPrimePowers(n As Uinteger, Byref cnt As Uinteger) As PrimePower Ptr
    Dim As PrimePower Ptr powers = Callocate(32 * Sizeof(PrimePower))
    cnt = 0
    Dim As Uinteger num = n
    
    For i As Uinteger = 2 To Sqr(n)
        If num Mod i = 0 Then
            powers[cnt].prime = i
            powers[cnt].power = 0
            While (num Mod i = 0)
                powers[cnt].power += 1
                num \= i
            Wend
            cnt += 1
        End If
    Next
    
    If num > 1 Then
        powers[cnt].prime = num
        powers[cnt].power = 1
        cnt += 1
    End If
    
    Return powers
End Function

Function phi(n As Uinteger) As Uinteger
    If n = 1 Then Return 1
    
    Dim As Uinteger cnt, result, i, p, k
    Dim As PrimePower Ptr powers = getPrimePowers(n, cnt)
    result = 1
    
    For i = 0 To cnt - 1
        p = powers[i].prime
        k = powers[i].power
        result *= (p - 1) * (p ^ (k - 1))
    Next
    
    Deallocate(powers)
    Return result
End Function

Function lambda(p As Uinteger, k As Uinteger) As Uinteger
    Dim As Uinteger res = phi(p ^ k)
    Return Iif(p = 2 Andalso k >= 3, res \ 2, res)
End Function

Function reducedTotient(n As Uinteger) As Uinteger
    If n = 1 Then Return 1
    
    Dim As Uinteger cnt, result, i
    Dim As PrimePower Ptr powers = getPrimePowers(n, cnt)
    result = 1
    
    For i = 0 To cnt - 1
        result = lcm(result, lambda(powers[i].prime, powers[i].power))
    Next
    
    Deallocate(powers)
    Return result
End Function

Function kIter(n As Uinteger) As Uinteger
    Return Iif(n <= 1, 0, 1 + kIter(reducedTotient(n)))
End Function

' Main program
Dim As Uinteger i, n
Print " n lambda(n) k(n)"
Print "-----------------"
For i = 1 To 25
    Print Using "##       ##   ##"; i; reducedTotient(i); kIter(i)
Next

Print !"\nIterations to 1     n     lambda(n)"
Print String(36, "=")

Dim As Double t0 = Timer
n = 1
For i = 0 To 17 '15
    While kIter(n) <> i
        n += 1
    Wend
    Print Using "##        ###,###,###  ###,###,###"; i; n; reducedTotient(n)
Next

Dim As Double t1 = Timer - t0
Print Using !"\nTime elapsed: # minutes ##.### seconds"; t1 \ 60; t1 Mod 60

Sleep

Output:

 n lambda(n) k(n)
-----------------
 1        1    0
 2        1    1
 3        2    2
 4        2    2
 5        4    3
 6        2    2
 7        6    3
 8        2    2
 9        6    3
10        4    3
11       10    4
12        2    2
13       12    3
14        6    3
15        4    3
16        4    3
17       16    4
18        6    3
19       18    4
20        4    3
21        6    3
22       10    4
23       22    5
24        2    2
25       20    4

Iterations to 1     n     lambda(n)
====================================
 0                  1            1
 1                  2            1
 2                  3            2
 3                  5            4
 4                 11           10
 5                 23           22
 6                 47           46
 7                283          282
 8                719          718
 9              1,439        1,438
10              2,879        2,878
11             34,549       34,548
12            138,197      138,196
13            531,441      354,294
14          1,594,323    1,062,882
15          4,782,969    3,188,646

Time elapsed: 2 minutes 9.137 seconds

stretch

16         14,348,907    9,565,938
17         43,046,721   28,697,814

Time elapsed:  47 minutes 54.000 seconds

jq

Adapted from Wren

Works with jq, the C implementation of jq (*)

Works with gojq, the Go implementation of jq (**)

Works with jaq, the Rust implementation of jq (***)

A point of interest in the following implementation is the way in which the filter `iteratedToOne/1` prevents contamination of the input object by encapsulating it using the pattern: `def F: {"global": .} |.... | .global;`

(*) The program seemingly made no progress after i == 12 in the second table and was terminated.

(**) gojq's space requirements for this task beyond i == 10 in the second table are prohibitive.

(***) The program shown below is compatible with jaq except that jaq does not allow arrays to be grown implicitly by assignment; a simple workaround would be to change .cache to be a JSON object. However the space requirements effectively prevent computation of the rows in the second table beyond the row with i == 11.

### Generic functions

def lpad($len): tostring | ($len - length) as $l | (" " * $l) + .;

def power($b): . as $in | reduce range(0;$b) as $i (1; . * $in);

# Emit an array of the prime factors of 'n' in order using a wheel with basis [2, 3, 5]
# e.g. 44 | primeFactors => [2,2,11]
# From https://rosettacode.org/wiki/Product_of_min_and_max_prime_factors#jq
def primeFactors:
  def out($i): until (.n % $i != 0; .factors += [$i] | .n = ((.n/$i)|floor) );
  if . < 2 then []
  else [4, 2, 4, 2, 4, 6, 2, 6] as $inc
    | { n: .,
        factors: [] }
    | out(2)
    | out(3)
    | out(5)
    | .k = 7
    | .i = 0
    | until(.k * .k > .n;
        if .n % .k == 0
        then .factors += [.k]
        | .n = ((.n/.k)|floor)
        else .k += $inc[.i]
        | .i = ((.i + 1) % 8)
        end)
    | if .n > 1 then .factors += [ .n ] else . end
  | .factors
  end;

# Define the helper function to take advantage of jq's tail-recursion optimization
def lcm($m; $n):
  def _lcm:
    # state is [m, n, i]
    if (.[2] % .[1]) == 0 then .[2] else (.[0:2] + [.[2] + $m]) | _lcm end;
  [$m, $n, $m] | _lcm;

def lcm(stream):
  reduce stream as $s (1; lcm(.; $s));

### Carmichael Lambda Function

def primePows:
  . as $n
  | { factPows: [] }
  | ($n | primeFactors) as $pf
  | .currFact = $pf[0]
  | .count = 1
  | reduce $pf[1:][] as $fact (.;
      if $fact != .currFact
      then .factPows += [[.currFact, .count]]
      | .currFact = $fact
      | .count = 1
      else .count += 1
      end)
  | .factPows + [[.currFact, .count]] ;

def phi($p; $r):
  ($p | power($r-1)) * ($p - 1);

# Initial cache values
def cache: [0, 1, 1, null, 2];   # [_, 1, phi(2; 1), _, phi(2; 2)]

# input: {cache}
# output: {cache, result}
def CarmichaelHelper($p; $r):
  ($p|power($r)) as $n
  | .cache[$n] as $cached
  | if $cached then .result=$cached
    else
      if $p > 2
      then .cache[$n] = phi($p; $r)
      else .cache[$n] = phi($p; $r - 1)
      end
      | .result = .cache[$n]
    end;

# input: {cache}
# output: {cache, result, a}
def CarmichaelLambda($n):
  .cache |= (. // cache)
  | .result = null
  | if $n < 1 then "CarmichaelLambda: argument must be a positive integer (vs\($n))" | error
    else
      if .cache[$n] then .result = .cache[$n]
      else ($n|primePows) as $pps
      | if ($pps|length) == 1
        then $pps[0][0] as $p
        | $pps[0][1] as $r
        | if $p > 2 then .cache[$n] = phi($p; $r)
          else .cache[$n] = phi($p; $r - 1)
          end
          | .result = .cache[$n]
        else .a = []
        | reduce $pps[] as $pp (.;
            CarmichaelHelper($pp[0]; $pp[1])
            | .a += [ .result ] )
        | .cache[$n] = lcm(.a[])
        | .result = .cache[$n]
        end
      end
    end ;

# input: {cache}
# output: {cache, result}
def iteratedToOne($j):
  {global: .,
   k: 0,
   $j}
  | until( .j <= 1;
        .j as $j
        | .global |= CarmichaelLambda($j)
        | .j = .global.result
        | .k +=  1 )
  | .global.result = .k
  | .global;

def task1($m):
  " n   λ   k",
  "----------",
  (foreach range(1; 1+$m) as $n ({};
     CarmichaelLambda($n)
     | .result as $lambda
     | iteratedToOne($n)
     | .result as $k
     | .emit = "\($n|lpad(2)) \($lambda|lpad(2)) \($k|lpad(2))" )
   | .emit ) ;

def task2($maxI; $maxN):
  "\nIterations to 1       i     lambda(i)",
  "=====================================",
  "   0                  1            1",
 ( . // {}
  | .cache |= (. // cache)
  | .found = [true, (range(0; $maxN)|false)]
  | .i=1
  | .n=null
  | while (.i <= $maxI and .n != $maxN;
      .emit = null
      | iteratedToOne(.i)
      | .n = .result
      | if (.found[.n]|not)
        then .found[.n] = true
        | CarmichaelLambda(.i)
        | .result as $lambda 
        | .emit = "\(.n|lpad(4)) \(.i|lpad(18)) \($lambda|lpad(12))"
        end
      | .i += 1 )
  | select(.emit).emit );

task1(25),
"",
task2(5*1e6; 15)

Output:

 n   λ   k
----------
 1  1  0
 2  1  1
 3  2  2
 4  2  2
 5  4  3
 6  2  2
 7  6  3
 8  2  2
 9  6  3
10  4  3
11 10  4
12  2  2
13 12  3
14  6  3
15  4  3
16  4  3
17 16  4
18  6  3
19 18  4
20  4  3
21  6  3
22 10  4
23 22  5
24  2  2
25 20  4

Iterations to 1       i     lambda(i)
=====================================
   0                  1            1
   1                  2            1
   2                  3            2
   3                  5            4
   4                 11           10
   5                 23           22
   6                 47           46
   7                283          282
   8                719          718
   9               1439         1438
  10               2879         2878
  11              34549        34548
  12             138197       138196
[...terminated...]

Julia

using DelimitedFiles
using Primes

const maxcached = typemax(Int32)
const carmicache = zeros(Int32, maxcached) # NB: 8 gigabytes memory used for this size cache
const rlock = ReentrantLock()

""" Carmichael reduced totient function lambda(n) """
function lambda(n::Integer)
    @assert n > 0
    n <= maxcached && carmicache[n] > 0 && return Int(carmicache[n])
    lam = 1
    for (p, e) in factor(n)
        if p == 2 && e > 2
            lam = lcm(lam, 2^(e - 2))
        else
            lam = lcm(lam, (p - 1) * p^(e - 1))
        end
    end
    if n <= maxcached
        lock(rlock)
        carmicache[n] = lam
        unlock(rlock)
    end
    return lam
end

""" 
Return k for the k-fold iterated lambda function, where k is count of iterations to 1
"""
function iterated_lambdas_to_one(i)
    k = 0
    while i > 1
        i = lambda(i)
        k += 1
    end
    return k
end

""" 
    Calculate values for the task. Due to a runtime of several days for a sequence length
    of 25, the function switches to multithreading for values over about 2 billion. Note
    that (empirically) adjacent pairs of values within the sequence above 100000 are less
    than a factor of 6 apart, so each next large value is sought within that range.
"""
function print_iterated_lambdas(upto = 25)
    println("Listing of (n, lambda(n), k for iteration to 1) for integers from 1 to 25:")
    firsts = zeros(Int, upto)
    if isfile("lambdas.txt")
        precalc = readdlm("lambdas.txt", Int64)
        if length(precalc) == upto
            firsts .= precalc # restore work done if resterted (since it takes days)
        end
    end
    for i = 1:25
        lam = lambda(i)
        k = iterated_lambdas_to_one(i)
        print("($i, $lam, $k)", i % 5 == 0 ? "\n" : "  ")
    end
    target = something(findlast(!iszero, firsts), 1)
    if firsts[target] < maxcached ÷ 7 # get first ones with single thread
        for i = 2:maxcached
            n = iterated_lambdas_to_one(i)
            if n <= upto && (firsts[n] == 0 || firsts[n] > i)
                firsts[n] = i
            end
        end
    end
    target = something(findlast(!iszero, firsts), 1) + 1
    while target <= upto # multithreaded for larger targets
        startval = firsts[target - 1] + 1
        for j in startval:startval:startval*6
            j += iseven(j) # make odd
            @Threads.threads for i in j:2:j+startval-1 # should not be even if i > 2
                n = iterated_lambdas_to_one(i)
                if n == target
                    if firsts[n] == 0 || firsts[n] > i
                        try
                            lock(rlock)
                            firsts[n] = i
                            writedlm("lambdas.txt", firsts) # save work done for restarts
                            unlock(rlock)
                        catch y
                            unlock(rlock)
                            @warn y
                            rethrow()
                        end
                    end
                    break
                end
                firsts[target] > 0 && firsts[target] < i && break
            end
        end
        target += 1
    end
    println("\nIterations to 1     n      lambda(n)\n=====================================")
    println("   0                1             1")
    for n = 1:upto
        println(lpad(n, 4) * lpad(firsts[n], 17) * lpad(lambda(firsts[n]), 14))
    end
end

print_iterated_lambdas()

Output:

Listing of (n, lambda(n), k for iteration to 1) for integers from 1 to 25:
(1, 1, 0)  (2, 1, 1)  (3, 2, 2)  (4, 2, 2)  (5, 4, 3)
(6, 2, 2)  (7, 6, 3)  (8, 2, 2)  (9, 6, 3)  (10, 4, 3)
(11, 10, 4)  (12, 2, 2)  (13, 12, 3)  (14, 6, 3)  (15, 4, 3)
(16, 4, 3)  (17, 16, 4)  (18, 6, 3)  (19, 18, 4)  (20, 4, 3)
(21, 6, 3)  (22, 10, 4)  (23, 22, 5)  (24, 2, 2)  (25, 20, 4)

Iterations to 1     n      lambda(n)
=====================================
   0                1             1
   1                2             1
   2                3             2
   3                5             4
   4               11            10
   5               23            22
   6               47            46
   7              283           282
   8              719           718
   9             1439          1438
  10             2879          2878
  11            34549         34548
  12           138197        138196
  13           531441        354294
  14          1594323       1062882
  15          4782969       3188646
  16         14348907       9565938
  17         43046721      28697814
  18         86093443      86093442
  19        258280327     258280326
  20        688747547     688747546
  21       3486784401    2324522934
  22      10460353203    6973568802
  23      31381059609   20920706406 
  24      94143178827   62762119218
  25     282429536481  188286357654

Mathematica |Wolfram Language

Translation of: Python

IteratedToOne[i_] := Module[{k = 0, iter = i}, 
  While[iter > 1, iter = CarmichaelLambda[iter];
   k++;];
  k]
Print["Listing of (n, lambda(n), k for iteration to 1) for integers from 1 to 25:"];
Do[lam = CarmichaelLambda[i];
 k = IteratedToOne[i];
 If[Mod[i, 5] == 0, Print[{i, lam, k}], Print[{i, lam, k}, "  "]], {i, 1, 25}]

upTo = 20;
maxToTest = 100000000000;
firsts = Table[0, {upTo + 1}];
firsts[[1]] = 1;

Print["\nIterations to 1     n   lambda(n)\n=================================="];
Print[StringForm["``     ``     ``", 0, 1, 1]];

Do[n = IteratedToOne[i];
 If[0 < n <= upTo && firsts[[n + 1]] == 0, firsts[[n + 1]] = i;
  j = If[n < 1, 0, CarmichaelLambda[i]];
  Print[StringForm["``     ``     ``", n, i, j]];
  If[n >= upTo, Break[]];], {i, 2, maxToTest}]

Output:

Listing of (n, lambda(n), k for iteration to 1) for integers from 1 to 25:
{1, 1, 0}  
{2, 1, 1}  
{3, 2, 2}  
{4, 2, 2}  
{5, 4, 3}
{6, 2, 2}  
{7, 6, 3}  
{8, 2, 2}  
{9, 6, 3}  
{10, 4, 3}
{11, 10, 4}  
{12, 2, 2}  
{13, 12, 3}  
{14, 6, 3}  
{15, 4, 3}
{16, 4, 3}  
{17, 16, 4}  
{18, 6, 3}  
{19, 18, 4}  
{20, 4, 3}
{21, 6, 3}  
{22, 10, 4}  
{23, 22, 5}  
{24, 2, 2}  
{25, 20, 4}
Iterations to 1     n   lambda(n)
==================================
0     1     1
1     2     1
2     3     2
3     5     4
4     11     10
5     23     22
6     47     46
7     283     282
8     719     718
9     1439     1438
10     2879     2878
11     34549     34548
12     138197     138196
13     531441     354294
14     1594323     1062882
15     4782969     3188646

Perl

# 20240925 Perl programming solution

use strict;
use warnings;
use ntheory qw(is_carmichael carmichael_lambda);

sub iterated_to_one {
   my ($n) = @_;
   for (my $k = 0;;) { ($n = carmichael_lambda($n)) > 1 ? $k++ : return ++$k }
}

print " n   λ   k\n----------\n";
for my $n (1..25) {
   printf "%2d  %2d  %2d\n", $n, carmichael_lambda($n), iterated_to_one($n)
}

print "\nIterations to 1       i     lambda(i)\n",'='x37,"\n";
print "   0                  1            1\n";

my ($max_n, $max_i, @found) = (15, 5_000_000, (1));
for my $i (1 .. $max_i) {
   unless ($found[ my $n = iterated_to_one($i) ]) {
      printf "%4d %18d %12d\n", $n, $i, carmichael_lambda($i);
      $n == $max_n ?  last  : ( $found[$n] = 1 ) 
   }
}

Output:

 n   λ   k
----------
 1   1   1
 2   1   1
 3   2   2
 4   2   2
 5   4   3
 6   2   2
 7   6   3
 8   2   2
 9   6   3
10   4   3
11  10   4
12   2   2
13  12   3
14   6   3
15   4   3
16   4   3
17  16   4
18   6   3
19  18   4
20   4   3
21   6   3
22  10   4
23  22   5
24   2   2
25  20   4

Iterations to 1       i     lambda(i)
=====================================
   0                  1            1
   1                  1            1
   2                  3            2
   3                  5            4
   4                 11           10
   5                 23           22
   6                 47           46
   7                283          282
   8                719          718
   9               1439         1438
  10               2879         2878
  11              34549        34548
  12             138197       138196
  13             531441       354294
  14            1594323      1062882
  15            4782969      3188646

Phix

with javascript_semantics
function lambda(sequence pk)
    integer {p,k} = pk, res = phi(power(p,k))
    return iff(p=2 and k>=3 ? res/2 : res)
end function

function reduced_totient(integer n)
    sequence p = prime_powers(n)
    for i,pk in p do
        p[i] = lambda(pk)
    end for
    return lcm(p)
-- alt: neater but ~2* slower
--  return lcm(apply(prime_powers(n),lambda))
end function

function k_iter(integer i)
    return iff(i>1?1+k_iter(reduced_totient(i)):0)
end function

function ident(integer i) return i end function

sequence n = tagset(25)
for i,d in {"n","eulers_totient","reduced_totient","k"} do
    sequence r = apply(n,{ident,phi,reduced_totient,k_iter}[i])
    printf(1,"%17s: %s\n",{d,join(r," ",fmt:="%2d")})
end for

atom t0 = time()
printf(1,"\n k           i   lambda(i)\n")
printf(1,  "==========================\n")
integer i = 1
for k=0 to iff(platform()=JS?12:15) do -- (keep JS < 3s)
    while k_iter(i)!=k do i += 1 end while
    printf(1,"%2d %,9d %,9d\n",{k,i,reduced_totient(i)})
end for
?elapsed(time()-t0)

Output:

                n:  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
   eulers_totient:  1  1  2  2  4  2  6  4  6  4 10  4 12  6  8  8 16  6 18  8 12 10 22  8 20
  reduced_totient:  1  1  2  2  4  2  6  2  6  4 10  2 12  6  4  4 16  6 18  4  6 10 22  2 20
                k:  0  1  2  2  3  2  3  2  3  3  4  2  3  3  3  3  4  3  4  3  3  4  5  2  4

 k           i   lambda(i)
==========================
 0           1           1
 1           2           1
 2           3           2
 3           5           4
 4          11          10
 5          23          22
 6          47          46
 7         283         282
 8         719         718
 9       1,439       1,438
10       2,879       2,878
11      34,549      34,548
12     138,197     138,196
13     531,441     354,294
14   1,594,323   1,062,882
15   4,782,969   3,188,646
"2 minutes and 27s"

stretch

16  14,348,907   9,565,938
17  43,046,721  28,697,814

I estimate that took ~47mins, at which point I killed it.
The Python code took 13mins 58 secs to get to 15, on the same box.

Python

Python has the Carmichael function in the SymPy library, where it is called reduced_totient.

from sympy import reduced_totient

def iterated_to_one(i):
    """ return k for the k-fold iterated lambda function where k is the first time iteration reaches 1 """
    k = 0
    while i > 1:
        i = reduced_totient(i)
        k += 1
    return k


if __name__ == "__main__":
    print("Listing of (n, lambda(n), k for iteration to 1) for integers from 1 to 25:")
    for i in range(1, 26):
        lam = reduced_totient(i)
        k = iterated_to_one(i)
        print(f'({i}, {lam}, {k})', end = '\n' if i % 5 == 0 else '  ')

    UP_TO = 20
    MAX_TO_TEST = 100_000_000_000
    FIRSTS = [0] * (UP_TO + 1)
    FIRSTS[0] = 1

    print('\nIterations to 1     n   lambda(n)\n==================================')
    print('   0                1          1')
    for i in range(2, MAX_TO_TEST):
        n = iterated_to_one(i)
        if 0 < n <= UP_TO and FIRSTS[n] == 0:
            FIRSTS[n] = i
            j = 0 if n < 1 else int(reduced_totient(i))
            print(f"{n:>4}{i:>17}{j:>11}")    
            if n >= UP_TO:
                break

Output:

Listing of (n, lambda(n), k for iteration to 1) for integers from 1 to 25:
(1, 1, 0)  (2, 1, 1)  (3, 2, 2)  (4, 2, 2)  (5, 4, 3)
(6, 2, 2)  (7, 6, 3)  (8, 2, 2)  (9, 6, 3)  (10, 4, 3)
(11, 10, 4)  (12, 2, 2)  (13, 12, 3)  (14, 6, 3)  (15, 4, 3)
(16, 4, 3)  (17, 16, 4)  (18, 6, 3)  (19, 18, 4)  (20, 4, 3)
(21, 6, 3)  (22, 10, 4)  (23, 22, 5)  (24, 2, 2)  (25, 20, 4)

Iterations to 1     n   lambda(n)
==================================
   0                1          1
   1                2          1
   2                3          2
   3                5          4
   4               11         10
   5               23         22
   6               47         46
   7              283        282
   8              719        718
   9             1439       1438
  10             2879       2878
  11            34549      34548
  12           138197     138196
  13           531441     354294
  14          1594323    1062882
  15          4782969    3188646
  16         14348907    9565938
  17         43046721   28697814
  18         86093443   86093442
  19        258280327  258280326
  20        688747547  688747546

Raku

Translation of: Wren

# 20240228 Raku programming solution

use Prime::Factor;

sub phi(Int $p, Int $r) { return $p**($r - 1) * ($p - 1) }

sub CarmichaelLambda(Int $n) {

   state %cache = 1 => 1, 2 => 1, 4 => 2;

   sub CarmichaelHelper(Int $p, Int $r) {
      my Int $n = $p ** $r;
      return %cache{$n} if %cache{$n}:exists;
      return %cache{$n} = $p > 2 ?? phi($p, $r) !! phi($p, $r - 1)
   }

   if $n < 1 { die "'n' must be a positive integer." }
   return %cache{$n} if %cache{$n}:exists;
   if ( my %pps = prime-factors($n).Bag ).elems == 1 { 
      my ($p, $r) = %pps.kv>>.Int;
      return %cache{$n} = $p > 2 ?? phi($p, $r) !! phi($p, $r - 1)
   }
   return [lcm] %pps.kv.map: -> $k, $v { CarmichaelHelper($k.Int, $v) } 
}

sub iteratedToOne($i is copy) {
   my $k = 0;
   while $i > 1 { $i = CarmichaelLambda($i) andthen $k++ } 
   return $k
}

say " n   λ   k";
say "----------";
for 1..25 -> $n {
   printf "%2d  %2d  %2d\n", $n, CarmichaelLambda($n), iteratedToOne($n)
}

say "\nIterations to 1       i     lambda(i)";
say "=====================================";
say "   0                  1            1";

my ($maxI, $maxN) = 5e6, 10; # for N=15, takes around 47 minutes with an i5-10500T
my @found = True, |( False xx $maxN );
for 1 .. $maxI -> $i {
   unless @found[ my $n = iteratedToOne($i) ] {
      printf "%4d %18d %12d\n", $n, $i, CarmichaelLambda($i);
      @found[$n] = True andthen ( last if $n == $maxN )
   }
}

Output:

Same as Wren example .. and it is probably 30X times slower 😅

REXX

Translation of: Julia

Libraries: How to use
Library: Functions
Library: Numbers
Library: Settings
Library: Abend
Library: Sequences

The used recurring algorithm was also proposed by ChatGPT. Shown here is the main program and the procedure to calculate Carmichael lambda function.
Procedure Efactors (in Sequences) returns the unique prime factors and their exponents.

include Settings

say version; say 'Carmichael''s lambda function'; say
numeric digits 10
arg n
if n = '' then
   n = 25
call Time 'r'
say '----------'
say ' n   l   k'
say '----------'
do i = 1 to n
   k = 1; a = Lambda2(i); b = a
   do while b > 1
      k = k+1; b = Lambda2(b)
   end
   say Right(i,2) Right(a,3) Right(k,3)
end
say '----------'
say Format(Time('e'),,3) 'seconds'
say

say '----------'
say ' k       n'
say '----------'
call Time 'r'
a = 1
do i = 1 to 15
   do j = a until k = i
      k = 0; b = j
      do until b = 1
         k = k+1; b = Lambda2(b)
      end
   end
   say Right(k,2) Right(j,7)
   a = j+1
end
say '----------'
say Format(Time('e'),,3) 'seconds'
exit

Lambda2:
/* Carmichael's lambda function = reduced totient function */
procedure expose efac.
arg x
/* Recurring calculation */
f = Efactors(x)
y = 1
do i = 1 to f
   p = efac.factor.i; e = efac.exponent.i
   if p = 2 & e > 2 then
      a = 2**(e-2)
   else
      a = (p-1)*p**(e-1)
   y = Lcm(y,a)
end
return y

include Numbers
include Functions
include Sequences
include Abend

Output:

REXX-ooRexx_5.0.0(MT)_64-bit 6.05 23 Dec 2022
Carmichael's lambda function

----------
 n   l   k
----------
 1   1   1
 2   1   1
 3   2   2
 4   2   2
 5   4   3
 6   2   2
 7   6   3
 8   2   2
 9   6   3
10   4   3
11  10   4
12   2   2
13  12   3
14   6   3
15   4   3
16   4   3
17  16   4
18   6   3
19  18   4
20   4   3
21   6   3
22  10   4
23  22   5
24   2   2
25  20   4
----------
  0.000 seconds

----------
 k       n
----------
 1       1
 2       3
 3       5
 4      11
 5      23
 6      47
 7     283
 8     719
 9    1439
10    2879
11   34549
12  138197
13  531441
14 1594323
15 4782969
----------
483.303 seconds

Because up to k = 15 this already takes 8 minutes, I didn't try the stretched task.

Sidef

Built-in as Number#carmichael_lambda:

func iteratedToOne(n) {
    var k = 0
    while (n > 1) {
        n.carmichael_lambda!
        ++k
    }
    return k
}

say " n   λ   k"
say "----------"

for n in (1..25) {
    printf("%2d  %2d  %2d\n", n, n.carmichael_lambda, iteratedToOne(n))
}

say "\nIterations to 1       i     lambda(i)"
say "====================================="

var table = []

1..Inf -> each {|k|
    var n = iteratedToOne(k)
    if (!table[n]) {
        table[n] = true
        printf("%4s %18s %12s\n", n, k, k.carmichael_lambda)
        break if (n == 15)
    }
}

Output:

 n   λ   k
----------
 1   1   0
 2   1   1
 3   2   2
 4   2   2
 5   4   3
 6   2   2
 7   6   3
 8   2   2
 9   6   3
10   4   3
11  10   4
12   2   2
13  12   3
14   6   3
15   4   3
16   4   3
17  16   4
18   6   3
19  18   4
20   4   3
21   6   3
22  10   4
23  22   5
24   2   2
25  20   4

Iterations to 1       i     lambda(i)
=====================================
   0                  1            1
   1                  2            1
   2                  3            2
   3                  5            4
   4                 11           10
   5                 23           22
   6                 47           46
   7                283          282
   8                719          718
   9               1439         1438
  10               2879         2878
  11              34549        34548
  12             138197       138196
  13             531441       354294
  14            1594323      1062882
  15            4782969      3188646

Wren

Library: Wren-math

Library: Wren-fmt

Takes about 88 seconds on my machine (Core i7) to get up to n = 15 for the third part of the task, so haven't attempted the stretch goal.

import "./math" for Int
import "./fmt" for Fmt

var primePows = Fn.new { |n|
    var factPows = []
    var pf = Int.primeFactors(n)
    var currFact = pf[0]
    var count = 1
    for (fact in pf.skip(1)) {
        if (fact != currFact) {
            factPows.add([currFact, count])
            currFact = fact
            count = 1
        } else {
            count = count + 1
        }
    }
    factPows.add([currFact, count])
    return factPows
}

var phi = Fn.new { |p, r| p.pow(r-1) * (p - 1) }

var cache = { 1: 1, 2: phi.call(2, 1), 4: phi.call(2, 2) }

var CarmichaelHelper = Fn.new { |p, r|
    var n = p.pow(r)
    if (cache.containsKey(n)) return cache[n]
    if (p > 2) return cache[n] = phi.call(p, r)
    return cache[n] = phi.call(p, r - 1)
}

var CarmichaelLambda = Fn.new { |n|
    if (n < 1) Fiber.abort("'n' must be a positive integer.")
    if (cache.containsKey(n)) return cache[n]
    var pps = primePows.call(n)
    if (pps.count == 1) {
        var p = pps[0][0]
        var r = pps[0][1]
        if (p > 2) return cache[n] = phi.call(p, r)
        return cache[n] = phi.call(p, r - 1)       
    }
    var a = []
    for (pp in pps) a.add(CarmichaelHelper.call(pp[0], pp[1]))
    return cache[n] = Int.lcm(a)
}

var iteratedToOne = Fn.new { |i|
    var k = 0
    while (i > 1) {
        i = CarmichaelLambda.call(i)
        k = k + 1
    }
    return k
}

System.print(" n   λ   k")
System.print("----------")
for (n in 1..25) {
    var lambda = CarmichaelLambda.call(n)
    var k = iteratedToOne.call(n)
    Fmt.print("$2d  $2d  $2d", n, lambda, k)
}

System.print("\nIterations to 1       i     lambda(i)")
System.print("=====================================")
System.print("   0                  1            1")
var maxI = 5 * 1e6
var maxN = 15
var found = List.filled(maxN + 1, false)
found[0] = true
var i = 1
while (i <= maxI) {
    var n = iteratedToOne.call(i)
    if (!found[n]) {
        found[n] = true
        var lambda = CarmichaelLambda.call(i)
        Fmt.print("$4d $,18d $,12d", n, i, lambda)
        if (n == maxN) break
    }
    i = i + 1
}

Output:

 n   λ   k
----------
 1   1   0
 2   1   1
 3   2   2
 4   2   2
 5   4   3
 6   2   2
 7   6   3
 8   2   2
 9   6   3
10   4   3
11  10   4
12   2   2
13  12   3
14   6   3
15   4   3
16   4   3
17  16   4
18   6   3
19  18   4
20   4   3
21   6   3
22  10   4
23  22   5
24   2   2
25  20   4

Iterations to 1       i     lambda(i)
=====================================
   0                  1            1
   1                  2            1
   2                  3            2
   3                  5            4
   4                 11           10
   5                 23           22
   6                 47           46
   7                283          282
   8                719          718
   9              1,439        1,438
  10              2,879        2,878
  11             34,549       34,548
  12            138,197      138,196
  13            531,441      354,294
  14          1,594,323    1,062,882
  15          4,782,969    3,188,646