-- |
-- Module:      Data.Poly.Internal.Dense.GcdDomain
-- Copyright:   (c) 2019 Andrew Lelechenko
-- Licence:     BSD3
-- Maintainer:  Andrew Lelechenko <andrew.lelechenko@gmail.com>
--
-- 'GcdDomain' instance with a 'GcdDomain' constraint on the coefficient type.
--

{-# LANGUAGE FlexibleInstances          #-}
{-# LANGUAGE MultiWayIf                 #-}
{-# LANGUAGE ScopedTypeVariables        #-}
{-# LANGUAGE TypeFamilies               #-}

{-# OPTIONS_GHC -fno-warn-orphans #-}

module Data.Poly.Internal.Dense.GcdDomain
  () where

import Prelude hiding (gcd, lcm, (^))
import Control.Exception
import Control.Monad
import Control.Monad.ST
import Data.Euclidean
import Data.Maybe
import Data.Semiring (Semiring(..), Ring(), isZero, minus)
import qualified Data.Vector.Generic as G
import qualified Data.Vector.Generic.Mutable as MG

import Data.Poly.Internal.Dense

-- | @since 0.3.0.0
instance (Eq a, Ring a, GcdDomain a, G.Vector v a) => GcdDomain (Poly v a) where
  divide :: Poly v a -> Poly v a -> Maybe (Poly v a)
divide (Poly v a
xs) (Poly v a
ys) =
    v a -> Poly v a
forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
v a -> Poly v a
toPoly' (v a -> Poly v a) -> Maybe (v a) -> Maybe (Poly v a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> v a -> v a -> Maybe (v a)
forall a (v :: * -> *).
(Eq a, Ring a, GcdDomain a, Vector v a) =>
v a -> v a -> Maybe (v a)
quotient v a
xs v a
ys
  {-# INLINABLE divide #-}

  gcd :: Poly v a -> Poly v a -> Poly v a
gcd (Poly v a
xs) (Poly v a
ys)
    | v a -> Bool
forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
xs = v a -> Poly v a
forall (v :: * -> *) a. v a -> Poly v a
Poly v a
ys
    | v a -> Bool
forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
ys = v a -> Poly v a
forall (v :: * -> *) a. v a -> Poly v a
Poly v a
xs
    | v a -> Int
forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
1 = v a -> Poly v a
forall (v :: * -> *) a. v a -> Poly v a
Poly (v a -> Poly v a) -> v a -> Poly v a
forall a b. (a -> b) -> a -> b
$ a -> v a
forall (v :: * -> *) a. Vector v a => a -> v a
G.singleton (a -> v a) -> a -> v a
forall a b. (a -> b) -> a -> b
$ (a -> a -> a) -> a -> v a -> a
forall (v :: * -> *) b a.
Vector v b =>
(a -> b -> a) -> a -> v b -> a
G.foldl' a -> a -> a
forall a. GcdDomain a => a -> a -> a
gcd (v a -> a
forall (v :: * -> *) a. Vector v a => v a -> a
G.unsafeHead v a
xs) v a
ys
    | v a -> Int
forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
ys Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
1 = v a -> Poly v a
forall (v :: * -> *) a. v a -> Poly v a
Poly (v a -> Poly v a) -> v a -> Poly v a
forall a b. (a -> b) -> a -> b
$ a -> v a
forall (v :: * -> *) a. Vector v a => a -> v a
G.singleton (a -> v a) -> a -> v a
forall a b. (a -> b) -> a -> b
$ (a -> a -> a) -> a -> v a -> a
forall (v :: * -> *) b a.
Vector v b =>
(a -> b -> a) -> a -> v b -> a
G.foldl' a -> a -> a
forall a. GcdDomain a => a -> a -> a
gcd (v a -> a
forall (v :: * -> *) a. Vector v a => v a -> a
G.unsafeHead v a
ys) v a
xs
    | Bool
otherwise = v a -> Poly v a
forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
v a -> Poly v a
toPoly' (v a -> Poly v a) -> v a -> Poly v a
forall a b. (a -> b) -> a -> b
$ v a -> v a -> v a
forall a (v :: * -> *).
(Eq a, Ring a, GcdDomain a, Vector v a) =>
v a -> v a -> v a
gcdNonEmpty v a
xs v a
ys
  {-# INLINE gcd #-}

  lcm :: Poly v a -> Poly v a -> Poly v a
lcm x :: Poly v a
x@(Poly v a
xs) y :: Poly v a
y@(Poly v a
ys)
    | v a -> Bool
forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
xs Bool -> Bool -> Bool
|| v a -> Bool
forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
ys = Poly v a
forall a. Semiring a => a
zero
    | Bool
otherwise = (Poly v a
x Poly v a -> Poly v a -> Poly v a
forall a. GcdDomain a => a -> a -> a
`divide'` Poly v a -> Poly v a -> Poly v a
forall a. GcdDomain a => a -> a -> a
gcd Poly v a
x Poly v a
y) Poly v a -> Poly v a -> Poly v a
forall a. Semiring a => a -> a -> a
`times` Poly v a
y
  {-# INLINABLE lcm #-}

  coprime :: Poly v a -> Poly v a -> Bool
coprime Poly v a
x Poly v a
y = Maybe (Poly v a) -> Bool
forall a. Maybe a -> Bool
isJust (Poly v a
forall a. Semiring a => a
one Poly v a -> Poly v a -> Maybe (Poly v a)
forall a. GcdDomain a => a -> a -> Maybe a
`divide` Poly v a -> Poly v a -> Poly v a
forall a. GcdDomain a => a -> a -> a
gcd Poly v a
x Poly v a
y)
  {-# INLINABLE coprime #-}

gcdNonEmpty
  :: (Eq a, Ring a, GcdDomain a, G.Vector v a)
  => v a
  -> v a
  -> v a
gcdNonEmpty :: forall a (v :: * -> *).
(Eq a, Ring a, GcdDomain a, Vector v a) =>
v a -> v a -> v a
gcdNonEmpty v a
xs v a
ys = (forall s. ST s (v a)) -> v a
forall a. (forall s. ST s a) -> a
runST ((forall s. ST s (v a)) -> v a) -> (forall s. ST s (v a)) -> v a
forall a b. (a -> b) -> a -> b
$ do
    let x :: a
x = (a -> a -> a) -> v a -> a
forall (v :: * -> *) a. Vector v a => (a -> a -> a) -> v a -> a
G.foldl1' a -> a -> a
forall a. GcdDomain a => a -> a -> a
gcd v a
xs
        y :: a
y = (a -> a -> a) -> v a -> a
forall (v :: * -> *) a. Vector v a => (a -> a -> a) -> v a -> a
G.foldl1' a -> a -> a
forall a. GcdDomain a => a -> a -> a
gcd v a
ys
        xy :: a
xy = a
x a -> a -> a
forall a. GcdDomain a => a -> a -> a
`gcd` a
y
    Mutable v s a
xs' <- v a -> ST s (Mutable v (PrimState (ST s)) a)
forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
v a -> m (Mutable v (PrimState m) a)
G.thaw v a
xs
    Mutable v s a
ys' <- v a -> ST s (Mutable v (PrimState (ST s)) a)
forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
v a -> m (Mutable v (PrimState m) a)
G.thaw v a
ys
    Mutable v s a
zs' <- Mutable v s a -> Mutable v s a -> ST s (Mutable v s a)
forall a (v :: * -> *) s.
(Eq a, Ring a, GcdDomain a, Vector v a) =>
Mutable v s a -> Mutable v s a -> ST s (Mutable v s a)
gcdM Mutable v s a
xs' Mutable v s a
ys'

    let lenZs :: Int
lenZs = Mutable v s a -> Int
forall (v :: * -> * -> *) a s. MVector v a => v s a -> Int
MG.length Mutable v s a
zs'
        go :: a -> Int -> ST s a
go a
acc Int
0 = a -> ST s a
forall a. a -> ST s a
forall (f :: * -> *) a. Applicative f => a -> f a
pure a
acc
        go a
acc Int
n = do
          a
t <- Mutable v (PrimState (ST s)) a -> Int -> ST s a
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> m a
MG.unsafeRead Mutable v s a
Mutable v (PrimState (ST s)) a
zs' (Int
n Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1)
          a -> Int -> ST s a
go (a
acc a -> a -> a
forall a. GcdDomain a => a -> a -> a
`gcd` a
t) (Int
n Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1)
    a
a <- Mutable v (PrimState (ST s)) a -> Int -> ST s a
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> m a
MG.unsafeRead Mutable v s a
Mutable v (PrimState (ST s)) a
zs' (Int
lenZs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1)
    a
z <- a -> Int -> ST s a
go a
a (Int
lenZs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1)

    [Int] -> (Int -> ST s ()) -> ST s ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenZs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1] ((Int -> ST s ()) -> ST s ()) -> (Int -> ST s ()) -> ST s ()
forall a b. (a -> b) -> a -> b
$ \Int
i ->
      Mutable v (PrimState (ST s)) a -> (a -> a) -> Int -> ST s ()
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> (a -> a) -> Int -> m ()
MG.unsafeModify Mutable v s a
Mutable v (PrimState (ST s)) a
zs'((a -> a -> a
forall a. Semiring a => a -> a -> a
`times` a
xy) (a -> a) -> (a -> a) -> a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. (a -> a -> a
forall a. GcdDomain a => a -> a -> a
`divide'` a
z)) Int
i

    Mutable v (PrimState (ST s)) a -> ST s (v a)
forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> m (v a)
G.unsafeFreeze Mutable v s a
Mutable v (PrimState (ST s)) a
zs'
{-# INLINABLE gcdNonEmpty #-}

gcdM
  :: (Eq a, Ring a, GcdDomain a, G.Vector v a)
  => G.Mutable v s a
  -> G.Mutable v s a
  -> ST s (G.Mutable v s a)
gcdM :: forall a (v :: * -> *) s.
(Eq a, Ring a, GcdDomain a, Vector v a) =>
Mutable v s a -> Mutable v s a -> ST s (Mutable v s a)
gcdM Mutable v s a
xs Mutable v s a
ys
  | Mutable v s a -> Bool
forall (v :: * -> * -> *) a s. MVector v a => v s a -> Bool
MG.null Mutable v s a
xs = Mutable v s a -> ST s (Mutable v s a)
forall a. a -> ST s a
forall (f :: * -> *) a. Applicative f => a -> f a
pure Mutable v s a
ys
  | Mutable v s a -> Bool
forall (v :: * -> * -> *) a s. MVector v a => v s a -> Bool
MG.null Mutable v s a
ys = Mutable v s a -> ST s (Mutable v s a)
forall a. a -> ST s a
forall (f :: * -> *) a. Applicative f => a -> f a
pure Mutable v s a
xs
  | Bool
otherwise = do
  let lenXs :: Int
lenXs = Mutable v s a -> Int
forall (v :: * -> * -> *) a s. MVector v a => v s a -> Int
MG.length Mutable v s a
xs
      lenYs :: Int
lenYs = Mutable v s a -> Int
forall (v :: * -> * -> *) a s. MVector v a => v s a -> Int
MG.length Mutable v s a
ys
  a
xLast <- Mutable v (PrimState (ST s)) a -> Int -> ST s a
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> m a
MG.unsafeRead Mutable v s a
Mutable v (PrimState (ST s)) a
xs (Int
lenXs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1)
  a
yLast <- Mutable v (PrimState (ST s)) a -> Int -> ST s a
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> m a
MG.unsafeRead Mutable v s a
Mutable v (PrimState (ST s)) a
ys (Int
lenYs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1)
  let z :: a
z  = a
xLast a -> a -> a
forall a. GcdDomain a => a -> a -> a
`lcm` a
yLast
      zx :: a
zx = a
z a -> a -> a
forall a. GcdDomain a => a -> a -> a
`divide'` a
xLast
      zy :: a
zy = a
z a -> a -> a
forall a. GcdDomain a => a -> a -> a
`divide'` a
yLast

  if
    | Int
lenYs Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
lenXs
    , Just a
xy <- a
xLast a -> a -> Maybe a
forall a. GcdDomain a => a -> a -> Maybe a
`divide` a
yLast -> do
      [Int] -> (Int -> ST s ()) -> ST s ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenYs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1] ((Int -> ST s ()) -> ST s ()) -> (Int -> ST s ()) -> ST s ()
forall a b. (a -> b) -> a -> b
$ \Int
i -> do
        a
y <- Mutable v (PrimState (ST s)) a -> Int -> ST s a
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> m a
MG.unsafeRead Mutable v s a
Mutable v (PrimState (ST s)) a
ys Int
i
        Bool -> ST s () -> ST s ()
forall (f :: * -> *). Applicative f => Bool -> f () -> f ()
when (a
y a -> a -> Bool
forall a. Eq a => a -> a -> Bool
/= a
forall a. Semiring a => a
zero) (ST s () -> ST s ()) -> ST s () -> ST s ()
forall a b. (a -> b) -> a -> b
$
          Mutable v (PrimState (ST s)) a -> (a -> a) -> Int -> ST s ()
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> (a -> a) -> Int -> m ()
MG.unsafeModify
            Mutable v s a
Mutable v (PrimState (ST s)) a
xs
            (\a
x -> a
x a -> a -> a
forall a. Ring a => a -> a -> a
`minus` a
y a -> a -> a
forall a. Semiring a => a -> a -> a
`times` a
xy)
            (Int
i Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
lenXs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
lenYs)
      Mutable v s a
xs' <- (a -> Bool) -> Mutable v s a -> ST s (Mutable v s a)
forall (v :: * -> *) a s.
Vector v a =>
(a -> Bool) -> Mutable v s a -> ST s (Mutable v s a)
dropWhileEndM a -> Bool
forall a. (Eq a, Semiring a) => a -> Bool
isZero Mutable v s a
xs
      Mutable v s a -> Mutable v s a -> ST s (Mutable v s a)
forall a (v :: * -> *) s.
(Eq a, Ring a, GcdDomain a, Vector v a) =>
Mutable v s a -> Mutable v s a -> ST s (Mutable v s a)
gcdM Mutable v s a
xs' Mutable v s a
ys
    | Int
lenXs Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
lenYs
    , Just a
yx <- a
yLast a -> a -> Maybe a
forall a. GcdDomain a => a -> a -> Maybe a
`divide` a
xLast -> do
      [Int] -> (Int -> ST s ()) -> ST s ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenXs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1] ((Int -> ST s ()) -> ST s ()) -> (Int -> ST s ()) -> ST s ()
forall a b. (a -> b) -> a -> b
$ \Int
i -> do
        a
x <- Mutable v (PrimState (ST s)) a -> Int -> ST s a
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> m a
MG.unsafeRead Mutable v s a
Mutable v (PrimState (ST s)) a
xs Int
i
        Bool -> ST s () -> ST s ()
forall (f :: * -> *). Applicative f => Bool -> f () -> f ()
when (a
x a -> a -> Bool
forall a. Eq a => a -> a -> Bool
/= a
forall a. Semiring a => a
zero) (ST s () -> ST s ()) -> ST s () -> ST s ()
forall a b. (a -> b) -> a -> b
$
          Mutable v (PrimState (ST s)) a -> (a -> a) -> Int -> ST s ()
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> (a -> a) -> Int -> m ()
MG.unsafeModify
            Mutable v s a
Mutable v (PrimState (ST s)) a
ys
            (\a
y -> a
y a -> a -> a
forall a. Ring a => a -> a -> a
`minus` a
x a -> a -> a
forall a. Semiring a => a -> a -> a
`times` a
yx)
            (Int
i Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
lenYs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
lenXs)
      Mutable v s a
ys' <- (a -> Bool) -> Mutable v s a -> ST s (Mutable v s a)
forall (v :: * -> *) a s.
Vector v a =>
(a -> Bool) -> Mutable v s a -> ST s (Mutable v s a)
dropWhileEndM a -> Bool
forall a. (Eq a, Semiring a) => a -> Bool
isZero Mutable v s a
ys
      Mutable v s a -> Mutable v s a -> ST s (Mutable v s a)
forall a (v :: * -> *) s.
(Eq a, Ring a, GcdDomain a, Vector v a) =>
Mutable v s a -> Mutable v s a -> ST s (Mutable v s a)
gcdM Mutable v s a
xs Mutable v s a
ys'
    | Int
lenYs Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
lenXs -> do
      [Int] -> (Int -> ST s ()) -> ST s ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenYs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1] ((Int -> ST s ()) -> ST s ()) -> (Int -> ST s ()) -> ST s ()
forall a b. (a -> b) -> a -> b
$ \Int
i -> do
        a
y <- Mutable v (PrimState (ST s)) a -> Int -> ST s a
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> m a
MG.unsafeRead Mutable v s a
Mutable v (PrimState (ST s)) a
ys Int
i
        Mutable v (PrimState (ST s)) a -> (a -> a) -> Int -> ST s ()
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> (a -> a) -> Int -> m ()
MG.unsafeModify
          Mutable v s a
Mutable v (PrimState (ST s)) a
xs
          (\a
x -> a
x a -> a -> a
forall a. Semiring a => a -> a -> a
`times` a
zx a -> a -> a
forall a. Ring a => a -> a -> a
`minus` a
y a -> a -> a
forall a. Semiring a => a -> a -> a
`times` a
zy)
          (Int
i Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
lenXs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
lenYs)
      [Int] -> (Int -> ST s ()) -> ST s ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenXs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
lenYs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1] ((Int -> ST s ()) -> ST s ()) -> (Int -> ST s ()) -> ST s ()
forall a b. (a -> b) -> a -> b
$
        Mutable v (PrimState (ST s)) a -> (a -> a) -> Int -> ST s ()
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> (a -> a) -> Int -> m ()
MG.unsafeModify Mutable v s a
Mutable v (PrimState (ST s)) a
xs (a -> a -> a
forall a. Semiring a => a -> a -> a
`times` a
zx)
      Mutable v s a
xs' <- (a -> Bool) -> Mutable v s a -> ST s (Mutable v s a)
forall (v :: * -> *) a s.
Vector v a =>
(a -> Bool) -> Mutable v s a -> ST s (Mutable v s a)
dropWhileEndM a -> Bool
forall a. (Eq a, Semiring a) => a -> Bool
isZero Mutable v s a
xs
      Mutable v s a -> Mutable v s a -> ST s (Mutable v s a)
forall a (v :: * -> *) s.
(Eq a, Ring a, GcdDomain a, Vector v a) =>
Mutable v s a -> Mutable v s a -> ST s (Mutable v s a)
gcdM Mutable v s a
xs' Mutable v s a
ys
    | Bool
otherwise -> do
      [Int] -> (Int -> ST s ()) -> ST s ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenXs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1] ((Int -> ST s ()) -> ST s ()) -> (Int -> ST s ()) -> ST s ()
