{-# LANGUAGE GeneralizedNewtypeDeriving #-}
-- | Proof that Z[sqrt(-5)] is a Prufer domain. This implies that it is 
-- possible to solve systems of equations over Z[sqrt(-5)].
module Algebra.ZSqrt5 (ZSqrt5(..)) where

import Test.QuickCheck

import Algebra.Structures.IntegralDomain
import Algebra.Structures.EuclideanDomain (quotient, genEuclidAlg)
import Algebra.Structures.BezoutDomain (toPrincipal)
import Algebra.Structures.PruferDomain
import Algebra.Structures.Coherent
import Algebra.Ideal
import Algebra.Z
import Algebra.Q

-- | Z[sqrt(-5)] is a pair such that (a,b) = a + b*sqrt(-5)
newtype ZSqrt5 = ZSqrt5 (Z,Z)
  deriving (Eq,Ord,Arbitrary)

instance Show ZSqrt5 where
  show (ZSqrt5 (a,b)) = show a ++ " + " ++ show b ++ " * sqrt(-5)"

-- Arithmetical properties
instance Ring ZSqrt5 where
  (ZSqrt5 (a,b)) <+> (ZSqrt5 (c,d)) = ZSqrt5 (a + c, b + d)
  (ZSqrt5 (a,b)) <*> (ZSqrt5 (c,d)) = ZSqrt5 (a*c - 5*b*d, a*d + b*c)
  neg (ZSqrt5 (a,b))                = ZSqrt5 (neg a, neg b)
  zero                              = ZSqrt5 (0,0)
  one                               = ZSqrt5 (1,0)
  
instance CommutativeRing ZSqrt5 where

instance IntegralDomain ZSqrt5 where

propIntDomZSqrt5 :: ZSqrt5 -> ZSqrt5 -> ZSqrt5 -> Property
propIntDomZSqrt5 = propIntegralDomain

--------------------------------------------------------------------------------
-- Useful auxiliary functions:


(+>) :: Z -> ZSqrt5 -> ZSqrt5
r +> (ZSqrt5 (a,b)) = ZSqrt5 (r+a,b)

(<+) :: ZSqrt5 -> Z -> ZSqrt5
(ZSqrt5 (a,b)) <+ r = ZSqrt5 (a+r,b)

infixl 6 +>, <+

--------------------------------------------------------------------------------

instance PruferDomain ZSqrt5 where
  -- Assume /= 0
  calcUVW (ZSqrt5 (a,b)) (ZSqrt5 (c,d)) = (u,v,w) 
    where
    -- Let s = (a+b*sqrt(-5))/(c+d*sqrt(-5))
    -- Compute p and q such that: s = p + q*sqrt(-5)
    p = toQ (a*c+5*b*d) / n
    q = toQ (b*c-a*d)   / n
    n = toQ (c^2+5*d^2)

    s :: (Q,Q)
    s = (p,q)

    -- Rewrite: 
    -- s = p + q*sqrt(-5)  <==>  s^2 - 2*sp + p^2 + 5*q^2 = 0
    --                     <==> a0*s^2 + a1*s + a2 = 0
    --
    -- Some computations give:
    a0' = toZ n
    a1' = -2 * toZ (p*n)
    a2' = a^2 + 5*b^2 
    
    -- Normalize:
    g  = genEuclidAlg [a0',a1',a2']
    a0 = a0' `quotient` g
    a1 = a1' `quotient` g
    a2 = a2' `quotient` g

    -- Compute m0, m1, m2, m3 such that:
    -- m0a0 + m1a0s + m2(a0s+a1) + m3(a0s+a1)s = 1
    a0s    = ZSqrt5 (toZ (a0 *> p), toZ (a0 *> q))
    a0sa1  = a0s <+ a1
    a0sa1s = ZSqrt5 (-a2,0)

    (Id [1],[n0,n1,n2],_) = toPrincipal (Id [a0,a1,a2]) 

    m0 = n0
    m1 = -n1
    m2 = n1
    m3 = -n2

    -- Finally we get u, v and w:
    u = m0 * a0 +> m2 *> a0sa1
    v = m0 *> a0s <+> m2 *> a0sa1s
    w = m1 * a0 +> m3 *> a0sa1

propPruferDomZSqrt5 :: ZSqrt5 -> ZSqrt5 -> ZSqrt5 -> Property
propPruferDomZSqrt5 x@(ZSqrt5 (a,b)) y@(ZSqrt5 (c,d)) z@(ZSqrt5 (e,f)) = 
  a /= 0 && b /= 0 && c /= 0 && d /= 0 && e /= 0 && f /= 0 ==> propPruferDomain x y z

instance Coherent ZSqrt5 where
  solve = solvePD