Amended on Mon 16 Aug 2010 01:53 AM (UTC) by Nick Gammon
Message
Given the above, we can also work out the LCM (least common multiple):
-- See: http://en.wikipedia.org/wiki/Least_common_multiple
function lcm (a, b)
a = bc.number (a)
b = bc.number (b)
return (a * b) / gcd (a, b)
end -- lcm
print (lcm (21, 6)) --> 42
print (lcm (4, 6)) --> 12
In arithmetic and number theory, the least common multiple or (LCM) of two rational numbers a and b is the smallest positive rational number that is an integer multiple of both a and b.
Amended on Mon 16 Aug 2010 01:52 AM (UTC) by Nick Gammon
Message
Another interesting algorithm - GCD (Greatest common divisor)
-- See: http://en.wikipedia.org/wiki/Greatest_common_divisor
-- and http://en.wikipedia.org/wiki/Euclidean_algorithm
function gcd (a, b)
local old_digits = bc.digits (0) -- we are using integers
local temp
a = bc.number (a)
b = bc.number (b)
while not bc.iszero (b) do
temp = b
b = bc.mod (a, b)
a = temp
end -- while
bc.digits (old_digits) -- put old decimal places back
return a
end -- gcd
print (gcd (8, 12)) --> 4
print (gcd (42, 56)) --> 14
print (gcd (206, 40)) --> 2
print (gcd (12070, 41820)) --> 170
print (gcd (10362, 12397)) --> 11
In mathematics, the greatest common divisor (GCD) of two or more non-zero integers, is the largest positive integer that divides the numbers without a remainder. For example, the GCD of 8 and 12 is 4.
Yes, well, I bought an interesting book recently (Coincidences, Chaos and All That Math Jazz), and it suggested this experiment:
In your favourite spreadsheet, put the number 0.5 in cell A1, and then make A2 be "=A1^2-2", and then "fill down" the same calculation over a few dozen rows. My first 20 look like this:
Quote: I think someone famous said "after the third decimal place, no-one gives a damn".
Hrmph, if only I could say the same. :-) In some of the code I am writing, I unfortunately care very much about some numbers after the third decimal... if I can make the inner loop a few tens of microseconds faster, that could mean speeding up the entire program considerably...
David Haley aka Ksilyan
Head Programmer,
Legends of the Darkstone
Amended on Fri 27 Apr 2007 06:51 AM (UTC) by Nick Gammon
Message
And now the inverse operation, taking the natural logarithm:
bc.digits (100)
function ln (x)
local x = bc.number (x) -- argument converted to bignum
-- if x <= 0.5 or >= 2 take the square root and double the result
if bc.compare (x, 0.5) ~= 1 or bc.compare (x, 2) == 1 then
return ln (bc.sqrt (x)) * 2
end -- if x <= 0.5 or >= 2
local a = bc.number ((x - 1) / (x + 1))
local e = bc.number (0)
local t = bc.number (a)
for i = 1, 200, 2 do
local E = e
e = e + (t / i)
t = t * a * a
-- break if no change
if e == E then
break
end -- if no change
end -- for
return e * 2
end -- function ln
print (ln (100))
Output
4.6051701859880913680359829093687284152022029772575459520666558019351452193547049604719944101791965216
Amended on Thu 26 Apr 2007 05:36 AM (UTC) by Nick Gammon
Message
And now to raise e to a power:
bc.digits (100)
function exp (x)
local x = bc.number (x) -- argument converted to bignum
-- if x > 1 halve it and square the result
if bc.compare (x, 1) == 1 then
return exp (x / 2) ^ 2
end -- if x > 1
local e = bc.number (1) -- result
local n = bc.number (1) -- denominator: i!
local t = bc.number (1) -- numerator: x^i
for i = 1, 200 do
local E = e
n = n * i -- n is i! (factorial of i)
t = t * x -- t is x^i (x to the power i)
e = e + (t / n) -- new value
-- break if no change
if e == E then
break
end -- if no change
end -- for
return e
end -- function exp
print (exp (5))
Output
148.4131591025766034211155800405522796234876675938789890467528451109120648209585760796884094598990189279
Amended on Thu 26 Apr 2007 03:38 AM (UTC) by Nick Gammon
Message
This algorithm calculates sin (x) to a specified precision - the precision is the count of times you want to go through the loop:
bc.digits (100)
function sine (x, precision)
local val = bc.number (1)
x = bc.number (x)
while precision > 0 do
val = 1 - val * x * x / (2 * precision) / (2 * precision + 1)
precision = precision - 1
end -- while
val = x * val
return val
end -- function sine
print (sine (1, 35))
Output
0.8414709848078965066525023216302989996225630607983710656727517099919104043912396689486397435430526959
Amended on Wed 25 Apr 2007 07:32 AM (UTC) by Nick Gammon
Message
Here's another interesting algorithm. This calculates e to 100 decimal places:
bc.digits (103)
n = bc.number (1)
e = bc.number (1)
for i = 1, 200 do
local E = e
n = n * i -- n is i! (factorial of i)
e = e + (1 / n)
if e == E then
break
end -- if no change
end -- for
print (e)
Output
2.7182818284590452353602874713526624977572470936999595749669676277240766303535475945713821785251664274229
Note that the last 3 digits are wrong - to get the first 100 right I had to calculate to 103 digits of precision.
A debugging display shows it stops on the 72nd iteration, which means it got it right on the 71st iteration.
The dates and times for posts above are shown in Universal Co-ordinated Time (UTC).
To show them in your local time you can join the forum, and then set the 'time correction' field in your profile to the number of hours difference between your location and UTC time.