Rosetta Code

Giuga numbers: Difference between revisions

Page Discussion

Giuga numbers (view source)

Revision as of 08:55, 24 April 2024

58,779 bytes added ,  29 days ago
→‎{{header|jq}}
 
(59 intermediate revisions by 24 users not shown)
Line 1:
{{draft task}}
 
;Definition
Line 26:
 
 
=={{header|FreeBASIC11l}}==
{{trans|Python}}
<lang freebasic>Function isGiuga(m As Uinteger) As Boolean
 
Dim As Uinteger n = m
<syntaxhighlight lang="11l">
Dim As Uinteger f = 2, l = Sqr(n)
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++
</syntaxhighlight>
 
{{out}}
<pre>
The first 4 Giuga numbers are:
30
858
1722
66198
</pre>
 
=={{header|ABC}}==
Based on the Algol 68 sample,
<syntaxhighlight lang="abc">
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
</syntaxhighlight>
{{out}}
<pre>
30 858 1722 66198
</pre>
 
=={{header|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.
<syntaxhighlight lang="algol68">
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
</syntaxhighlight>
{{out}}
<pre>
30 858 1722 66198
</pre>
 
=={{header|ALGOL W}}==
{{Trans|ALGOL 68}}
<syntaxhighlight lang="algolw">
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.
</syntaxhighlight>
{{out}}
Same as Algol 68
 
=={{header|Arturo}}==
<syntaxhighlight lang="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?</syntaxhighlight>
 
{{out}}
 
<pre>30
858
1722
66198</pre>
 
=={{header|AWK}}==
<syntaxhighlight lang="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) }
}
}
}
</syntaxhighlight>
{{out}}
<pre>
Giuga numbers 1-4: 30 858 1722 66198
</pre>
 
=={{header|BASIC}}==
==={{header|BASIC256}}===
<syntaxhighlight lang="vb">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</syntaxhighlight>
{{out}}
<pre>Same as FreeBASIC entry.</pre>
 
==={{header|FreeBASIC}}===
<syntaxhighlight lang="vbnet">Function isGiuga(m As Uinteger) As Boolean
Dim As Uinteger n = m, f = 2, l = Sqr(n)
Do
If n Mod f = 0 Then
Line 47 ⟶ 284:
If isGiuga(n) Then c += 1: Print n; " ";
n += 1
Loop Until c = limit</langsyntaxhighlight>
{{out}}
<pre>The first 4 Giuga numbers are: 30 858 1722 66198</pre>
 
==={{header|Gambas}}===
<syntaxhighlight lang="vbnet">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</syntaxhighlight>
{{out}}
<pre>Same as FreeBASIC entry.</pre>
 
==={{header|PureBasic}}===
<syntaxhighlight lang="vb">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()</syntaxhighlight>
{{out}}
<pre>Same as FreeBASIC entry.</pre>
 
=={{header|C++}}==
===Brute force===
Based on the Go solution. Takes 26 minutes on my system (Intel Core i5 3.2GHz).
<syntaxhighlight lang="cpp">#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;
}
}
}</syntaxhighlight>
 
{{out}}
<pre>
First 5 Giuga numbers:
30
858
1722
66198
2214408306
</pre>
 
===Faster version===
{{libheader|Boost}}
{{trans|Wren}}
<syntaxhighlight lang="cpp">#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);
}</syntaxhighlight>
 
{{out}}
<pre>
n = 3: 30
n = 4: 1722 858
n = 5: 66198
n = 6: 24423128562 2214408306
</pre>
 
=={{header|Delphi}}==
{{works with|Delphi|6.0}}
{{libheader|SysUtils,StdCtrls}}
 
 
<syntaxhighlight lang="Delphi">
 
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;
 
 
 
</syntaxhighlight>
{{out}}
<pre>
30
858
1722
66198
 
Elapsed Time: 4.991 Sec.
 