forall a b. (a -> b) -> a -> b
$ \Int
i -> do
        a
x <- Mutable v (PrimState (ST s)) a -> Int -> ST s a
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> m a
MG.unsafeRead Mutable v s a
Mutable v (PrimState (ST s)) a
xs Int
i
        Mutable v (PrimState (ST s)) a -> (a -> a) -> Int -> ST s ()
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> (a -> a) -> Int -> m ()
MG.unsafeModify
          Mutable v s a
Mutable v (PrimState (ST s)) a
ys
          (\a
y -> a
y a -> a -> a
forall a. Semiring a => a -> a -> a
`times` a
zy a -> a -> a
forall a. Ring a => a -> a -> a
`minus` a
x a -> a -> a
forall a. Semiring a => a -> a -> a
`times` a
zx)
          (Int
i Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
lenYs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
lenXs)
      [Int] -> (Int -> ST s ()) -> ST s ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenYs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
lenXs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1] ((Int -> ST s ()) -> ST s ()) -> (Int -> ST s ()) -> ST s ()
forall a b. (a -> b) -> a -> b
$
        Mutable v (PrimState (ST s)) a -> (a -> a) -> Int -> ST s ()
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> (a -> a) -> Int -> m ()
MG.unsafeModify Mutable v s a
Mutable v (PrimState (ST s)) a
ys (a -> a -> a
forall a. Semiring a => a -> a -> a
`times` a
zy)
      Mutable v s a
