soluciones :: Integer -> [(Integer, Integer)]

   λ> take 3 (soluciones 7)
   [(0,0),(14,14),(-42,6)]
   λ> sum [x+y | (x,y) <- take 6 (soluciones (primes !! (10^6)))]
   185830404

import Data.List (sort)
import Data.Numbers.Primes (primes)

-- 1ª solución
-- ===========

soluciones :: Integer -> [(Integer, Integer)]
soluciones p =
  [(x,y) | (x,y) <- pares,
           p * (x + y) == x * y]

-- pares es la lista de los pares de números enteros ordenados por la
-- suma de los valores absolutos de sus componentes. Por ejemplo,
--    λ> take 38 pares
--    [(0,0),(-1,0),(0,-1),(0,1),(1,0),(-2,0),(-1,-1),(-1,1),(0,-2),
--     (0,2),(1,-1),(1,1),(2,0),(-3,0),(-2,-1),(-2,1),(-1,-2),(-1,2),
--     (0,-3),(0,3),(1,-2),(1,2),(2,-1),(2,1),(3,0),(-4,0),(-3,-1),
--     (-3,1),(-2,-2),(-2,2),(-1,-3),(-1,3),(0,-4),(0,4),(1,-3),(1,3),
--     (2,-2),(2,2)]
pares :: [(Integer, Integer)]
pares =
  concat [sort (concat [aux (x,n-x) | x <- [0..n]])
         | n <- [0..]]
  where
    aux (0,0) = [(0,0)]
    aux (0,y) = [(0,y), (0,-y)]
    aux (x,0) = [(x,0), (-x,0)]
    aux (x,y) = [(x,y), (-x,y), (x,-y), (-x,-y)]

-- 2ª solución
-- ===========

-- Observando los siguientes cálculos:
--    take 6 (soluciones 2)  == [(0,0),(-2,1),(1,-2),(4,4),(3,6),(6,3)]
--    take 6 (soluciones 3)  == [(0,0),(-6,2),(2,-6),(6,6),(4,12),(12,4)]
--    take 6 (soluciones 5)  == [(0,0),(10,10),(-20,4),(4,-20),(6,30),(30,6)]
--    take 6 (soluciones 7)  == [(0,0),(14,14),(-42,6),(6,-42),(8,56),(56,8)]
--    take 6 (soluciones 11) == [(0,0),(22,22),(-110,10),(10,-110),(12,132),(132,12)]

soluciones2 :: Integer -> [(Integer, Integer)]
soluciones2 2 = [(0,0),(-2,1),(1,-2),(4,4),(3,6),(6,3)]
soluciones2 3 = [(0,0),(-6,2),(2,-6),(6,6),(4,12),(12,4)]
soluciones2 p =
  [(0,0), (2*p,2*p), (p*(1-p),p-1), (p-1,p*(1-p)), (p+1,p*(p+1)), (p*(p+1),p+1)]

-- Comparación de equivalencia
-- ===========================

-- La propiedad es
prop_equivalencia :: Bool
prop_equivalencia =
  and [ take 6 (soluciones p) == soluciones2 p
      | p <- takeWhile (<= 37) primes]

-- La comprobación es
--    λ> prop_equivalencia
--    True

-- Comparación de eficiencia
-- =========================

-- La comparación es
--    λ> take 6 (soluciones 37)
--    [(0,0),(74,74),(-1332,36),(36,-1332),(38,1406),(1406,38)]
--    (4.68 secs, 3,362,009,344 bytes)
--    λ> take 6 (soluciones2 37)
--    [(0,0),(74,74),(-1332,36),(36,-1332),(38,1406),(1406,38)]
--    (0.01 secs, 144,192 bytes)