</pre>
 
 
=={{header|EasyLang}}==
{{trans|Python}}
<syntaxhighlight lang="easylang">
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
.
</syntaxhighlight>
 
=={{header|Euler}}==
Uses the C style for loop procedure from the [[Sieve of Eratosthenes]] task
'''begin'''
'''new''' for; '''new''' n; '''new''' gCount;
for &lt;- ` '''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'''
&apos;
;
gCount &lt;- 0;
for( ` n &lt;- 2 &apos;, ` gCount &lt; 4 &apos;, ` n &lt;- n + 4 &apos;
, ` '''begin'''
'''new''' v; '''new''' f; '''new''' isGiuga; '''new''' fCount;
v &lt;- n % 2;
isGiuga &lt;- '''true''';
fCount &lt;- 1;
for( ` f &lt;- 3 &apos;, ` f &lt;= v '''and''' isGiuga &apos;, ` f &lt;- f + 2 &apos;
, ` '''if''' v '''mod''' f = 0 '''then''' '''begin'''
fCount &lt;- fCount + 1;
isGiuga &lt;- [ [ n % f ] - 1 ] '''mod''' f = 0;
v &lt;- v % f
'''end''' '''else''' 0
&apos;
);
'''if''' isGiuga '''then''' '''begin'''
'''if''' fCount &gt; 1 '''then''' '''begin'''
gCount &lt;- gCount + 1;
'''out''' n
'''end''' '''else''' 0
'''end''' '''else''' 0
'''end'''
&apos;
)
'''end''' $
 
=={{header|Go}}==
{{trans|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.
{{libheader|Go-rcu}}
<syntaxhighlight lang="go">package main
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 around 2 hours 40 minutes.
<lang go>package main
 
import ("fmt"
 
"fmt"
var factors []int
"rcu"
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
for// n :=can't 4;be len(giuga)2 <or limit;divisible nby += 2 {4
for n := 6; len(giuga) < limit; n += 4 {
factors := rcu.PrimeFactors(n)
primeFactors(n)
// can't be prime or semi-prime
if len(factors) > 2 {
isSquareFreeisGiuga := true
for i_, f := 1;range i < len(factors); i++ {
if factors[i](n/f-1)%f =!= factors[i-1]0 {
isSquareFreeisGiuga = false
break
}
}
if isSquareFreeisGiuga {
isGiugagiuga := trueappend(giuga, n)
for _, f := range factors {
if (n/f-1)%f != 0 {
isGiuga = false
break
}
}
if isGiuga {
giuga = append(giuga, n)
}
}
}
Line 91 ⟶ 693:
fmt.Println("The first", limit, "Giuga numbers are:")
fmt.Println(giuga)
}</langsyntaxhighlight>
 
{{out}}
Line 97 ⟶ 699:
The first 5 Giuga numbers are:
[30 858 1722 66198 2214408306]
</pre>
 
=={{header|Haskell}}==
<syntaxhighlight lang="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..]</syntaxhighlight>
{{out}}
<pre>
[30,858,1722,66198]
</pre>
 
Line 103 ⟶ 732:
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):
 
<langsyntaxhighlight Jlang="j">giguaP=: {{ (1<y)*(-.1 p:y)**/(=<.) y ((_1+%)%]) q: y }}"0</langsyntaxhighlight>
 
<langsyntaxhighlight Jlang="j"> 1+I.giguaP 1+i.1e5
30 858 1722 66198</langsyntaxhighlight>
 
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:
 
<syntaxhighlight lang="j">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
}}</syntaxhighlight>
 
Example use:<syntaxhighlight lang="j"> giuga 1
 
giuga 2
 
giuga 3
30
giuga 4
1722 858
giuga 5
66198
giuga 6
24423128562 2214408306</syntaxhighlight>
 
=={{header|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.
 
<syntaxhighlight lang="java">
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>();
 
}
</syntaxhighlight>
{{ out }}
<pre>
Found Giuga numbers: [30, 858, 1722, 66198]
</pre>
 
=={{header|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.
<syntaxhighlight lang="jq">
# `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))
</syntaxhighlight>
{{output}}
<pre>
30
858
1722
66198
</pre>
 
