# Giuga numbers

Giuga numbers
You are encouraged to solve this task according to the task description, using any language you may know.
Definition

A Giuga number is a composite number n which is such that each of its distinct prime factors f divide (n/f - 1) exactly.

All known Giuga numbers are even though it is not known for certain that there are no odd examples.

Example

30 is a Giuga number because its distinct prime factors are 2, 3 and 5 and:

• 30/2 - 1 = 14 is divisible by 2
• 30/3 - 1 = 9 is divisible by 3
• 30/5 - 1 = 5 is divisible by 5

Task

Determine and show here the first four Giuga numbers.

Stretch

Determine the fifth Giuga number and any more you have the patience for.

References

## 11l

Translation of: Python
```F isGiuga(m)
V n = m
V f = 2
V l = sqrt(n)
L
I n % f == 0
I ((m I/ f) - 1) % f != 0
R 0B
n I/= f
I f > n
R 1B
E
f++
I f > l
R 0B

V n = 3
V c = 0
print(‘The first 4 Giuga numbers are:’)
L c < 4
I isGiuga(n)
c++
print(n)
n++```
Output:
```The first 4 Giuga numbers are:
30
858
1722
66198
```

## ABC

Based on the Algol 68 sample,

```HOW TO REPORT is.giuga n:
\ each prime factor must appear only once, e.g.: for 2:
\ [ ( n / 2 ) - 1 ] mod 2 = 0  => n / 2 is odd => n isn't divisible by 4
\ similarly for other primes
PUT  1, 3, 1, floor( n / 2 ) IN f.count, f, giuga, v
WHILE f <= v AND giuga = 1:
IF v mod f = 0:
PUT f.count + 1 IN f.count
IF ( ( floor( n / f ) ) - 1 ) mod f <> 0: PUT 0 IN giuga
PUT floor( v / f ) IN v
PUT f + 2 IN f
IF giuga = 1:             \ n is still a candidate, check it is not prime
IF f.count = 1: FAIL       \ only 1 factor - it is prime so not giuga
REPORT giuga = 1

PUT 0, -2 IN g.count, n
WHILE g.count < 4:
PUT n + 4 IN n                          \ assume the numbers are all even
IF is.giuga n:
PUT g.count + 1 IN g.count
WRITE n```
Output:
```30 858 1722 66198
```

## ALGOL 68

As with the Wren and other samples, assumes the first four Giuga numbers are even and also uses the fact that no prime squared can be a divisor of a Giuga number.

```BEGIN # find some Giuga numbers, composites n such that all their distinct   #
#                        prime factors f exactly divide ( n / f ) - 1  #

# find the first four Giuga numbers                                      #
# each prime factor must appear only once, e.g.: for 2:                  #
# [ ( n / 2 ) - 1 ] mod 2 = 0  => n / 2 is odd => n isn't divisible by 4 #
# similarly for other primes                                             #
INT g count := 0;
FOR n FROM 2 BY 4 WHILE g count < 4 DO # assume the numbers are all even #
INT  v        := n OVER 2;
BOOL is giuga := TRUE;
INT  f count  := 1;
FOR f FROM 3 BY 2 WHILE f <= v AND is giuga DO
IF v MOD f = 0 THEN
# have a prime factor                                        #
f count +:= 1;
is giuga := ( ( n OVER f ) - 1 ) MOD f = 0;
v OVERAB f
FI
OD;
IF is giuga THEN
# n is still a candidate, check it is not prime                  #
IF f count > 1 THEN
g count +:= 1;
print( ( " ", whole( n, 0 ) ) )
FI
FI
OD
END```
Output:
``` 30 858 1722 66198
```

## ALGOL W

Translation of: ALGOL 68
```begin % find some Giuga numbers, composites n such that all their distinct   %
%                        prime factors f exactly divide ( n / f ) - 1  %

% find the first four Giuga numbers                                      %
% each prime factor must appear only once, e.g.: for 2:                  %
% [ ( n / 2 ) - 1 ] mod 2 = 0  => n / 2 is odd => n isn't divisible by 4 %
% similarly for other primes                                             %

integer gCount, n;
gCount :=  0;
n      := -2;
while begin n := n + 4;
gCount < 4
end
do begin    % assume the numbers are all even                            %
integer v, f, fCount;
logical isGiuga;
v       := n div 2;
isGiuga := TRUE;
fCount  := 1;
f       := 1;
while begin f := f + 2;
f <= v and isGiuga
end
do begin
if v rem f = 0 then begin
% have a prime factor                                        %
fCount  := fCount + 1;
isGiuga := ( ( n div f ) - 1 ) rem f = 0;
v := v div f
end if_v_rem_f_eq_0
end while_f_le_v_and_isGiuga ;
if isGiuga then begin
% n is still a candidate, check it is not prime                  %
if fCount > 1 then begin
gCount := gCount + 1;
writeon( i_w := 1, s_w := 0, " ", n )
end if_fCount_gt_1
end if_isGiuga
end while_gCount_lt_4
end.```
Output:

Same as Algol 68

## Arturo

```giuga?: function [n]->
and? -> not? prime? n
-> every? factors.prime n 'f
-> zero? (dec n/f) % f

print.lines select.first:4 1..∞ => giuga?
```
Output:
```30
858
1722
66198```

## AWK

```# syntax: GAWK -f GIUGA_NUMBER.AWK
BEGIN {
n = 3
stop = 4
printf("Giuga numbers 1-%d:",stop)
while (count < stop) {
if (is_giuga(n)) {
count++
printf(" %d",n)
}
n++
}
printf("\n")
exit(0)
}
function is_giuga(m,  f,l,n) {
n = m
f = 2
l = sqrt(n)
while (1) {
if (n % f == 0) {
if (((m / f) - 1) % f != 0) { return(0) }
n /= f
if (f > n) { return(1) }
}
else {
if (++f > l) { return(0) }
}
}
}
```
Output:
```Giuga numbers 1-4: 30 858 1722 66198
```

## BASIC

### BASIC256

```n = 3
c = 0
limit = 4
print "The first"; limit; " Giuga numbers are: ";
do
if isGiuga(N) then c += 1: print n; "  ";
n += 1
until c = limit
end

function isGiuga(m)
n = m
f = 2
l = sqr(n)
while True
if n mod f = 0 then
if ((m / f) - 1) mod f <> 0 then return false
n /= f
if f > n then return true
else
f += 1
if f > l then return false
end if
end while
end function```
Output:
`Same as FreeBASIC entry.`

### FreeBASIC

```Function isGiuga(m As Uinteger) As Boolean
Dim As Uinteger n = m, f = 2, l = Sqr(n)
Do
If n Mod f = 0 Then
If ((m / f) - 1) Mod f <> 0 Then Return False
n /= f
If f > n Then Return True
Else
f += 1
If f > l Then Return False
End If
Loop
End Function

Dim As Uinteger n = 3, c = 0, limit = 4
Print "The first "; limit; " Giuga numbers are: ";
Do
If isGiuga(n) Then c += 1: Print n; "  ";
n += 1
Loop Until c = limit
```
Output:
`The first 4 Giuga numbers are: 30  858  1722  66198`

### Gambas

