Continued fraction/Arithmetic/G(matrix ng, continued fraction n1, continued fraction n2): Difference between revisions
Continued fraction/Arithmetic/G(matrix ng, continued fraction n1, continued fraction n2) (view source)
Revision as of 17:20, 20 November 2023
, 5 months ago→{{header|Wren}}: Changed to Wren S/H
m (→{{header|Wren}}: Changed to Wren S/H) |
|||
| (13 intermediate revisions by 3 users not shown) | |||
Line 1:
{{
This task performs the basic mathematical functions on 2 continued fractions. This requires the full version of matrix NG:
: <math>\begin{bmatrix}
Line 1,823:
===Using multiple precision numbers===
{{trans|Standard ML}}
{{libheader|ats2-xprelude}}
{{libheader|GMP}}
For this program you need [https://sourceforge.net/p/chemoelectric/ats2-xprelude ats2-xprelude].
Line 6,124 ⟶ 6,126:
484/49 / 22/7 is {0 1 0 0 0 0 1 0} [9; 1, 7, 6] [3; 7] gives [3; 7]
√2 * √2 = [1; 0, 1]
</pre>
=={{header|Haskell}}==
{{trans|Mercury}}
This Haskell follows the Mercury, in using infinitely long lazy lists to represent continued fractions. There are two kinds of terms: "infinite" and "finite integer".
<syntaxhighlight lang="Haskell">
----------------------------------------------------------------------
data Term = InfiniteTerm | IntegerTerm Integer
type ContinuedFraction = [Term] -- The list should be infinitely long.
type NG8 = (Integer, Integer, Integer, Integer,
Integer, Integer, Integer, Integer)
----------------------------------------------------------------------
cf2string (cf :: ContinuedFraction) =
loop 0 "[" cf
where loop i s lst =
case lst of {
(InfiniteTerm : _) -> s ++ "]" ;
(IntegerTerm value : tail) ->
(if i == 20 then
s ++ ",...]"
else
let {
sepStr =
case i of {
0 -> "";
1 -> ";";
_ -> ","
};
termStr = show value;
s1 = s ++ sepStr ++ termStr
}
in loop (i + 1) s1 tail)
}
----------------------------------------------------------------------
repeatingTerm (term :: Term) =
term : repeatingTerm term
infiniteContinuedFraction = repeatingTerm InfiniteTerm
i2cf (i :: Integer) =
-- Continued fraction representing an integer.
IntegerTerm i : infiniteContinuedFraction
r2cf (n :: Integer) (d :: Integer) =
-- Continued fraction representing a rational number.
let (q, r) = divMod n d in
(if r == 0 then
(IntegerTerm q : infiniteContinuedFraction)
else
(IntegerTerm q : r2cf d r))
----------------------------------------------------------------------
add_cf = apply_ng8 (0, 1, 1, 0, 0, 0, 0, 1)
sub_cf = apply_ng8 (0, 1, -1, 0, 0, 0, 0, 1)
mul_cf = apply_ng8 (1, 0, 0, 0, 0, 0, 0, 1)
div_cf = apply_ng8 (0, 1, 0, 0, 0, 0, 1, 0)
apply_ng8
(ng :: NG8)
(x :: ContinuedFraction)
(y :: ContinuedFraction) =
--
let (a12, a1, a2, a, b12, b1, b2, b) = ng in
if iseqz [b12, b1, b2, b] then
infiniteContinuedFraction -- No more finite terms to output.
else if iseqz [b2, b] then
let (ng1, x1, y1) = absorb_x_term ng x y in
apply_ng8 ng1 x1 y1
else if atLeastOne_iseqz [b2, b] then
let (ng1, x1, y1) = absorb_y_term ng x y in
apply_ng8 ng1 x1 y1
else if iseqz [b1] then
let (ng1, x1, y1) = absorb_x_term ng x y in
apply_ng8 ng1 x1 y1
else
let {
(q12, r12) = maybeDivide a12 b12;
(q1, r1) = maybeDivide a1 b1;
(q2, r2) = maybeDivide a2 b2;
(q, r) = maybeDivide a b
}
in
if not (iseqz [b12]) && q == q12 && q == q1 && q == q2 then
-- Output a term.
(if integerExceedsInfinitizingThreshold q then
infiniteContinuedFraction
else
let new_ng = (b12, b1, b2, b, r12, r1, r2, r) in
(IntegerTerm q : apply_ng8 new_ng x y))
else
-- Put a1, a2, and a over a common denominator and compare
-- some magnitudes.
let {
n1 = a1 * b2 * b;
n2 = a2 * b1 * b;
n = a * b1 * b2
}
in
(if abs (n1 - n) > abs (n2 - n) then
let (ng1, x1, y1) = absorb_x_term ng x y in
apply_ng8 ng1 x1 y1
else
let (ng1, x1, y1) = absorb_y_term ng x y in
apply_ng8 ng1 x1 y1)
absorb_x_term
((a12, a1, a2, a, b12, b1, b2, b) :: NG8)
(x :: ContinuedFraction)
(y :: ContinuedFraction) =
--
case x of {
(IntegerTerm n : xtail) -> (
let new_ng = (a2 + (a12 * n), a + (a1 * n), a12, a1,
b2 + (b12 * n), b + (b1 * n), b12, b1) in
if (ng8ExceedsProcessingThreshold new_ng) then
-- Pretend we have reached an infinite term.
((a12, a1, a12, a1, b12, b1, b12, b1),
infiniteContinuedFraction, y)
else
(new_ng, xtail, y)
);
(InfiniteTerm : _) ->
((a12, a1, a12, a1, b12, b1, b12, b1),
infiniteContinuedFraction, y)
}
absorb_y_term
((a12, a1, a2, a, b12, b1, b2, b) :: NG8)
(x :: ContinuedFraction)
(y :: ContinuedFraction) =
--
case y of {
(IntegerTerm n : ytail) -> (
let new_ng = (a1 + (a12 * n), a12, a + (a2 * n), a2,
b1 + (b12 * n), b12, b + (b2 * n), b2) in
if (ng8ExceedsProcessingThreshold new_ng) then
-- Pretend we have reached an infinite term.
((a12, a12, a2, a2, b12, b12, b2, b2),
x, infiniteContinuedFraction)
else
(new_ng, x, ytail)
);
(InfiniteTerm : _) ->
((a12, a12, a2, a2, b12, b12, b2, b2),
x, infiniteContinuedFraction)
}
ng8ExceedsProcessingThreshold (a12, a1, a2, a,
b12, b1, b2, b) =
(integerExceedsProcessingThreshold a12 ||
integerExceedsProcessingThreshold a1 ||
integerExceedsProcessingThreshold a2 ||
integerExceedsProcessingThreshold a ||
integerExceedsProcessingThreshold b12 ||
integerExceedsProcessingThreshold b1 ||
integerExceedsProcessingThreshold b2 ||
integerExceedsProcessingThreshold b)
integerExceedsProcessingThreshold i =
abs i >= 2 ^ 512
integerExceedsInfinitizingThreshold i =
abs i >= 2 ^ 64
maybeDivide a b =
if b == 0
then (0, 0)
else divMod a b
iseqz [] = True
iseqz (head : tail) = head == 0 && iseqz tail
atLeastOne_iseqz [] = False
atLeastOne_iseqz (head : tail) = head == 0 || atLeastOne_iseqz tail
----------------------------------------------------------------------
zero = i2cf 0
one = i2cf 1
two = i2cf 2
three = i2cf 3
four = i2cf 4
one_fourth = r2cf 1 4
one_third = r2cf 1 3
one_half = r2cf 1 2
two_thirds = r2cf 2 3
three_fourths = r2cf 3 4
goldenRatio = repeatingTerm (IntegerTerm 1)
silverRatio = repeatingTerm (IntegerTerm 2)
sqrt2 = IntegerTerm 1 : silverRatio
sqrt5 = IntegerTerm 2 : repeatingTerm (IntegerTerm 4)
----------------------------------------------------------------------
padLeft n s
| length s < n = replicate (n - length s) ' ' ++ s
| otherwise = s
padRight n s
| length s < n = s ++ replicate (n - length s) ' '
| otherwise = s
show_cf (expression, cf, note) =
let exprStr = padLeft 19 expression in
do { putStr exprStr;
putStr " => ";
if note == "" then
putStrLn (cf2string cf)
else
let cfStr = padRight 48 (cf2string cf) in
do { putStr cfStr;
putStrLn note }
}
thirteen_elevenths = r2cf 13 11
twentytwo_sevenths = r2cf 22 7
main = do {
show_cf ("golden ratio", goldenRatio, "(1 + sqrt(5))/2");
show_cf ("silver ratio", silverRatio, "(1 + sqrt(2))");
show_cf ("sqrt(2)", sqrt2, "from the module");
show_cf ("sqrt(2)", silverRatio `sub_cf` one,
"from the silver ratio");
show_cf ("sqrt(5)", sqrt5, "from the module");
show_cf ("sqrt(5)", (two `mul_cf` goldenRatio) `sub_cf` one,
"from the golden ratio");
show_cf ("13/11", thirteen_elevenths, "");
show_cf ("22/7", twentytwo_sevenths, "approximately pi");
show_cf ("13/11 + 1/2", thirteen_elevenths `add_cf` one_half, "");
show_cf ("22/7 + 1/2", twentytwo_sevenths `add_cf` one_half, "");
show_cf ("(22/7) * 1/2", twentytwo_sevenths `mul_cf` one_half, "");
show_cf ("(22/7) / 2", twentytwo_sevenths `div_cf` two, "");
show_cf ("sqrt(2) + sqrt(2)", sqrt2 `add_cf` sqrt2, "");
show_cf ("sqrt(2) - sqrt(2)", sqrt2 `sub_cf` sqrt2, "");
show_cf ("sqrt(2) * sqrt(2)", sqrt2 `mul_cf` sqrt2, "");
show_cf ("sqrt(2) / sqrt(2)", sqrt2 `div_cf` sqrt2, "");
return ()
}
----------------------------------------------------------------------
</syntaxhighlight>
{{out}}
<pre>$ ghc continued_fraction_task.hs && ./continued_fraction_task
[1 of 1] Compiling Main ( continued_fraction_task.hs, continued_fraction_task.o )
Linking continued_fraction_task ...
golden ratio => [1;1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,...] (1 + sqrt(5))/2
silver ratio => [2;2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,...] (1 + sqrt(2))
sqrt(2) => [1;2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,...] from the module
sqrt(2) => [1;2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,...] from the silver ratio
sqrt(5) => [2;4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,...] from the module
sqrt(5) => [2;4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,...] from the golden ratio
13/11 => [1;5,2]
22/7 => [3;7] approximately pi
13/11 + 1/2 => [1;1,2,7]
22/7 + 1/2 => [3;1,1,1,4]
(22/7) * 1/2 => [1;1,1,3]
(22/7) / 2 => [1;1,1,3]
sqrt(2) + sqrt(2) => [2;1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,...]
sqrt(2) - sqrt(2) => [0]
sqrt(2) * sqrt(2) => [2]
sqrt(2) / sqrt(2) => [1]
</pre>
Line 6,476 ⟶ 6,752:
sqrt(2) * sqrt(2) => [2]
sqrt(2) / sqrt(2) => [1]
</pre>
=={{header|Java}}==
<syntaxhighlight lang="java">
import java.util.ArrayList;
import java.util.List;
public final class ContinuedFractionArithmeticG2 {
public static void main(String[] aArgs) {
test("[3; 7] + [0; 2]", new NG( new NG8(0, 1, 1, 0, 0, 0, 0, 1), new R2cf(1, 2), new R2cf(22, 7) ),
new NG( new NG4(2, 1, 0, 2), new R2cf(22, 7) ));
test("[1; 5, 2] * [3; 7]", new NG( new NG8(1, 0, 0, 0, 0, 0, 0, 1), new R2cf(13, 11), new R2cf(22, 7) ),
new R2cf(286, 77) );
test("[1; 5, 2] - [3; 7]", new NG( new NG8(0, 1, -1, 0, 0, 0, 0, 1), new R2cf(13, 11), new R2cf(22, 7) ),
new R2cf(-151, 77) );
test("Divide [] by [3; 7]",
new NG( new NG8(0, 1, 0, 0, 0, 0, 1, 0), new R2cf(22 * 22, 7 * 7), new R2cf(22,7)) );
test("([0; 3, 2] + [1; 5, 2]) * ([0; 3, 2] - [1; 5, 2])",
new NG( new NG8(1, 0, 0, 0, 0, 0, 0, 1),
new NG( new NG8(0, 1, 1, 0, 0, 0, 0, 1),
new R2cf(2, 7), new R2cf(13, 11)),
new NG( new NG8(0, 1, -1, 0, 0, 0, 0, 1), new R2cf(2, 7), new R2cf(13, 11) ) ),
new R2cf(-7797, 5929) );
}
private static void test(String aDescription, ContinuedFraction... aFractions) {
System.out.println("Testing: " + aDescription);
for ( ContinuedFraction fraction : aFractions ) {
while ( fraction.hasMoreTerms() ) {
System.out.print(fraction.nextTerm() + " ");
}
System.out.println();
}
System.out.println();
}
private static abstract class MatrixNG {
protected abstract void consumeTerm();
protected abstract void consumeTerm(int aN);
protected abstract boolean needsTerm();
protected int configuration = 0;
protected int currentTerm = 0;
protected boolean hasTerm = false;
}
private static class NG4 extends MatrixNG {
public NG4(int aA1, int aA, int aB1, int aB) {
a1 = aA1; a = aA; b1 = aB1; b = aB;
}
public void consumeTerm() {
a = a1;
b = b1;
}
public void consumeTerm(int aN) {
int temp = a; a = a1; a1 = temp + a1 * aN;
temp = b; b = b1; b1 = temp + b1 * aN;
}
public boolean needsTerm() {
if ( b1 == 0 && b == 0 ) {
return false;
}
if ( b1 == 0 || b == 0 ) {
return true;
}
currentTerm = a / b;
if ( currentTerm == a1 / b1 ) {
int temp = a; a = b; b = temp - b * currentTerm;
temp = a1; a1 = b1; b1 = temp - b1 * currentTerm;
hasTerm = true;
return false;
}
return true;
}
private int a1, a, b1, b;
}
private static class NG8 extends MatrixNG {
public NG8(int aA12, int aA1, int aA2, int aA, int aB12, int aB1, int aB2, int aB) {
a12 = aA12; a1 = aA1; a2 = aA2; a = aA; b12 = aB12; b1 = aB1; b2 = aB2; b = aB;
}
public void consumeTerm() {
if ( configuration == 0 ) {
a = a1; a2 = a12;
b = b1; b2 = b12;
} else {
a = a2; a1 = a12;
b = b2; b1 = b12;
}
}
public void consumeTerm(int aN) {
if ( configuration == 0 ) {
int temp = a; a = a1; a1 = temp + a1 * aN;
temp = a2; a2 = a12; a12 = temp + a12 * aN;
temp = b; b = b1; b1 = temp + b1 * aN;
temp = b2; b2 = b12; b12 = temp + b12 * aN;
} else {
int temp = a; a = a2; a2 = temp + a2 * aN;
temp = a1; a1 = a12; a12 = temp + a12 * aN;
temp = b; b = b2; b2 = temp + b2 * aN;
temp = b1; b1 = b12; b12 = temp + b12 * aN;
}
}
public boolean needsTerm() {
if ( b1 == 0 && b == 0 && b2 == 0 && b12 == 0 ) {
return false;
}
if ( b == 0 ) {
configuration = ( b2 == 0 ) ? 0 : 1;
return true;
}
ab = (double) a / b;
if ( b2 == 0 ) {
configuration = 1;
return true;
}
a2b2 = (double) a2 / b2;
if ( b1 == 0 ) {
configuration = 0;
return true;
}
a1b1 = (double) a1 / b1;
if ( b12 == 0 ) {
configuration = setConfiguration();
return true;
}
a12b12 = (double) a12 / b12;
currentTerm = (int) ab;
if ( currentTerm == (int) a1b1 && currentTerm == (int) a2b2 && currentTerm == (int) a12b12 ) {
int temp = a; a = b; b = temp - b * currentTerm;
temp = a1; a1 = b1; b1 = temp - b1 * currentTerm;
temp = a2; a2 = b2; b2 = temp - b2 * currentTerm;
temp = a12; a12 = b12; b12 = temp - b12 * currentTerm;
hasTerm = true;
return false;
}
configuration = setConfiguration();
return true;
}
private int setConfiguration() {
return ( Math.abs(a1b1 - ab) > Math.abs(a2b2 - ab) ) ? 0 : 1;
}
private int a12, a1, a2, a, b12, b1, b2, b;
private double ab, a1b1, a2b2, a12b12;
}
private static interface ContinuedFraction {
public boolean hasMoreTerms();
public int nextTerm();
}
private static class R2cf implements ContinuedFraction {
public R2cf(int aN1, int aN2) {
n1 = aN1; n2 = aN2;
}
public boolean hasMoreTerms() {
return Math.abs(n2) > 0;
}
public int nextTerm() {
final int term = n1 / n2;
final int temp = n2;
n2 = n1 - term * n2;
n1 = temp;
return term;
}
private int n1, n2;
}
private static class NG implements ContinuedFraction {
public NG(NG4 aNG, ContinuedFraction aCF) {
matrixNG = aNG;
cf.add(aCF);
}
public NG(NG8 aNG, ContinuedFraction aCF1, ContinuedFraction aCF2) {
matrixNG = aNG;
cf.add(aCF1); cf.add(aCF2);
}
public boolean hasMoreTerms() {
while ( matrixNG.needsTerm() ) {
if ( cf.get(matrixNG.configuration).hasMoreTerms() ) {
matrixNG.consumeTerm(cf.get(matrixNG.configuration).nextTerm());
} else {
matrixNG.consumeTerm();
}
}
return matrixNG.hasTerm;
}
public int nextTerm() {
matrixNG.hasTerm = false;
return matrixNG.currentTerm;
}
private MatrixNG matrixNG;
private List<ContinuedFraction> cf = new ArrayList<ContinuedFraction>();
}
}
</syntaxhighlight>
{{ out }}
<pre>
Testing: [3; 7] + [0; 2]
3 1 1 1 4
3 1 1 1 4
Testing: [1; 5, 2] * [3; 7]
3 1 2 2
3 1 2 2
Testing: [1; 5, 2] - [3; 7]
-1 -1 -24 -1 -2
-1 -1 -24 -1 -2
Testing: Divide [] by [3; 7]
3 7
Testing: ([0; 3, 2] + [1; 5, 2]) * ([0; 3, 2] - [1; 5, 2])
-1 -3 -5 -1 -2 -1 -26 -3
-1 -3 -5 -1 -2 -1 -26 -3
</pre>
Line 6,922 ⟶ 7,456:
-1 -3 -5 -1 -2 -1 -26 -3
-1 -3 -5 -1 -2 -1 -26 -3
</pre>
=={{header|Mercury}}==
{{works with|Mercury|22.01.1}}
{{trans|Standard ML}}
This program was written with reference to the Standard ML, but really is a different kind of implementation: the continued fractions are represented as lazy lists.
The program is not very fast, but this might be due mainly to the <code>integer</code> type in the Mercury standard library not being very fast. I do not know. If so, an interface to the GNU Multiple Precision Arithmetic Library might speed things up quite a bit.
The program comes in two source files. The main program goes in file <code>continued_fraction_task.m</code>:
<syntaxhighlight lang="Mercury">
%%% -*- mode: mercury; prolog-indent-width: 2; -*-
%%%
%%% A program in two files:
%%% continued_fraction_task.m (this file)
%%% continued_fraction.m (the continued_fraction module)
%%%
%%% Compile with:
%%% mmc --make --use-subdirs continued_fraction_task
%%%
:- module continued_fraction_task.
:- interface.
:- import_module io.
:- pred main(io::di, io::uo) is det.
:- implementation.
:- import_module continued_fraction.
:- import_module integer.
:- import_module rational.
:- import_module string.
:- pred show(string::in, continued_fraction::in, string::in,
io::di, io::uo) is det.
:- pred show(string::in, continued_fraction::in,
io::di, io::uo) is det.
show(Expression, CF, Note, !IO) :-
pad_left(Expression, (' '), 19, Expr1),
print(Expr1, !IO),
print(" => ", !IO),
(if (Note = "")
then (print(to_string(CF), !IO),
nl(!IO))
else (pad_right(to_string(CF), (' '), 48, CF1_Str),
print(CF1_Str, !IO),
print(Note, !IO),
nl(!IO))).
show(Expression, CF, !IO) :- show(Expression, CF, "", !IO).
:- func thirteen_elevenths = continued_fraction.
thirteen_elevenths = from_rational(rational(13, 11)).
:- func twentytwo_sevenths = continued_fraction.
twentytwo_sevenths = from_rational(rational(22, 7)).
main(!IO) :-
show("golden ratio", golden_ratio, "(1 + sqrt(5))/2", !IO),
show("silver ratio", silver_ratio, "(1 + sqrt(2))", !IO),
show("sqrt(2)", sqrt2, "from the module", !IO),
show("sqrt(2)", silver_ratio - one, "from the silver ratio", !IO),
show("sqrt(5)", sqrt5, "from the module", !IO),
show("sqrt(5)", (two * golden_ratio) - one, "from the golden ratio", !IO),
show("13/11", thirteen_elevenths, !IO),
show("22/7", twentytwo_sevenths, "approximately pi", !IO),
show("13/11 + 1/2", thirteen_elevenths + one_half, !IO),
show("22/7 + 1/2", twentytwo_sevenths + one_half, !IO),
show("(22/7) * 1/2", twentytwo_sevenths * one_half, !IO),
show("(22/7) / 2", twentytwo_sevenths / two, !IO),
show("sqrt(2) + sqrt(2)", sqrt2 + sqrt2, !IO),
show("sqrt(2) - sqrt(2)", sqrt2 - sqrt2, !IO),
show("sqrt(2) * sqrt(2)", sqrt2 * sqrt2, !IO),
show("sqrt(2) / sqrt(2)", sqrt2 / sqrt2, !IO),
true.
:- end_module continued_fraction_task.
</syntaxhighlight>
The <code>continued_fraction</code> module source goes in file <code>continued_fraction.m</code>:
<syntaxhighlight lang="Mercury">
%%% -*- mode: mercury; prolog-indent-width: 2; -*-
:- module continued_fraction.
:- interface.
:- import_module int.
:- import_module integer.
:- import_module lazy.
:- import_module rational.
%% A continued fraction is a kind of lazy list. The list is always
%% infinitely long. However, you need not consider terms that come
%% after an infinite term.
:- type continued_fraction
---> continued_fraction(lazy(node)).
:- type node
---> cons(term, continued_fraction).
:- type term
---> infinite_term
; integer_term(integer).
%% ng8 = {A12, A1, A2, A, B12, B1, B2, B}.
:- type ng8 == {integer, integer, integer, integer,
integer, integer, integer, integer}.
%% Get a human-readable string. The second form takes a "MaxTerms"
%% argument. The first form uses a default value for MaxTerms.
:- func to_string(continued_fraction) = string.
:- func to_string(int, continued_fraction) = string.
%% Make a term from a regular int.
:- func int_term(int) = term.
%% A "continued fraction" with only infinite terms.
:- func infinite_continued_fraction = continued_fraction.
%% A continued fraction whose term repeats infinitely.
:- func repeating_term(term) = continued_fraction.
%% A continued fraction representing an integer.
:- func from_integer(integer) = continued_fraction.
:- func from_int(int) = continued_fraction.
%% A continued fraction representing a rational number.
:- func from_rational(rational) = continued_fraction.
%% A continued fraction that is a bihomographic function of two other
%% continued fractions.
:- func apply_ng8(ng8, continued_fraction,
continued_fraction) = continued_fraction.
:- func '+'(continued_fraction,
continued_fraction) = continued_fraction.
:- func '-'(continued_fraction,
continued_fraction) = continued_fraction.
:- func '*'(continued_fraction,
continued_fraction) = continued_fraction.
:- func '/'(continued_fraction,
continued_fraction) = continued_fraction.
%% Miscellaneous continued fractions.
:- func zero = continued_fraction.
:- func one = continued_fraction.
:- func two = continued_fraction.
:- func three = continued_fraction.
:- func four = continued_fraction.
%%
:- func one_fourth = continued_fraction.
:- func one_third = continued_fraction.
:- func one_half = continued_fraction.
:- func two_thirds = continued_fraction.
:- func three_fourths = continued_fraction.
%%
:- func golden_ratio = continued_fraction. % (1 + sqrt(5))/2
:- func silver_ratio = continued_fraction. % (1 + sqrt(2))
:- func sqrt2 = continued_fraction. % The square root of two.
:- func sqrt5 = continued_fraction. % The square root of five.
:- implementation.
:- import_module string.
%%--------------------------------------------------------------------
to_string(CF) = to_string(20, CF).
to_string(MaxTerms, CF) = S :-
to_string(MaxTerms, CF, 0, "[", S).
:- pred to_string(int::in, continued_fraction::in, int::in,
string::in, string::out) is det.
to_string(MaxTerms, CF, I, S0, S) :-
CF = continued_fraction(Node),
force(Node) = cons(Term, CF1),
(if (Term = integer_term(N))
then (if (I = MaxTerms) then (S = S0 ++ ",...]")
else (TermStr = (integer.to_string(N)),
Separator = (if (I = 0) then ""
else if (I = 1) then ";"
else ","),
to_string(MaxTerms, CF1, I + 1,
S0 ++ Separator ++ TermStr, S)))
else (S = S0 ++ "]")).
%%--------------------------------------------------------------------
int_term(I) = integer_term(integer(I)).
infinite_continued_fraction = CF :-
CF = repeating_term(infinite_term).
repeating_term(T) = CF :-
CF = continued_fraction(Node),
Node = delay((func) = cons(T, repeating_term(T))).
from_integer(I) = CF :-
CF = continued_fraction(Node),
Node = delay((func) = cons(integer_term(I), infinite_continued_fraction)).
from_int(I) = from_integer(integer(I)).
from_rational(R) = CF :-
N = numer(R),
D = denom(R),
CF = from_rational_integers(N, D).
:- func from_rational_integers(integer, integer) = continued_fraction.
from_rational_integers(N, D) = CF :-
if (D = zero) then (CF = infinite_continued_fraction)
else (divide_with_rem(N, D, Q, R),
CF = continued_fraction(
delay((func) = cons(integer_term(Q),
from_rational_integers(D, R)))
)).
%%--------------------------------------------------------------------
(X : continued_fraction) + (Y : continued_fraction) = (
apply_ng8({zero, one, one, zero, zero, zero, zero, one}, X, Y)
).
(X : continued_fraction) - (Y : continued_fraction) = (
apply_ng8({zero, one, negative_one, zero, zero, zero, zero, one}, X, Y)
).
(X : continued_fraction) * (Y : continued_fraction) = (
apply_ng8({one, zero, zero, zero, zero, zero, zero, one}, X, Y)
).
(X : continued_fraction) / (Y : continued_fraction) = (
apply_ng8({zero, one, zero, zero, zero, zero, one, zero}, X, Y)
).
apply_ng8(NG, X, Y) = CF :-
NG = {A12, A1, A2, A, B12, B1, B2, B},
(if iseqz(B12, B1, B2, B) then (
%% There are no more finite terms to output.
CF = infinite_continued_fraction
)
else if iseqz(B2, B) then (
absorb_x_term(NG, NG1, X, X1, Y, Y1),
CF = apply_ng8(NG1, X1, Y1)
)
else if at_least_one_iseqz(B2, B) then (
absorb_y_term(NG, NG1, X, X1, Y, Y1),
CF = apply_ng8(NG1, X1, Y1)
)
else if iseqz(B1) then (
absorb_x_term(NG, NG1, X, X1, Y, Y1),
CF = apply_ng8(NG1, X1, Y1)
)
else (
maybe_divide(A12, B12, Q12, R12),
maybe_divide(A1, B1, Q1, R1),
maybe_divide(A2, B2, Q2, R2),
maybe_divide(A, B, Q, R),
(if (not (iseqz(B12)), Q = Q12, Q = Q1, Q = Q2)
then (
%% Output a term.
if (integer_exceeds_infinitizing_threshold(Q))
then (CF = infinite_continued_fraction)
else (NG1 = {B12, B1, B2, B, R12, R1, R2, R},
CF = continued_fraction(
delay((func) = cons(integer_term(Q),
apply_ng8(NG1, X, Y)))
))
)
else (
%% Put A1, A2, and A over a common denominator and compare some
%% magnitudes.
N1 = A1 * B2 * B,
N2 = A2 * B1 * B,
N = A * B1 * B2,
(if (abs(N1 - N) > abs(N2 - N))
then (absorb_x_term(NG, NG1, X, X1, Y, Y1),
CF = apply_ng8(NG1, X1, Y1))
else (absorb_y_term(NG, NG1, X, X1, Y, Y1),
CF = apply_ng8(NG1, X1, Y1)))
))
)).
:- pred absorb_x_term(ng8::in, ng8::out,
continued_fraction::in, continued_fraction::out,
continued_fraction::in, continued_fraction::out)
is det.
absorb_x_term(!NG, !X, !Y) :-
(!.NG) = {A12, A1, A2, A, B12, B1, B2, B},
(!.X) = continued_fraction(XNode),
force(XNode) = cons(XTerm, X1),
(if (XTerm = integer_term(N))
then (New_NG = {A2 + (A12 * N), A + (A1 * N), A12, A1,
B2 + (B12 * N), B + (B1 * N), B12, B1},
(if (ng8_exceeds_processing_threshold(New_NG))
then (
%% Pretend we have reached an infinite term.
!:NG = {A12, A1, A12, A1, B12, B1, B12, B1},
!:X = infinite_continued_fraction
)
else (!:NG = New_NG, !:X = X1)))
else (!:NG = {A12, A1, A12, A1, B12, B1, B12, B1},
!:X = infinite_continued_fraction)).
:- pred absorb_y_term(ng8::in, ng8::out,
continued_fraction::in, continued_fraction::out,
continued_fraction::in, continued_fraction::out)
is det.
absorb_y_term(!NG, !X, !Y) :-
(!.NG) = {A12, A1, A2, A, B12, B1, B2, B},
(!.Y) = continued_fraction(YNode),
force(YNode) = cons(YTerm, Y1),
(if (YTerm = integer_term(N))
then (New_NG = {A1 + (A12 * N), A12, A + (A2 * N), A2,
B1 + (B12 * N), B12, B + (B2 * N), B2},
(if (ng8_exceeds_processing_threshold(New_NG))
then (
%% Pretend we have reached an infinite term.
!:NG = {A12, A12, A2, A2, B12, B12, B2, B2},
!:Y = infinite_continued_fraction
)
else (!:NG = New_NG, !:Y = Y1)))
else (!:NG = {A12, A12, A2, A2, B12, B12, B2, B2},
!:Y = infinite_continued_fraction)).
:- pred ng8_exceeds_processing_threshold(ng8::in) is semidet.
:- pred integer_exceeds_processing_threshold(integer::in) is semidet.
ng8_exceeds_processing_threshold({A12, A1, A2, A,
B12, B1, B2, B}) :-
(integer_exceeds_processing_threshold(A12) ;
integer_exceeds_processing_threshold(A1) ;
integer_exceeds_processing_threshold(A2) ;
integer_exceeds_processing_threshold(A) ;
integer_exceeds_processing_threshold(B12) ;
integer_exceeds_processing_threshold(B1) ;
integer_exceeds_processing_threshold(B2) ;
integer_exceeds_processing_threshold(B)).
integer_exceeds_processing_threshold(Integer) :-
abs(Integer) >= pow(two, integer(512)).
:- pred integer_exceeds_infinitizing_threshold(integer::in) is semidet.
integer_exceeds_infinitizing_threshold(Integer) :-
abs(Integer) >= pow(two, integer(64)).
:- pred maybe_divide(integer::in, integer::in,
integer::out, integer::out) is det.
maybe_divide(N, D, Q, R) :-
if iseqz(D) then (Q = zero, R = zero)
else divide_with_rem(N, D, Q, R).
:- pred iseqz(integer::in) is semidet.
:- pred iseqz(integer::in, integer::in) is semidet.
:- pred iseqz(integer::in, integer::in,
integer::in, integer::in) is semidet.
iseqz(Integer) :- Integer = zero.
iseqz(A, B) :- iseqz(A), iseqz(B).
iseqz(A, B, C, D) :- iseqz(A), iseqz(B), iseqz(C), iseqz(D).
:- pred at_least_one_iseqz(integer::in, integer::in) is semidet.
at_least_one_iseqz(A, B) :- (A = zero; B = zero).
%%--------------------------------------------------------------------
:- func two_plus_sqrt5 = continued_fraction.
two_plus_sqrt5 = repeating_term(int_term(4)).
zero = from_int(0).
one = from_int(1).
two = from_int(2).
three = from_int(3).
four = from_int(4).
one_fourth = from_rational(rational(1, 4)).
one_third = from_rational(rational(1, 3)).
one_half = from_rational(rational(1, 2)).
two_thirds = from_rational(rational(2, 3)).
three_fourths = from_rational(rational(3, 4)).
golden_ratio = repeating_term(int_term(1)).
silver_ratio = repeating_term(int_term(2)).
sqrt2 = continued_fraction(delay((func) = cons(int_term(1), silver_ratio))).
sqrt5 = continued_fraction(delay((func) = cons(int_term(2), two_plus_sqrt5))).
%%--------------------------------------------------------------------
:- end_module continued_fraction.
</syntaxhighlight>
{{out}}
<pre>$ mmc --make --use-subdirs continued_fraction_task && ./continued_fraction_task
Making Mercury/int3s/continued_fraction_task.int3
Making Mercury/int3s/continued_fraction.int3
Making Mercury/ints/continued_fraction.int
Making Mercury/ints/continued_fraction_task.int
Making Mercury/cs/continued_fraction.c
Making Mercury/cs/continued_fraction_task.c
Making Mercury/os/continued_fraction.o
Making Mercury/os/continued_fraction_task.o
Making continued_fraction_task
golden ratio => [1;1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,...] (1 + sqrt(5))/2
silver ratio => [2;2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,...] (1 + sqrt(2))
sqrt(2) => [1;2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,...] from the module
sqrt(2) => [1;2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,...] from the silver ratio
sqrt(5) => [2;4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,...] from the module
sqrt(5) => [2;4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,...] from the golden ratio
13/11 => [1;5,2]
22/7 => [3;7] approximately pi
13/11 + 1/2 => [1;1,2,7]
22/7 + 1/2 => [3;1,1,1,4]
(22/7) * 1/2 => [1;1,1,3]
(22/7) / 2 => [1;1,1,3]
sqrt(2) + sqrt(2) => [2;1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,...]
sqrt(2) - sqrt(2) => [0]
sqrt(2) * sqrt(2) => [2]
sqrt(2) / sqrt(2) => [1]
</pre>
Line 7,540 ⟶ 8,489:
=={{header|OCaml}}==
{{trans|Standard ML}}
{{libheader|Zarith}}
To facilitate comparison of
You will need the Zarith module, which is a commonly used interface to GNU Multiple Precision.
<syntaxhighlight lang="ocaml">
Line 8,041 ⟶ 8,993:
sqrt2 / sqrt2 => [1]
</pre>
=={{header|Owl Lisp}}==
{{trans|Mercury}}
{{works with|Owl Lisp|0.2.1}}
<syntaxhighlight lang="scheme">
;;--------------------------------------------------------------------
;;
;; A continued fraction will be represented as an infinite lazy list
;; of terms, where a term is either an integer or the symbol 'infinity
;;
;;--------------------------------------------------------------------
(define (cf2maxterms-string maxterms cf)
(let repeat ((p cf)
(i 0)
(s "["))
(let ((term (lcar p))
(rest (lcdr p)))
(if (eq? term 'infinity)
(string-append s "]")
(if (>= i maxterms)
(string-append s ",...]")
(let ((separator (case i
((0) "")
((1) ";")
(else ",")))
(term-str (number->string term)))
(repeat rest (+ i 1)
(string-append s separator term-str))))))))
(define (cf2string cf)
(cf2maxterms-string 20 cf))
;;--------------------------------------------------------------------
(define (repeated-term-cf term)
(lunfold (lambda (x) (values x x))
term
(lambda (x) #false)))
(define cf-end (repeated-term-cf 'infinity))
(define (i2cf term) (pair term cf-end))
(define (r2cf fraction)
;; The entire finite-length continued fraction is constructed in
;; reverse order. The list is then rebuilt in the correct order, and
;; given an infinite number of 'infinity terms as its tail.
(let repeat ((n (numerator fraction))
(d (denominator fraction))
(revlst #null))
(if (zero? d)
(lappend (reverse revlst) cf-end)
(let-values (((q r) (truncate/ n d)))
(repeat d r (cons q revlst))))))
(define (->cf x)
(cond ((integer? x) (i2cf x))
((rational? x) (r2cf x))
(else x)))
;;--------------------------------------------------------------------
(define (maybe-divide a b)
(if (zero? b)
(values #false #false)
(truncate/ a b)))
(define integer-exceeds-processing-threshold?
(let ((threshold+1 (expt 2 512)))
(lambda (i) (>= (abs i) threshold+1))))
(define (any-exceed-processing-threshold? lst)
(any integer-exceeds-processing-threshold? lst))
(define integer-exceeds-infinitizing-threshold?
(let ((threshold+1 (expt 2 64)))
(lambda (i) (>= (abs i) threshold+1))))
(define (ng8-values ng)
(values (list-ref ng 0)
(list-ref ng 1)
(list-ref ng 2)
(list-ref ng 3)
(list-ref ng 4)
(list-ref ng 5)
(list-ref ng 6)
(list-ref ng 7)))
(define (absorb-x-term ng x y)
(let-values (((a12 a1 a2 a b12 b1 b2 b) (ng8-values ng)))
(define (no-more-x)
(values (list a12 a1 a12 a1 b12 b1 b12 b1) cf-end y))
(let ((term (lcar x)))
(if (eq? term 'infinity)
(no-more-x)
(let ((new-ng (list (+ a2 (* a12 term))
(+ a (* a1 term)) a12 a1
(+ b2 (* b12 term))
(+ b (* b1 term)) b12 b1)))
(cond ((any-exceed-processing-threshold? new-ng)
;; Pretend we have reached the end of x.
(no-more-x))
(else (values new-ng (lcdr x) y))))))))
(define (absorb-y-term ng x y)
(let-values (((a12 a1 a2 a b12 b1 b2 b) (ng8-values ng)))
(define (no-more-y)
(values (list a12 a12 a2 a2 b12 b12 b2 b2) x cf-end))
(let ((term (lcar y)))
(if (eq? term 'infinity)
(no-more-y)
(let ((new-ng (list (+ a1 (* a12 term)) a12
(+ a (* a2 term)) a2
(+ b1 (* b12 term)) b12
(+ b (* b2 term)) b2)))
(cond ((any-exceed-processing-threshold? new-ng)
;; Pretend we have reached the end of y.
(no-more-y))
(else (values new-ng x (lcdr y)))))))))
(define (apply-ng8 ng x y)
(let repeat ((ng ng)
(x (->cf x))
(y (->cf y)))
(let-values (((a12 a1 a2 a b12 b1 b2 b) (ng8-values ng)))
(cond ((every zero? (list b12 b1 b2 b)) cf-end)
((every zero? (list b2 b))
(let-values (((ng x y) (absorb-x-term ng x y)))
(repeat ng x y)))
((any zero? (list b2 b))
(let-values (((ng x y) (absorb-y-term ng x y)))
(repeat ng x y)))
((zero? b1)
(let-values (((ng x y) (absorb-x-term ng x y)))
(repeat ng x y)))
(else
(let-values (((q12 r12) (maybe-divide a12 b12))
((q1 r1) (maybe-divide a1 b1))
((q2 r2) (maybe-divide a2 b2))
((q r) (maybe-divide a b)))
(cond
((and (not (zero? b12))
(= q q12) (= q q1) (= q q2))
;; Output a term.
(if (integer-exceeds-infinitizing-threshold? q)
cf-end ; Pretend the term is infinite.
(let ((new-ng (list b12 b1 b2 b r12 r1 r2 r)))
(pair q (repeat new-ng x y)))))
(else
;; Put a1, a2, and a over a common denominator
;; and compare some magnitudes.
(let ((n1 (* a1 b2 b))
(n2 (* a2 b1 b))
(n (* a b1 b2)))
(let ((absorb-term
(if (> (abs (- n1 n)) (abs (- n2 n)))
absorb-x-term
absorb-y-term)))
(let-values (((ng x y) (absorb-term ng x y)))
(repeat ng x y))))))))))))
;;--------------------------------------------------------------------
(define (make-ng8-operation ng) (lambda (x y) (apply-ng8 ng x y)))
(define cf+ (make-ng8-operation '(0 1 1 0 0 0 0 1)))
(define cf- (make-ng8-operation '(0 1 -1 0 0 0 0 1)))
(define cf* (make-ng8-operation '(1 0 0 0 0 0 0 1)))
(define cf/ (make-ng8-operation '(0 1 0 0 0 0 1 0)))
;;--------------------------------------------------------------------
(define golden-ratio (repeated-term-cf 1))
(define silver-ratio (repeated-term-cf 2))
(define sqrt2 (pair 1 silver-ratio))
(define sqrt5 (pair 2 (repeated-term-cf 4)))
;;--------------------------------------------------------------------
(define (show expression cf note)
(let* ((cf (cf2string cf))
(expr-len (string-length expression))
(expr-pad-len (max 0 (- 18 expr-len)))
(expr-pad (make-string expr-pad-len #\space)))
(display expr-pad)
(display expression)
(display " => ")
(display cf)
(unless (string=? note "")
(let* ((cf-len (string-length cf))
(cf-pad-len (max 0 (- 48 cf-len)))
(cf-pad (make-string cf-pad-len #\space)))
(display cf-pad)
(display note)))
(newline)))
(show "13/11 + 1/2" (cf+ 13/11 1/2) "(cf+ 13/11 1/2)")
(show "22/7 + 1/2" (cf+ 22/7 1/2) "(cf+ 22/7 1/2)")
(show "22/7 * 1/2" (cf* 22/7 1/2) "(cf* 22/7 1/2)")
(show "golden ratio" golden-ratio "(repeated-term-cf 1)")
(show "silver ratio" silver-ratio "(repeated-term-cf 2)")
(show "sqrt(2)" sqrt2 "(pair 1 silver-ratio)")
(show "sqrt(2)" (cf- silver-ratio 1) "(cf- silver-ratio 1)")
(show "sqrt(5)" sqrt5 "(pair 2 (repeated-term-cf 4)")
(show "sqrt(5)" (cf- (cf* 2 golden-ratio) 1)
"(cf- (cf* 2 golden-ratio) 1)")
(show "sqrt(2) + sqrt(2)" (cf+ sqrt2 sqrt2) "(cf+ sqrt2 sqrt2)")
(show "sqrt(2) - sqrt(2)" (cf- sqrt2 sqrt2) "(cf- sqrt2 sqrt2)")
(show "sqrt(2) * sqrt(2)" (cf* sqrt2 sqrt2) "(cf* sqrt2 sqrt2)")
(show "sqrt(2) / sqrt(2)" (cf/ sqrt2 sqrt2) "(cf/ sqrt2 sqrt2)")
(show "(1 + 1/sqrt(2))/2" (cf/ (cf+ 1 (cf/ 1 sqrt2)) 2)
"(cf/ (cf+ 1 (cf/ 1 sqrt2)) 2)")
(show "(1 + 1/sqrt(2))/2" (apply-ng8 '(0 1 0 0 0 0 2 0)
silver-ratio sqrt2)
"(apply-ng8 '(0 1 0 0 0 0 2 0) sqrt2 sqrt2)")
(show "(1 + 1/sqrt(2))/2" (apply-ng8 '(1 0 0 1 0 0 0 8)
silver-ratio silver-ratio)
"(apply-ng8 '(1 0 0 1 0 0 0 8) silver-ratio silver-ratio)")
;;--------------------------------------------------------------------
</syntaxhighlight>
{{out}}
<pre>$ ol continued-fraction-task-Owl.scm
13/11 + 1/2 => [1;1,2,7] (cf+ 13/11 1/2)
22/7 + 1/2 => [3;1,1,1,4] (cf+ 22/7 1/2)
22/7 * 1/2 => [1;1,1,3] (cf* 22/7 1/2)
golden ratio => [1;1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,...] (repeated-term-cf 1)
silver ratio => [2;2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,...] (repeated-term-cf 2)
sqrt(2) => [1;2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,...] (pair 1 silver-ratio)
sqrt(2) => [1;2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,...] (cf- silver-ratio 1)
sqrt(5) => [2;4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,...] (pair 2 (repeated-term-cf 4)
sqrt(5) => [2;4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,...] (cf- (cf* 2 golden-ratio) 1)
sqrt(2) + sqrt(2) => [2;1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,...] (cf+ sqrt2 sqrt2)
sqrt(2) - sqrt(2) => [0] (cf- sqrt2 sqrt2)
sqrt(2) * sqrt(2) => [2] (cf* sqrt2 sqrt2)
sqrt(2) / sqrt(2) => [1] (cf/ sqrt2 sqrt2)
(1 + 1/sqrt(2))/2 => [0;1,5,1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,...] (cf/ (cf+ 1 (cf/ 1 sqrt2)) 2)
(1 + 1/sqrt(2))/2 => [0;1,5,1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,...] (apply-ng8 '(0 1 0 0 0 0 2 0) sqrt2 sqrt2)
(1 + 1/sqrt(2))/2 => [0;1,5,1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,4,1,...] (apply-ng8 '(1 0 0 1 0 0 0 8) silver-ratio silver-ratio)</pre>
=={{header|Phix}}==
Line 9,691 ⟶ 10,885:
=={{header|Wren}}==
{{trans|Kotlin}}
<syntaxhighlight lang="
construct new() {
_cfn = 0
| |||