ys' <- (a -> Bool) -> Mutable v s a -> ST s (Mutable v s a)
forall (v :: * -> *) a s.
Vector v a =>
(a -> Bool) -> Mutable v s a -> ST s (Mutable v s a)
dropWhileEndM a -> Bool
forall a. (Eq a, Semiring a) => a -> Bool
isZero Mutable v s a
ys
      Mutable v s a -> Mutable v s a -> ST s (Mutable v s a)
forall a (v :: * -> *) s.
(Eq a, Ring a, GcdDomain a, Vector v a) =>
Mutable v s a -> Mutable v s a -> ST s (Mutable v s a)
gcdM Mutable v s a
xs Mutable v s a
ys'
{-# INLINABLE gcdM #-}

divide' :: GcdDomain a => a -> a -> a
divide' :: forall a. GcdDomain a => a -> a -> a
divide' = (a -> Maybe a -> a
forall a. a -> Maybe a -> a
fromMaybe ([Char] -> a
forall a. HasCallStack => [Char] -> a
error [Char]
"gcd: violated internal invariant") (Maybe a -> a) -> (a -> Maybe a) -> a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.) ((a -> Maybe a) -> a -> a) -> (a -> a -> Maybe a) -> a -> a -> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. a -> a -> Maybe a
forall a. GcdDomain a => a -> a -> Maybe a
divide