```Public Sub Main()

Dim n As Integer = 3, c As Integer = 0, limit As Integer = 4

Print "The first "; limit; " Giuga numbers are: ";
Do
If isGiuga(n) Then
c += 1
Print n; "  ";
Endif
n += 1
Loop Until c = limit

End

Function isGiuga(m As Integer) As Boolean

Dim n As Integer = m, f As Integer = 2, l As Integer = Sqr(n)

Do
If n Mod f = 0 Then
Dim q As Integer = (m / f)
If (q - 1) Mod f <> 0 Then Return False
n /= f
If f > n Then Return True
Else
f += 1
If f > l Then Return False
End If
Loop

End Function
```
Output:
`Same as FreeBASIC entry.`

### PureBasic

```Procedure.b isGiuga(m.i)
Define.i n = m,	f = 2, l = Sqr(n)
While #True
If Mod(n, f) = 0:
If Mod(((m / f) - 1), f) <> 0: ProcedureReturn #False: EndIf
n = n / f
If f > n: ProcedureReturn #True: EndIf
Else
f + 1
If f > l: ProcedureReturn #False: EndIf
EndIf
Wend
EndProcedure

OpenConsole()
Define.i n = 3, c = 0, limit = 4
Print("The first " + Str(limit) + " Giuga numbers are: ")
Repeat
If isGiuga(N):
c + 1
Print(Str(n) + "  ")
EndIf
n + 1
Until c = limit

Input()
CloseConsole()```
Output:
`Same as FreeBASIC entry.`

## C++

### Brute force

Based on the Go solution. Takes 26 minutes on my system (Intel Core i5 3.2GHz).

```#include <iostream>

// Assumes n is even with exactly one factor of 2.
bool is_giuga(unsigned int n) {
unsigned int m = n / 2;
auto test_factor = [&m, n](unsigned int p) -> bool {
if (m % p != 0)
return true;
m /= p;
return m % p != 0 && (n / p - 1) % p == 0;
};
if (!test_factor(3) || !test_factor(5))
return false;
static constexpr unsigned int wheel[] = {4, 2, 4, 2, 4, 6, 2, 6};
for (unsigned int p = 7, i = 0; p * p <= m; ++i) {
if (!test_factor(p))
return false;
p += wheel[i & 7];
}
return m == 1 || (n / m - 1) % m == 0;
}

int main() {
std::cout << "First 5 Giuga numbers:\n";
// n can't be 2 or divisible by 4
for (unsigned int i = 0, n = 6; i < 5; n += 4) {
if (is_giuga(n)) {
std::cout << n << '\n';
++i;
}
}
}
```
Output:
```First 5 Giuga numbers:
30
858
1722
66198
2214408306
```

### Faster version

Library: Boost
Translation of: Wren
```#include <boost/rational.hpp>

#include <algorithm>
#include <cstdint>
#include <iostream>
#include <vector>

using rational = boost::rational<uint64_t>;

bool is_prime(uint64_t n) {
if (n < 2)
return false;
if (n % 2 == 0)
return n == 2;
if (n % 3 == 0)
return n == 3;
for (uint64_t p = 5; p * p <= n; p += 4) {
if (n % p == 0)
return false;
p += 2;
if (n % p == 0)
return false;
}
return true;
}

uint64_t next_prime(uint64_t n) {
while (!is_prime(n))
++n;
return n;
}

std::vector<uint64_t> divisors(uint64_t n) {
std::vector<uint64_t> div{1};
for (uint64_t i = 2; i * i <= n; ++i) {
if (n % i == 0) {
div.push_back(i);
if (i * i != n)
div.push_back(n / i);
}
}
div.push_back(n);
sort(div.begin(), div.end());
return div;
}

void giuga_numbers(uint64_t n) {
std::cout << "n = " << n << ":";
std::vector<uint64_t> p(n, 0);
std::vector<rational> s(n, 0);
p[2] = 2;
p[1] = 2;
s[1] = rational(1, 2);
for (uint64_t t = 2; t > 1;) {
p[t] = next_prime(p[t] + 1);
s[t] = s[t - 1] + rational(1, p[t]);
if (s[t] == 1 || s[t] + rational(n - t, p[t]) <= rational(1)) {
--t;
} else if (t < n - 2) {
++t;
uint64_t c = s[t - 1].numerator();
uint64_t d = s[t - 1].denominator();
p[t] = std::max(p[t - 1], c / (d - c));
} else {
uint64_t c = s[n - 2].numerator();
uint64_t d = s[n - 2].denominator();
uint64_t k = d * d + c - d;
auto div = divisors(k);
uint64_t count = (div.size() + 1) / 2;
for (uint64_t i = 0; i < count; ++i) {
uint64_t h = div[i];
if ((h + d) % (d - c) == 0 && (k / h + d) % (d - c) == 0) {
uint64_t r1 = (h + d) / (d - c);
uint64_t r2 = (k / h + d) / (d - c);
if (r1 > p[n - 2] && r2 > p[n - 2] && r1 != r2 &&
is_prime(r1) && is_prime(r2)) {
std::cout << ' ' << d * r1 * r2;
}
}
}
}
}
std::cout << '\n';
}

int main() {
for (uint64_t n = 3; n < 7; ++n)
giuga_numbers(n);
}
```
Output:
```n = 3: 30
n = 4: 1722 858
n = 5: 66198
n = 6: 24423128562 2214408306
```

## Delphi

Works with: Delphi version 6.0

```function IsGiugaNumber(N: integer): boolean;
var IA: TIntegerDynArray;
var I,V: integer;
begin
Result:=False;
if IsPrime(N) then exit;
GetPrimeFactors(N,IA);
for I:=0 to High(IA) do
begin
V:=N div IA[I] - 1;
if (V mod IA[I])<>0 then exit;
end;
Result:=True;
end;

procedure ShowGiugaNumbers(Memo: TMemo);
var I,Cnt: integer;
begin
Cnt:=0;
for I:=4 to High(integer) do
if IsGiugaNumber(I) then
begin
Inc(Cnt);
Memo.Lines.Add(IntToStr(I));
if Cnt>=4 then break;
end;
end;
```
Output:
```30
858
1722
66198

Elapsed Time: 4.991 Sec.

```

## EasyLang

Translation of: Python
```func giuga m .
n = m
for f = 2 to sqrt n
while n mod f = 0
if (m div f - 1) mod f <> 0
return 0
.
n = n div f
if f > n
return 1
.
.
.
.
n = 3
while cnt < 4
if giuga n = 1
cnt += 1
print n
.
n += 1
.```

## Euler

Uses the C style for loop procedure from the Sieve of Eratosthenes task

```begin
new for; new n; new gCount;
for   <- ` formal init; formal test; formal incr; formal body;
begin
label again;
init;
again:          if test then begin body; incr; goto again end else 0
end
'
;
gCount <- 0;
for( ` n <- 2 ', ` gCount < 4 ', ` n <- n + 4 '
, ` begin
new v; new f; new isGiuga; new fCount;
v       <- n % 2;
isGiuga <- true;
fCount  <- 1;
for( ` f <- 3 ', ` f <= v and isGiuga ', ` f <- f + 2 '
, ` if v mod f = 0 then begin
fCount  <- fCount + 1;
isGiuga <- [ [ n % f ] - 1 ] mod f = 0;
v <- v % f
end else 0
'
);
if isGiuga then begin
if fCount > 1 then begin
gCount <- gCount + 1;
out n
end else 0
end else 0
end
'
)
end \$
```

## Go

Translation of: Wren

I thought I'd see how long it would take to 'brute force' the fifth Giuga number and the answer (without using parallelization, Core i7) is about 1 hour 38 minutes.

```package main

