Reliable cube root in Haskell
I am doing question 62 at project euler and came up with the following to test whether a number is cubic:
isInt x = x == fromInteger (round x)
isCube x= isInt $ x**(1/3)
But due to floating point error, it returns incorrect results:
*Main> isCube (384^3)
False
Is there a way to implement a more reliable cube test?
On a side-note, here is the rest of my solution, which doesn't work because of a type interface error on filter (isCube) (perms n)
:
cubes = [n^3|n<-[1..]]
perms n = map read $ permutations $ show n :: [Integer]
answer = he开发者_开发百科ad [n|n<-cubes,(length $ filter (isCube) (perms n)) == 5]
What do I need to do to fix the error?
No instances for (Floating Integer, RealFrac Integer)
arising from a use of `isCube' at prob62.hs:10:44-49
Any optimisations are also welcome ;-)
Try to avoid using floating point numbers as much as possible, especially when you have a problem which concerns integer values. Floating point numbers have problems with rounding and that certain values (like 1/3) cannot be represented exactly. So it's no surprise that you get mysterious answers.
First of all, in order to fix your type error you have to redefine isCube
. If you check it's type signature it looks like this:
isCube :: (RealFrac a, Floating a) => a -> Bool
Note that it expects something that is of class Floating
as its first argument. Your problem is that you want to use this function on integer values and integers are not an instance of Floating
. You can redefine isCube
like this to make the function type check.
isCube x = isInt $ (fromIntegral x) ** (1/3)
However, that will not make your program correct.
One way to make your program more correct is to do what Henrik suggested. It would look like this:
isCube x = (round (fromIntegral x ** (1/3))) ^ 3 == x
Good luck!
Don't know much about Haskell, but I would take the cube root, round to the nearerst integer, take the cube, and compare to the original value.
For another approach useful for Integer
values have a look at the integerCubeRoot
function in the arithmoi package.
Example:
ghci> import Math.NumberTheory.Powers.Cube
ghci> let x = 12345^3333
ghci> length $ show x
13637
ghci> isCube x
True
ghci> isCube (x+1)
False
ghci> length $ show $ integerCubeRoot x
4546
perms
has the type [Integer]
. isCube
has the type (RealFrac a, Floating a) => a -> Bool
(as you can check in GHCI). The RealFrac
constraint comes from round x
, the Floating
constraint comes from x**(1/3)
. Since Integer
is neither RealFrac
nor Floating
, isCube
can't be used as Integer -> Bool
. So filter isCube (perms n)
doesn't make sense.
So you need to fix isCube
to work properly on Integer
s:
isCube x = isInt $ (fromInteger x)**(1/3)
In fact, the reason isCube (384^3)
even compiles is that it "really" means isCube ((fromInteger 384)^(fromInteger 3))
.
Of course, this will still work badly due to floating point errors. Basically, checking floating numbers for equality, as you do in isInt
, is almost always a bad idea. See other answers for explanation how to make a better test.
精彩评论