Arithmetic/Rational
The objective of this task is to create a reasonably complete implementation of rational arithmetic in the particular language using the idioms of the language.
You are encouraged to solve this task according to the task description, using any language you may know.
For example: Define a new type called frac with binary operator "//" of two integers that returns a structure made up of the numerator and the denominator (as per a rational number).
Further define the appropriate rational unary operators abs and '-', with the binary operators for addition '+', subtraction '-', multiplication '×', division '/', integer division '÷', modulo division, the comparison operators (e.g. '<', '≤', '>', & '≥') and equality operators (e.g. '=' & '≠').
Define standard coercion operators for casting int to frac etc.
If space allows, define standard increment and decrement operators (e.g. '+:=' & '-:=' etc.).
Finally test the operators: Use the new type frac to find all perfect numbers less than 219 by summing the reciprocal of the factors.
See also
Ada
The generic package specification:
generic
type Number is range <>;
package Generic_Rational is
type Rational is private;
function "abs" (A : Rational) return Rational;
function "+" (A : Rational) return Rational;
function "-" (A : Rational) return Rational;
function Inverse (A : Rational) return Rational;
function "+" (A : Rational; B : Rational) return Rational;
function "+" (A : Rational; B : Number ) return Rational;
function "+" (A : Number; B : Rational) return Rational;
function "-" (A : Rational; B : Rational) return Rational;
function "-" (A : Rational; B : Number ) return Rational;
function "-" (A : Number; B : Rational) return Rational;
function "*" (A : Rational; B : Rational) return Rational;
function "*" (A : Rational; B : Number ) return Rational;
function "*" (A : Number; B : Rational) return Rational;
function "/" (A : Rational; B : Rational) return Rational;
function "/" (A : Rational; B : Number ) return Rational;
function "/" (A : Number; B : Rational) return Rational;
function "/" (A : Number; B : Number) return Rational;
function ">" (A : Rational; B : Rational) return Boolean;
function ">" (A : Number; B : Rational) return Boolean;
function ">" (A : Rational; B : Number) return Boolean;
function "<" (A : Rational; B : Rational) return Boolean;
function "<" (A : Number; B : Rational) return Boolean;
function "<" (A : Rational; B : Number) return Boolean;
function ">=" (A : Rational; B : Rational) return Boolean;
function ">=" (A : Number; B : Rational) return Boolean;
function ">=" (A : Rational; B : Number) return Boolean;
function "<=" (A : Rational; B : Rational) return Boolean;
function "<=" (A : Number; B : Rational) return Boolean;
function "<=" (A : Rational; B : Number) return Boolean;
function "=" (A : Number; B : Rational) return Boolean;
function "=" (A : Rational; B : Number) return Boolean;
function Numerator (A : Rational) return Number;
function Denominator (A : Rational) return Number;
Zero : constant Rational;
One : constant Rational;
private
type Rational is record
Numerator : Number;
Denominator : Number;
end record;
Zero : constant Rational := (0, 1);
One : constant Rational := (1, 1);
end Generic_Rational;
The package can be instantiated with any integer type. It provides rational numbers represented by a numerator and denominator cleaned from the common divisors. Mixed arithmetic of the base integer type and the rational type is supported. Division to zero raises Constraint_Error. The implementation of the specification above is as follows:
package body Generic_Rational is
function GCD (A, B : Number) return Number is
begin
if A = 0 then
return B;
end if;
if B = 0 then
return A;
end if;
if A > B then
return GCD (B, A mod B);
else
return GCD (A, B mod A);
end if;
end GCD;
function Inverse (A : Rational) return Rational is
begin
if A.Numerator > 0 then
return (A.Denominator, A.Numerator);
elsif A.Numerator < 0 then
return (-A.Denominator, -A.Numerator);
else
raise Constraint_Error;
end if;
end Inverse;
function "abs" (A : Rational) return Rational is
begin
return (abs A.Numerator, A.Denominator);
end "abs";
function "+" (A : Rational) return Rational is
begin
return A;
end "+";
function "-" (A : Rational) return Rational is
begin
return (-A.Numerator, A.Denominator);
end "-";
function "+" (A : Rational; B : Rational) return Rational is
Common : constant Number := GCD (A.Denominator, B.Denominator);
A_Denominator : constant Number := A.Denominator / Common;
B_Denominator : constant Number := B.Denominator / Common;
begin
return (A.Numerator * B_Denominator + B.Numerator * A_Denominator) /
(A_Denominator * B.Denominator);
end "+";
function "+" (A : Rational; B : Number) return Rational is
begin
return (A.Numerator + B * A.Denominator) / A.Denominator;
end "+";
function "+" (A : Number; B : Rational) return Rational is
begin
return B + A;
end "+";
function "-" (A : Rational; B : Rational) return Rational is
begin
return A + (-B);
end "-";
function "-" (A : Rational; B : Number) return Rational is
begin
return A + (-B);
end "-";
function "-" (A : Number; B : Rational) return Rational is
begin
return A + (-B);
end "-";
function "*" (A : Rational; B : Rational) return Rational is
begin
return (A.Numerator * B.Numerator) / (A.Denominator * B.Denominator);
end "*";
function "*" (A : Rational; B : Number) return Rational is
Common : constant Number := GCD (A.Denominator, abs B);
begin
return (A.Numerator * B / Common, A.Denominator / Common);
end "*";
function "*" (A : Number; B : Rational) return Rational is
begin
return B * A;
end "*";
function "/" (A : Rational; B : Rational) return Rational is
begin
return A * Inverse (B);
end "/";
function "/" (A : Rational; B : Number) return Rational is
Common : constant Number := GCD (abs A.Numerator, abs B);
begin
if B > 0 then
return (A.Numerator / Common, A.Denominator * (B / Common));
else
return ((-A.Numerator) / Common, A.Denominator * ((-B) / Common));
end if;
end "/";
function "/" (A : Number; B : Rational) return Rational is
begin
return Inverse (B) * A;
end "/";
function "/" (A : Number; B : Number) return Rational is
Common : constant Number := GCD (abs A, abs B);
begin
if B = 0 then
raise Constraint_Error;
elsif A = 0 then
return (0, 1);
elsif A > 0 xor B > 0 then
return (-(abs A / Common), abs B / Common);
else
return (abs A / Common, abs B / Common);
end if;
end "/";
function ">" (A, B : Rational) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator > 0;
end ">";
function ">" (A : Number; B : Rational) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator > 0;
end ">";
function ">" (A : Rational; B : Number) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator > 0;
end ">";
function "<" (A, B : Rational) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator < 0;
end "<";
function "<" (A : Number; B : Rational) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator < 0;
end "<";
function "<" (A : Rational; B : Number) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator < 0;
end "<";
function ">=" (A, B : Rational) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator >= 0;
end ">=";
function ">=" (A : Number; B : Rational) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator >= 0;
end ">=";
function ">=" (A : Rational; B : Number) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator >= 0;
end ">=";
function "<=" (A, B : Rational) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator <= 0;
end "<=";
function "<=" (A : Number; B : Rational) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator <= 0;
end "<=";
function "<=" (A : Rational; B : Number) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator <= 0;
end "<=";
function "=" (A : Number; B : Rational) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator = 0;
end "=";
function "=" (A : Rational; B : Number) return Boolean is
Diff : constant Rational := A - B;
begin
return Diff.Numerator = 0;
end "=";
function Numerator (A : Rational) return Number is
begin
return A.Numerator;
end Numerator;
function Denominator (A : Rational) return Number is
begin
return A.Denominator;
end Denominator;
end Generic_Rational;
The implementation uses solution of the greatest common divisor task. Here is the implementation of the test:
with Ada.Numerics.Elementary_Functions; use Ada.Numerics.Elementary_Functions;
with Ada.Text_IO; use Ada.Text_IO;
with Generic_Rational;
procedure Test_Rational is
package Integer_Rational is new Generic_Rational (Integer);
use Integer_Rational;
begin
for Candidate in 2..2**15 loop
declare
Sum : Rational := 1 / Candidate;
begin
for Divisor in 2..Integer (Sqrt (Float (Candidate))) loop
if Candidate mod Divisor = 0 then -- Factor is a divisor of Candidate
Sum := Sum + One / Divisor + Rational'(Divisor / Candidate);
end if;
end loop;
if Sum = 1 then
Put_Line (Integer'Image (Candidate) & " is perfect");
end if;
end;
end loop;
end Test_Rational;
The perfect numbers are searched by summing of the reciprocal of each of the divisors of a candidate except 1. This sum must be 1 for a perfect number.
- Output:
6 is perfect 28 is perfect 496 is perfect 8128 is perfect
ALGOL 68
<lang algol68> MODE FRAC = STRUCT( INT num #erator#, den #ominator#);
FORMAT frac repr = $g(-0)"//"g(-0)$;
PROC gcd = (INT a, b) INT: # greatest common divisor #
(a = 0 | b |: b = 0 | a |: ABS a > ABS b | gcd(b, a MOD b) | gcd(a, b MOD a));
PROC lcm = (INT a, b)INT: # least common multiple #
a OVER gcd(a, b) * b;
PROC raise not implemented error = ([]STRING args)VOID: (
put(stand error, ("Not implemented error: ",args, newline));
stop
);
PRIO // = 9; # higher then the ** operator #
OP // = (INT num, den)FRAC: ( # initialise and normalise #
INT common = gcd(num, den);
IF den < 0 THEN
( -num OVER common, -den OVER common)
ELSE
( num OVER common, den OVER common)
FI
);
OP + = (FRAC a, b)FRAC: (
INT common = lcm(den OF a, den OF b);
FRAC result := ( common OVER den OF a * num OF a + common OVER den OF b * num OF b, common );
num OF result//den OF result
);
OP - = (FRAC a, b)FRAC: a + -b,
* = (FRAC a, b)FRAC: (
INT num = num OF a * num OF b,
den = den OF a * den OF b;
INT common = gcd(num, den);
(num OVER common) // (den OVER common)
);
OP / = (FRAC a, b)FRAC: a * FRAC(den OF b, num OF b),# real division #
% = (FRAC a, b)INT: ENTIER (a / b), # integer divison #
%* = (FRAC a, b)FRAC: a/b - FRACINIT ENTIER (a/b), # modulo division #
** = (FRAC a, INT exponent)FRAC:
IF exponent >= 0 THEN
(num OF a ** exponent, den OF a ** exponent )
ELSE
(den OF a ** exponent, num OF a ** exponent )
FI;
OP REALINIT = (FRAC frac)REAL: num OF frac / den OF frac,
FRACINIT = (INT num)FRAC: num // 1,
FRACINIT = (REAL num)FRAC: (
# express real number as a fraction # # a future execise! #
raise not implemented error(("Convert a REAL to a FRAC","!"));
SKIP
);
OP < = (FRAC a, b)BOOL: num OF (a - b) < 0,
> = (FRAC a, b)BOOL: num OF (a - b) > 0,
<= = (FRAC a, b)BOOL: NOT ( a > b ),
>= = (FRAC a, b)BOOL: NOT ( a < b ),
= = (FRAC a, b)BOOL: (num OF a, den OF a) = (num OF b, den OF b),
/= = (FRAC a, b)BOOL: (num OF a, den OF a) /= (num OF b, den OF b);
# Unary operators #
OP - = (FRAC frac)FRAC: (-num OF frac, den OF frac),
ABS = (FRAC frac)FRAC: (ABS num OF frac, ABS den OF frac),
ENTIER = (FRAC frac)INT: (num OF frac OVER den OF frac) * den OF frac;
COMMENT Operators for extended characters set, and increment/decrement:
OP +:= = (REF FRAC a, FRAC b)REF FRAC: ( a := a + b ),
+=: = (FRAC a, REF FRAC b)REF FRAC: ( b := a + b ),
-:= = (REF FRAC a, FRAC b)REF FRAC: ( a := a - b ),
*:= = (REF FRAC a, FRAC b)REF FRAC: ( a := a * b ),
/:= = (REF FRAC a, FRAC b)REF FRAC: ( a := a / b ),
%:= = (REF FRAC a, FRAC b)REF FRAC: ( a := FRACINIT (a % b) ),
%*:= = (REF FRAC a, FRAC b)REF FRAC: ( a := a %* b );
# OP aliases for extended character sets (eg: Unicode, APL, ALCOR and GOST 10859) #
OP × = (FRAC a, b)FRAC: a * b,
÷ = (FRAC a, b)INT: a OVER b,
÷× = (FRAC a, b)FRAC: a MOD b,
÷* = (FRAC a, b)FRAC: a MOD b,
%× = (FRAC a, b)FRAC: a MOD b,
≤ = (FRAC a, b)FRAC: a <= b,
≥ = (FRAC a, b)FRAC: a >= b,
≠ = (FRAC a, b)BOOL: a /= b,
↑ = (FRAC frac, INT exponent)FRAC: frac ** exponent,
÷×:= = (REF FRAC a, FRAC b)REF FRAC: ( a := a MOD b ),
%×:= = (REF FRAC a, FRAC b)REF FRAC: ( a := a MOD b ),
÷*:= = (REF FRAC a, FRAC b)REF FRAC: ( a := a MOD b );
# BOLD aliases for CPU that only support uppercase for 6-bit bytes - wrist watches #
OP OVER = (FRAC a, b)INT: a % b,
MOD = (FRAC a, b)FRAC: a %*b,
LT = (FRAC a, b)BOOL: a < b,
GT = (FRAC a, b)BOOL: a > b,
LE = (FRAC a, b)BOOL: a <= b,
GE = (FRAC a, b)BOOL: a >= b,
EQ = (FRAC a, b)BOOL: a = b,
NE = (FRAC a, b)BOOL: a /= b,
UP = (FRAC frac, INT exponent)FRAC: frac**exponent;
# the required standard assignment operators #
OP PLUSAB = (REF FRAC a, FRAC b)REF FRAC: ( a +:= b ), # PLUS #
PLUSTO = (FRAC a, REF FRAC b)REF FRAC: ( a +=: b ), # PRUS #
MINUSAB = (REF FRAC a, FRAC b)REF FRAC: ( a *:= b ),
DIVAB = (REF FRAC a, FRAC b)REF FRAC: ( a /:= b ),
OVERAB = (REF FRAC a, FRAC b)REF FRAC: ( a %:= b ),
MODAB = (REF FRAC a, FRAC b)REF FRAC: ( a %*:= b );
END COMMENT Example: searching for Perfect Numbers.
FRAC sum:= FRACINIT 0;
FORMAT perfect = $b(" perfect!","")$;
FOR i FROM 2 TO 2**19 DO
INT candidate := i;
FRAC sum := 1 // candidate;
REAL real sum := 1 / candidate;
FOR factor FROM 2 TO ENTIER sqrt(candidate) DO
IF candidate MOD factor = 0 THEN
sum := sum + 1 // factor + 1 // ( candidate OVER factor);
real sum +:= 1 / factor + 1 / ( candidate OVER factor)
FI
OD;
IF den OF sum = 1 THEN
printf(($"Sum of reciprocal factors of "g(-0)" = "g(-0)" exactly, about "g(0,real width) f(perfect)l$,
candidate, ENTIER sum, real sum, ENTIER sum = 1))
FI
OD</lang>
- Output:
Sum of reciprocal factors of 6 = 1 exactly, about 1.0000000000000000000000000001 perfect! Sum of reciprocal factors of 28 = 1 exactly, about 1.0000000000000000000000000001 perfect! Sum of reciprocal factors of 120 = 2 exactly, about 2.0000000000000000000000000002 Sum of reciprocal factors of 496 = 1 exactly, about 1.0000000000000000000000000001 perfect! Sum of reciprocal factors of 672 = 2 exactly, about 2.0000000000000000000000000001 Sum of reciprocal factors of 8128 = 1 exactly, about 1.0000000000000000000000000001 perfect! Sum of reciprocal factors of 30240 = 3 exactly, about 3.0000000000000000000000000002 Sum of reciprocal factors of 32760 = 3 exactly, about 3.0000000000000000000000000003 Sum of reciprocal factors of 523776 = 2 exactly, about 2.0000000000000000000000000005
BBC BASIC
<lang bbcbasic> *FLOAT64
DIM frac{num, den}
DIM Sum{} = frac{}, Kf{} = frac{}, One{} = frac{}
One.num = 1 : One.den = 1
FOR n% = 2 TO 2^19-1
Sum.num = 1 : Sum.den = n%
FOR k% = 2 TO SQR(n%)
IF (n% MOD k%) = 0 THEN
Kf.num = 1 : Kf.den = k%
PROCadd(Sum{}, Kf{})
PROCnormalise(Sum{})
Kf.den = n% DIV k%
PROCadd(Sum{}, Kf{})
PROCnormalise(Sum{})
ENDIF
NEXT
IF FNeq(Sum{}, One{}) PRINT n% " is perfect"
NEXT n%
END
DEF PROCabs(a{}) : a.num = ABS(a.num) : ENDPROC
DEF PROCneg(a{}) : a.num = -a.num : ENDPROC
DEF PROCadd(a{}, b{})
LOCAL t : t = a.den * b.den
a.num = a.num * b.den + b.num * a.den
a.den = t
ENDPROC
DEF PROCsub(a{}, b{})
LOCAL t : t = a.den * b.den
a.num = a.num * b.den - b.num * a.den
a.den = t
ENDPROC
DEF PROCmul(a{}, b{})
a.num *= b.num : a.den *= b.den
ENDPROC
DEF PROCdiv(a{}, b{})
a.num *= b.den : a.den *= b.num
ENDPROC
DEF FNeq(a{}, b{}) = a.num * b.den = b.num * a.den
DEF FNlt(a{}, b{}) = a.num * b.den < b.num * a.den
DEF FNgt(a{}, b{}) = a.num * b.den > b.num * a.den
DEF FNne(a{}, b{}) = a.num * b.den <> b.num * a.den
DEF FNle(a{}, b{}) = a.num * b.den <= b.num * a.den
DEF FNge(a{}, b{}) = a.num * b.den >= b.num * a.den
DEF PROCnormalise(a{})
LOCAL a, b, t
a = a.num : b = a.den
WHILE b <> 0
t = a
a = b
b = t - b * INT(t / b)
ENDWHILE
a.num /= a : a.den /= a
IF a.den < 0 a.num *= -1 : a.den *= -1
ENDPROC</lang>
Output:
6 is perfect
28 is perfect
496 is perfect
8128 is perfect
C
C does not have overloadable operators. The following implementation does not define all operations so as to keep the example short. Note that the code passes around struct values instead of pointers to keep it simple, a practice normally avoided for efficiency reasons. <lang c>#include <stdio.h>
- include <stdlib.h>
- define FMT "%lld"
typedef long long int fr_int_t; typedef struct { fr_int_t num, den; } frac;
fr_int_t gcd(fr_int_t m, fr_int_t n) { fr_int_t t; while (n) { t = n; n = m % n; m = t; } return m; }
frac frac_new(fr_int_t num, fr_int_t den) { frac a; if (!den) { printf("divide by zero: "FMT"/"FMT"\n", num, den); abort(); }
int g = gcd(num, den);
if (g) { num /= g; den /= g; } else { num = 0; den = 1; }
if (den < 0) { den = -den; num = -num; } a.num = num; a.den = den; return a; }
- define BINOP(op, n, d) frac frac_##op(frac a, frac b) { return frac_new(n,d); }
BINOP(add, a.num * b.den + b.num * a.den, a.den * b.den); BINOP(sub, a.num * b.den - b.num + a.den, a.den * b.den); BINOP(mul, a.num * b.num, a.den * b.den); BINOP(div, a.num * b.den, a.den * b.num);
int frac_cmp(frac a, frac b) { int l = a.num * b.den, r = a.den * b.num; return l < r ? -1 : l > r; }
- define frac_cmp_int(a, b) frac_cmp(a, frac_new(b, 1))
int frtoi(frac a) { return a.den / a.num; } double frtod(frac a) { return (double)a.den / a.num; }
int main() { int n, k; frac sum, kf;
for (n = 2; n < 1<<19; n++) { sum = frac_new(1, n);
for (k = 2; k * k < n; k++) { if (n % k) continue; kf = frac_new(1, k); sum = frac_add(sum, kf);
kf = frac_new(1, n / k); sum = frac_add(sum, kf); } if (frac_cmp_int(sum, 1) == 0) printf("%d\n", n); }
return 0; }</lang> See Rational Arithmetic/C
C#
using System;
struct Fraction : IEquatable<Fraction>, IComparable<Fraction>
{
public readonly long Num;
public readonly long Denom;
public Fraction(long num, long denom)
{
if (num == 0)
{
denom = 1;
}
else if (denom == 0)
{
throw new ArgumentException("Denominator may not be zero", "denom");
}
else if (denom < 0)
{
num = -num;
denom = -denom;
}
long d = GCD(num, denom);
this.Num = num / d;
this.Denom = denom / d;
}
private static long GCD(long x, long y)
{
return y == 0 ? x : GCD(y, x % y);
}
private static long LCM(long x, long y)
{
return x / GCD(x, y) * y;
}
public Fraction Abs()
{
return new Fraction(Math.Abs(Num), Denom);
}
public Fraction Reciprocal()
{
return new Fraction(Denom, Num);
}
#region Conversion Operators
public static implicit operator Fraction(long i)
{
return new Fraction(i, 1);
}
public static explicit operator double(Fraction f)
{
return f.Num == 0 ? 0 : (double)f.Num / f.Denom;
}
#endregion
#region Arithmetic Operators
public static Fraction operator -(Fraction f)
{
return new Fraction(-f.Num, f.Denom);
}
public static Fraction operator +(Fraction a, Fraction b)
{
long m = LCM(a.Denom, b.Denom);
long na = a.Num * m / a.Denom;
long nb = b.Num * m / b.Denom;
return new Fraction(na + nb, m);
}
public static Fraction operator -(Fraction a, Fraction b)
{
return a + (-b);
}
public static Fraction operator *(Fraction a, Fraction b)
{
return new Fraction(a.Num * b.Num, a.Denom * b.Denom);
}
public static Fraction operator /(Fraction a, Fraction b)
{
return a * b.Reciprocal();
}
public static Fraction operator %(Fraction a, Fraction b)
{
long l = a.Num * b.Denom, r = a.Denom * b.Num;
long n = l / r;
return new Fraction(l - n * r, a.Denom * b.Denom);
}
#endregion
#region Comparison Operators
public static bool operator ==(Fraction a, Fraction b)
{
return a.Num == b.Num && a.Denom == b.Denom;
}
public static bool operator !=(Fraction a, Fraction b)
{
return a.Num != b.Num || a.Denom != b.Denom;
}
public static bool operator <(Fraction a, Fraction b)
{
return (a.Num * b.Denom) < (a.Denom * b.Num);
}
public static bool operator >(Fraction a, Fraction b)
{
return (a.Num * b.Denom) > (a.Denom * b.Num);
}
public static bool operator <=(Fraction a, Fraction b)
{
return !(a > b);
}
public static bool operator >=(Fraction a, Fraction b)
{
return !(a < b);
}
#endregion
#region Object Members
public override bool Equals(object obj)
{
if (obj is Fraction)
return ((Fraction)obj) == this;
else
return false;
}
public override int GetHashCode()
{
return Num.GetHashCode() ^ Denom.GetHashCode();
}
public override string ToString()
{
return Num.ToString() + "/" + Denom.ToString();
}
#endregion
#region IEquatable<Fraction> Members
public bool Equals(Fraction other)
{
return other == this;
}
#endregion
#region IComparable<Fraction> Members
public int CompareTo(Fraction other)
{
return (this.Num * other.Denom).CompareTo(this.Denom * other.Num);
}
#endregion
}
Test program:
using System;
static class Program
{
static void Main(string[] args)
{
int max = 1 << 19;
for (int candidate = 2; candidate < max; candidate++)
{
Fraction sum = new Fraction(1, candidate);
int max2 = (int)Math.Sqrt(candidate);
for (int factor = 2; factor <= max2; factor++)
{
if (candidate % factor == 0)
{
sum += new Fraction(1, factor);
sum += new Fraction(1, candidate / factor);
}
}
if (sum == 1)
Console.WriteLine("{0} is perfect", candidate);
}
}
}
- Output:
6 is perfect 28 is perfect 496 is perfect 8128 is perfect
C++
Boost provides a rational number template. <lang cpp>#include <iostream>
- include "math.h"
- include "boost/rational.hpp"
typedef boost::rational<int> frac;
bool is_perfect(int c) {
frac sum(1, c);
for (int f = 2;f < sqrt(static_cast<float>(c)); ++f){
if (c % f == 0) sum += frac(1,f) + frac(1, c/f);
}
if (sum.denominator() == 1){
return (sum == 1);
}
return false;
}
int main() {
for (int candidate = 2; candidate < 0x80000; ++candidate){
if (is_perfect(candidate))
std::cout << candidate << " is perfect" << std::endl;
} return 0;
}</lang>
Clojure
Ratios are built in to Clojure and support math operations already. They automatically reduce and become Integers if possible. <lang Clojure>user> 22/7 22/7 user> 34/2 17 user> (+ 37/5 42/9) 181/15</lang>
Common Lisp
Common Lisp has rational numbers built-in and integrated with all other number types. Common Lisp's number system is not extensible so reimplementing rational arithmetic would require all-new operator names. <lang lisp>(loop for candidate from 2 below (expt 2 19)
for sum = (+ (/ candidate)
(loop for factor from 2 to (isqrt candidate)
when (zerop (mod candidate factor))
sum (+ (/ factor) (/ (floor candidate factor)))))
when (= sum 1)
collect candidate)</lang>
D
Rational implementation based on BigInt. <lang d>import std.bigint, std.traits;
T gcd(T)(/*in*/ T a, /*in*/ T b) /*pure nothrow*/ {
// std.numeric.gcd doesn't work with BigInt return (b != 0) ? gcd(b, a % b) : (a < 0) ? -a : a;
}
T lcm(T)(/*in*/ T a, /*in*/ T b) {
return a / gcd(a, b) * b;
}
BigInt toBig(T : BigInt)(/*const*/ ref T n) pure nothrow { return n; }
BigInt toBig(T)(in ref T n) pure nothrow if (isIntegral!T) {
return BigInt(n);
}
struct Rational {
/*const*/ private BigInt num, den; // numerator & denominator
private enum Type { NegINF = -2,
NegDEN = -1,
NaRAT = 0,
NORMAL = 1,
PosINF = 2 };
this(U : Rational)(U n) pure nothrow {
num = n.num;
den = n.den;
}
this(U)(in U n) pure nothrow if (isIntegral!U) {
num = toBig(n);
den = 1UL;
}
this(U, V)(/*in*/ U n, /*in*/ V d) /*pure nothrow*/ {
num = toBig(n);
den = toBig(d);
/*const*/ BigInt common = gcd(num, den);
if (common != 0) {
num /= common;
den /= common;
} else { // infinite or NOT a Number
num = (num == 0) ? 0 : (num < 0) ? -1 : 1;
den = 0;
}
if (den < 0) { // assure den is non-negative
num = -num;
den = -den;
}
}
BigInt nomerator() /*const*/ pure nothrow @property {
return num;
}
BigInt denominator() /*const*/ pure nothrow @property {
return den;
}
string toString() /*const*/ {
if (den == 0) {
if (num == 0)
return "NaRat";
else
return ((num < 0) ? "-" : "+") ~ "infRat";
}
return toDecimalString(num) ~
(den == 1 ? "" : ("/" ~ toDecimalString(den)));
}
Rational opBinary(string op)(/*in*/ Rational r)
/*const pure nothrow*/ if (op == "+" || op == "-") {
BigInt common = lcm(den, r.den);
BigInt n = mixin("common / den * num" ~ op ~
"common / r.den * r.num" );
return Rational(n, common);
}
Rational opBinary(string op)(/*in*/ Rational r)
/*const pure nothrow*/ if (op == "*") {
return Rational(num * r.num, den * r.den);
}
Rational opBinary(string op)(/*in*/ Rational r)
/*const pure nothrow*/ if (op == "/") {
return Rational(num * r.den, den * r.num);
}
Rational opBinary(string op, T)(in T r)
/*const pure nothrow*/ if (isIntegral!T && (op == "+" ||
op == "-" || op == "*" || op == "/")) {
return opBinary!op(Rational(r));
}
Rational opBinary(string op)(in size_t p)
/*const pure nothrow*/ if (op == "^^") {
return Rational(num ^^ p, den ^^ p);
}
Rational opBinaryRight(string op, T)(in T l)
/*const pure nothrow*/ if (isIntegral!T) {
return Rational(l).opBinary!op(Rational(num, den));
}
Rational opUnary(string op)()
/*const pure nothrow*/ if (op == "+" || op == "-") {
return Rational(mixin(op ~ "num"), den);
}
int opCmp(T)(/*in*/ T r) /*const pure nothrow*/ {
Rational rhs = Rational(r);
if (type() == Type.NaRAT || rhs.type() == Type.NaRAT)
throw new Exception("Compare invlove an NaRAT.");
if (type() != Type.NORMAL ||
rhs.type() != Type.NORMAL) // for infinite
return (type() == rhs.type()) ? 0 :
((type() < rhs.type()) ? -1 : 1);
BigInt diff = num * rhs.den - den * rhs.num;
return (diff == 0) ? 0 : ((diff < 0) ? -1 : 1);
}
int opEquals(T)(/*in*/ T r) /*const pure nothrow*/ {
Rational rhs = Rational(r);
if (type() == Type.NaRAT || rhs.type() == Type.NaRAT)
return false;
return num == rhs.num && den == rhs.den;
}
Type type() /*const pure nothrow*/ {
if (den > 0) return Type.NORMAL;
if (den < 0) return Type.NegDEN;
if (num > 0) return Type.PosINF;
if (num < 0) return Type.NegINF;
return Type.NaRAT;
}
}
version (arithmetic_rational_main) { // test part void main() {
import std.stdio, std.math;
foreach (p; 2 .. 2 ^^ 19) {
auto sum = Rational(1, p);
immutable limit = 1 + cast(uint)sqrt(cast(real)p);
foreach (factor; 2 .. limit)
if (p % factor == 0)
sum = sum + Rational(1, factor) + Rational(factor, p);
if (sum.denominator == 1)
writefln("Sum of recipr. factors of %6s = %s exactly%s",
p, sum, (sum == 1) ? ", perfect." : ".");
}
}
}</lang>
Use the -version=rational_arithmetic_main compiler switch to run the test code.
- Output:
Sum of recipr. factors of 6 = 1 exactly, perfect. Sum of recipr. factors of 28 = 1 exactly, perfect. Sum of recipr. factors of 120 = 2 exactly. Sum of recipr. factors of 496 = 1 exactly, perfect. Sum of recipr. factors of 672 = 2 exactly. Sum of recipr. factors of 8128 = 1 exactly, perfect. Sum of recipr. factors of 30240 = 3 exactly. Sum of recipr. factors of 32760 = 3 exactly.
Run-time is about 9.5 seconds. It's quite slow because in DMD v.2.060 BigInts have no memory optimizations.
Output using p up to 2^^19, as requested by the task:
Sum of recipr. factors of 6 = 1 exactly, perfect. Sum of recipr. factors of 28 = 1 exactly, perfect. Sum of recipr. factors of 120 = 2 exactly. Sum of recipr. factors of 496 = 1 exactly, perfect. Sum of recipr. factors of 672 = 2 exactly. Sum of recipr. factors of 8128 = 1 exactly, perfect. Sum of recipr. factors of 30240 = 3 exactly. Sum of recipr. factors of 32760 = 3 exactly. Sum of recipr. factors of 523776 = 2 exactly.
Elisa
<lang Elisa>component RationalNumbers;
type Rational;
Rational(Numerator = integer, Denominater = integer) -> Rational;
Rational + Rational -> Rational;
Rational - Rational -> Rational;
Rational * Rational -> Rational;
Rational / Rational -> Rational;
Rational == Rational -> boolean;
Rational <> Rational -> boolean;
Rational >= Rational -> boolean;
Rational <= Rational -> boolean;
Rational > Rational -> boolean;
Rational < Rational -> boolean;
+ Rational -> Rational;
- Rational -> Rational;
abs(Rational) -> Rational;
Rational(integer) -> Rational;
Numerator(Rational) -> integer;
Denominator(Rational) -> integer;
begin
Rational(A,B) = Rational:[A;B];
R1 + R2 = Normalize( R1.A * R2.B + R1.B * R2.A, R1.B * R2.B);
R1 - R2 = Normalize( R1.A * R2.B - R1.B * R2.A, R1.B * R2.B);
R1 * R2 = Normalize( R1.A * R2.A, R1.B * R2.B);
R1 / R2 = Normalize( R1.A * R2.B, R1.B * R2.A);
R1 == R2 = [ R = (R1 - R2); R.A == 0];
R1 <> R2 = [ R = (R1 - R2); R.A <> 0];
R1 >= R2 = [ R = (R1 - R2); R.A >= 0];
R1 <= R2 = [ R = (R1 - R2); R.A <= 0];
R1 > R2 = [ R = (R1 - R2); R.A > 0];
R1 < R2 = [ R = (R1 - R2); R.A < 0];
+ R = R;
- R = Rational(-R.A, R.B);
abs(R) = Rational(abs(R.A), abs(R.B));
Rational(I) = Rational (I, 1);
Numerator(R) = R.A;
Denominator(R) = R.B;
<< internal definitions >>
Normalize (A = integer, B = integer) -> Rational;
Normalize (A, B) = [ exception( B == 0, "Illegal Rational Number");
Common = GCD(abs(A), abs(B)); if B < 0 then Rational(-A / Common, -B / Common) else Rational( A / Common, B / Common) ];
GCD (A = integer, B = integer) -> integer;
GCD (A, B) = [ if A == 0 then return(B);
if B == 0 then return(A); if A > B then GCD (B, mod(A,B))
else GCD (A, mod(B,A)) ];
end component RationalNumbers;</lang> Tests <lang Elisa>use RationalNumbers;
PerfectNumbers( Limit = integer) -> multi(integer); PerfectNumbers( Limit) =
[ Candidate = 2 .. Limit;
Sum:= Rational(1,Candidate); [ Divisor = 2 .. integer(sqrt(real(Candidate))); if mod(Candidate, Divisor) == 0 then Sum := Sum + Rational(1, Divisor) + Rational(Divisor, Candidate); ]; if Sum == Rational(1,1) then Candidate
];
PerfectNumbers(10000)?</lang>
- Output:
6 28 496 8128
Forth
<lang forth>\ Rationals can use any double cell operations: 2!, 2@, 2dup, 2swap, etc. \ Uses the stack convention of the built-in "*/" for int * frac -> int
- numerator drop ;
- denominator nip ;
- s>rat 1 ; \ integer to rational (n/1)
- rat>s / ; \ integer
- rat>frac mod ; \ fractional part
- rat>float swap s>f s>f f/ ;
- rat. swap 1 .r [char] / emit . ;
\ normalize: factors out gcd and puts sign into numerator
- gcd ( a b -- gcd ) begin ?dup while tuck mod repeat ;
- rat-normalize ( rat -- rat ) 2dup gcd tuck / >r / r> ;
- rat-abs swap abs swap ;
- rat-negate swap negate swap ;
- 1/rat over 0< if negate swap negate else swap then ;
- rat+ ( a b c d -- ad+bc bd )
rot 2dup * >r rot * >r * r> + r> rat-normalize ;
- rat- rat-negate rat+ ;
- rat* ( a b c d -- ac bd )
rot * >r * r> rat-normalize ;
- rat/ swap rat* ;
- rat-equal d= ;
- rat-less ( a b c d -- ad<bc )
-rot * >r * r> < ;
- rat-more 2swap rat-less ;
- rat-inc tuck + swap ;
- rat-dec tuck - swap ;</lang>
Fortran
<lang fortran>module module_rational
implicit none private public :: rational public :: rational_simplify public :: assignment (=) public :: operator (//) public :: operator (+) public :: operator (-) public :: operator (*) public :: operator (/) public :: operator (<) public :: operator (<=) public :: operator (>) public :: operator (>=) public :: operator (==) public :: operator (/=) public :: abs public :: int public :: modulo type rational integer :: numerator integer :: denominator end type rational interface assignment (=) module procedure assign_rational_int, assign_rational_real end interface interface operator (//) module procedure make_rational end interface interface operator (+) module procedure rational_add end interface interface operator (-) module procedure rational_minus, rational_subtract end interface interface operator (*) module procedure rational_multiply end interface interface operator (/) module procedure rational_divide end interface interface operator (<) module procedure rational_lt end interface interface operator (<=) module procedure rational_le end interface interface operator (>) module procedure rational_gt end interface interface operator (>=) module procedure rational_ge end interface interface operator (==) module procedure rational_eq end interface interface operator (/=) module procedure rational_ne end interface interface abs module procedure rational_abs end interface interface int module procedure rational_int end interface interface modulo module procedure rational_modulo end interface
contains
recursive function gcd (i, j) result (res)
integer, intent (in) :: i
integer, intent (in) :: j
integer :: res
if (j == 0) then
res = i
else
res = gcd (j, modulo (i, j))
end if
end function gcd
function rational_simplify (r) result (res) type (rational), intent (in) :: r type (rational) :: res integer :: g g = gcd (r % numerator, r % denominator) res = r % numerator / g // r % denominator / g end function rational_simplify
function make_rational (numerator, denominator) result (res) integer, intent (in) :: numerator integer, intent (in) :: denominator type (rational) :: res res = rational (numerator, denominator) end function make_rational
subroutine assign_rational_int (res, i) type (rational), intent (out), volatile :: res integer, intent (in) :: i res = i // 1 end subroutine assign_rational_int
subroutine assign_rational_real (res, x)
type (rational), intent(out), volatile :: res
real, intent (in) :: x
integer :: x_floor
real :: x_frac
x_floor = floor (x)
x_frac = x - x_floor
if (x_frac == 0) then
res = x_floor // 1
else
res = (x_floor // 1) + (1 // floor (1 / x_frac))
end if
end subroutine assign_rational_real
function rational_add (r, s) result (res)
type (rational), intent (in) :: r
type (rational), intent (in) :: s
type (rational) :: res
res = r % numerator * s % denominator + r % denominator * s % numerator // &
& r % denominator * s % denominator
end function rational_add
function rational_minus (r) result (res) type (rational), intent (in) :: r type (rational) :: res res = - r % numerator // r % denominator end function rational_minus
function rational_subtract (r, s) result (res)
type (rational), intent (in) :: r
type (rational), intent (in) :: s
type (rational) :: res
res = r % numerator * s % denominator - r % denominator * s % numerator // &
& r % denominator * s % denominator
end function rational_subtract
function rational_multiply (r, s) result (res) type (rational), intent (in) :: r type (rational), intent (in) :: s type (rational) :: res res = r % numerator * s % numerator // r % denominator * s % denominator end function rational_multiply
function rational_divide (r, s) result (res) type (rational), intent (in) :: r type (rational), intent (in) :: s type (rational) :: res res = r % numerator * s % denominator // r % denominator * s % numerator end function rational_divide
function rational_lt (r, s) result (res)
type (rational), intent (in) :: r
type (rational), intent (in) :: s
type (rational) :: r_simple
type (rational) :: s_simple
logical :: res
r_simple = rational_simplify (r)
s_simple = rational_simplify (s)
res = r_simple % numerator * s_simple % denominator < &
& s_simple % numerator * r_simple % denominator
end function rational_lt
function rational_le (r, s) result (res)
type (rational), intent (in) :: r
type (rational), intent (in) :: s
type (rational) :: r_simple
type (rational) :: s_simple
logical :: res
r_simple = rational_simplify (r)
s_simple = rational_simplify (s)
res = r_simple % numerator * s_simple % denominator <= &
& s_simple % numerator * r_simple % denominator
end function rational_le
function rational_gt (r, s) result (res)
type (rational), intent (in) :: r
type (rational), intent (in) :: s
type (rational) :: r_simple
type (rational) :: s_simple
logical :: res
r_simple = rational_simplify (r)
s_simple = rational_simplify (s)
res = r_simple % numerator * s_simple % denominator > &
& s_simple % numerator * r_simple % denominator
end function rational_gt
function rational_ge (r, s) result (res)
type (rational), intent (in) :: r
type (rational), intent (in) :: s
type (rational) :: r_simple
type (rational) :: s_simple
logical :: res
r_simple = rational_simplify (r)
s_simple = rational_simplify (s)
res = r_simple % numerator * s_simple % denominator >= &
& s_simple % numerator * r_simple % denominator
end function rational_ge
function rational_eq (r, s) result (res) type (rational), intent (in) :: r type (rational), intent (in) :: s logical :: res res = r % numerator * s % denominator == s % numerator * r % denominator end function rational_eq
function rational_ne (r, s) result (res) type (rational), intent (in) :: r type (rational), intent (in) :: s logical :: res res = r % numerator * s % denominator /= s % numerator * r % denominator end function rational_ne
function rational_abs (r) result (res) type (rational), intent (in) :: r type (rational) :: res res = sign (r % numerator, r % denominator) // r % denominator end function rational_abs
function rational_int (r) result (res) type (rational), intent (in) :: r integer :: res res = r % numerator / r % denominator end function rational_int
function rational_modulo (r) result (res) type (rational), intent (in) :: r integer :: res res = modulo (r % numerator, r % denominator) end function rational_modulo
end module module_rational</lang> Example: <lang fortran>program perfect_numbers
use module_rational implicit none integer, parameter :: n_min = 2 integer, parameter :: n_max = 2 ** 19 - 1 integer :: n integer :: factor type (rational) :: sum
do n = n_min, n_max
sum = 1 // n
factor = 2
do
if (factor * factor >= n) then
exit
end if
if (modulo (n, factor) == 0) then
sum = rational_simplify (sum + (1 // factor) + (factor // n))
end if
factor = factor + 1
end do
if (sum % numerator == 1 .and. sum % denominator == 1) then
write (*, '(i0)') n
end if
end do
end program perfect_numbers</lang>
- Output:
6 28 496 8128
GAP
Rational numbers are built-in. <lang gap>2/3 in Rationals;
- true
2/3 + 3/4;
- 17/12</lang>
Go
Go does not have user defined operators. For this reason the task as written cannot be done in Go.
Go does however have arbitrary-precision integers and rational numbers in the big package of the standard library. Code here implements the perfect number test described in the task using the standard library.
<lang go>package main
import (
"fmt" "math" "math/big"
)
func main() {
var recip big.Rat
max := int64(1 << 19)
for candidate := int64(2); candidate < max; candidate++ {
sum := big.NewRat(1, candidate)
max2 := int64(math.Sqrt(float64(candidate)))
for factor := int64(2); factor <= max2; factor++ {
if candidate%factor == 0 {
sum.Add(sum, recip.SetFrac64(1, factor))
if f2 := candidate / factor; f2 != factor {
sum.Add(sum, recip.SetFrac64(1, f2))
}
}
}
if sum.Denom().Int64() == 1 {
perfectstring := ""
if sum.Num().Int64() == 1 {
perfectstring = "perfect!"
}
fmt.Printf("Sum of recipr. factors of %d = %d exactly %s\n",
candidate, sum.Num().Int64(), perfectstring)
}
}
}</lang>
- Output:
Sum of recipr. factors of 6 = 1 exactly perfect! Sum of recipr. factors of 28 = 1 exactly perfect! Sum of recipr. factors of 120 = 2 exactly Sum of recipr. factors of 496 = 1 exactly perfect! Sum of recipr. factors of 672 = 2 exactly Sum of recipr. factors of 8128 = 1 exactly perfect! Sum of recipr. factors of 30240 = 3 exactly Sum of recipr. factors of 32760 = 3 exactly Sum of recipr. factors of 523776 = 2 exactly
Groovy
Groovy does not provide any built-in facility for rational arithmetic. However, it does support arithmetic operator overloading. Thus it is not too hard to build a fairly robust, complete, and intuitive rational number class, such as the following: <lang groovy>class Rational implements Comparable {
final BigInteger numerator, denominator
static final Rational ONE = new Rational(1, 1)
static final Rational ZERO = new Rational(0, 1)
Rational(BigInteger whole) { this(whole, 1) }
Rational(BigDecimal decimal) {
this(
decimal.scale() < 0 ? decimal.unscaledValue()*10**(-decimal.scale()) : decimal.unscaledValue(),
decimal.scale() < 0 ? 1 : 10**(decimal.scale())
)
}
Rational(num, denom) {
assert denom != 0 : "Denominator must not be 0"
def values = denom > 0 ? [num, denom] : [-num, -denom] //reduce(num, denom)
numerator = values[0]
denominator = values[1]
}
private List reduce(BigInteger num, BigInteger denom) {
BigInteger sign = ((num < 0) != (denom < 0)) ? -1 : 1
num = num.abs()
denom = denom.abs()
BigInteger commonFactor = gcd(num, denom)
[num.intdiv(commonFactor) * sign, denom.intdiv(commonFactor)]
}
public Rational toLeastTerms() {
def reduced = reduce(numerator, denominator)
new Rational(reduced[0], reduced[1])
}
private BigInteger gcd(BigInteger n, BigInteger m) { n == 0 ? m : { while(m%n != 0) { def t=n; n=m%n; m=t }; n }() }
Rational plus (Rational r) { new Rational(numerator*r.denominator + r.numerator*denominator, denominator*r.denominator) }
Rational plus (BigInteger n) { new Rational(numerator + n*denominator, denominator) }
Rational next () { new Rational(numerator + denominator, denominator) }
Rational minus (Rational r) { new Rational(numerator*r.denominator - r.numerator*denominator, denominator*r.denominator) }
Rational minus (BigInteger n) { new Rational(numerator - n*denominator, denominator) }
Rational previous () { new Rational(numerator - denominator, denominator) }
Rational multiply (Rational r) { new Rational(numerator*r.numerator, denominator*r.denominator) }
Rational multiply (BigInteger n) { new Rational(numerator*n, denominator) }
Rational div (Rational r) { new Rational(numerator*r.denominator, denominator*r.numerator) }
Rational div (BigInteger n) { new Rational(numerator, denominator*n) }
BigInteger intdiv (BigInteger n) { numerator.intdiv(denominator*n) }
Rational negative () { new Rational(-numerator, denominator) }
Rational abs () { new Rational(numerator.abs(), denominator) }
Rational reciprocal() { new Rational(denominator, numerator) }
Rational power(BigInteger n) { new Rational(numerator ** n, denominator ** n) }
boolean asBoolean() { numerator != 0 }
BigDecimal toBigDecimal() { (numerator as BigDecimal)/(denominator as BigDecimal) }
BigInteger toBigInteger() { numerator.intdiv(denominator) }
Double toDouble() { toBigDecimal().toDouble() }
double doubleValue() { toDouble() as double }
Float toFloat() { toBigDecimal().toFloat() }
float floatValue() { toFloat() as float }
Integer toInteger() { toBigInteger().toInteger() }
int intValue() { toInteger() as int }
Long toLong() { toBigInteger().toLong() }
long longValue() { toLong() as long }
Object asType(Class type) {
switch (type) {
case this.getClass(): return this
case Boolean.class: return asBoolean()
case BigDecimal.class: return toBigDecimal()
case BigInteger.class: return toBigInteger()
case Double.class: return toDouble()
case Float.class: return toFloat()
case Integer.class: return toInteger()
case Long.class: return toLong()
case String.class: return toString()
default: throw new ClassCastException("Cannot convert from type Rational to type " + type)
}
}
boolean equals(o) {
compareTo(o) == 0
}
int compareTo(o) {
o instanceof Rational \
? compareTo(o as Rational) \
: o instanceof Number \
? compareTo(o as Number)\
: (Double.NaN as int)
}
int compareTo(Rational r) { numerator*r.denominator <=> denominator*r.numerator }
int compareTo(Number n) { numerator <=> denominator*(n as BigInteger) }
int hashCode() { [numerator, denominator].hashCode() }
String toString() {
def reduced = reduce(numerator, denominator)
"${reduced[0]}//${reduced[1]}"
}
}</lang> The following script tests some of this class's features: <lang groovy>def x = new Rational(5, 20) def y = new Rational(9, 12) def z = new Rational(0, 10000)
println x println y println z println (x <=> y) println ((x as Rational).compareTo(y)) assert x*3 == y assert (z + 1) <= y*4 assert x != y
println "x + y == ${x} + ${y} == ${x + y}" println "x + z == ${x} + ${z} == ${x + z}" println "x - y == ${x} - ${y} == ${x - y}" println "x - z == ${x} - ${z} == ${x - z}" println "x * y == ${x} * ${y} == ${x * y}" println "y ** 3 == ${y} ** 3 == ${y ** 3}" println "x * z == ${x} * ${z} == ${x * z}" println "x / y == ${x} / ${y} == ${x / y}" try { print "x / z == ${x} / ${z} == "; println "${x / z}" } catch (Throwable t) { println t.message }
println "-x == -${x} == ${-x}" println "-y == -${y} == ${-y}" println "-z == -${z} == ${-z}"
print "x as int == ${x} as int == "; println x.intValue() print "x as double == ${x} as double == "; println x.doubleValue() print "1 / x as int == 1 / ${x} as int == "; println x.reciprocal().intValue() print "1.0 / x == 1.0 / ${x} == "; println x.reciprocal().doubleValue() print "y as int == ${y} as int == "; println y.intValue() print "y as double == ${y} as double == "; println y.doubleValue() print "1 / y as int == 1 / ${y} as int == "; println y.reciprocal().intValue() print "1.0 / y == 1.0 / ${y} == "; println y.reciprocal().doubleValue() print "z as int == ${z} as int == "; println z.intValue() print "z as double == ${z} as double == "; println z.doubleValue() try { print "1 / z as int == 1 / ${z} as int == "; println z.reciprocal().intValue() } catch (Throwable t) { println t.message } try { print "1.0 / z == 1.0 / ${z} == "; println z.reciprocal().doubleValue() } catch (Throwable t) { println t.message }
println "++x == ++ ${x} == ${++x}" println "++y == ++ ${y} == ${++y}" println "++z == ++ ${z} == ${++z}" println "-- --x == -- -- ${x} == ${-- (--x)}" println "-- --y == -- -- ${y} == ${-- (--y)}" println "-- --z == -- -- ${z} == ${-- (--z)}" println x println y println z
println (x <=> y) assert x*3 == y assert (z + 1) <= y*4 assert (x < y)
println (new Rational(25)) println (new Rational(25.0)) println (new Rational(0.25))
println Math.PI println (new Rational(Math.PI)) println ((new Rational(Math.PI)).toBigDecimal()) println ((new Rational(Math.PI)) as BigDecimal) println ((new Rational(Math.PI)) as Double) println ((new Rational(Math.PI)) as double) println ((new Rational(Math.PI)) as boolean) println (z as boolean) try { println ((new Rational(Math.PI)) as Date) } catch (Throwable t) { println t.message } try { println ((new Rational(Math.PI)) as char) } catch (Throwable t) { println t.message }</lang>
- Output:
1//4 3//4 0//1 -1 -1 x + y == 1//4 + 3//4 == 1//1 x + z == 1//4 + 0//1 == 1//4 x - y == 1//4 - 3//4 == -1//2 x - z == 1//4 - 0//1 == 1//4 x * y == 1//4 * 3//4 == 3//16 y ** 3 == 3//4 ** 3 == 27//64 x * z == 1//4 * 0//1 == 0//1 x / y == 1//4 / 3//4 == 1//3 x / z == 1//4 / 0//1 == Denominator must not be 0. Expression: (denom != 0). Values: denom = 0 -x == -1//4 == -1//4 -y == -3//4 == -3//4 -z == -0//1 == 0//1 x as int == 1//4 as int == 0 x as double == 1//4 as double == 0.25 1 / x as int == 1 / 1//4 as int == 4 1.0 / x == 1.0 / 1//4 == 4.0 y as int == 3//4 as int == 0 y as double == 3//4 as double == 0.75 1 / y as int == 1 / 3//4 as int == 1 1.0 / y == 1.0 / 3//4 == 1.3333333333 z as int == 0//1 as int == 0 z as double == 0//1 as double == 0.0 1 / z as int == 1 / 0//1 as int == Denominator must not be 0. Expression: (denom != 0). Values: denom = 0 1.0 / z == 1.0 / 0//1 == Denominator must not be 0. Expression: (denom != 0). Values: denom = 0 ++x == ++ 1//4 == 5//4 ++y == ++ 3//4 == 7//4 ++z == ++ 0//1 == 1//1 -- --x == -- -- 5//4 == -3//4 -- --y == -- -- 7//4 == -1//4 -- --z == -- -- 1//1 == -1//1 1//4 3//4 0//1 -1 25//1 25//1 1//4 3.141592653589793 884279719003555//281474976710656 3.141592653589793115997963468544185161590576171875 3.141592653589793115997963468544185161590576171875 3.141592653589793 3.141592653589793 true false Cannot convert from type Rational to type class java.util.Date Cannot convert from type Rational to type class java.lang.Character
The following uses the Rational class to find all perfect numbers less than 219: <lang groovy>def factorize = { target ->
if (target == 1L) {
return [1L]
} else if ([2L, 3L].contains(target)) {
return [1L, target]
}
def targetSqrt = Math.ceil(Math.sqrt(target)) as long
def lowfactors = (2L..(targetSqrt)).findAll { (target % it) == 0 }
if (lowfactors.isEmpty()) {
return [1L, target]
}
def nhalf = lowfactors.size() - ((lowfactors[-1] == targetSqrt) ? 1 : 0)
return ([1L] + lowfactors + ((nhalf-1)..0).collect { target.intdiv(lowfactors[it]) } + [target]).unique()
}
1.upto(2**19) {
if ((it % 100000) == 0) { println "HT" }
else if ((it % 1000) == 0) { print "." }
def factors = factorize(it)
def isPerfect = factors.collect{ factor -> new Rational( factor ).reciprocal() }.sum() == new Rational(2)
if (isPerfect) { println() ; println ([perfect: it, factors: factors]) }
}</lang>
- Output:
[perfect:6, factors:[1, 2, 3, 6]] [perfect:28, factors:[1, 2, 4, 7, 14, 28]] [perfect:496, factors:[1, 2, 4, 8, 16, 31, 62, 124, 248, 496]] ........ [perfect:8128, factors:[1, 2, 4, 8, 16, 32, 64, 127, 254, 508, 1016, 2032, 4064, 8128]] ...........................................................................................HT ...................................................................................................HT ...................................................................................................HT ...................................................................................................HT ...................................................................................................HT ........................
Haskell
Haskell provides a Rational type, which is really an alias for Ratio Integer (Ratio being a polymorphic type implementing rational numbers for any Integral type of numerators and denominators). The fraction is constructed using the % operator.
<lang haskell>import Data.Ratio
-- simply prints all the perfect numbers main = mapM_ print [candidate
| candidate <- [2 .. 2^19],
getSum candidate == 1]
where getSum candidate = 1 % candidate +
sum [1 % factor + 1 % (candidate `div` factor)
| factor <- [2 .. floor(sqrt(fromIntegral(candidate)))],
candidate `mod` factor == 0]</lang>
For a sample implementation of Ratio, see the Haskell 98 Report.
Icon and Unicon
The IPL provides support for rational arithmetic
- The data type is called 'rational' not 'frac'.
- Use the record constructor 'rational' to create a rational. Sign must be 1 or -1.
- Neither Icon nor Unicon supports operator overloading. Augmented assignments make little sense w/o this.
- Procedures include 'negrat' (unary -), 'addrat' (+), 'subrat' (-), 'mpyrat' (*), 'divrat' (modulo /).
Additional procedures are implemented here to complete the task:
- 'makerat' (make), 'absrat' (abs), 'eqrat' (=), 'nerat' (~=), 'ltrat' (<), 'lerat' (<=), 'gerat' (>=), 'gtrat' (>)
<lang Icon>procedure main()
limit := 2^19
write("Perfect numbers up to ",limit," (using rational arithmetic):")
every write(is_perfect(c := 2 to limit))
write("End of perfect numbers")
# verify the rest of the implementation
zero := makerat(0) # from integer
half := makerat(0.5) # from real
qtr := makerat("1/4") # from strings ...
one := makerat("1")
mone := makerat("-1")
verifyrat("eqrat",zero,zero)
verifyrat("ltrat",zero,half)
verifyrat("ltrat",half,zero)
verifyrat("gtrat",zero,half)
verifyrat("gtrat",half,zero)
verifyrat("nerat",zero,half)
verifyrat("nerat",zero,zero)
verifyrat("absrat",mone,)
end
procedure is_perfect(c) #: test for perfect numbers using rational arithmetic
rsum := rational(1, c, 1)
every f := 2 to sqrt(c) do
if 0 = c % f then
rsum := addrat(rsum,addrat(rational(1,f,1),rational(1,integer(c/f),1)))
if rsum.numer = rsum.denom = 1 then
return c
end</lang>
- Output:
Perfect numbers up to 524288 (using rational arithmetic): 6 28 496 8128 End of perfect numbers Testing eqrat( (0/1), (0/1) ) ==> returned (0/1) Testing ltrat( (0/1), (1/2) ) ==> returned (1/2) Testing ltrat( (1/2), (0/1) ) ==> failed Testing gtrat( (0/1), (1/2) ) ==> failed Testing gtrat( (1/2), (0/1) ) ==> returned (0/1) Testing nerat( (0/1), (1/2) ) ==> returned (1/2) Testing nerat( (0/1), (0/1) ) ==> failed Testing absrat( (-1/1), ) ==> returned (1/1)
The following task functions are missing from the IPL: <lang Icon>procedure verifyrat(p,r1,r2) #: verification tests for rational procedures return write("Testing ",p,"( ",rat2str(r1),", ",rat2str(\r2) | &null," ) ==> ","returned " || rat2str(p(r1,r2)) | "failed") end
procedure makerat(x) #: make rational (from integer, real, or strings) local n,d static c initial c := &digits++'+-'
return case type(x) of {
"real" : real2rat(x)
"integer" : ratred(rational(x,1,1))
"string" : if x ? ( n := integer(tab(many(c))), ="/", d := integer(tab(many(c))), pos(0)) then
ratred(rational(n,d,1))
else
makerat(numeric(x))
}
end
procedure absrat(r1) #: abs(rational)
r1 := ratred(r1) r1.sign := 1 return r1
end
invocable all # for string invocation
procedure xoprat(op,r1,r2) #: support procedure for binary operations that cross denominators
local numer, denom, div
r1 := ratred(r1) r2 := ratred(r2)
return if op(r1.numer * r2.denom,r2.numer * r1.denom) then r2 # return right argument on success
end
procedure eqrat(r1,r2) #: rational r1 = r2 return xoprat("=",r1,r2) end
procedure nerat(r1,r2) #: rational r1 ~= r2 return xoprat("~=",r1,r2) end
procedure ltrat(r1,r2) #: rational r1 < r2 return xoprat("<",r1,r2) end
procedure lerat(r1,r2) #: rational r1 <= r2 return xoprat("<=",r1,r2) end
procedure gerat(r1,r2) #: rational r1 >= r2 return xoprat(">=",r1,r2) end
procedure gtrat(r1,r2) #: rational r1 > r2 return xoprat(">",r1,r2) end
link rational</lang>
The
provides rational and gcd in numbers. Record definition and usage is shown below:
<lang Icon> record rational(numer, denom, sign) # rational type
addrat(r1,r2) # Add rational numbers r1 and r2. divrat(r1,r2) # Divide rational numbers r1 and r2. medrat(r1,r2) # Form mediant of r1 and r2. mpyrat(r1,r2) # Multiply rational numbers r1 and r2. negrat(r) # Produce negative of rational number r. rat2real(r) # Produce floating-point approximation of r rat2str(r) # Convert the rational number r to its string representation. real2rat(v,p) # Convert real to rational with precision p (default 1e-10). Warning: excessive p gives ugly fractions reciprat(r) # Produce the reciprocal of rational number r. str2rat(s) # Convert the string representation (such as "3/2") to a rational number subrat(r1,r2) # Subtract rational numbers r1 and r2.
gcd(i, j) # returns greatest common divisor of i and j</lang>
J
J implements rational numbers: <lang j> 3r4*2r5 3r10</lang> That said, note that J's floating point numbers work just fine for the stated problem: <lang j> is_perfect_rational=: 2 = (1 + i.) +/@:%@([ #~ 0 = |) ]</lang> Faster version (but the problem, as stated, is still tremendously inefficient): <lang j>factors=: */&>@{@((^ i.@>:)&.>/)@q:~&__ is_perfect_rational=: 2= +/@:%@,@factors</lang> Exhaustive testing would take forever: <lang j> I.is_perfect_rational@"0 i.2^19 6 28 496 8128
I.is_perfect_rational@x:@"0 i.2^19x
6 28 496 8128</lang> More limited testing takes reasonable amounts of time: <lang j> (#~ is_perfect_rational"0) (* <:@+:) 2^i.10x 6 28 496 8128</lang>
Java
Uses BigRational class: Arithmetic/Rational/Java <lang java>class BigRationalFindPerfectNumbers {
public static void main(String[] args) {
System.out.println("Running BigRational built-in tests");
if (BigRational.testFeatures()) {
int MAX_NUM = (1 << 19);
System.out.println();
System.out.println("Searching for perfect numbers in the range [1, " + (MAX_NUM - 1) + "]");
BigRational TWO = BigRational.valueOf(2);
for (int i = 1; i < MAX_NUM; i++) {
BigRational reciprocalSum = BigRational.ONE;
if (i > 1)
reciprocalSum = reciprocalSum.add(BigRational.valueOf(i).reciprocal());
int maxDivisor = (int)Math.sqrt(i);
if (maxDivisor >= i)
maxDivisor--;
for (int divisor = 2; divisor <= maxDivisor; divisor++) {
if ((i % divisor) == 0) {
reciprocalSum = reciprocalSum.add(BigRational.valueOf(divisor).reciprocal());
int dividend = i / divisor;
if (divisor != dividend)
reciprocalSum = reciprocalSum.add(BigRational.valueOf(dividend).reciprocal());
}
}
if (reciprocalSum.equals(TWO))
System.out.println(String.valueOf(i) + " is a perfect number");
}
}
return;
}
}</lang>
- Output:
Running BigRational built-in tests PASS: BaseConstructor-1 PASS: BaseConstructor-2 PASS: BaseConstructor-3 PASS: BaseConstructor-4 PASS: Inequality-1 PASS: Inequality-2 PASS: IntegerConstructor-1 PASS: IntegerConstructor-2 ...(omitted for brevity)... PASS: Reciprocal-4 PASS: Signum-1 PASS: Signum-2 PASS: Signum-3 PASS: Numerator-1 PASS: Numerator-2 PASS: Denominator-1 PASS: Denominator-2 Passed all tests Searching for perfect numbers in the range [1, 524287] 6 is a perfect number 28 is a perfect number 496 is a perfect number 8128 is a perfect number
JavaScript
- The core of the Rational class
<lang javascript>// the constructor function Rational(numerator, denominator) {
if (denominator === undefined)
denominator = 1;
else if (denominator == 0)
throw "divide by zero";
this.numer = numerator;
if (this.numer == 0)
this.denom = 1;
else
this.denom = denominator;
this.normalize();
}
// getter methods Rational.prototype.numerator = function() {return this.numer}; Rational.prototype.denominator = function() {return this.denom};
// clone a rational Rational.prototype.dup = function() {
return new Rational(this.numerator(), this.denominator());
};
// conversion methods Rational.prototype.toString = function() {
if (this.denominator() == 1) {
return this.numerator().toString();
} else {
// implicit conversion of numbers to strings
return this.numerator() + '/' + this.denominator()
}
}; Rational.prototype.toFloat = function() {return eval(this.toString())} Rational.prototype.toInt = function() {return Math.floor(this.toFloat())};
// reduce Rational.prototype.normalize = function() {
// greatest common divisor
var a=Math.abs(this.numerator()), b=Math.abs(this.denominator())
while (b != 0) {
var tmp = a;
a = b;
b = tmp % b;
}
// a is the gcd
this.numer /= a;
this.denom /= a;
if (this.denom < 0) {
this.numer *= -1;
this.denom *= -1;
}
return this;
}
// absolute value // returns a new rational Rational.prototype.abs = function() {
return new Rational(Math.abs(this.numerator()), this.denominator());
};
// inverse // returns a new rational Rational.prototype.inv = function() {
return new Rational(this.denominator(), this.numerator());
};
// // arithmetic methods
// variadic, modifies receiver Rational.prototype.add = function() {
for (var i = 0; i < arguments.length; i++) {
this.numer = this.numer * arguments[i].denominator() + this.denom * arguments[i].numerator();
this.denom = this.denom * arguments[i].denominator();
}
return this.normalize();
};
// variadic, modifies receiver Rational.prototype.subtract = function() {
for (var i = 0; i < arguments.length; i++) {
this.numer = this.numer * arguments[i].denominator() - this.denom * arguments[i].numerator();
this.denom = this.denom * arguments[i].denominator();
}
return this.normalize();
};
// unary "-" operator // returns a new rational Rational.prototype.neg = function() {
return (new Rational(0)).subtract(this);
};
// variadic, modifies receiver Rational.prototype.multiply = function() {
for (var i = 0; i < arguments.length; i++) {
this.numer *= arguments[i].numerator();
this.denom *= arguments[i].denominator();
}
return this.normalize();
};
// modifies receiver Rational.prototype.divide = function(rat) {
return this.multiply(rat.inv());
}
// increment
// modifies receiver
Rational.prototype.inc = function() {
this.numer += this.denominator(); return this.normalize();
}
// decrement // modifies receiver Rational.prototype.dec = function() {
this.numer -= this.denominator(); return this.normalize();
}
// // comparison methods
Rational.prototype.isZero = function() {
return (this.numerator() == 0);
} Rational.prototype.isPositive = function() {
return (this.numerator() > 0);
} Rational.prototype.isNegative = function() {
return (this.numerator() < 0);
}
Rational.prototype.eq = function(rat) {
return this.dup().subtract(rat).isZero();
} Rational.prototype.ne = function(rat) {
return !(this.eq(rat));
} Rational.prototype.lt = function(rat) {
return this.dup().subtract(rat).isNegative();
} Rational.prototype.gt = function(rat) {
return this.dup().subtract(rat).isPositive();
} Rational.prototype.le = function(rat) {
return !(this.gt(rat));
} Rational.prototype.ge = function(rat) {
return !(this.lt(rat));
}</lang>
- Testing
<lang javascript>function assert(cond, msg) { if (!cond) throw msg; }
print('testing') var a, b, c, d, e, f;
//test creation a = new Rational(0); assert(a.toString() == "0", "Rational(0).toString() == '0'") a = new Rational(2); assert(a.toString() == "2", "Rational(2).toString() == '2'") a = new Rational(1,2); assert(a.toString() == "1/2", "Rational(1,2).toString() == '1/2'") b = new Rational(2,-12); assert(b.toString() == "-1/6", "Rational(1,6).toString() == '1/6'") f = new Rational(0,9)
a = new Rational(1,3) b = new Rational(1,2) c = new Rational(1,3)
assert(!(a.eq(b)), "1/3 == 1/2") assert(a.eq(c), "1/3 == 1/3") assert(a.ne(b), "1/3 != 1/2") assert(!(a.ne(c)), "1/3 != 1/3") assert(a.lt(b), "1/3 < 1/2") assert(!(b.lt(a)), "1/2 < 1/3") assert(!(a.lt(c)), "1/3 < 1/3") assert(!(a.gt(b)), "1/3 > 1/2") assert(b.gt(a), "1/2 > 1/3") assert(!(a.gt(c)), "1/3 > 1/3")
assert(a.le(b), "1/3 <= 1/2") assert(!(b.le(a)), "1/2 <= 1/3") assert(a.le(c), "1/3 <= 1/3") assert(!(a.ge(b)), "1/3 >= 1/2") assert(b.ge(a), "1/2 >= 1/3") assert(a.ge(c), "1/3 >= 1/3")
a = new Rational(1,2) b = new Rational(1,6) a.add(b); assert(a.eq(new Rational(2,3)), "1/2 + 1/6 == 2/3") c = a.neg(); assert(a.eq(new Rational(2,3)), "neg(1/2) == -1/2")
assert(c.eq(new Rational(2,-3)), "neg(1/2) == -1/2")
d = c.abs(); assert(c.eq(new Rational(-2,3)), "abs(neg(1/2)) == 1/2")
assert(d.eq(new Rational(2,3)), "abs(neg(1/2)) == 1/2")
b.subtract(a); assert(b.eq(new Rational(-1,2)), "1/6 - 1/2 == -1/3")
c = a.neg().abs(); assert(c.eq(a), "abs(neg(1/2)) == 1/2") c = (new Rational(-1,3)).inv(); assert(c.toString() == '-3', "inv(1/6 - 1/2) == -3") try {
e = f.inv(); throw "should have been an error: " +f + '.inv() = ' + e
} catch (e) {
assert(e == "divide by zero", "0.inv() === error")
}
b = new Rational(1,6) b.add(new Rational(2,3), new Rational(4,2)); assert(b.toString() == "17/6", "1/6+2/3+4/2 == 17/6");
a = new Rational(1,3); b = new Rational(1,6) c = new Rational(5,6); d = new Rational(1/5); e = new Rational(2); f = new Rational(0,9);
assert(c.dup().multiply(d).eq(b), "5/6 * 1/5 = 1/6")
assert(c.dup().multiply(d,e).eq(a), "5/6 * 1/5 *2 = 1/3")
assert(c.dup().multiply(d,e,f).eq(f), "5/6 * 1/5 *2*0 = 0")
c.divide(new Rational(5)); assert(c.eq(b), "5/6 / 5 = 1/6b")
try {
e = c.divide(f) throw "should have been an error: " + c + "/" + f + '= ' + e
} catch (e) {
assert(e == "divide by zero", "0.inv() === error")
}
print('all tests passed');</lang>
- Finding perfect numbers
<lang javascript>function factors(num) {
var factors = new Array();
var sqrt = Math.floor(Math.sqrt(num));
for (var i = 1; i <= sqrt; i++) {
if (num % i == 0) {
factors.push(i);
if (num / i != i)
factors.push(num / i);
}
}
factors.sort(function(a,b){return a-b}); // numeric sort
return factors;
}
function isPerfect(n) {
var sum = new Rational(0);
var fctrs = factors(n);
for (var i = 0; i < fctrs.length; i++)
sum.add(new Rational(1, fctrs[i]));
// note, fctrs includes 1, so sum should be 2 return sum.toFloat() == 2.0;
}
// find perfect numbers less than 2^19 for (var n = 2; n < Math.pow(2,19); n++)
if (isPerfect(n))
print("perfect: " + n);
// test 5th perfect number var n = Math.pow(2,12) * (Math.pow(2,13) - 1); if (isPerfect(n))
print("perfect: " + n);</lang>
- Output:
perfect: 6 perfect: 28 perfect: 496 perfect: 8128 perfect: 33550336
Lua
<lang lua>function gcd(a,b) return a == 0 and b or gcd(b % a, a) end
do
local function coerce(a, b)
if type(a) == "number" then return rational(a, 1), b end
if type(b) == "number" then return a, rational(b, 1) end
return a, b
end
rational = setmetatable({
__add = function(a, b)
local a, b = coerce(a, b)
return rational(a.num * b.den + a.den * b.num, a.den * b.den)
end,
__sub = function(a, b)
local a, b = coerce(a, b)
return rational(a.num * b.den - a.den * b.num, a.den * b.den)
end,
__mul = function(a, b)
local a, b = coerce(a, b)
return rational(a.num * b.num, a.den * b.den)
end,
__div = function(a, b)
local a, b = coerce(a, b)
return rational(a.num * b.den, a.den * b.num)
end,
__pow = function(a, b)
if type(a) == "number" then return a ^ (b.num / b.den) end
return rational(a.num ^ b, a.den ^ b) --runs into a problem if these aren't integers
end,
__concat = function(a, b)
if getmetatable(a) == rational then return a.num .. "/" .. a.den .. b end
return a .. b.num .. "/" .. b.den
end,
__unm = function(a) return rational(-a.num, -a.den) end}, {
__call = function(z, a, b) return setmetatable({num = a / gcd(a, b),den = b / gcd(a, b)}, z) end} )
end
print(rational(2, 3) + rational(3, 5) - rational(1, 10) .. "") --> 7/6 print((rational(4, 5) * rational(5, 9)) ^ rational(1, 2) .. "") --> 2/3 print(rational(45, 60) / rational(5, 2) .. "") --> 3/10 print(5 + rational(1, 3) .. "") --> 16/3
function findperfs(n)
local ret = {}
for i = 1, n do
sum = rational(1, i)
for fac = 2, i^.5 do
if i % fac == 0 then
sum = sum + rational(1, fac) + rational(fac, i)
end
end
if sum.den == sum.num then
ret[#ret + 1] = i
end
end
return table.concat(ret, '\n')
end print(findperfs(2^19))</lang>
Maple
Maple has full built-in support for arithmetic with fractions (rational numbers). Fractions are treated like any other number in Maple. <lang Maple> > a := 3 / 5;
a := 3/5
> numer( a );
3
> denom( a );
5
</lang> However, while you can enter a fraction such as "4/6", it will automatically be reduced so that the numerator and denominator have no common factor: <lang Maple> > b := 4 / 6;
b := 2/3
</lang> All the standard arithmetic operators work with rational numbers. It is not necessary to call any special routines. <lang Maple> > a + b;
19
--
15
> a * b;
2/5
> a / b;
9/10
> a - b;
-1
--
15
> a + 1;
8/5
> a - 1;
-2/5
</lang> Notice that fractions are treated as exact quantities; they are not converted to floats. However, you can get a floating point approximation to any desired accuracy by applying the function evalf to a fraction. <lang Maple> > evalf( 22 / 7 ); # default is 10 digits
3.142857143
> evalf[100]( 22 / 7 ); # 100 digits 3.142857142857142857142857142857142857142857142857142857142857142857\
142857142857142857142857142857143
</lang>
Mathematica
Mathematica has full support for fractions built-in. If one divides two exact numbers it will be left as a fraction if it can't be simplified. Comparison, addition, division, product et cetera are built-in: <lang Mathematica>4/16 3/8 8/4 4Pi/2 16!/10! Sqrt[9/16] Sqrt[3/4] (23/12)^5 2 + 1/(1 + 1/(3 + 1/4))
1/2+1/3+1/5 8/Pi+Pi/8 //Together 13/17 + 7/31 Sum[1/n,{n,1,100}] (*summation of 1/1 + 1/2 + 1/3 + 1/4+ .........+ 1/99 + 1/100*)
1/2-1/3 a=1/3;a+=1/7
1/4==2/8 1/4>3/8 Pi/E >23/20 1/3!=123/370 Sin[3]/Sin[2]>3/20
Numerator[6/9] Denominator[6/9]</lang> gives back:
1/4 3/8 2 2 Pi 5765760 3/4 Sqrt[3]/2 6436343 / 248832 47/17 31/30 (64+Pi^2) / (8 Pi) 522 / 527 14466636279520351160221518043104131447711 / 2788815009188499086581352357412492142272 1/6 10/21 True False True True True 2 3
As you can see, Mathematica automatically handles fraction as exact things, it doesn't evaluate the fractions to a float. It only does this when either the numerator or the denominator is not exact. I only showed integers above, but Mathematica can handle symbolic fraction in the same and complete way: <lang Mathematica>c/(2 c) (b^2 - c^2)/(b - c) // Cancel 1/2 + b/c // Together</lang> gives back: <lang Mathematica>1/2 b+c (2 b+c) / (2 c)</lang> Moreover it does simplification like Sin[x]/Cos[x] => Tan[x]. Division, addition, subtraction, powering and multiplication of a list (of any dimension) is automatically threaded over the elements: <lang Mathematica>1+2*{1,2,3}^3</lang> gives back: <lang Mathematica>{3, 17, 55}</lang> To check for perfect numbers in the range 1 to 2^25 we can use: <lang Mathematica>found={}; CheckPerfect[num_Integer]:=If[Total[1/Divisors[num]]==2,AppendTo[found,num]]; Do[CheckPerfect[i],{i,1,2^25}]; found</lang> gives back: <lang Mathematica>{6, 28, 496, 8128, 33550336}</lang> Final note; approximations of fractions to any precision can be found using the function N.
Maxima
<lang maxima>/* Rational numbers are builtin */ a: 3 / 11; 3/11
b: 117 / 17; 117/17
a + b; 1338/187
a - b; -1236/187
a * b; 351/187
a / b; 17/429
a^5; 243/161051
num(a); 3
denom(a); 11
ratnump(a); true</lang>
Modula-2
This is incomplete as the Perfect Numbers task has not been addressed.
- Definition Module
<lang modula2>DEFINITION MODULE Rational;
TYPE RAT = RECORD
numerator : INTEGER;
denominator : INTEGER;
END;
PROCEDURE IGCD( i : INTEGER; j : INTEGER ) : INTEGER; PROCEDURE ILCM( i : INTEGER; j : INTEGER ) : INTEGER; PROCEDURE IABS( i : INTEGER ) : INTEGER;
PROCEDURE RNormalize( i : RAT ) : RAT; PROCEDURE RCreate( num : INTEGER; dem : INTEGER ) : RAT; PROCEDURE RAdd( i : RAT; j : RAT ) : RAT; PROCEDURE RSubtract( i : RAT; j : RAT ) : RAT; PROCEDURE RMultiply( i : RAT; j : RAT ) : RAT; PROCEDURE RDivide( i : RAT; j : RAT ) : RAT; PROCEDURE RAbs( i : RAT ) : RAT; PROCEDURE RInv( i : RAT ) : RAT; PROCEDURE RNeg( i : RAT ) : RAT;
PROCEDURE RInc( i : RAT ) : RAT; PROCEDURE RDec( i : RAT ) : RAT; PROCEDURE REQ( i : RAT; j : RAT ) : BOOLEAN; PROCEDURE RNE( i : RAT; j : RAT ) : BOOLEAN; PROCEDURE RLT( i : RAT; j : RAT ) : BOOLEAN; PROCEDURE RLE( i : RAT; j : RAT ) : BOOLEAN; PROCEDURE RGT( i : RAT; j : RAT ) : BOOLEAN; PROCEDURE RGE( i : RAT; j : RAT ) : BOOLEAN;
PROCEDURE RIsZero( i : RAT ) : BOOLEAN; PROCEDURE RIsNegative( i : RAT ) : BOOLEAN; PROCEDURE RIsPositive( i : RAT ) : BOOLEAN;
PROCEDURE RToString( i : RAT; VAR S : ARRAY OF CHAR ); PROCEDURE RToRational( s : ARRAY OF CHAR ) : RAT;
PROCEDURE WriteRational( i : RAT );
END Rational.</lang>
- Implementation Module
<lang modula2>IMPLEMENTATION MODULE Rational;
FROM Strings IMPORT Assign, Append, Pos, Copy, Length; FROM NumberConversion IMPORT IntToString, StringToInt;
FROM InOut IMPORT WriteString (*, WriteCard,WriteLine, WriteInt, WriteLn *);
PROCEDURE IGCD( i : INTEGER; j : INTEGER ) : INTEGER;
VAR
res : INTEGER;
BEGIN
IF j = 0 THEN
res := i;
ELSE
res := IGCD( j, i MOD j );
END;
RETURN res; END IGCD;
PROCEDURE ILCM( i : INTEGER; j : INTEGER ) : INTEGER;
VAR
res : INTEGER;
BEGIN
res := (i DIV IGCD( i, j ) ) * j;
RETURN res;
END ILCM;
PROCEDURE IABS( i : INTEGER ) : INTEGER;
VAR
res : INTEGER;
BEGIN
IF i < 0 THEN
res := i * (-1);
ELSE
res := i;
END;
RETURN res;
END IABS;
PROCEDURE RNormalize( i : RAT ) : RAT;
VAR
gcd : INTEGER;
res : RAT;
BEGIN
gcd := IGCD( ABS( i.numerator ), ABS( i.denominator ) );
IF gcd <> 0 THEN
res.numerator := i.numerator DIV gcd;
res.denominator := i.denominator DIV gcd;
IF ( res.denominator < 0 ) THEN
res.numerator := res.numerator * (-1);
res.denominator := res.denominator * (-1);
END;
ELSE
WITH res DO
numerator := 0;
denominator := 0;
END;
END;
RETURN res;
END RNormalize;
PROCEDURE RCreate( num : INTEGER; dem : INTEGER ) : RAT;
VAR
rat : RAT;
BEGIN
WITH rat DO
numerator := num;
denominator := dem;
END;
RETURN RNormalize(rat);
END RCreate;
PROCEDURE RAdd( i : RAT; j : RAT ) : RAT;
BEGIN
RETURN RCreate( i.numerator * j.denominator + j.numerator * i.denominator, i.denominator * j.denominator );
END RAdd;
PROCEDURE RSubtract( i : RAT; j : RAT ) : RAT;
BEGIN
RETURN RCreate( i.numerator * j.denominator - j.numerator * i.denominator, i.denominator * j.denominator );
END RSubtract;
PROCEDURE RMultiply( i : RAT; j : RAT ) : RAT;
BEGIN
RETURN RCreate( i.numerator * j.numerator, i.denominator * j.denominator );
END RMultiply;
PROCEDURE RDivide( i : RAT; j : RAT ) : RAT;
BEGIN
RETURN RCreate( i.numerator * j.denominator, i.denominator * j.numerator );
END RDivide;
PROCEDURE RAbs( i : RAT ) : RAT;
BEGIN
RETURN RCreate( IABS( i.numerator ), i.denominator );
END RAbs;
PROCEDURE RInv( i : RAT ) : RAT;
BEGIN
RETURN RCreate( i.denominator, i.numerator );
END RInv;
PROCEDURE RNeg( i : RAT ) : RAT;
BEGIN
RETURN RCreate( i.numerator * (-1), i.denominator );
END RNeg;
PROCEDURE RInc( i : RAT ) : RAT;
BEGIN
RETURN RCreate( i.numerator + i.denominator, i.denominator );
END RInc;
PROCEDURE RDec( i : RAT ) : RAT;
BEGIN
RETURN RCreate( i.numerator - i.denominator, i.denominator );
END RDec;
PROCEDURE REQ( i : RAT; j : RAT ) : BOOLEAN;
VAR
ii : RAT;
jj : RAT;
BEGIN
ii := RNormalize( i );
jj := RNormalize( j );
RETURN ( ( ii.numerator = jj.numerator ) AND ( ii.denominator = jj.denominator ) );
END REQ;
PROCEDURE RNE( i : RAT; j : RAT ) : BOOLEAN;
BEGIN
RETURN NOT REQ( i, j );
END RNE;
PROCEDURE RLT( i : RAT; j : RAT ) : BOOLEAN;
BEGIN
RETURN RIsNegative( RSubtract( i, j ) );
END RLT;
PROCEDURE RLE( i : RAT; j : RAT ) : BOOLEAN;
BEGIN
RETURN NOT RGT( i, j );
END RLE;
PROCEDURE RGT( i : RAT; j : RAT ) : BOOLEAN;
BEGIN
RETURN RIsPositive( RSubtract( i, j ) );
END RGT;
PROCEDURE RGE( i : RAT; j : RAT ) : BOOLEAN;
BEGIN
RETURN NOT RLT( i, j );
END RGE;
PROCEDURE RIsZero( i : RAT ) : BOOLEAN;
BEGIN
RETURN i.numerator = 0;
END RIsZero;
PROCEDURE RIsNegative( i : RAT ) : BOOLEAN;
BEGIN
RETURN i.numerator < 0;
END RIsNegative;
PROCEDURE RIsPositive( i : RAT ) : BOOLEAN;
BEGIN
RETURN i.numerator > 0;
END RIsPositive;
PROCEDURE RToString( i : RAT; VAR S : ARRAY OF CHAR );
VAR
num : ARRAY [1..15] OF CHAR;
den : ARRAY [1..15] OF CHAR;
BEGIN
IF RIsZero( i ) THEN
Assign("0", S );
ELSE
IntToString( i.numerator, num, 1 );
Assign( num, S );
IF ( i.denominator <> 1 ) THEN
IntToString( i.denominator, den, 1 );
Append( S, "/" );
Append( S, den );
END;
END;
END RToString;
PROCEDURE RToRational( s : ARRAY OF CHAR ) : RAT;
VAR
n : CARDINAL;
numer : INTEGER;
denom : INTEGER;
LHS, RHS : ARRAY [ 1..20 ] OF CHAR;
Flag : BOOLEAN;
BEGIN
numer := 0;
denom := 0;
n := Pos( "/", s );
IF n > HIGH( s ) THEN
StringToInt( s, numer, Flag );
IF Flag THEN
denom := 1;
END;
ELSE
Copy( s, 0, n, LHS );
Copy( s, n+1, Length( s ), RHS );
StringToInt( LHS, numer, Flag );
IF Flag THEN
StringToInt( RHS, denom, Flag );
END;
END;
RETURN RCreate( numer, denom );
END RToRational;
PROCEDURE WriteRational( i : RAT );
VAR
res : ARRAY [0 .. 80] OF CHAR;
BEGIN
RToString( i, res );
WriteString( res );
END WriteRational;
END Rational.</lang>
- Test Program
<lang modula2>MODULE TestRat;
FROM InOut IMPORT WriteString, WriteLine;
FROM Terminal IMPORT KeyPressed;
FROM Strings IMPORT CompareStr;
FROM Rational IMPORT RAT, IGCD, RCreate, RToString, RIsZero, RNormalize,
RToRational, REQ, RNE, RLT, RLE, RGT, RGE, WriteRational,
RAdd, RSubtract, RMultiply, RDivide, RAbs, RNeg, RInv;
VAR
res : INTEGER; a, b, c, d, e, f : RAT; ans : ARRAY [1..100] OF CHAR;
PROCEDURE Assert( F : BOOLEAN; S : ARRAY OF CHAR ); BEGIN
IF ( NOT F) THEN
WriteLine( S );
END;
END Assert;
BEGIN
a := RCreate( 0, 0 ); Assert( RIsZero( a ), "RIsZero( a )");
a := RToRational("2");
RToString( a, ans );
res := CompareStr( ans, "2" );
Assert( (res = 0), "CompareStr( RToString( a ), '2' ) = 0");
a := RToRational("1/2");
RToString( a, ans );
res := CompareStr( ans, "1/2");
Assert( res = 0, "CompareStr( RToString( a, ans ), '1/2') = 0");
b := RToRational( "2/-12" ); RToString( b, ans ); res := CompareStr( ans, "-1/6"); Assert( res = 0, "CompareStr( RToString( b, ans ), '-1/6') = 0");
f := RCreate( 0, 9 ); (* rationalizes internally to zero *)
a := RToRational("1/3");
b := RToRational("1/2");
c := RCreate( 1, 3 );
Assert( NOT REQ( a, b ), "1/3 == 1/2" ); Assert( REQ( a, c ), "1/3 == 1/3" ); Assert( RNE( a, b ), "1/3 != 1/2" ); Assert( RLT( a, b ), "1/3 < 1/2" ); Assert( NOT RLT(b,a), "1/2 < 1/3" ); Assert( NOT RLT(a,c), "1/3 < 1/3" ); Assert( NOT RGT(a,b), "1/3 > 1/2" ); Assert( RGT(b,a), "1/2 > 1/3" ); Assert( NOT RGT(a,c), "1/3 > 1/3" );
Assert( RLE( a, b ), "1/3 <= 1/2" ); Assert( NOT RLE( b, a ), "1/2 <= 1/3" ); Assert( RLE( a, c ), "1/3 <= 1/3" ); Assert( NOT RGE(a,b), "1/3 >= 1/2" ); Assert( RGE(b,a), "1/2 >= 1/3" ); Assert( RGE( a,c ), "1/3 >= 1/3" );
a := RCreate(1,2);
b := RCreate(1,6);
a := RAdd( a, b );
Assert( REQ( a, RToRational("2/3")), "1/2 + 1/6 == 2/3" );
c := RNeg( a ); Assert( REQ( a, RCreate( 2,3)), "2/3 == 2/3" ); Assert( REQ( c, RCreate( 2,-3)), "Neg 1/2 == -1/2" ); a := RCreate( 2,-3);
d := RAbs( c ); Assert( REQ( d, RCreate( 2,3 ) ), "abs(neg(1/2))==1/2" );
a := RToRational( "1/2"); b := RSubtract( b, a );
Assert( REQ( b, RCreate(-1,3) ), "1/6 - 1/2 == -1/3" );
c := RInv(b); RToString( c, ans ); res := CompareStr( ans, "-3" ); Assert( res = 0, "inv(1/6 - 1/2) == -3" );
f := RInv( f ); (* as f normalized to zero, the reciprocal is still zero *)
b := RCreate( 1, 6); b := RAdd( b, RAdd( RCreate( 2,3), RCreate( 4,2 ))); RToString( b, ans ); res := CompareStr( ans, "17/6" ); Assert( res = 0, "1/6 + 2/3 + 4/2 = 17/6" );
a := RCreate(1,3);
b := RCreate(1,6);
c := RCreate(5,6);
d := RToRational("1/5");
e := RToRational("2");
f := RToRational("0/9");
Assert( REQ( RMultiply( c, d ), b ), "5/6 * 1/5 = 1/6" );
Assert( REQ( RMultiply( c, RMultiply( d, e ) ), a ), "5/6 * 1/5 * 2 = 1/3" );
Assert( REQ( RMultiply( c, RMultiply( d, RMultiply( e, f ) ) ), f ), "5/6 * 1/5 * 2 * 0" );
Assert( REQ( b, RDivide( c, RToRational("5" ) ) ), "5/6 / 5 = 1/6" );
e := RDivide( c, f ); (* RDivide multiplies so no divide by zero *)
WriteString("Press any key..."); WHILE NOT KeyPressed() DO END;
END TestRat.</lang>
Objective-C
File frac.h
#import <Foundation/Foundation.h>
@interface RCRationalNumber : NSObject
{
@private
int numerator;
int denominator;
BOOL autoSimplify;
BOOL withSign;
}
+(instancetype)valueWithNumerator:(int)num andDenominator: (int)den;
+(instancetype)valueWithDouble: (double)fnum;
+(instancetype)valueWithInteger: (int)inum;
+(instancetype)valueWithRational: (RCRationalNumber *)rnum;
-(instancetype)initWithNumerator: (int)num andDenominator: (int)den;
-(instancetype)initWithDouble: (double)fnum precision: (int)prec;
-(instancetype)initWithInteger: (int)inum;
-(instancetype)initWithRational: (RCRationalNumber *)rnum;
-(NSComparisonResult)compare: (RCRationalNumber *)rnum;
-(id)simplify: (BOOL)act;
-(void)setAutoSimplify: (BOOL)v;
-(void)setWithSign: (BOOL)v;
-(BOOL)autoSimplify;
-(BOOL)withSign;
-(NSString *)description;
// ops
-(id)multiply: (RCRationalNumber *)rnum;
-(id)divide: (RCRationalNumber *)rnum;
-(id)add: (RCRationalNumber *)rnum;
-(id)sub: (RCRationalNumber *)rnum;
-(id)abs;
-(id)neg;
-(id)mod: (RCRationalNumber *)rnum;
-(int)sign;
-(BOOL)isNegative;
-(id)reciprocal;
// getter
-(int)numerator;
-(int)denominator;
//setter
-(void)setNumerator: (int)num;
-(void)setDenominator: (int)num;
// defraction
-(double)number;
-(int)integer;
@end
- File frac.m
#import <Foundation/Foundation.h>
#import <math.h>
#import "frac.h"
// gcd: [[Greatest common divisor#Recursive_Euclid_algorithm]]
// if built in as "private" function, add static.
static int lcm(int a, int b)
{
return a / gcd(a,b) * b;
}
@implementation RCRationalNumber
// initializers
-(instancetype)init
{
NSLog(@"initialized to unity");
return [self initWithInteger: 1];
}
-(instancetype)initWithNumerator: (int)num andDenominator: (int)den
{
if ((self = [super init]) != nil) {
if (den == 0) {
NSLog(@"denominator is zero");
return nil;
}
[self setNumerator: num];
[self setDenominator: den];
[self setWithSign: YES];
[self setAutoSimplify: YES];
[self simplify: YES];
}
return self;
}
-(instancetype)initWithInteger:(int)inum
{
return [self initWithNumerator: inum andDenominator: 1];
}
-(instancetype)initWithDouble: (double)fnum precision: (int)prec
{
if ( prec > 9 ) prec = 9;
double p = pow(10.0, (double)prec);
int nd = (int)(fnum * p);
return [self initWithNumerator: nd andDenominator: (int)p ];
}
-(instancetype)initWithRational: (RCRationalNumber *)rnum
{
return [self initWithNumerator: [rnum numerator] andDenominator: [rnum denominator]];
}
// comparing
-(NSComparisonResult)compare: (RCRationalNumber *)rnum
{
if ( [self number] > [rnum number] ) return NSOrderedDescending;
if ( [self number] < [rnum number] ) return NSOrderedAscending;
return NSOrderedSame;
}
// string rapresentation of the Q
-(NSString *)description
{
[self simplify: [self autoSimplify]];
return [NSString stringWithFormat: @"%@%d/%d", [self isNegative] ? @"-" :
( [self withSign] ? @"+" : @"" ),
abs([self numerator]), [self denominator]];
}
// setter options
-(void)setAutoSimplify: (BOOL)v
{
autoSimplify = v;
[self simplify: v];
}
-(void)setWithSign: (BOOL)v
{
withSign = v;
}
// getter for options
-(BOOL)autoSimplify
{
return autoSimplify;
}
-(BOOL)withSign
{
return withSign;
}
// "simplify" the fraction ...
-(id)simplify: (BOOL)act
{
if ( act ) {
int common = gcd([self numerator], [self denominator]);
[self setNumerator: [self numerator]/common];
[self setDenominator: [self denominator]/common];
}
return self;
}
// diadic operators
-(id)multiply: (RCRationalNumber *)rnum
{
int newnum = [self numerator] * [rnum numerator];
int newden = [self denominator] * [rnum denominator];
return [RCRationalNumber valueWithNumerator: newnum
andDenominator: newden];
}
-(id)divide: (RCRationalNumber *)rnum
{
return [self multiply: [rnum reciprocal]];
}
-(id)add: (RCRationalNumber *)rnum
{
int common = lcm([self denominator], [rnum denominator]);
int resnum = common / [self denominator] * [self numerator] +
common / [rnum denominator] * [rnum numerator];
return [RCRationalNumber valueWithNumerator: resnum andDenominator: common];
}
-(id)sub: (RCRationalNumber *)rnum
{
return [self add: [rnum neg]];
}
-(id)mod: (RCRationalNumber *)rnum
{
return [[self divide: rnum]
sub: [RCRationalNumber valueWithInteger: [[self divide: rnum] integer]]];
}
// unary operators
-(id)neg
{
return [RCRationalNumber valueWithNumerator: -1*[self numerator]
andDenominator: [self denominator]];
}
-(id)abs
{
return [RCRationalNumber valueWithNumerator: abs([self numerator])
andDenominator: [self denominator]];
}
-(id)reciprocal
{
return [RCRationalNumber valueWithNumerator: [self denominator]
andDenominator: [self numerator]];
}
// get the sign
-(int)sign
{
return ([self numerator] < 0) ? -1 : 1;
}
// or just test if negative
-(BOOL)isNegative
{
return [self numerator] < 0;
}
// Q as real floating point
-(double)number
{
return (double)[self numerator] / (double)[self denominator];
}
// Q as (truncated) integer
-(int)integer
{
return [self numerator] / [self denominator];
}
// set num and den indipendently, fixing sign accordingly
-(void)setNumerator: (int)num
{
numerator = num;
}
-(void)setDenominator: (int)num
{
if ( num < 0 ) numerator = -numerator;
denominator = abs(num);
}
// getter
-(int)numerator
{
return numerator;
}
-(int)denominator
{
return denominator;
}
// class method
+(instancetype)valueWithNumerator:(int)num andDenominator: (int)den
{
return [[self alloc] initWithNumerator: num andDenominator: den];
}
+(instancetype)valueWithDouble: (double)fnum
{
return [[self alloc] initWithDouble: fnum];
}
+(instancetype)valueWithInteger: (int)inum
{
return [[self alloc] initWithInteger: inum];
}
+(instancetype)valueWithRational: (RCRationalNumber *)rnum
{
return [[self alloc] initWithRational: rnum];
}
@end
- Testing
#import <Foundation/Foundation.h>
#import "frac.h"
#import <math.h>
int main()
{
@autoreleasepool {
int i;
for(i=2; i < 0x80000; i++) {
int candidate = i;
RCRationalNumber *sum = [RCRationalNumber valueWithNumerator: 1
andDenominator: candidate];
int factor;
for(factor=2; factor < sqrt((double)candidate); factor++) {
if ( (candidate % factor) == 0 ) {
sum = [[sum add: [RCRationalNumber valueWithNumerator: 1
andDenominator: factor]]
add: [RCRationalNumber valueWithNumerator: 1
andDenominator: (candidate/factor)]];
}
}
if ( [sum denominator] == 1 ) {
printf("Sum of recipr. factors of %d = %d exactly %s\n",
candidate, [sum integer], ([sum integer]==1) ? "perfect!" : "");
}
}
}
return 0;
}
OCaml
OCaml's Num library implements arbitrary-precision rational numbers: <lang ocaml>#load "nums.cma";; open Num;;
for candidate = 2 to 1 lsl 19 do
let sum = ref (num_of_int 1 // num_of_int candidate) in
for factor = 2 to truncate (sqrt (float candidate)) do
if candidate mod factor = 0 then
sum := !sum +/ num_of_int 1 // num_of_int factor
+/ num_of_int 1 // num_of_int (candidate / factor)
done;
if is_integer_num !sum then
Printf.printf "Sum of recipr. factors of %d = %d exactly %s\n%!"
candidate (int_of_num !sum) (if int_of_num !sum = 1 then "perfect!" else "")
done;;</lang> Delimited overloading can be used to make the arithmetic expressions more readable: <lang ocaml>let () =
for candidate = 2 to 1 lsl 19 do
let sum = ref Num.(1 / of_int candidate) in
for factor = 2 to truncate (sqrt (float candidate)) do
if candidate mod factor = 0 then
sum := Num.(!sum + 1 / of_int factor + of_int factor / of_int candidate)
done;
if Num.is_integer_num !sum then
Printf.printf "Sum of recipr. factors of %d = %d exactly %s\n%!"
candidate Num.(to_int !sum) (if Num.(!sum = 1) then "perfect!" else "")
done</lang>
It might be implemented like this:
[insert implementation here]
ooRexx
<lang ooRexx> loop candidate = 6 to 2**19
sum = .fraction~new(1, candidate) max2 = rxcalcsqrt(candidate)~trunc
loop factor = 2 to max2
if candidate // factor == 0 then do
sum += .fraction~new(1, factor)
sum += .fraction~new(1, candidate / factor)
end
end
if sum == 1 then say candidate "is a perfect number"
end
- class fraction inherit orderable
- method init
expose numerator denominator use strict arg numerator, denominator = 1
if numerator == 0 then denominator = 0
else if denominator == 0 then raise syntax 98.900 array("Fraction denominator cannot be zero")
-- if the denominator is negative, make the numerator carry the sign
if denominator < 0 then do
numerator = -numerator
denominator = - denominator
end
-- find the greatest common denominator and reduce to -- the simplest form gcd = self~gcd(numerator~abs, denominator~abs)
numerator /= gcd denominator /= gcd
-- fraction instances are immutable, so these are -- read only attributes
- attribute numerator GET
- attribute denominator GET
-- calculate the greatest common denominator of a numerator/denominator pair
- method gcd private
use arg x, y
loop while y \= 0
-- check if they divide evenly
temp = x // y
x = y
y = temp
end
return x
-- calculate the least common multiple of a numerator/denominator pair
- method lcm private
use arg x, y return x / self~gcd(x, y) * y
- method abs
expose numerator denominator -- the denominator is always forced to be positive return self~class~new(numerator~abs, denominator)
- method reciprocal
expose numerator denominator return self~class~new(denominator, numerator)
-- convert a fraction to regular Rexx number
- method toNumber
expose numerator denominator
if numerator == 0 then return 0 return numerator/denominator
- method negative
expose numerator denominator return self~class~new(-numerator, denominator)
- method add
expose numerator denominator use strict arg other -- convert to a fraction if a regular number if \other~isa(.fraction) then other = self~class~new(other, 1)
multiple = self~lcm(denominator, other~denominator) newa = numerator * multiple / denominator newb = other~numerator * multiple / other~denominator return self~class~new(newa + newb, multiple)
- method subtract
use strict arg other return self + (-other)
- method times
expose numerator denominator use strict arg other -- convert to a fraction if a regular number if \other~isa(.fraction) then other = self~class~new(other, 1) return self~class~new(numerator * other~numerator, denominator * other~denominator)
- method divide
use strict arg other -- convert to a fraction if a regular number if \other~isa(.fraction) then other = self~class~new(other, 1) -- and multiply by the reciprocal return self * other~reciprocal
-- compareTo method used by the orderable interface to implement -- the operator methods
- method compareTo
expose numerator denominator -- convert to a fraction if a regular number if \other~isa(.fraction) then other = self~class~new(other, 1)
return (numerator * other~denominator - denominator * other~numerator)~sign
-- we still override "==" and "\==" because we want to bypass the -- checks for not being an instance of the class
- method "=="
expose numerator denominator use strict arg other
-- convert to a fraction if a regular number if \other~isa(.fraction) then other = self~class~new(other, 1) -- Note: these are numeric comparisons, so we're using the "=" -- method so those are handled correctly return numerator = other~numerator & denominator = other~denominator
- method "\=="
use strict arg other return \self~"\=="(other)
-- some operator overrides -- these only work if the left-hand-side of the -- subexpression is a quaternion
- method "*"
forward message("TIMES")
- method "/"
forward message("DIVIDE")
- method "-"
-- need to check if this is a prefix minus or a subtract
if arg() == 0 then
forward message("NEGATIVE")
else
forward message("SUBTRACT")
- method "+"
-- need to check if this is a prefix plus or an addition
if arg() == 0 then
return self -- we can return this copy since it is imutable
else
forward message("ADD")
- method string
expose numerator denominator if denominator == 1 then return numerator return numerator"/"denominator
-- override hashcode for collection class hash uses
- method hashCode
expose numerator denominator return numerator~hashcode~bitxor(numerator~hashcode)
- requires rxmath library
</lang> Output:
6 is a perfect number 28 is a perfect number 496 is a perfect number 8128 is a perfect number
PARI/GP
Pari handles rational arithmetic natively. <lang parigp>for(n=2,1<<19,
s=0; fordiv(n,d,s+=1/d); if(s==2,print(n))
)</lang>
Perl
Perl's Math::BigRat core module implements arbitrary-precision rational numbers. The bigrat pragma can be used to turn on transparent BigRat support:
<lang perl>use bigrat;
foreach my $candidate (2 .. 2**19) {
my $sum = 1 / $candidate;
foreach my $factor (2 .. sqrt($candidate)+1) {
if ($candidate % $factor == 0) {
$sum += 1 / $factor + 1 / ($candidate / $factor);
}
}
if ($sum->denominator() == 1) {
print "Sum of recipr. factors of $candidate = $sum exactly ", ($sum == 1 ? "perfect!" : ""), "\n";
}
}</lang> It might be implemented like this:
[insert implementation here]
Perl 6
Perl 6 supports rational arithmetic natively. <lang perl6>for 2..2**19 -> $candidate {
my $sum = 1 / $candidate;
for 2 .. ceiling(sqrt($candidate)) -> $factor {
if $candidate %% $factor {
$sum += 1 / $factor + 1 / ($candidate / $factor);
}
}
if $sum.denominator == 1 {
say "Sum of reciprocal factors of $candidate = $sum exactly", ($sum == 1 ?? ", perfect!" !! ".");
}
}</lang> Note also that ordinary decimal literals are stored as Rats, so the following loop always stops exactly on 10 despite 0.1 not being exactly representable in floating point: <lang perl6>for 1.0, 1.1, 1.2 ... 10 { .say }</lang> The arithmetic is all done in rationals, which are converted to floating-point just before display so that people don't have to puzzle out what 53/10 means.
PicoLisp
<lang PicoLisp>(load "@lib/frac.l")
(for (N 2 (> (** 2 19) N) (inc N))
(let (Sum (frac 1 N) Lim (sqrt N))
(for (F 2 (>= Lim F) (inc F))
(when (=0 (% N F))
(setq Sum
(f+ Sum
(f+ (frac 1 F) (frac 1 (/ N F))) ) ) ) )
(when (= 1 (cdr Sum))
(prinl
"Perfect " N
", sum is " (car Sum)
(and (= 1 (car Sum)) ": perfect") ) ) ) )</lang>
- Output:
Perfect 6, sum is 1: perfect Perfect 28, sum is 1: perfect Perfect 120, sum is 2 Perfect 496, sum is 1: perfect Perfect 672, sum is 2 Perfect 8128, sum is 1: perfect Perfect 30240, sum is 3 Perfect 32760, sum is 3 Perfect 523776, sum is 2
Python
Python 3's standard library already implements a Fraction class: <lang python>from fractions import Fraction
for candidate in range(2, 2**19):
sum = Fraction(1, candidate)
for factor in range(2, int(candidate**0.5)+1):
if candidate % factor == 0:
sum += Fraction(1, factor) + Fraction(1, candidate // factor)
if sum.denominator == 1:
print("Sum of recipr. factors of %d = %d exactly %s" %
(candidate, int(sum), "perfect!" if sum == 1 else ""))</lang>
It might be implemented like this: <lang python>def lcm(a, b):
return a // gcd(a,b) * b
def gcd(u, v):
return gcd(v, u%v) if v else abs(u)
class Fraction:
def __init__(self, numerator, denominator):
common = gcd(numerator, denominator)
self.numerator = numerator//common
self.denominator = denominator//common
def __add__(self, frac):
common = lcm(self.denominator, frac.denominator)
n = common // self.denominator * self.numerator + common // frac.denominator * frac.numerator
return Fraction(n, common)
def __sub__(self, frac):
return self.__add__(-frac)
def __neg__(self):
return Fraction(-self.numerator, self.denominator)
def __abs__(self):
return Fraction(abs(self.numerator), abs(self.denominator))
def __mul__(self, frac):
return Fraction(self.numerator * frac.numerator, self.denominator * frac.denominator)
def __div__(self, frac):
return self.__mul__(frac.reciprocal())
def reciprocal(self):
return Fraction(self.denominator, self.numerator)
def __cmp__(self, n):
return int(float(self) - float(n))
def __float__(self):
return float(self.numerator / self.denominator)
def __int__(self):
return (self.numerator // self.denominator)</lang>
Ruby
Ruby's standard library already implements a Rational class: <lang ruby>require 'rational'
for candidate in 2 .. 2**19:
sum = Rational(1, candidate)
for factor in 2 ... candidate**0.5
if candidate % factor == 0
sum += Rational(1, factor) + Rational(1, candidate / factor)
end
end
if sum.denominator == 1
puts "Sum of recipr. factors of %d = %d exactly %s" %
[candidate, sum.to_i, sum == 1 ? "perfect!" : ""]
end
end</lang> It might be implemented like this:
[insert implementation here]
Scheme
Scheme has native rational numbers.
<lang scheme>; simply prints all the perfect numbers (do ((candidate 2 (+ candidate 1))) ((>= candidate (expt 2 19)))
(let ((sum (/ 1 candidate)))
(do ((factor 2 (+ factor 1))) ((>= factor (sqrt candidate)))
(if (= 0 (modulo candidate factor))
(set! sum (+ sum (/ 1 factor) (/ factor candidate)))))
(if (= 1 (denominator sum))
(begin (display candidate) (newline)))))</lang>
It might be implemented like this:
[insert implementation here]
Scala
<lang scala>class Rational(n: Long, d:Long) extends Ordered[Rational] {
require(d!=0) private val g:Long = gcd(n, d) val numerator:Long = n/g val denominator:Long = d/g
def this(n:Long)=this(n,1)
def +(that:Rational):Rational=new Rational(
numerator*that.denominator + that.numerator*denominator,
denominator*that.denominator)
def -(that:Rational):Rational=new Rational(
numerator*that.denominator - that.numerator*denominator,
denominator*that.denominator)
def *(that:Rational):Rational=
new Rational(numerator*that.numerator, denominator*that.denominator)
def /(that:Rational):Rational=
new Rational(numerator*that.denominator, that.numerator*denominator)
def unary_~ :Rational=new Rational(denominator, numerator)
def unary_- :Rational=new Rational(-numerator, denominator)
def abs :Rational=new Rational(Math.abs(numerator), Math.abs(denominator))
override def compare(that:Rational):Int=
(this.numerator*that.denominator-that.numerator*this.denominator).toInt
override def toString()=numerator+"/"+denominator
private def gcd(x:Long, y:Long):Long=
if(y==0) x else gcd(y, x%y)
}
object Rational {
def apply(n: Long, d:Long)=new Rational(n,d) def apply(n:Long)=new Rational(n) implicit def longToRational(i:Long)=new Rational(i)
}</lang>
<lang scala>def find_perfects():Unit= {
for (candidate <- 2 until 1<<19)
{
var sum= ~Rational(candidate)
for (factor <- 2 until (Math.sqrt(candidate)+1).toInt)
{
if (candidate%factor==0)
sum+= ~Rational(factor)+ ~Rational(candidate/factor)
}
if (sum.denominator==1 && sum.numerator==1)
printf("Perfect number %d sum is %s\n", candidate, sum)
}
}</lang>
Slate
Slate uses infinite-precision fractions transparently. <lang slate>54 / 7. 20 reciprocal. (5 / 6) reciprocal. (5 / 6) as: Float.</lang>
Smalltalk
Smalltalk uses naturally and transparently fractions (through the class Fraction):
st> 54/7 54/7 st> 54/7 + 1 61/7 st> 54/7 < 50 true st> 20 reciprocal 1/20 st> (5/6) reciprocal 6/5 st> (5/6) asFloat 0.8333333333333334
<lang smalltalk>| sum | 2 to: (2 raisedTo: 19) do: [ :candidate |
sum := candidate reciprocal.
2 to: (candidate sqrt) do: [ :factor |
( (candidate \\ factor) = 0 )
ifTrue: [
sum := sum + (factor reciprocal) + ((candidate / factor) reciprocal)
]
].
( (sum denominator) = 1 )
ifTrue: [
('Sum of recipr. factors of %1 = %2 exactly %3' %
{ candidate printString .
(sum asInteger) printString .
( sum = 1 ) ifTrue: [ 'perfect!' ]
ifFalse: [ ' ' ] }) displayNl
]
].</lang>
Tcl
Code to find factors of a number not shown:
namespace eval rat {}
proc rat::new {args} {
if {[llength $args] == 0} {
set args {0}
}
lassign [split {*}$args] n d
if {$d == 0} {
error "divide by zero"
}
if {$d < 0} {
set n [expr {-1 * $n}]
set d [expr {abs($d)}]
}
return [normalize $n $d]
}
proc rat::split {args} {
if {[llength $args] == 1} {
lassign [::split $args /] n d
if {$d eq ""} {
set d 1
}
} else {
lassign $args n d
}
return [list $n $d]
}
proc rat::join {rat} {
lassign $rat n d
if {$n == 0} {
return 0
} elseif {$d == 1} {
return $n
} else {
return $n/$d
}
}
proc rat::normalize {n d} {
set gcd [gcd $n $d]
return [join [list [expr {$n/$gcd}] [expr {$d/$gcd}]]]
}
proc rat::gcd {a b} {
while {$b != 0} {
lassign [list $b [expr {$a % $b}]] a b
}
return $a
}
proc rat::abs {rat} {
lassign [split $rat] n d
return [join [list [expr {abs($n)}] $d]]
}
proc rat::inv {rat} {
lassign [split $rat] n d
return [normalize $d $n]
}
proc rat::+ {args} {
set n 0
set d 1
foreach arg $args {
lassign [split $arg] an ad
set n [expr {$n*$ad + $an*$d}]
set d [expr {$d * $ad}]
}
return [normalize $n $d]
}
proc rat::- {args} {
lassign [split [lindex $args 0]] n d
if {[llength $args] == 1} {
return [join [list [expr {-1 * $n}] $d]]
}
foreach arg [lrange $args 1 end] {
lassign [split $arg] an ad
set n [expr {$n*$ad - $an*$d}]
set d [expr {$d * $ad}]
}
return [normalize $n $d]
}
proc rat::* {args} {
set n 1
set d 1
foreach arg $args {
lassign [split $arg] an ad
set n [expr {$n * $an}]
set d [expr {$d * $ad}]
}
return [normalize $n $d]
}
proc rat::/ {a b} {
set r [* $a [inv $b]]
if {[string match */0 $r]} {
error "divide by zero"
}
return $r
}
proc rat::== {a b} {
return [expr {[- $a $b] == 0}]
}
proc rat::!= {a b} {
return [expr { ! [== $a $b]}]
}
proc rat::< {a b} {
lassign [split [- $a $b]] n d
return [expr {$n < 0}]
}
proc rat::> {a b} {
lassign [split [- $a $b]] n d
return [expr {$n > 0}]
}
proc rat::<= {a b} {
return [expr { ! [> $a $b]}]
}
proc rat::>= {a b} {
return [expr { ! [< $a $b]}]
}
################################################
proc is_perfect {num} {
set sum [rat::new 0]
foreach factor [all_factors $num] {
set sum [rat::+ $sum [rat::new 1/$factor]]
}
# note, all_factors includes 1, so sum should be 2
return [rat::== $sum 2]
}
proc get_perfect_numbers {} {
set t [clock seconds]
set limit [expr 2**19]
for {set num 2} {$num < $limit} {incr num} {
if {[is_perfect $num]} {
puts "perfect: $num"
}
}
puts "elapsed: [expr {[clock seconds] - $t}] seconds"
set num [expr {2**12 * (2**13 - 1)}] ;# 5th perfect number
if {[is_perfect $num]} {
puts "perfect: $num"
}
}
source primes.tcl
get_perfect_numbers
- Output:
perfect: 6 perfect: 28 perfect: 496 perfect: 8128 elapsed: 477 seconds perfect: 33550336
TI-89 BASIC
While TI-89 BASIC has built-in rational and symbolic arithmetic, it does not have user-defined data types.