import "fmt"

var factors []int
var inc = []int{4, 2, 4, 2, 4, 6, 2, 6}

// Assumes n is even with exactly one factor of 2.
// Empties 'factors' if any other prime factor is repeated.
func primeFactors(n int) {
factors = factors[:0]
factors = append(factors, 2)
last := 2
n /= 2
for n%3 == 0 {
if last == 3 {
factors = factors[:0]
return
}
last = 3
factors = append(factors, 3)
n /= 3
}
for n%5 == 0 {
if last == 5 {
factors = factors[:0]
return
}
last = 5
factors = append(factors, 5)
n /= 5
}
for k, i := 7, 0; k*k <= n; {
if n%k == 0 {
if last == k {
factors = factors[:0]
return
}
last = k
factors = append(factors, k)
n /= k
} else {
k += inc[i]
i = (i + 1) % 8
}
}
if n > 1 {
factors = append(factors, n)
}
}

func main() {
const limit = 5
var giuga []int
// n can't be 2 or divisible by 4
for n := 6; len(giuga) < limit; n += 4 {
primeFactors(n)
// can't be prime or semi-prime
if len(factors) > 2 {
isGiuga := true
for _, f := range factors {
if (n/f-1)%f != 0 {
isGiuga = false
break
}
}
if isGiuga {
giuga = append(giuga, n)
}
}
}
fmt.Println("The first", limit, "Giuga numbers are:")
fmt.Println(giuga)
}
```
Output:
```The first 5 Giuga numbers are:
[30 858 1722 66198 2214408306]
```

## Haskell

```--for obvious theoretical reasons the smallest divisor of a number bare 1
--must be prime
primeFactors :: Int -> [Int]
primeFactors n = snd \$ until ( (== 1) . fst ) step (n , [] )
where
step :: (Int , [Int] ) -> (Int , [Int] )
step (n , li) = ( div n h , li ++ [h] )
where
h :: Int
h = head \$ tail \$ divisors n --leave out 1

divisors :: Int -> [Int]
divisors n = [d | d <- [1 .. n] , mod n d == 0]

isGiuga :: Int -> Bool
isGiuga n = (divisors n /= [1,n]) && all (\i -> mod ( div n i - 1 ) i == 0 )
(primeFactors n)