=={{header|Julia}}==
<langsyntaxhighlight rubylang="julia">using Primes
 
isGiuga(n) = all(f -> f != n && rem(n ÷ f - 1, f) == 0, factor(Vector, n))
Line 126 ⟶ 914:
 
getGiuga(4)
</langsyntaxhighlight>{{out}}
<pre>
30
Line 134 ⟶ 922:
</pre>
=== Ad hoc faster version ===
<langsyntaxhighlight rubylang="julia">using Primes
 
isGiuga(n) = all(f -> f != n && rem(n ÷ f - 1, f) == 0, factor(Vector, n))
 
function getgiugas(numberwanted, verbose = true)
Line 144 ⟶ 930:
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, ")")
Line 158 ⟶ 944:
@time getgiugas(2, false)
@time getgiugas(6)
</langsyntaxhighlight>{{out}}
<pre>
30 (elapsed: 0.0)
Line 168 ⟶ 954:
432.066249 seconds (235 allocations: 12.523 KiB)
</pre>
 
=={{header|Lua}}==
{{Trans|ALGOL 68}}
<syntaxhighlight lang="lua">
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</syntaxhighlight>
{{out}}
<pre>
30 858 1722 66198
</pre>
 
=={{header|Mathematica}}/{{header|Wolfram Language}}==
{{trans|Julia}}
<syntaxhighlight lang="Mathematica">
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]]
</syntaxhighlight>
{{out}}
<pre>
{30, 858, 1722, 66198}
 
</pre>
 
 
=={{header|Maxima}}==
<syntaxhighlight lang="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);
</syntaxhighlight>
{{out}}
<pre>
[30,858,1722,66198]
</pre>
 
=={{header|Nim}}==
{{trans|FreeBASIC}}
<syntaxhighlight lang="Nim">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()
</syntaxhighlight>
 
{{out}}
<pre>The first 4 Giuga numbers are: 30 858 1722 66198
</pre>
 
=={{header|PARI/GP}}==
{{trans|Julia}}
<syntaxhighlight lang="PARI/GP">
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))
</syntaxhighlight>
{{out}}
<pre>
[30, 858, 1722, 66198]
 
