Jump to content

A/B = c. Given c, derive A & B knowing (only) that they are both integers.


Browseruk

Recommended Posts

You correctly state that there may not be a unique answer, since A and B may be multiplied by any non-zero number. That begs  the question: What is your point? Contrary to your claim, your thread title does not completely describe the problem. As a rough guess of what you meant: There are numbers that cannot be represented as fractions of natural numbers. The most prominent cases are pi (as in "ratio of circle circumference to radius") and "e" (as in the exponential function).

Link to comment
Share on other sites

Sorry, I thought I made my point clear when I wrote: 

Quote
1 hour ago, Browseruk said:

I'm after the smallest integers that produce c.

 

Given an improper fraction (as a decimal value)(eg.1.05) how can I derive the smallest pairing of integers A / B (eg. 21/20 ) that produce that value?

Edited by Browseruk
clarification
Link to comment
Share on other sites

Why not just write "1.05" as "105/100" and simplify the fraction?
I'm not sure if that's what you are asking, but if it is, from there do a prime factorization of the numbers and  eliminate common factors 

Link to comment
Share on other sites

2 hours ago, moth said:

Why not just write "1.05" as "105/100" and simplify the fraction?
 

But what if the value I get is 1.2105263157894736842105263157895? (23/19)

1 minute ago, Sensei said:

..there is infinite quantity of answers..

102/100 = 1.02

204/200 = 1.02

306/300 = 1.02

[...]

1020/1000 = 1.02

etc. etc.

 

Yes. That' s what I was indicating with "multiple [sic]"; but there can only be one that cannot be reduced by a common factor.

Link to comment
Share on other sites

19 minutes ago, Browseruk said:

But what if the value I get is 1.2105263157894736842105263157895? (23/19)

Yes. That' s what I was indicating with "multiple [sic]"; but there can only be one that cannot be reduced by a common factor.

then hopefully you can do prime factorization quickly 


1020/1000
2 * 510 / 2 * 500
  2 * 255 / 2 * 250
   5 *51 / 2 * 125
              / 5 * 25
              / 5 * 5
51/50=1.02

1.2105263157894736842105263157895 / 100000000000000000000000000000000 =?

Edited by moth
Link to comment
Share on other sites

I made little brute-force algorithm for you. Executable is in Release folder.

// RationalSeeker v1.0.1 (c) 2019 by Sensei : main project file.
#include "stdafx.h"
using namespace System;
int main(array<System::String ^> ^args)
{
    if( args->Length == 1 ) {
        try {
            String ^arg = args[ 0 ];
            Decimal rational = Decimal::Parse( arg );
            Console::WriteLine( L"Rational " + rational );
            Decimal temp_rational = rational;
            Int64 temp = 1;
            while( Math::Floor( temp_rational ) != temp_rational ) {
                temp_rational *= 10;
                temp *= 10;
            }
            Int64 integer = (Int64) temp_rational;
            Console::WriteLine( L"Integer " + integer );
	            Int64 half = integer / 2;
            for( Int64 i = 1; i < half; i++ ) {
                if( ( ( integer % i ) == 0 ) && ( ( temp % i ) == 0 ) ) {
                    Int64 numerator = integer / i;
                    Int64 denominator = temp / i;
                    Console::WriteLine( numerator+ L"/" + denominator + L"=" + rational );
                    //break;
                }
            }
        } catch( Exception ^e ) {
            Console::Error->WriteLine( e->ToString() );
        }
    } else {
        Console::Error->WriteLine( L"Required parameter missing!" );
    }
    return 0;
}

It first get rid of fractional part by multiplication by 10 in the loop.

Then goes in the second loop through the all numbers from 1 to half of integer (it could go through the all primes including 1 instead, but first we would have to have array of the all primes).

Edit: fixed variable names, and decreased archive size considerably..

 

 

RationalSeeker.zip