solution :: [Int]
solution = take 4 \$ filter isGiuga [2..]
```
Output:
```[30,858,1722,66198]
```

## J

We can brute force this task building a test for giuga numbers and checking the first hundred thousand integers (which takes a small fraction of a second):

```giguaP=: {{ (1<y)*(-.1 p:y)**/(=<.) y ((_1+%)%]) q: y }}"0
```
```   1+I.giguaP 1+i.1e5
30 858 1722 66198
```

These numbers have some interesting properties but there's an issue with guaranteeing correctness of more sophisticated approaches. That said, here's a translation of the pari/gp implementation on the talk page:

```divisors=: [: /:~@, */&>@{@((^ i.@>:)&.>/)@q:~&__

giuga=: {{
r=. i.0
p=. (2) 0 1} s=. 1r2,}.(2>.y-1+t=.1)\$0
while. t do.
p=. p t}~ 4 p:t{p
s=. s t}~ (s{~t-1)+1%t{p
if. (1=t{s) +. 1 >: (t{s)+(y-t+1)%t{p do.
t=. t-1
elseif. t<y-3 do.
p=. p (t+1)}~ (p{~t) >. (%-.)s{~t
t=. t+1
else.
'c d'=. 2 x: s{~y-3
dc=. d-c
k=. (d^2)-dc
for_h. ({.~ <.@-:@>:@#) f=. divisors k do.
if. 0=dc|h+d do.
if. 0=dc|dkh=. d+k%h do.
py3=. p{~y-3
if. py3 < r1=. (h+d)%dc do.
if. py3 < r2=. dkh%dc do.
if. r1~:r2 do.
if. 1 p: r1 do.
if. 1 p: r2 do.
r=. r, d*r1*r2
end.
end.
end.
end.
end.
end.
end.
end.
end.
end.
r
}}
```

Example use:

```   giuga 1

giuga 2

giuga 3
30
giuga 4
1722 858
giuga 5
66198
giuga 6
24423128562 2214408306
```

## Java

Algorithm uses the mathematical facts that a Giuga number must be square-free and cannot be a semi-prime.

It does not assume that a Giuga number is even, and exhaustively tests all possible candidates up to approximately 147,000 in around 20 milliseconds.

```import java.util.ArrayList;
import java.util.Collections;
import java.util.List;

public final class GiugaNumbers {

public static void main(String[] aArgs) {
primes = List.of( 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59 );

List<Integer> primeCounts = List.of( 3, 4, 5 );
for ( int primeCount : primeCounts ) {
primeFactors = new ArrayList<Integer>(Collections.nCopies(primeCount, 0));
combinations(primeCount, 0, 0);
}

Collections.sort(results);
System.out.println("Found Giuga numbers: " + results);
}

private static void checkIfGiugaNumber(List<Integer> aPrimeFactors) {
final int product = aPrimeFactors.stream().reduce(1, Math::multiplyExact);

for ( int factor : aPrimeFactors ) {
final int divisor = factor * factor;
if ( ( product - factor ) % divisor != 0 ) {
return;
}
}

results.add(product);
}

private static void combinations(int aPrimeCount, int aIndex, int aLevel) {
if ( aLevel == aPrimeCount ) {
checkIfGiugaNumber(primeFactors);
return;
}

for ( int i = aIndex; i < primes.size(); i++ ) {
primeFactors.set(aLevel, primes.get(i));
combinations(aPrimeCount, i + 1, aLevel + 1);
}
}

private static List<Integer> primes;
private static List<Integer> primeFactors;
private static List<Integer> results = new ArrayList<Integer>();

}
```
Output:
```Found Giuga numbers: [30, 858, 1722, 66198]
```

## jq

Works with jq, the C implementation of jq

Works with gojq, the Go implementation of jq

Works with jaq, the Rust implementation of jq

The algorithm used here is naive but suffices to find the first four Giuga numbers in a reasonable amount of time whether using the C, Go, or Rust implementations. On my machine, gojq is fastest at about 12s.

```# `div/1` is defined firstly to take advantage of gojq's infinite-precision
# integer arithmetic, and secondly to ensure jaq returns an integer.
def div(\$j):
(. - (. % \$j)) / \$j | round;  # round is for the sake of jaq

# For convenience
def div(\$i; \$j): \$i|div(\$j);

def is_giuga:
. as \$m
| sqrt as \$limit
| {n: \$m, f: 2}
| until(.ans;
if (.n % .f) == 0
then if ((div(\$m; .f) - 1) % .f) != 0 then .ans = 0
else .n = div(.n; .f)
| if .f > .n then .ans = 1 end
end
else .f += 1
| if .f > \$limit then .ans = 0 end
end)
| .ans == 1 ;

limit(4; range(4; infinite) | select(is_giuga))```
Output:
```30
858
1722
66198
```

## Julia

```using Primes

isGiuga(n) = all(f -> f != n && rem(n ÷ f - 1, f) == 0, factor(Vector, n))

function getGiuga(N)
gcount = 0
for i in 4:typemax(Int)
if isGiuga(i)
println(i)
(gcount += 1) >= N && break
end
end
end

getGiuga(4)
```
Output:
```30
858
1722
66198
```

### Ad hoc faster version

```using Primes

function getgiugas(numberwanted, verbose = true)
n, found, nfound = 6, Int[], 0
starttime = time()
while nfound < numberwanted
if n % 5 == 0 || n % 7 == 0 || n % 11 == 0
for (p, e) in eachfactor(n)
(e != 1 || rem(n ÷ p - 1, p) != 0) && @goto nextnumber
end
verbose && println(n, "  (elapsed: ", time() - starttime, ")")
push!(found, n)
nfound += 1
end
@label nextnumber
n += 6 # all mult of 6
end
return found
end

@time getgiugas(2, false)
@time getgiugas(6)
```
Output:
```30  (elapsed: 0.0)
858  (elapsed: 0.0)
1722  (elapsed: 0.0)
66198  (elapsed: 0.0009999275207519531)
2214408306  (elapsed: 18.97099995613098)
24423128562  (elapsed: 432.06500005722046)
432.066249 seconds (235 allocations: 12.523 KiB)
```

## Lua

Translation of: ALGOL 68
```do --[[ find the first 4 Giuga numbers, composites n such that all their
distinct prime factors f exactly divide ( n / f ) - 1

each prime factor must appear only once, e.g.: for 2:
[ ( n / 2 ) - 1 ] mod 2 = 0  => n / 2 is odd => n isn't divisible by 4
similarly for other primes
]]

local gCount, n = 0, 2
while gCount < 4 do
n = n + 4                           -- assume the numbers are all even
local fCount, f, isGiuga, v = 1, 1, true, math.floor( n / 2 )
while f <= ( v - 2 ) and isGiuga do
f = f + 2
if v % f == 0 then                          -- have a prime factor
fCount  = fCount + 1
isGiuga = ( math.floor( n / f ) - 1 ) % f == 0
v       = math.floor( v / f )
end
end
if isGiuga then       -- n is still a candidate, check it is not prime
if fCount > 1 then
gCount = gCount + 1
io.write( " ", n )
end
end
end
end
```
Output:
``` 30 858 1722 66198
```

## Mathematica/Wolfram Language

Translation of: Julia
```IsGiuga[n_Integer] := Module[{factors},
factors = FactorInteger[n];
AllTrue[factors, Function[{f},
Mod[Quotient[n, f[[1]]] - 1, f[[1]]] == 0 && f[[1]] != n]]
]

GetGiuga[N_Integer] := Module[{giugaNumbers = {}, i = 4},
While[Length[giugaNumbers] < N,
If[IsGiuga[i], AppendTo[giugaNumbers, i]];
i++
];
giugaNumbers
]

Print[GetGiuga[4]]
```
Output:
```{30, 858, 1722, 66198}

```

## Maxima

```/* Predicate function that checks wether an integer is a Giuga number or not */
giugap(n):=if not primep(n) then block(ifactors(n),map(lambda([x],mod((n/x)-1,x)=0),map(first,%%)),
if length(unique(%%))=1 and apply(lhs,unique(%%))=0 then true)\$

/* Function that returns a list of the first len Giuga integers */
giuga_count(len):=block(
[i:1,count:0,result:[]],
while count<len do (if giugap(i) then (result:endcons(i,result),count:count+1),i:i+1),
result)\$

/* Test case */
giuga_count(4);
```
Output:
```[30,858,1722,66198]
```

## Nim

Translation of: FreeBASIC
```import std/math

func isGiuga(m: Natural): bool =
var n = m
var f = 2
var l = int(sqrt(n.toFloat))
while true:
if n mod f == 0:
if (m div f - 1) mod f != 0:
return false
n = n div f
if f > n:
return true
else:
inc f
if f > l:
return false

var n = 3
var c = 0
const Limit = 4
stdout.write "The first ", Limit, " Giuga numbers are: "
while true:
if n.isGiuga:
inc c
stdout.write n, " "
if c == Limit: break
inc n
echo()
```
Output:
```The first 4 Giuga numbers are: 30 858 1722 66198
```

## PARI/GP

Translation of: Julia
```is_giuga(n) = {
my(factors = factor(n));
for (i = 1, #factors[,1],
if (factors[i,1] == n, return(0));
if (Mod(n \ factors[i,1] - 1, factors[i,1]) != 0, return(0));
);
return(1);
}

get_giuga(N) = {
my(giugaNumbers = [], i = 4);
while(#giugaNumbers < N,
if (is_giuga(i), giugaNumbers = concat(giugaNumbers, i));
i++;
);
giugaNumbers
}

print(get_giuga(4))```
Output:
```[30, 858, 1722, 66198]

```

## Pascal

### Free Pascal

OK.Cheating to find square free numbers like julia in distance 6
That means always factors 2,3 and minimum one of 5,7,11.

```program Giuga;

{\$IFDEF FPC}
{\$MODE DELPHI}  {\$OPTIMIZATION ON,ALL}  {\$COPERATORS ON}
{\$ELSE}
{\$APPTYPE CONSOLE}
{\$ENDIF}
uses
sysutils
{\$IFDEF WINDOWS},Windows{\$ENDIF}
;
//######################################################################
//prime decomposition only squarefree and multiple of 6

type
tprimeFac = packed record
pfpotPrimIdx : array[0..9] of Uint64;
pfMaxIdx  : Uint32;
end;
tpPrimeFac = ^tprimeFac;
tPrimes = array[0..65535] of Uint32;

var
{\$ALIGN 8}
SmallPrimes: tPrimes;
{\$ALIGN 32}

procedure InitSmallPrimes;
//get primes. #0..65535.Sieving only odd numbers
const
MAXLIMIT = (821641-1) shr 1;
var
pr : array[0..MAXLIMIT] of byte;
p,j,d,flipflop :NativeUInt;
Begin
SmallPrimes[0] := 2;
fillchar(pr[0],SizeOf(pr),#0);
p := 0;
repeat
repeat
p +=1
until pr[p]= 0;
j := (p+1)*p*2;
if j>MAXLIMIT then
BREAK;
d := 2*p+1;
repeat
pr[j] := 1;
j += d;
until j>MAXLIMIT;
until false;

SmallPrimes[1] := 3;
SmallPrimes[2] := 5;
j := 3;
d := 7;
flipflop := (2+1)-1;//7+2*2,11+2*1,13,17,19,23
p := 3;
repeat
if pr[p] = 0 then
begin
SmallPrimes[j] := d;
inc(j);
end;
d += 2*flipflop;
p+=flipflop;
flipflop := 3-flipflop;
until (p > MAXLIMIT) OR (j>High(SmallPrimes));
end;

function OutPots(pD:tpPrimeFac;n:NativeInt):Ansistring;
var
s: String[31];
chk,p: NativeInt;
Begin
str(n,s);
result := s+' :';
with pd^ do
begin
chk := 1;
For n := 0 to pfMaxIdx-1 do
Begin
if n>0 then
result += '*';
p := pfpotPrimIdx[n];
chk *= p;
str(p,s);
result += s;
end;
str(chk,s);
result += '_chk_'+s+'<';
end;
end;

function IsSquarefreeDecomp6(var res:tPrimeFac;n:Uint64):boolean;inline;
//factorize only not prime/semiprime and squarefree n= n div 6
var
pr,i,q,idx :NativeUInt;
Begin
with res do
Begin
Idx := 2;

q := n DIV 5;
if n = 5*q then
Begin
pfpotPrimIdx[2] := 5;
n := q;
q := q div 5;
if q*5=n then
EXIT(false);
inc(Idx);
end;

q := n DIV 7;
if n = 7*q then
Begin
pfpotPrimIdx[Idx] := 7;
n := q;
q := q div 7;
if q*7=n then
EXIT(false);
inc(Idx);
end;

q := n DIV 11;
if n = 11*q then
Begin
pfpotPrimIdx[Idx] := 11;
n := q;
q := q div 11;
if q*11=n then
EXIT(false);
inc(Idx);
end;

if Idx < 3 then
Exit(false);

i := 5;
while i < High(SmallPrimes) do
begin
pr := SmallPrimes[i];
q := n DIV pr;
//if n < pr*pr
if pr > q then
BREAK;
if n = pr*q then
Begin
pfpotPrimIdx[Idx] := pr;
n := q;
q := n div pr;
if pr*q = n then
EXIT(false);
inc(Idx);
end;
inc(i);
end;
if n <> 1 then
begin
pfpotPrimIdx[Idx] := n;
inc(Idx);
end;
pfMaxIdx := idx;
end;
exit(true);
end;

function ChkGiuga(n:Uint64;pPrimeDecomp :tpPrimeFac):boolean;inline;
var
p : Uint64;
idx: NativeInt;
begin
with pPrimeDecomp^ do
Begin
idx := pfMaxIdx-1;
repeat
p := pfpotPrimIdx[idx];
result := (((n DIV p)-1)MOD p) = 0;
if not(result) then
EXIT;
dec(idx);
until idx<0;
end;
end;

const
LMT = 24423128562;//2214408306;//
var
PrimeDecomp :tPrimeFac;
T0:Int64;
n,n6 : UInt64;
cnt:Uint32;
Begin
InitSmallPrimes;

T0 := GetTickCount64;
with PrimeDecomp do
begin
pfpotPrimIdx[0]:= 2;
pfpotPrimIdx[1]:= 3;
end;
n := 0;
n6 := 0;
cnt := 0;
repeat
//only multibles of 6
inc(n,6);
inc(n6);
//no square factor of 2
if n AND 3 = 0 then
continue;
//no square factor of 3
if n MOD 9 = 0 then
continue;

if IsSquarefreeDecomp6(PrimeDecomp,n6)then
if ChkGiuga(n,@PrimeDecomp) then
begin
inc(cnt);
writeln(cnt:3,'..',OutPots(@PrimeDecomp,n),'  ',(GettickCount64-T0)/1000:6:3,' s');
end;
until n >= LMT;
T0 := GetTickCount64-T0;
writeln('Found ',cnt);
writeln('Tested til ',n,' runtime ',T0/1000:0:3,' s');
writeln;
writeln(OutPots(@PrimeDecomp,n));
end.
```
@home AMD 5600G ( 4,4 Ghz ) fpc3.2.2 -O4 -Xs:
```  1..30 :2*3*5_chk_30<   0.000 s
2..858 :2*3*11*13_chk_858<   0.000 s
3..1722 :2*3*7*41_chk_1722<   0.000 s
4..66198 :2*3*11*17*59_chk_66198<   0.000 s
5..2214408306 :2*3*11*23*31*47057_chk_2214408306<  17.120 s
6..24423128562 :2*3*7*43*3041*4447_chk_24423128562<  450.180 s
Found 6
Tested til 24423128562 runtime 450.180 s

24423128562 :2*3*7*43*3041*4447_chk_24423128562
TIO.RUN (~2.3 Ghz )takes ~4x runtime ? ( 2214408306 DIV 2 ) in 36 secs :-(
```

#### alternative version

Generating recursive squarefree numbers of ascending primes and check those numbers.
2*3 are set fixed. 2*3*5 followed 2*3*7 than 2*3*11. Therefor the results are unsorted.

```program Giuga;
{
30 = 2 * 3 * 5.
858 = 2 * 3 * 11 * 13.
1722 = 2 * 3 * 7 * 41.
66198 = 2 * 3 * 11 * 17 * 59.
2214408306 = 2 * 3 * 11 * 23 * 31 * 47057.
24423128562 = 2 * 3 * 7 * 43 * 3041 * 4447.
432749205173838 = 2 * 3 * 7 * 59 * 163 * 1381 * 775807.
14737133470010574 = 2 * 3 * 7 * 71 * 103 * 67213 * 713863.
550843391309130318 = 2 * 3 * 7 * 71 * 103 * 61559 * 29133437.
244197000982499715087866346 = 2 * 3 * 11 * 23 * 31 * 47137 * 28282147 * 3892535183.
554079914617070801288578559178 = 2 * 3 * 11 * 23 * 31 * 47059 * 2259696349 * 110725121051.
1910667181420507984555759916338506 = 2 * 3 * 7 * 43 * 1831 * 138683 * 2861051 * 1456230512169437. }
{\$IFDEF FPC}
{\$MODE DELPHI}  {\$OPTIMIZATION ON,ALL}  {\$COPERATORS ON}
{\$ELSE}
{\$APPTYPE CONSOLE}
{\$ENDIF}
uses
sysutils
{\$IFDEF WINDOWS},Windows{\$ENDIF}
;
const
LMT =14737133470010574;// 432749205173838;//24423128562;//2214408306;
type
tFac = packed record
fMul     :Uint64;
fPrime,
fPrimIdx,
fprimMaxIdx,dummy :Uint32;
dummy2: Uint64;
end;
tFacs  = array[0..15] of tFac;
tPrimes = array[0..62157] of Uint32;//775807 see factor of 432749205173838
//  tPrimes = array[0..4875{14379}] of Uint32;//sqrt 24423128562
//  tPrimes = array[0..1807414] of Uint32;//29133437
//  tPrimes = array[0..50847533] of Uint32;// 1e9
//  tPrimes = array[0..5761454] of Uint32;//1e8
var
SmallPrimes: tPrimes;
T0 : Int64;
cnt:Uint32;

procedure InitSmallPrimes;
//get primes. #0..65535.Sieving only odd numbers
const
//MAXLIMIT = (trunc(sqrt(LMT)+1)-1) shr 1+4;
MAXLIMIT = 775807 DIV 2+1;//(trunc(sqrt(LMT)+1)-1) shr 1+4;
var
pr : array of byte;
pPr :pByte;
p,j,d,flipflop :NativeUInt;
Begin
SmallPrimes[0] := 2;
setlength(pr,MAXLIMIT);
pPr := @pr[0];
p := 0;
repeat
repeat
p +=1
until pPr[p]= 0;
j := (p+1)*p*2;
if j>MAXLIMIT then
BREAK;
d := 2*p+1;
repeat
pPr[j] := 1;
j += d;
until j>MAXLIMIT;
until false;

SmallPrimes[1] := 3;
SmallPrimes[2] := 5;
j := 3;
d := 7;
flipflop := (2+1)-1;//7+2*2,11+2*1,13,17,19,23
p := 3;
repeat
if pPr[p] = 0 then
begin
SmallPrimes[j] := d;
inc(j);
end;
d += 2*flipflop;
p+=flipflop;
flipflop := 3-flipflop;
until (p > MAXLIMIT) OR (j>High(SmallPrimes));
setlength(pr,0);
end;

procedure OutFac(var F:tFacs;maxIdx:Uint32);
var
i : integer;
begin
write(cnt:3,'  ');
For i := 0 to maxIdx do
write(F[i].fPrime,'*');
write(#8,' = ',F[maxIdx].fMul);
writeln('      ',(GetTickCount64-T0)/1000:10:3,' s');
//readln;
end;

function ChkGiuga(var F:tFacs;MaxIdx:Uint32):boolean;inline;
var
n : Uint64;
idx: NativeInt;
p : Uint32;
begin
n := F[MaxIdx].fMul;
idx := MaxIdx;
repeat
p := F[idx].fPrime;
result := (((n DIV p)-1)MOD p) = 0;
if not(result) then
EXIT;
dec(idx);
until idx<0;
inc(cnt);
end;

procedure InsertNextPrimeFac(var F:tFacs;idx:Uint32);
var
Mul : Uint64;
i,p : uint32;
begin
with F[idx-1] do
begin
Mul := fMul;
i := fPrimIdx;
end;

while i<High(SmallPrimes) do
begin
inc(i);
with F[idx] do
begin
if i >fprimMaxIdx then
BREAK;
p := SmallPrimes[i];
if p*Mul>LMT then
BREAK;
fMul := p*Mul;
fPrime := p;
fPrimIdx := i;
IF (Mul-1) MOD p = 0 then
IF ChkGiuga(F,idx) then
OutFac(F,idx);
end;
// max 6 factors //for lmt 24e9 need 7 factors
if idx <5 then
InsertNextPrimeFac(F,idx+1);
end;
end;

var
{\$ALIGN 32}
Facs : tFacs;
i : integer;
Begin
InitSmallPrimes;

T0 := GetTickCount64;
with Facs[0] do
begin
fMul := 2;fPrime := 2;fPrimIdx:= 0;
end;
with Facs[1] do
begin
fMul := 2*3;fPrime := 3;fPrimIdx:= 1;
end;
i := 2;
//search the indices of mx found factor
while SmallPrimes[i] < 11 do inc(i); Facs[2].fprimMaxIdx := i;
while SmallPrimes[i] < 71 do inc(i); Facs[3].fprimMaxIdx := i;
while SmallPrimes[i] < 3041 do inc(i); Facs[4].fprimMaxIdx := i;
while SmallPrimes[i] < 67213 do inc(i); Facs[5].fprimMaxIdx := i;
while SmallPrimes[i] < 775807 do inc(i); Facs[6].fprimMaxIdx := i;
{
writeln('Found ',cnt,' in ',(GetTickCount64-T0)/1000:10:3,' s');
//start with 2*3*7
with Facs[2] do
begin
fMul := 2*3*7;fPrime := 7;fPrimIdx:= 3;
end;
InsertNextPrimeFac(Facs,3);
//start with 2*3*11
writeln('Found ',cnt,' in ',(GetTickCount64-T0)/1000:10:3,' s');
with Facs[2] do
begin
fMul := 2*3*11;fPrime := 11;fPrimIdx:= 4;
end;
InsertNextPrimeFac(Facs,3);
}
InsertNextPrimeFac(Facs,2);
writeln('Found ',cnt,' in ',(GetTickCount64-T0)/1000:10:3,' s');
writeln;
end.
```
@TIO.RUN:
```  1  2*3*5 = 30           0.000 s
2  2*3*7*41 = 1722           0.810 s
3  2*3*7*43*3041*4447 = 24423128562           0.871 s
4  2*3*11*13 = 858           1.057 s
5  2*3*11*17*59 = 66198           1.089 s
6  2*3*11*23*31*47057 = 2214408306           1.152 s
Found 6 in      1.526 s
Real time: 1.682 s CPU share: 99.42 %

@home:
Limit:432749205173838
start 2*3*7
Found 0 in      0.000 s
1  2*3*7*41 = 1722         147.206 s
2  2*3*7*43*3041*4447 = 24423128562         163.765 s
3  2*3*7*59*163*1381*775807 = 432749205173838         179.124 s
Found 3 in    197.002 s
start 2*3*11
4  2*3*11*13 = 858         197.002 s
5  2*3*11*17*59 = 66198         219.166 s
6  2*3*11*23*31*47057 = 2214408306         244.468 s

real  5m10,271s

Limit :14737133470010574
start 2*3*7
Found 0 in      0.000 s
1  2*3*7*41 = 1722        1330.819 s
2  2*3*7*43*3041*4447 = 24423128562        1567.028 s
3  2*3*7*59*163*1381*775807 = 432749205173838        1788.203 s
4  2*3*7*71*103*67213*713863 = 14737133470010574        2051.552 s
Found 4 in   2129.801 s
start 2*3*11
5  2*3*11*13 = 858        2129.801 s
6  2*3*11*17*59 = 66198        2305.752 s
7  2*3*11*23*31*47057 = 2214408306        2591.984 s
Found 7 in   3654.610 s

real  60m54,612s
```

## Perl

```#!/usr/bin/perl

use strict; # https://rosettacode.org/wiki/Giuga_numbers
use warnings;
use ntheory qw( factor forcomposites );
use List::Util qw( all );

forcomposites
{
my \$n = \$_;
all { (\$n / \$_ - 1) % \$_ == 0 } factor \$n and print "\$n\n";
} 4, 67000;
```
Output:
```30
858
1722
66198
```

## Phix

```with javascript_semantics
constant limit = 4
sequence giuga = {}
integer n = 4
while length(giuga)<limit do
sequence pf = prime_factors(n)
for f in pf do
if remainder(n/f-1,f) then pf={} exit end if
end for
if length(pf) then giuga &= n end if
n += 2
end while
printf(1,"The first %d Giuga numbers are: %v\n",{limit,giuga})
```
Output:
```The first 4 Giuga numbers are: {30,858,1722,66198}
```

### Faster version

Translation of: Wren
```--
-- demo\rosetta\Giuga_number.exw
-- =============================
--
with javascript_semantics
requires("1.0.2") -- (is_prime2 tweak)

procedure giuga(integer n)
printf(1,"n = %d:",n)
sequence p = repeat(0,n),
s = repeat(0,n)
p[1..2] = 1
s[1] = {1,2}
integer t = 2
while t>1 do
integer pt = p[t]+1
p[t] = pt
pt = get_prime(pt)
integer {c,d} = s[t-1]
{c,d} = {c*pt+d,d*pt}
s[t] = {c,d}
if c=d
or c*pt+(n-t)*d<=d*pt then
t -= 1
elsif t < (n - 2) then
t += 1
{c,d} = s[t-1]
p[t] = max(p[t-1], is_prime2(floor(c/(d-c)),true))
else
{c,d} = s[n-2]
integer dmc = d-c,
k = d*d-dmc
sequence f = factors(k,1)
for i=1 to floor(length(f)/2) do
integer h = f[i]
if remainder(h+d,dmc) == 0
and remainder(k/h+d,dmc) == 0 then
integer r1 = (h + d) / dmc,
r2 = (k/h + d) / dmc,
pn2 = get_prime(p[n-2])
if  r1 > pn2
and r2 > pn2
and r1 != r2
and is_prime(r1)
and is_prime(r2) then
printf(1," %d",d * r1 * r2)
end if
end if
end for
end if
end while
printf(1,"\n")
end procedure

papply({3,4,5,6},giuga)
```
Output:
```n = 3: 30
n = 4: 1722 858
n = 5: 66198
n = 6: 24423128562 2214408306
```

You can (almost) push things a little further on 64-bit:

```if machine_bits()=64 then giuga(7) end if
```

and get

```n = 7: 432749205173838 550843391309130318 14737133470010574
```

It took about 3 minutes for those to show, but then carried on in a doomed quest for an hour before I killed it.

## PL/M

Works with: 8080 PL/M Compiler

... under CP/M (or an emulator)

Only finds 3 Giuga numbers as the fourth is larger than 65535, the largest integer supported by the 8080 PL/M compiler.

```100H: /* FIND SOME GIUGA NUMBERS, COMPOSITES N SUCH THAT ALL THEIR DISTINCT */
/*                       PRIME FACTORS F EXACTLY DIVIDE ( N / F ) - 1 */

/* CP/M BDOS SYSTEM CALL AND I/O ROUTINES                                */
BDOS: PROCEDURE( FN, ARG ); DECLARE FN BYTE, ARG ADDRESS; GOTO 5; END;
PR\$CHAR:   PROCEDURE( C ); DECLARE C BYTE;    CALL BDOS( 2, C );  END;
PR\$STRING: PROCEDURE( S ); DECLARE S ADDRESS; CALL BDOS( 9, S );  END;
PR\$NL:     PROCEDURE;   CALL PR\$CHAR( 0DH ); CALL PR\$CHAR( 0AH ); END;
PR\$NUMBER: PROCEDURE( N ); /* PRINTS A NUMBER IN THE MINIMUN FIELD WIDTH  */
DECLARE N ADDRESS;
DECLARE V ADDRESS, N\$STR ( 6 )BYTE, W BYTE;
V = N;
W = LAST( N\$STR );
N\$STR( W ) = '\$';
N\$STR( W := W - 1 ) = '0' + ( V MOD 10 );
DO WHILE( ( V := V / 10 ) > 0 );
N\$STR( W := W - 1 ) = '0' + ( V MOD 10 );
END;
CALL PR\$STRING( .N\$STR( W ) );
END PR\$NUMBER;

/* TASK                                                                  */

/* FIND THE FIRST THREE GIUGA NUMBERS (THE FOURTH IS > 65535)            */
/* EACH PRIME FACTOR CAN ONLY APPEAR ONCE, E.G.,: FOR 2:                 */
/* (( N / 2 ) - 1) MOD 2 = 0  => N / 2 IS ODD => N NOT DIVISIBLE BY 4    */
/* SIMILARLY FOR OTHER PRIMES                                            */
DECLARE ( N, V, FCOUNT, F ) ADDRESS;
DECLARE IS\$GIUGA            BYTE;
N     = 2;
DO WHILE N < 65000;                   /* ASSUME THE NUMEBRS ARE ALL EVEN */
V        = N / 2;
IS\$GIUGA = 1;
FCOUNT   = 1;
F        = 1;
DO WHILE ( F := F + 2 ) <= V AND IS\$GIUGA;
IF V MOD F = 0 THEN DO;
/* HAVE A PRIME FACTOR                                          */
FCOUNT   = FCOUNT + 1;
IS\$GIUGA = ( ( N / F ) - 1 ) MOD F = 0;
V = V / F;
END;
END;
IF IS\$GIUGA THEN DO;
IF FCOUNT > 1 THEN DO;
/* N IS NOT PRIME, SO IS GIUGA                                  */
CALL PR\$CHAR( ' ' );CALL PR\$NUMBER( N );
END;
END;
N = N + 4;
END;

EOF```
Output:
``` 30 858 1722
```

## Python

Translation of: FreeBASIC
```#!/usr/bin/python

from math import sqrt

def isGiuga(m):
n = m
f = 2
l = sqrt(n)
while True:
if n % f == 0:
if ((m / f) - 1) % f != 0:
return False
n /= f
if f > n:
return True
else:
f += 1
if f > l:
return False

if __name__ == '__main__':
n = 3
c = 0
print("The first 4 Giuga numbers are: ")
while c < 4:
if isGiuga(n):
c += 1
print(n)
n += 1
```

## Quackery

`dpfs` is Distinct Prime Factors.

``` [ [] swap
dup times
[ dup i^ 2 + /mod
0 = if
[ nip dip
[ i^ 2 + join ]
[ dup i^ 2 + /mod
0 = iff
nip again ] ]
drop
dup 1 = if conclude ]
drop ]                    is dpfs  ( n --> [ )

[ dup dpfs
dup size 2 < iff
[ 2drop false ]
done
true unrot
witheach
[ 2dup / 1 -
swap mod 0 != if
[ dip not
conclude ] ]
drop ]                    is giuga ( n --> b )

[] 0
[ 1+ dup giuga if
[ tuck join swap ]
over size 4 = until ]
drop
echo```
Output:
`[ 30 858 1722 66198 ]`

## Raku

```my @primes = (3..60).grep: &is-prime;

print 'First four Giuga numbers: ';

put sort flat (2..4).map: -> \$c {
@primes.combinations(\$c).map: {
my \$n = [×] 2,|\$_;
\$n if all .map: { (\$n / \$_ - 1) %% \$_ };
}
}
```
Output:
`First 4 Giuga numbers: 30 858 1722 66198`

## Ring

```see "working..." + nl
see "The first 4 Giuga numbers are:" + nl
load "stdlibcore.ring"

Comp = []
num = 0
n = 1
while true
n++
if not isPrime(n)
Comp = []
for p = 1 to n
if isPrime(p) AND (n % p = 0)
add(Comp,p)
ok
next
flag = 1
for ind = 1 to len(Comp)
f = Comp[ind]
res = (n/f)- 1
if res % f != 0
flag = 0
exit
ok
next
if flag = 1
see "" + n + " "
num++
ok
if num = 4
exit
ok
ok
end
see nl + "done..." + nl```
Output:
```working...
The first 4 Giuga numbers are:
30 858 1722 66198
done...
```

## RPL

Works with: HP version 49
```≪ DUP → n
≪ FACTORS 1 SF
1 OVER SIZE FOR j
DUP j GET
n OVER / 1 -
IF SWAP MOD THEN 1 CF END
2 STEP DROP
1 FS?
≫ ≫ 'GIUGA?' STO

≪ { } 4
WHILE OVER SIZE 4 < REPEAT
IF DUP GIUGA? THEN LOOK SWAP OVER + SWAP END
2 +
END DROP
≫ 'TASK' STO
```
Output:
```1: {30 858 1722 66198}
```

## Ruby

```require 'prime'

giuga = (1..).lazy.select do |n|
pd = n.prime_division
pd.sum{|_, d| d} > 1 &&  #composite
pd.all?{|f, _| (n/f - 1) % f == 0}
end

p giuga.take(4).to_a
```
Output:
```[30, 858, 1722, 66198]
```

## Rust

```use prime_tools ;

fn prime_decomposition( mut number : u32) -> Vec<u32> {
let mut divisors : Vec<u32> = Vec::new( ) ;
let mut divisor : u32 = 2 ;
while number != 1 {
if number % divisor == 0 {
divisors.push( divisor ) ;
number /= divisor ;
}
else {
divisor += 1 ;
}
}
divisors
}

fn is_giuga( num : u32 ) -> bool {
let prime_factors : Vec<u32> = prime_decomposition( num ) ;
! prime_tools::is_u32_prime( num ) &&
prime_factors.into_iter( ).all( |n : u32| (num/n -1) % n == 0 )
}

fn main() {
let mut giuga_numbers : Vec<u32> = Vec::new( ) ;
let mut num : u32 = 2 ;
while giuga_numbers.len( ) != 4 {
if is_giuga( num ) {
giuga_numbers.push( num ) ;
}
num += 1 ;
}
println!("{:?}" , giuga_numbers ) ;
}
```
Output:
```[30, 858, 1722, 66198]
```

## Scala

Translation of: Java
```object GiugaNumbers {

private var results: List[Int] = List()
def main(args: Array[String]): Unit = {
val primes: List[Int] = List(2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59)

val primeCounts: List[Int] = List(3, 4, 5)
for (primeCount <- primeCounts) {
var primeFactors: List[Int] = List.fill(primeCount)(0)
combinations(primes, primeCount, 0, 0, primeFactors)
}

val sortedResults = results.sorted
println(s"Found Giuga numbers: \$sortedResults")
}

private def checkIfGiugaNumber(primeFactors: List[Int]): Unit = {
val product: Int = primeFactors.reduce(Math.multiplyExact)
for (factor <- primeFactors) {
val divisor: Int = factor * factor
if ((product - factor) % divisor != 0) {
return
}
}
results :+= product
}

private def combinations(primes: List[Int], primeCount: Int, index: Int, level: Int, primeFactors: List[Int]): Unit = {
if (level == primeCount) {
checkIfGiugaNumber(primeFactors)
return
}

for (i <- index until primes.length) {
val updatedPrimeFactors = primeFactors.updated(level, primes(i))
combinations(primes, primeCount, i + 1, level + 1, updatedPrimeFactors)
}
}
}
```
Output:
```Found Giuga numbers: List(30, 858, 1722, 66198)

```

## Sidef

```4..Inf `by` 2 -> lazy.grep {|n|
n.factor.all {|p| p `divides` (n/p - 1) }
}.first(4).say
```
Output:
```[30, 858, 1722, 66198]
```

## Wren

### Version 1 (Brute force)

Simple brute force but assumes all Giuga numbers will be even, must be square-free and can't be semi-prime.

Takes only about 0.05 seconds to find the first four Giuga numbers but finding the fifth would take many hours using this approach, so I haven't bothered.

```var factors = []
var inc = [4, 2, 4, 2, 4, 6, 2, 6]

// Assumes n is even with exactly one factor of 2.
// Empties 'factors' if any other prime factor is repeated.
var primeFactors = Fn.new { |n|
factors.clear()
var last = 2
factors.add(2)
n = (n/2).truncate
while (n%3 == 0) {
if (last == 3) {
factors.clear()
return
}
last = 3
factors.add(3)
n = (n/3).truncate
}
while (n%5 == 0) {
if (last == 5) {
factors.clear()
return
}
last = 5
factors.add(5)
n = (n/5).truncate
}
var k = 7
var i = 0
while (k * k <= n) {
if (n%k == 0) {
if (last == k) {
factors.clear()
return
}
last = k
factors.add(k)
n = (n/k).truncate
} else {
k = k + inc[i]
i = (i + 1) % 8
}
}
if (n > 1) factors.add(n)
}

var limit = 4
var giuga = []
var n = 6 // can't be 2 or 4
while (giuga.count < limit) {
primeFactors.call(n)
// can't be prime or semi-prime
if (factors.count > 2 && factors.all { |f| (n/f - 1) % f == 0 }) {
giuga.add(n)
}
n = n + 4 // can't be divisible by 4
}
System.print("The first %(limit) Giuga numbers are:")
System.print(giuga)
```
Output:
```The first 4 Giuga numbers are:
[30, 858, 1722, 66198]
```

### Version 2 (Pari-GP translation)

Library: Wren-math
Library: Wren-rat

This is a translation of the very fast Pari-GP code in the talk page. Only takes 0.015 seconds to find the first six Giuga numbers.

```import "./math" for Math, Int
import "./rat" for Rat

var giuga = Fn.new { |n|
System.print("n = %(n):")
var p = List.filled(n, 0)
var s = List.filled(n, null)
for (i in 0..n-2) s[i] = Rat.zero
p[2] = 2
p[1] = 2
var t = 2
s[1] = Rat.half
while (t > 1) {
p[t] = Int.isPrime(p[t] + 1) ? p[t] + 1 : Int.nextPrime(p[t] + 1)
s[t] = s[t-1] + Rat.new(1, p[t])
if (s[t] == Rat.one || s[t] + Rat.new(n - t, p[t]) <= Rat.one) {
t = t - 1
} else if (t < n - 2) {
t = t + 1
p[t] = Math.max(p[t-1], (s[t-1] / (Rat.one - s[t-1])).toFloat).floor
} else {
var c = s[n-2].num
var d = s[n-2].den
var k = d * d + c - d
var f = Int.divisors(k)
for (i in 0...((f.count + 1)/2).floor) {
var h = f[i]
if ((h + d) % (d-c) == 0 && (k/h + d) % (d - c) == 0) {
var r1 = (h + d) / (d - c)
var r2 = (k/h + d) / (d - c)
if (r1 > p[n-2] && r2 > p[n-2] && r1 != r2 && Int.isPrime(r1) && Int.isPrime(r2)) {
var w = d * r1 * r2
System.print(w)
}
}
}
}
}
}

for (n in 3..6) {
giuga.call(n)
System.print()
}
```
Output:
```n = 3:
30

n = 4:
1722
858

n = 5:
66198

n = 6:
24423128562
2214408306
```

## XPL0

```func Giuga(N0);         \Return 'true' if Giuga number
int  N0;
int  N, F, Q1, Q2, L;
[N:= N0;  F:= 2;  L:= sqrt(N);
loop    [Q1:= N/F;
if rem(0) = 0 then      \found a prime factor
[Q2:= N0/F;
if rem((Q2-1)/F) # 0 then return false;
N:= Q1;
if F>N then quit;
]
else    [F:= F+1;
if F>L then return false;
];
];
return true;
];

int N, C;
[N:= 3;  C:= 0;
loop    [if Giuga(N) then
[IntOut(0, N);  ChOut(0, ^ );
C:= C+1;
if C >= 4 then quit;
];
N:= N+1;
];
]```
Output:
```30 858 1722 66198
```