</pre>
 
 
=={{header|Pascal}}==
==={{header|Free Pascal}}===
OK.Cheating to find square free numbers like julia in distance 6<BR>
changing main routine at the end of [[Factors_of_an_integer#using_Prime_decomposition]] <BR>
That means always factors 2,3 and minimum one of 5,7,11.<br>
Mostly, about 70%, the highest primefactor remains in pfRemain.<br>
<br>
<lang pascal>
<syntaxhighlight lang="pascal">program Giuga;
function ChkGuiga(n:Uint64;pPrimeDecomp :tpPrimeFac):boolean;inline;
 
{$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}
pFr : Uint64;
SmallPrimes: tPrimes;
idx: NativeInt;
{$ALIGN 32}
p : Uint32;
 
begin
procedure InitSmallPrimes;
with pPrimeDecomp^ do
//get primes. #0..65535.Sieving only odd numbers
Begin
const
pfR := pfRemain;
MAXLIMIT = (821641-1) shr 1;
idx := pfMaxIdx;//always>0 for pfDivcnt>4
var
//check squarefree. i primefactors -> count of diviors = 2^i
ppr := 1array[0..MAXLIMIT] shlof idxbyte;
p,j,d,flipflop :NativeUInt;
if pfR <> 1 then
Begin
begin
SmallPrimes[0] := 2;
if (p+p <> pfdivcnt) then
fillchar(pr[0],SizeOf(pr),#0);
EXIT(false);
p if sqr(pfR)>:=n then0;
repeat
EXIT(false);
//squarefreerepeat
p +=1
result := (((n DIV pfR)-1)MOD pfR) = 0;
ifuntil resultpr[p]= then0;
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
idxSmallPrimes[j] -:= 1d;
repeatinc(j);
p := SmallPrimes[pfpotPrimIdx[idx]];
result := (((n DIV p)-1)MOD p) = 0;
if Not(result) then
EXIT(result);
dec(idx);
until idx <0;
end;
d += 2*flipflop;
end
p+=flipflop;
else
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
ifchk (p:= <> pfdivcnt) then1;
For n := 0 to pfMaxIdx-1 do
EXIT(false);
result := true;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
dec(p := pfpotPrimIdx[idx)];
p := SmallPrimes[pfpotPrimIdx[idx]];
result := (((n DIV p)-1)MOD p) = 0;
if not(result) then
EXIT(false);
until dec(idx<=0);
until idx<0;
end;
end;
end;
 
const
LMT = 24423128562;//24423128562;//2214408306;//
var
PrimeDecomp :tPrimeFac;
pPrimeDecomp :tpPrimeFac;
T0:Int64;
n,n6 : NativeUIntUInt64;
cnt:Uint32;
Begin
Line 233 ⟶ 1,306:
 
T0 := GetTickCount64;
with PrimeDecomp do
Init_Sieve(1);
begin
pfpotPrimIdx[0]:= 2;
pfpotPrimIdx[1]:= 3;
end;
n := 0;
n6 := 0;
cnt := 0;
repeat
//only evenmultibles numbersof 6
inc(n,26);
GetNextPrimeDecompinc(n6);
//no square factor of 2
pPrimeDecomp:= GetNextPrimeDecomp;
//if notn prime/semiprimeAND 3 = 0 then
continue;
if pPrimeDecomp^.pfDivCnt >4 then
//no square factor of 3
if ChkGuiga(n,pPrimeDecomp) then
if n MOD 9 = 0 then
continue;
 
if IsSquarefreeDecomp6(PrimeDecomp,n6)then
if ChkGiuga(n,@PrimeDecomp) then
begin
inc(cnt);
writeln(cnt:3,'..',OutPots(pPrimeDecomp@PrimeDecomp,n),' ',(GettickCount64-T0)/1000:6:3,' s');
end;
until n >= LMT;
Line 253 ⟶ 1,336:
writeln('Tested til ',n,' runtime ',T0/1000:0:3,' s');
writeln;
writeln(OutPots(pPrimeDecomp@PrimeDecomp,n));
end.</langsyntaxhighlight>
{{out|@home AMD 5600G ( 4,4 Ghz ) fpc3.2.2 -O4 -Xs}}
<pre>
1..30 : 8 : 2*3*5_chk_30<_SoD_72< 0.001000 s
2..858 : 16 : 2*3*11*13_chk_858<_SoD_2016< 0.001000 s
3..1722 : 16 : 2*3*7*41_chk_1722<_SoD_4032< 0.001000 s
4..66198 : 32 : 2*3*11*17*59_chk_66198<_SoD_155520< 0.002000 s
5..2214408306 : 64 : 2*3*11*23*31*47057_chk_2214408306<_SoD_5204238336< 4117.656120 s
6..24423128562 : 64 : 2*3*7*43*3041*4447_chk_24423128562<_SoD_57154166784< 533450.047180 s
Found 6
Tested til 24423128562 runtime 533450.047180 s</pre>
 
24423128562 :2*3*7*43*3041*4447_chk_24423128562
TIO.RUN (~2.3 Ghz )takes ~4x runtime ? ( 2214408306 DIV 2 ) in 36 secs :-(
</pre>
====alternative version====
Generating recursive squarefree numbers of ascending primes and check those numbers.<BR>
2*3 are set fixed. 2*3*5 followed 2*3*7 than 2*3*11. Therefor the results are unsorted.
<syntaxhighlight lang="pascal">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.
</syntaxhighlight>
{{out|@TIO.RUN}}
<pre>
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
</pre>
 
=={{header|Perl}}==
<langsyntaxhighlight lang="perl">#!/usr/bin/perl
 
use strict; # https://rosettacode.org/wiki/Giuga_numbers
Line 278 ⟶ 1,609:
my $n = $_;
all { ($n / $_ - 1) % $_ == 0 } factor $n and print "$n\n";
} 4, 67000;</langsyntaxhighlight>
{{out}}
<pre>
Line 288 ⟶ 1,619:
 
=={{header|Phix}}==
<!--<langsyntaxhighlight Phixlang="phix">(phixonline)-->
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<span style="color: #008080;">constant</span> <span style="color: #000000;">limit</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">4</span>
Line 302 ⟶ 1,633:
<span style="color: #008080;">end</span> <span style="color: #008080;">while</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"The first %d Giuga numbers are: %v\n"</span><span style="color: #0000FF;">,{</span><span style="color: #000000;">limit</span><span style="color: #0000FF;">,</span><span style="color: #000000;">giuga</span><span style="color: #0000FF;">})</span>
<!--</langsyntaxhighlight>-->
{{out}}
<pre>
The first 4 Giuga numbers are: {30,858,1722,66198}
</pre>
===Faster version===
{{trans|Wren}}
<!--<syntaxhighlight lang="phix">(phixonline)-->
<span style="color: #000080;font-style:italic;">--
-- demo\rosetta\Giuga_number.exw
-- =============================
--</span>
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<span style="color: #7060A8;">requires</span><span style="color: #0000FF;">(</span><span style="color: #008000;">"1.0.2"</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">-- (is_prime2 tweak)</span>
<span style="color: #008080;">procedure</span> <span style="color: #000000;">giuga</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"n = %d:"</span><span style="color: #0000FF;">,</span><span style="color: #000000;">n</span><span style="color: #0000FF;">)</span>
<span style="color: #004080;">sequence</span> <span style="color: #000000;">p</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span><span style="color: #000000;">n</span><span style="color: #0000FF;">),</span>
<span style="color: #000000;">s</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span><span style="color: #000000;">n</span><span style="color: #0000FF;">)</span>
<span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">1</span><span style="color: #0000FF;">..</span><span style="color: #000000;">2</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1</span>
<span style="color: #000000;">s</span><span style="color: #0000FF;">[</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #000000;">2</span><span style="color: #0000FF;">}</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">t</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">2</span>
<span style="color: #008080;">while</span> <span style="color: #000000;">t</span><span style="color: #0000FF;">></span><span style="color: #000000;">1</span> <span style="color: #008080;">do</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">pt</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">t</span><span style="color: #0000FF;">]+</span><span style="color: #000000;">1</span>
<span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">t</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">pt</span>
<span style="color: #000000;">pt</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">get_prime</span><span style="color: #0000FF;">(</span><span style="color: #000000;">pt</span><span style="color: #0000FF;">)</span>
<span style="color: #004080;">integer</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">c</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">s</span><span style="color: #0000FF;">[</span><span style="color: #000000;">t</span><span style="color: #0000FF;">-</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]</span>
<span style="color: #0000FF;">{</span><span style="color: #000000;">c</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">c</span><span style="color: #0000FF;">*</span><span style="color: #000000;">pt</span><span style="color: #0000FF;">+</span><span style="color: #000000;">d</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">*</span><span style="color: #000000;">pt</span><span style="color: #0000FF;">}</span>
<span style="color: #000000;">s</span><span style="color: #0000FF;">[</span><span style="color: #000000;">t</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">c</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">}</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">c</span><span style="color: #0000FF;">=</span><span style="color: #000000;">d</span>
<span style="color: #008080;">or</span> <span style="color: #000000;">c</span><span style="color: #0000FF;">*</span><span style="color: #000000;">pt</span><span style="color: #0000FF;">+(</span><span style="color: #000000;">n</span><span style="color: #0000FF;">-</span><span style="color: #000000;">t</span><span style="color: #0000FF;">)*</span><span style="color: #000000;">d</span><span style="color: #0000FF;"><=</span><span style="color: #000000;">d</span><span style="color: #0000FF;">*</span><span style="color: #000000;">pt</span> <span style="color: #008080;">then</span>
<span style="color: #000000;">t</span> <span style="color: #0000FF;">-=</span> <span style="color: #000000;">1</span>
<span style="color: #008080;">elsif</span> <span style="color: #000000;">t</span> <span style="color: #0000FF;"><</span> <span style="color: #0000FF;">(</span><span style="color: #000000;">n</span> <span style="color: #0000FF;">-</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">then</span>
<span style="color: #000000;">t</span> <span style="color: #0000FF;">+=</span> <span style="color: #000000;">1</span>
<span style="color: #0000FF;">{</span><span style="color: #000000;">c</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">s</span><span style="color: #0000FF;">[</span><span style="color: #000000;">t</span><span style="color: #0000FF;">-</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]</span>
<span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">t</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">max</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">t</span><span style="color: #0000FF;">-</span><span style="color: #000000;">1</span><span style="color: #0000FF;">],</span> <span style="color: #7060A8;">is_prime2</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">floor</span><span style="color: #0000FF;">(</span><span style="color: #000000;">c</span><span style="color: #0000FF;">/(</span><span style="color: #000000;">d</span><span style="color: #0000FF;">-</span><span style="color: #000000;">c</span><span style="color: #0000FF;">)),</span><span style="color: #004600;">true</span><span style="color: #0000FF;">))</span>
<span style="color: #008080;">else</span>
<span style="color: #0000FF;">{</span><span style="color: #000000;">c</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">s</span><span style="color: #0000FF;">[</span><span style="color: #000000;">n</span><span style="color: #0000FF;">-</span><span style="color: #000000;">2</span><span style="color: #0000FF;">]</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">dmc</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">d</span><span style="color: #0000FF;">-</span><span style="color: #000000;">c</span><span style="color: #0000FF;">,</span>
<span style="color: #000000;">k</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">d</span><span style="color: #0000FF;">*</span><span style="color: #000000;">d</span><span style="color: #0000FF;">-</span><span style="color: #000000;">dmc</span>
<span style="color: #004080;">sequence</span> <span style="color: #000000;">f</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">factors</span><span style="color: #0000FF;">(</span><span style="color: #000000;">k</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">floor</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">f</span><span style="color: #0000FF;">)/</span><span style="color: #000000;">2</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">h</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">f</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]</span>
<span style="color: #008080;">if</span> <span style="color: #7060A8;">remainder</span><span style="color: #0000FF;">(</span><span style="color: #000000;">h</span><span style="color: #0000FF;">+</span><span style="color: #000000;">d</span><span style="color: #0000FF;">,</span><span style="color: #000000;">dmc</span><span style="color: #0000FF;">)</span> <span style="color: #0000FF;">==</span> <span style="color: #000000;">0</span>
<span style="color: #008080;">and</span> <span style="color: #7060A8;">remainder</span><span style="color: #0000FF;">(</span><span style="color: #000000;">k</span><span style="color: #0000FF;">/</span><span style="color: #000000;">h</span><span style="color: #0000FF;">+</span><span style="color: #000000;">d</span><span style="color: #0000FF;">,</span><span style="color: #000000;">dmc</span><span style="color: #0000FF;">)</span> <span style="color: #0000FF;">==</span> <span style="color: #000000;">0</span> <span style="color: #008080;">then</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">r1</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">(</span><span style="color: #000000;">h</span> <span style="color: #0000FF;">+</span> <span style="color: #000000;">d</span><span style="color: #0000FF;">)</span> <span style="color: #0000FF;">/</span> <span style="color: #000000;">dmc</span><span style="color: #0000FF;">,</span>
<span style="color: #000000;">r2</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">(</span><span style="color: #000000;">k</span><span style="color: #0000FF;">/</span><span style="color: #000000;">h</span> <span style="color: #0000FF;">+</span> <span style="color: #000000;">d</span><span style="color: #0000FF;">)</span> <span style="color: #0000FF;">/</span> <span style="color: #000000;">dmc</span><span style="color: #0000FF;">,</span>
<span style="color: #000000;">pn2</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">get_prime</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">n</span><span style="color: #0000FF;">-</span><span style="color: #000000;">2</span><span style="color: #0000FF;">])</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">r1</span> <span style="color: #0000FF;">></span> <span style="color: #000000;">pn2</span>
<span style="color: #008080;">and</span> <span style="color: #000000;">r2</span> <span style="color: #0000FF;">></span> <span style="color: #000000;">pn2</span>
<span style="color: #008080;">and</span> <span style="color: #000000;">r1</span> <span style="color: #0000FF;">!=</span> <span style="color: #000000;">r2</span>
<span style="color: #008080;">and</span> <span style="color: #7060A8;">is_prime</span><span style="color: #0000FF;">(</span><span style="color: #000000;">r1</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">and</span> <span style="color: #7060A8;">is_prime</span><span style="color: #0000FF;">(</span><span style="color: #000000;">r2</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">then</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">" %d"</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span> <span style="color: #0000FF;">*</span> <span style="color: #000000;">r1</span> <span style="color: #0000FF;">*</span> <span style="color: #000000;">r2</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">while</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"\n"</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">procedure</span>
<span style="color: #7060A8;">papply</span><span style="color: #0000FF;">({</span><span style="color: #000000;">3</span><span style="color: #0000FF;">,</span><span style="color: #000000;">4</span><span style="color: #0000FF;">,</span><span style="color: #000000;">5</span><span style="color: #0000FF;">,</span><span style="color: #000000;">6</span><span style="color: #0000FF;">},</span><span style="color: #000000;">giuga</span><span style="color: #0000FF;">)</span>
<!--</syntaxhighlight>-->
{{out}}
<pre>
n = 3: 30
n = 4: 1722 858
n = 5: 66198
n = 6: 24423128562 2214408306
</pre>
You can (almost) push things a little further on 64-bit:
<!--<syntaxhighlight lang="phix">-->
<span style="color: #008080;">if</span> <span style="color: #7060A8;">machine_bits</span><span style="color: #0000FF;">()=</span><span style="color: #000000;">64</span> <span style="color: #008080;">then</span> <span style="color: #000000;">giuga</span><span style="color: #0000FF;">(</span><span style="color: #000000;">7</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<!--</syntaxhighlight>-->
and get
<pre>
n = 7: 432749205173838 550843391309130318 14737133470010574
</pre>
It took about 3 minutes for those to show, but then carried on in a doomed quest for an hour before I killed it.
 
=={{header|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.
<syntaxhighlight lang="plm">
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
</syntaxhighlight>
{{out}}
<pre>
30 858 1722
</pre>
 
=={{header|Python}}==
{{trans|FreeBASIC}}
<langsyntaxhighlight lang="python">#!/usr/bin/python
 
from math import sqrt
Line 339 ⟶ 1,808:
c += 1
print(n)
n += 1</langsyntaxhighlight>
 
=={{header|Quackery}}==
 
<code>dpfs</code> is Distinct Prime Factors.
 
<syntaxhighlight lang="Quackery"> [ [] 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</syntaxhighlight>
 
{{out}}
 
<pre>[ 30 858 1722 66198 ]</pre>
 
=={{header|Raku}}==
 
<syntaxhighlight lang="raku" perl6line>my @primes = (3..60).grep: &is-prime;
 
print 'First four Giuga numbers: ';
Line 353 ⟶ 1,861:
$n if all .map: { ($n / $_ - 1) %% $_ };
}
}</langsyntaxhighlight>
{{out}}
<pre>First 4 Giuga numbers: 30 858 1722 66198</pre>
 
=={{header|Ring}}==
<syntaxhighlight lang="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
</syntaxhighlight>
{{out}}
<pre>
working...
The first 4 Giuga numbers are:
30 858 1722 66198
done...
</pre>
 
=={{header|RPL}}==
{{works with|HP|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?
≫ ≫ '<span style="color:blue">GIUGA?</span>' STO
≪ { } 4
'''WHILE''' OVER SIZE 4 < '''REPEAT'''
'''IF''' DUP <span style="color:blue">GIUGA?</span> '''THEN''' LOOK SWAP OVER + SWAP '''END'''
2 +
'''END''' DROP
≫ '<span style="color:blue">TASK</span>' STO
 
{{out}}
<pre>
1: {30 858 1722 66198}
</pre>
 
=={{header|Ruby}}==
 
<syntaxhighlight lang="ruby" line>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
</syntaxhighlight>
{{out}}
<pre>[30, 858, 1722, 66198]
</pre>
 
=={{header|Rust}}==
<syntaxhighlight lang="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 ) ;
}</syntaxhighlight>
{{out}}
<pre>
[30, 858, 1722, 66198]
</pre>
 
=={{header|Scala}}==
{{trans|Java}}
<syntaxhighlight lang="Scala">
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)
}
}
}
</syntaxhighlight>
{{out}}
<pre>
Found Giuga numbers: List(30, 858, 1722, 66198)
 
</pre>
 
=={{header|Sidef}}==
<syntaxhighlight lang="ruby">4..Inf `by` 2 -> lazy.grep {|n|
n.factor.all {|p| p `divides` (n/p - 1) }
}.first(4).say</syntaxhighlight>
{{out}}
<pre>
[30, 858, 1722, 66198]
</pre>
 
=={{header|Wren}}==
===Version 1 (Brute force)===
{{libheader|Wren-math}}
{{libheader|Wren-seq}}
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.105 seconds to find the first four Giuga numbers but finding the fifth would take many hours using this approach, so I haven't bothered. Hopefully, someone can come up with a better method.
<syntaxhighlight lang="wren">var factors = []
<lang ecmascript>import "./math" for Int
var inc = [4, 2, 4, 2, 4, 6, 2, 6]
import "./seq" for Lst
 
// 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) {
var factors = Int.primeFactors.call(n)
// can't be prime or semi-prime
var c = factors.count
if (cfactors.count > 2 && factors.all { |f| (n/f - 1) % f == 0 }) {
var factors2 = Lstgiuga.pruneadd(factorsn)
if (c == factors2.count && factors2.all { |f| (n/f - 1) % f == 0 }) {
giuga.add(n)
}
}
n = n + 24 // can't be divisible by 4
}
System.print("The first %(limit) Giuga numbers are:")
System.print(giuga)</langsyntaxhighlight>
 
{{out}}
Line 387 ⟶ 2,120:
The first 4 Giuga numbers are:
[30, 858, 1722, 66198]
</pre>
<br>
===Version 2 (Pari-GP translation)===
{{libheader|Wren-math}}
{{libheader|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.
<syntaxhighlight lang="wren">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()
}</syntaxhighlight>
 
{{out}}
<pre>
n = 3:
30
 
n = 4:
1722
858
 
n = 5:
66198
 
n = 6:
24423128562
2214408306
</pre>
 
=={{header|XPL0}}==
<langsyntaxhighlight XPL0lang="xpl0">func Giuga(N0); \Return 'true' if Giuga number
int N0;
int N, F, Q1, Q2, L;
Line 417 ⟶ 2,216:
N:= N+1;
];
]</langsyntaxhighlight>
 
{{out}}
Peak
2,442

edits

Retrieved from "https://rosettacode.org/wiki/Special:MobileDiff/321418...363916"