ps. If you want to find it faster, replace:

            for( Int64 i = 1; i < half; i++ ) {

by

            for( Int64 i = half; i > 0; i-- ) {
Edited by Sensei
Link to comment
Share on other sites

2 hours ago, Browseruk said:

But what if the value I get is 1.2105263157894736842105263157895? (23/19)

You cannot enter numbers like this to computer.. Decimal number parser (and entire decimal number C++ class) would have to be modified to support parenthesis, and then keep number as a p/q instead of regular number.

It has (infinitely) repeating decimal period.

1.(210526315789473684)

https://en.wikipedia.org/wiki/Repeating_decimal

Other examples:

1/3 = 0.33333(3)

4/3 = 1.33333(3)

Check section "Converting repeating decimals to fractions" in the above article.

Edited by Sensei
Link to comment
Share on other sites

59 minutes ago, uncool said:

What, exactly, are you "given"? Just a partial decimal expansion of the fraction?

The values I receive are the output from  a magneto-static analysis. I get them in the form of IEEE754 64-bit floating point values. they are all between 1 and 2 (not inclusive).

Visually inspecting the data, I notice some of the values corresponded to the division of two relatively small integers; and that leads me to suspect that possibly all the values are such. (harmonics); but I need to prove it.

So I knocked up an iterative solution:

sub improper(@) { 
    my( $n, $x, $i, $f ) = ( 1, shift, 0, 1e9 ); 
    do{ 
        ++$n; 
        $i = int( $n * $x ); 
        $f = $n * $x - $i; 
    } until abs( $f ) < 1e-9; 
    return $i, $n 
}

And used it to produce a table of all the 3 decimal place fractions  1<x<2; diacarding those where the divisor was 1000; and came up with a table of just over 600 pairings:

1.001:  33033 33000
1.002:  501 500
1.003:  33099 33000
1.004:  251 250
1.005:  4221 4200
1.006:  503 500
1.007:  17119 17000
1.008:  126 125
1.009:  17153 17000
1.010:  101 100
1.011:  17187 17000
1.012:  253 250
1.013:  9117 9000
1.014:  507 500
1.015:  2233 2200
1.016:  127 125
1.017:  9153 9000
1.018:  509 500
1.019:  9171 9000
1.020:   51 50
...
1.975:  79 40
1.976:  247 125
1.978:  989 500
1.980:  99 50
1.982:  991 500
1.984:  248 125
1.985:  397 200
1.986:  993 500
1.988:  497 250
1.990:  199 100
1.992:  249 125
1.994:  997 500
1.995:  399 200
1.996:  499 250
1.998:  999 500

But that takes 0.5 second for 3 decimal digit producing 607 integer pairings; and 45 seconds for 4 decimal digits producing 6230 integer pairings. ~450 seconds for 5 dps and ~60,000 solutions. 

The trend is fairly clear. IEEE supports ~15 dps, so producing a table would take a lifetime and more storage than exists (guess! :), so I was hoping for a non-iterative, or possibly a less iterative. Knowing that gcd() can be programmed to converge very quickly, I looked to see if that could help and came across the formula's section of wikipedia's lcm() page:wikipedia's lcm() page and thought that there was a possibility in there somewhere, but since that math is over my head, I asked here.

 

Link to comment
Share on other sites

3 minutes ago, Browseruk said:

But that takes 0.5 second for 3 decimal digit producing 607 integer pairings; and 45 seconds for 4 decimal digits producing 6230 integer pairings. ~450 seconds for 5 dps and ~60,000 solutions. 

You're going, like I in the first source code, from 1, and increasing loop counter. Try going from half of range, and then decreasing by 1, till it's bigger than 0.

BTW, we can see your code gave wrong answers, e.g.

1.001:  33033 33000

It should be 1001/1000, isn't?

 

Link to comment
Share on other sites

11 hours ago, Browseruk said:

The values I receive are the output from  a magneto-static analysis. I get them in the form of IEEE754 64-bit floating point values. they are all between 1 and 2 (not inclusive).

So in short, what you are given is a number between 1 and 2, up to about 16 decimal places. You suspect that it is a rational number with relatively small numerator and denominator, and want to find those quickly. My suggestion would be to use continued fractions with a "cut-off" of a large number. For example:

Let's say the number we want to find is 13/9 = 1.4444444...., but we have some source of error and instead get x0 = 1.4444387260. Choose a cutoff of 100. Continued fractions would work as follows:

a0 = 1, x1 = 1/(x0 - 1) = 2.25002895

a1 = 2, x2 = 1/(x1 - 2) = 3.99953685363

a2 = 3, x3 = 1/(x2 - 3) = 1.00046336097

a3 = 1, x4 = 1/(x3 - 1) = 2158.1446534

2158 is larger than the cutoff, so we stop and assume our number is: 1 + 1/(2 + 1/(3 + 1/1)) = 1 + 1/(2 + 1/4) = 1 + 4/9 = 13/9.

 

This all assumes that your "noise" is relatively small (to be precise, if I'm not mistaken it should be a small fraction of the square of the reciprocal of the largest denominator you allow). If the noise is larger than that, then there will be multiple fractions that "fit". 

Edited by uncool
Link to comment
Share on other sites

18 hours ago, Browseruk said:

My (overlong?) title describes the problem completely; my question is how?

I realise that there can be multiple [sic] answers; I'm after the smallest integers that produce c.

eg. Given 1.05; A/B = 21/20.

 

You want Euclid's Algorithm: en.wikipedia.org/wiki/Euclidean_algorithm

E.g. you have 1.05, which is 105/100. Euclid gives that 5 is the greatest common divisor of 105 and 100. Therefore 105/100 = (105/5)/(100/5) = 21/20. 

Link to comment
Share on other sites

21 minutes ago, swansont said:

The part to the right of the decimal is what matters. You have to focus on that. Invert it and that’s your denominator. Call it y. Multiply by y/y

1.05

1/ 0.05 = 20

multiply the number by 20/20, and you get 21/20

Inverting it only works if the number is a natural number plus the reciprocal of a natural number. Take the 13/9 example I gave above. 

The generalization that works (see if you have something close to an integer; if not, remove the integer part and try again) is continued fractions, as I described. 

Link to comment
Share on other sites

49 minutes ago, uncool said:

Inverting it only works if the number is a natural number plus the reciprocal of a natural number. Take the 13/9 example I gave above. 

The generalization that works (see if you have something close to an integer; if not, remove the integer part and try again) is continued fractions, as I described. 

In that case you haven’t been given c, as described in the title. 

But I realize now that inverting the decimal might not yield an integer, so you have to multiply by another factor to get it.

 

Link to comment
Share on other sites

23 hours ago, Sensei said:

Try going from half of range, and then decreasing by 1, till it's bigger than 0.

That helps a little, but not enough.

 

23 hours ago, Sensei said:

BTW, we can see your code gave wrong answers,

Yeah. Early version; but enough to get a feel for the performance.

12 hours ago, uncool said:

My suggestion would be to use continued fractions with a "cut-off" of a large number.

Thanks. That's the cluebat I was hoping for. This still barely tested code seems to work reasonably well:

sub rational{
    my @a = int( my( $x ) = shift );
    {
        my $e = ( $x - int( $x ) );
        last unless $e > 1e-3;
        push @a, int( $x = 1.0 / $e );
        redo;
    }
    my( $n, $d ) = ( 1, pop @a );
    $n += $d * $_, ( $n, $d ) = ( $d, $n ) for reverse @a;
    ( $n, $d ) = ( $d, $n );
    return $n, $d;
}

A few examples@

C:\test>rational -N=21 -D=20
21 / 20

C:\test>rational -N=23 -D=19
23 / 19

C:\test>rational -N=101 -D=97
101 / 97

C:\test>rational -N=101 -D=97.1
1010 / 971

C:\test>rational -N=101 -D=97.01
10100 / 9701

C:\test>rational -N=101 -D=97.000001
101 / 97

C:\test>rational -N=101 -D=97.001
101000 / 97001

C:\test>rational -N=101 -D=97.0001
1010000 / 970001

C:\test>rational -N=101 -D=97.00001
7.04919182283353e+185 / 6.77001660699773e+185

C:\test>rational -N=101 -D=97.000001
101 / 97

C:\test>rational -F=1.70834275
6833371 / 4000000

C:\test>rational -F=1.60759143
3.15561131237261e+181 / 1.96294360213939e+181

Thanks a lot.

Link to comment
Share on other sites

Hrm. That's less able to handle small deviations than I expected - 101/97 differs from 101/97.00001 by about 0.0000001, or about 1/10 the precision I thought it would be able to handle. Also, if you want your algorithm to finish quickly, it might be a good idea to stop it when the denominator gets too large (so you don't get these absurd 10^185 answers). The numerator and denominator can be calculated at each step (as described at https://en.wikipedia.org/wiki/Continued_fraction#Infinite_continued_fractions_and_convergents , storing 2 convergents at a time), which keeps the algorithm quick. 

Edited by uncool
Link to comment
Share on other sites

5 hours ago, uncool said:

The numerator and denominator can be calculated at each step (as described at https://en.wikipedia.org/wiki/Continued_fraction#Infinite_continued_fractions_and_convergents 

Really. Is that what that lot says?  I never would have guessed. Heck, even now I don't have to guess, I'm still having trouble deriving that meaning from those words :) 

5 hours ago, uncool said:

Hrm. That's less able to handle small deviations than I expected - 101/97 differs from 101/97.00001 by about 0.0000001, or about 1/10 the precision I thought it would be able to handle.

Note: I'm using a limit of 1e-3 rather than the 1e-2 you advised.  Primarily because it gives an exact answer for one of my real values

C:\test>rational -F=1.70834275
6833371 / 4000000

rather  than

C:\test>rational -F=1.70834275
41 / 24

Having just tried 3e-3 also get the exact answer. Maybe the thing to do is pass the limit as an argument.

I could also split the routine into two, returning the terms from the first and passing them to the second, so that the calling code can decide how many terms to calculate. I need to play with some more real values and see what falls out.

Link to comment
Share on other sites

2 hours ago, Browseruk said:

Really. Is that what that lot says?  I never would have guessed. Heck, even now I don't have to guess, I'm still having trouble deriving that meaning from those words :) 

Using the example I gave (1.4444387260):

Always start with p0 = 0, q0 = 1, p1 = 1, q1 = 0.

a1 = 1, p2 = a1*p1 + p0 = 1, q2 = a1*q1 + q0 = 1, convergent = 1/1

a2 = 2, p3 = a2*p2 + p1 = 3, q3 = a2*q2 + q1 = 2, convergent = 3/2

a3 = 3, p4 = a3*p3 + p2 = 10, q4 = a3*q3 + q2 = 7, convergent = 10/7

a4 = 1, p5 = a4*p4 + p3 = 13, q5 = a4*q4 + q3 = 9, convergent = 13/9.

a5 = 2158, cutoff reached. 

Alternatively, at the end, it could be:

a5 = 2158, p6 = a5*p5 + p4 = 28064, q6 = a5*q5 + q4 = 19429. 19429 is too large, use last convergent, return 13/9. 

Edited by uncool
Link to comment
Share on other sites

Create an account or sign in to comment

You need to be a member in order to leave a comment

Create an account

Sign up for a new account in our community. It's easy!

Register a new account

Sign in

Already have an account? Sign in here.

Sign In Now
×
×
  • Create New...

Important Information

We have placed cookies on your device to help make this website better. You can adjust your cookie settings, otherwise we'll assume you're okay to continue.