isZeroM
  :: (Eq a, Semiring a, G.Vector v a)
  => G.Mutable v s a
  -> ST s Bool
isZeroM :: forall a (v :: * -> *) s.
(Eq a, Semiring a, Vector v a) =>
Mutable v s a -> ST s Bool
isZeroM Mutable v s a
xs = Int -> ST s Bool
go (Mutable v s a -> Int
forall (v :: * -> * -> *) a s. MVector v a => v s a -> Int
MG.length Mutable v s a
xs)
  where
    go :: Int -> ST s Bool
go Int
0 = Bool -> ST s Bool
forall a. a -> ST s a
forall (f :: * -> *) a. Applicative f => a -> f a
pure Bool
True
    go Int
n = do
      a
x <- Mutable v (PrimState (ST s)) a -> Int -> ST s a
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> m a
MG.unsafeRead Mutable v s a
Mutable v (PrimState (ST s)) a
xs (Int
n Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1)
      if a
x a -> a -> Bool
forall a. Eq a => a -> a -> Bool
== a
forall a. Semiring a => a
zero then Int -> ST s Bool
go (Int
n Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1) else Bool -> ST s Bool
forall a. a -> ST s a
forall (f :: * -> *) a. Applicative f => a -> f a
pure Bool
False
{-# INLINE isZeroM #-}

quotient
  :: (Eq a, Ring a, GcdDomain a, G.Vector v a)
  => v a
  -> v a
  -> Maybe (v a)
quotient :: forall a (v :: * -> *).
(Eq a, Ring a, GcdDomain a, Vector v a) =>
v a -> v a -> Maybe (v a)
quotient v a
xs v a
ys
  | v a -> Bool
forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
ys = ArithException -> Maybe (v a)
forall a e. Exception e => e -> a
throw ArithException
DivideByZero
  | v a -> Bool
forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
xs = v a -> Maybe (v a)
forall a. a -> Maybe a
Just v a
xs
  | v a -> Int
forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< v a -> Int
forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
ys = Maybe (v a)
forall a. Maybe a
Nothing
  | Bool
otherwise = (forall s. ST s (Maybe (v a))) -> Maybe (v a)
forall a. (forall s. ST s a) -> a
runST ((forall s. ST s (Maybe (v a))) -> Maybe (v a))
-> (forall s. ST s (Maybe (v a))) -> Maybe (v a)
forall a b. (a -> b) -> a -> b
$ do
    let lenXs :: Int
lenXs = v a -> Int
forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs
        lenYs :: Int
lenYs = v a -> Int
forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
ys
        lenQs :: Int
lenQs = Int
lenXs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
lenYs Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
1
    Mutable v s a
qs <- Int -> ST s (Mutable v (PrimState (ST s)) a)
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
Int -> m (v (PrimState m) a)
MG.unsafeNew Int
lenQs
    Mutable v s a
rs <- Int -> ST s (Mutable v (PrimState (ST s)) a)
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
Int -> m (v (PrimState m) a)
MG.unsafeNew Int
lenXs
    Mutable v (PrimState (ST s)) a -> v a -> ST s ()
forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> v a -> m ()
G.unsafeCopy Mutable v s a
Mutable v (PrimState (ST s)) a
rs v a
xs

    let go :: Int -> ST s (Maybe (v a))
go Int
i
          | Int
i Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
0 = do
            Bool
b <- Mutable v s a -> ST s Bool
forall a (v :: * -> *) s.
(Eq a, Semiring a, Vector v a) =>
Mutable v s a -> ST s Bool
isZeroM Mutable v s a
rs
            if Bool
b
              then v a -> Maybe (v a)
forall a. a -> Maybe a
Just (v a -> Maybe (v a)) -> ST s (v a) -> ST s (Maybe (v a))
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
<$> Mutable v (PrimState (ST s)) a -> ST s (v a)
forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> m (v a)
G.unsafeFreeze Mutable v s a
Mutable v (PrimState (ST s)) a
qs
              else Maybe (v a) -> ST s (Maybe (v a))
forall a. a -> ST s a
forall (f :: * -> *) a. Applicative f => a -> f a
pure Maybe (v a)
forall a. Maybe a
Nothing
          | Bool
otherwise = do
            a
r <- Mutable v (PrimState (ST s)) a -> Int -> ST s a
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> m a
MG.unsafeRead Mutable v s a
Mutable v (PrimState (ST s)) a
rs (Int
lenYs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1 Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
i)
            case a
r a -> a -> Maybe a
forall a. GcdDomain a => a -> a -> Maybe a
`divide` v a -> a
forall (v :: * -> *) a. Vector v a => v a -> a
G.unsafeLast v a
ys of
              Maybe a
Nothing -> Maybe (v a) -> ST s (Maybe (v a))
forall a. a -> ST s a
forall (f :: * -> *) a. Applicative f => a -> f a
pure Maybe (v a)
forall a. Maybe a
Nothing
              Just a
q -> do
                Mutable v (PrimState (ST s)) a -> Int -> a -> ST s ()
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
MG.unsafeWrite Mutable v s a
Mutable v (PrimState (ST s)) a
qs Int
i a
q
                [Int] -> (Int -> ST s ()) -> ST s ()
forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenYs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1] ((Int -> ST s ()) -> ST s ()) -> (Int -> ST s ()) -> ST s ()
forall a b. (a -> b) -> a -> b
$ \Int
k ->
                  Mutable v (PrimState (ST s)) a -> (a -> a) -> Int -> ST s ()
forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> (a -> a) -> Int -> m ()
MG.unsafeModify
                    Mutable v s a
Mutable v (PrimState (ST s)) a
rs
                    (\a
c -> a
c a -> a -> a
forall a. Ring a => a -> a -> a
`minus` a
q a -> a -> a
forall a. Semiring a => a -> a -> a
`times` v a -> Int -> a
forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
ys Int
k)
                    (Int
i Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
k)
                Int -> ST s (Maybe (v a))
go (Int
i Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1)

    Int -> ST s (Maybe (v a))
go (Int
lenQs Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int
1)
{-# INLINABLE quotient #-}