This is is just a bunch of algebraic relationships. One could use brute force, but it would take several hours.

Instead, reassign these variables:

a = x + y

b = x – y

c = x + z

d = x – z

e = y + z

f = y – z

Then playing around with them, we get:

x = (a+b)/2 = (c+d)/2 because b = c – e and e = a – d

y = (e+f)/2 = (2a-d-c)/2 because e = a – d and f = a – c

z = (c-d)/2

x + y + z = (2a-d+c)/2

So now we just iterate on a,c,d over a list of squares testing for:

a > c > d;

c-d is even;

b,e,f are squares;

and x > y > z > 1

Routine runs in about 120 seconds on a single core 3Ghz system.

```
```import Data.List
isSqr :: Int -> Bool
isSqr n = sq * sq == n
where sq = floor $ sqrt $ (fromIntegral n::Double)
slist :: [Int]
slist = [x | x <- [2..1000000], isSqr x]
test :: Int -> Int -> Int -> Bool
test a c d =
a > c && c > d && even (c-d) && isSqr(c-a+d) && isSqr(a-d) && isSqr(a-c) && x > y && y > z && z > 1
where
x = (c+d) `div` (fromIntegral 2)
y = ((2*a)-d-c) `div` (fromIntegral 2)
z = (c-d) `div` (fromIntegral 2)
main = print $ head $ sort [((2*a)-d+c) `div` (fromIntegral 2)| a<-slist,c<-slist,d<-slist,test a c d]