-- |
-- Module:      Data.Poly.Internal.Dense
-- Copyright:   (c) 2019 Andrew Lelechenko
-- Licence:     BSD3
-- Maintainer:  Andrew Lelechenko <andrew.lelechenko@gmail.com>
--
-- Dense polynomials of one variable.
--

{-# LANGUAGE FlexibleInstances          #-}
{-# LANGUAGE GeneralizedNewtypeDeriving #-}
{-# LANGUAGE PatternSynonyms            #-}
{-# LANGUAGE ScopedTypeVariables        #-}
{-# LANGUAGE TypeApplications           #-}
{-# LANGUAGE TypeFamilies               #-}
{-# LANGUAGE ViewPatterns               #-}

module Data.Poly.Internal.Dense
  ( Poly(..)
  , VPoly
  , UPoly
  , leading
  , dropWhileEndM
  -- * Num interface
  , toPoly
  , monomial
  , scale
  , pattern X
  , eval
  , subst
  , deriv
  , integral
  -- * Semiring interface
  , toPoly'
  , monomial'
  , scale'
  , pattern X'
  , eval'
  , subst'
  , substitute'
  , deriv'
  , unscale'
  , integral'
  , timesRing
  ) where

import Prelude hiding (quotRem, quot, rem, gcd, lcm)
import Control.DeepSeq (NFData)
import Control.Monad
import Control.Monad.ST
import Data.Bits
import Data.Euclidean (Euclidean, Field, quot)
import Data.Kind
import Data.List (foldl', intersperse)
import Data.Semiring (Semiring(..), Ring(), minus)
import qualified Data.Semiring as Semiring
import qualified Data.Vector as V
import qualified Data.Vector.Generic as G
import qualified Data.Vector.Generic.Mutable as MG
import qualified Data.Vector.Unboxed as U
import GHC.Exts

-- | Polynomials of one variable with coefficients from @a@,
-- backed by a 'G.Vector' @v@ (boxed, unboxed, storable, etc.).
--
-- Use the pattern 'X' for construction:
--
-- >>> (X + 1) + (X - 1) :: VPoly Integer
-- 2 * X + 0
-- >>> (X + 1) * (X - 1) :: UPoly Int
-- 1 * X^2 + 0 * X + (-1)
--
-- Polynomials are stored normalized, without leading
-- zero coefficients, so 0 * 'X' + 1 equals to 1.
--
-- The 'Ord' instance does not make much sense mathematically,
-- it is defined only for the sake of 'Data.Set.Set', 'Data.Map.Map', etc.
--
-- Due to being polymorphic by multiple axis, the performance of `Poly` crucially
-- depends on specialisation of instances. Clients are strongly recommended
-- to compile with @ghc-options:@ @-fspecialise-aggressively@ and suggested to enable @-O2@.
--
-- @since 0.1.0.0
newtype Poly (v :: Type -> Type) (a :: Type) = Poly
  { forall (v :: * -> *) a. Poly v a -> v a
unPoly :: v a
  -- ^ Convert a 'Poly' to a vector of coefficients
  -- (first element corresponds to the constant term).
  --
  -- @since 0.1.0.0
  }
  deriving
  ( Poly v a -> Poly v a -> Bool
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
forall (v :: * -> *) a. Eq (v a) => Poly v a -> Poly v a -> Bool
/= :: Poly v a -> Poly v a -> Bool
$c/= :: forall (v :: * -> *) a. Eq (v a) => Poly v a -> Poly v a -> Bool
== :: Poly v a -> Poly v a -> Bool
$c== :: forall (v :: * -> *) a. Eq (v a) => Poly v a -> Poly v a -> Bool
Eq, Poly v a -> Poly v a -> Bool
Poly v a -> Poly v a -> Ordering
Poly v a -> Poly v a -> Poly v a
forall a.
Eq a
-> (a -> a -> Ordering)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> a)
-> (a -> a -> a)
-> Ord a
forall {v :: * -> *} {a}. Ord (v a) => Eq (Poly v a)
forall (v :: * -> *) a. Ord (v a) => Poly v a -> Poly v a -> Bool
forall (v :: * -> *) a.
Ord (v a) =>
Poly v a -> Poly v a -> Ordering
forall (v :: * -> *) a.
Ord (v a) =>
Poly v a -> Poly v a -> Poly v a
min :: Poly v a -> Poly v a -> Poly v a
$cmin :: forall (v :: * -> *) a.
Ord (v a) =>
Poly v a -> Poly v a -> Poly v a
max :: Poly v a -> Poly v a -> Poly v a
$cmax :: forall (v :: * -> *) a.
Ord (v a) =>
Poly v a -> Poly v a -> Poly v a
>= :: Poly v a -> Poly v a -> Bool
$c>= :: forall (v :: * -> *) a. Ord (v a) => Poly v a -> Poly v a -> Bool
> :: Poly v a -> Poly v a -> Bool
$c> :: forall (v :: * -> *) a. Ord (v a) => Poly v a -> Poly v a -> Bool
<= :: Poly v a -> Poly v a -> Bool
$c<= :: forall (v :: * -> *) a. Ord (v a) => Poly v a -> Poly v a -> Bool
< :: Poly v a -> Poly v a -> Bool
$c< :: forall (v :: * -> *) a. Ord (v a) => Poly v a -> Poly v a -> Bool
compare :: Poly v a -> Poly v a -> Ordering
$ccompare :: forall (v :: * -> *) a.
Ord (v a) =>
Poly v a -> Poly v a -> Ordering
Ord
  , Poly v a -> ()
forall a. (a -> ()) -> NFData a
forall (v :: * -> *) a. NFData (v a) => Poly v a -> ()
rnf :: Poly v a -> ()
$crnf :: forall (v :: * -> *) a. NFData (v a) => Poly v a -> ()
NFData -- ^ @since 0.3.2.0
  )

-- | @since 0.3.1.0
instance (Eq a, Semiring a, G.Vector v a) => IsList (Poly v a) where
  type Item (Poly v a) = a
  fromList :: [Item (Poly v a)] -> Poly v a
fromList = forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
v a -> Poly v a
toPoly' forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall (v :: * -> *) a. Vector v a => [a] -> v a
G.fromList
  fromListN :: Int -> [Item (Poly v a)] -> Poly v a
fromListN = (forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
v a -> Poly v a
toPoly' forall b c a. (b -> c) -> (a -> b) -> a -> c
.) forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall (v :: * -> *) a. Vector v a => Int -> [a] -> v a
G.fromListN
  toList :: Poly v a -> [Item (Poly v a)]
toList = forall (v :: * -> *) a. Vector v a => v a -> [a]
G.toList forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall (v :: * -> *) a. Poly v a -> v a
unPoly

instance (Show a, G.Vector v a) => Show (Poly v a) where
  showsPrec :: Int -> Poly v a -> ShowS
showsPrec Int
d (Poly v a
xs)
    | forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
xs
      = String -> ShowS
showString String
"0"
    | forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs forall a. Eq a => a -> a -> Bool
== Int
1
      = forall a. Show a => Int -> a -> ShowS
showsPrec Int
d (forall (v :: * -> *) a. Vector v a => v a -> a
G.head v a
xs)
    | Bool
otherwise
      = Bool -> ShowS -> ShowS
showParen (Int
d forall a. Ord a => a -> a -> Bool
> Int
0)
      forall a b. (a -> b) -> a -> b
$ forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl forall b c a. (b -> c) -> (a -> b) -> a -> c
(.) forall a. a -> a
id
      forall a b. (a -> b) -> a -> b
$ forall a. a -> [a] -> [a]
intersperse (String -> ShowS
showString String
" + ")
      forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) b a.
Vector v b =>
(a -> Int -> b -> a) -> a -> v b -> a
G.ifoldl (\[ShowS]
acc Int
i a
c -> Int -> a -> ShowS
showCoeff Int
i a
c forall a. a -> [a] -> [a]
: [ShowS]
acc) [] v a
xs
    where
      -- Powers are guaranteed to be non-negative
      showCoeff :: Int -> a -> String -> String
      showCoeff :: Int -> a -> ShowS
showCoeff Int
0 a
c = forall a. Show a => Int -> a -> ShowS
showsPrec Int
7 a
c
      showCoeff Int
1 a
c = forall a. Show a => Int -> a -> ShowS
showsPrec Int
7 a
c forall b c a. (b -> c) -> (a -> b) -> a -> c
. String -> ShowS
showString String
" * X"
      showCoeff Int
i a
c = forall a. Show a => Int -> a -> ShowS
showsPrec Int
7 a
c forall b c a. (b -> c) -> (a -> b) -> a -> c
. String -> ShowS
showString (String
" * X^" forall a. [a] -> [a] -> [a]
++ forall a. Show a => a -> String
show Int
i)

-- | Polynomials backed by boxed vectors.
--
-- @since 0.2.0.0
type VPoly = Poly V.Vector

-- | Polynomials backed by unboxed vectors.
--
-- @since 0.2.0.0
type UPoly = Poly U.Vector

-- | Make a 'Poly' from a list of coefficients
-- (first element corresponds to the constant term).
--
-- >>> :set -XOverloadedLists
-- >>> toPoly [1,2,3] :: VPoly Integer
-- 3 * X^2 + 2 * X + 1
-- >>> toPoly [0,0,0] :: UPoly Int
-- 0
--
-- @since 0.1.0.0
toPoly :: (Eq a, Num a, G.Vector v a) => v a -> Poly v a
toPoly :: forall a (v :: * -> *).
(Eq a, Num a, Vector v a) =>
v a -> Poly v a
toPoly = forall (v :: * -> *) a. v a -> Poly v a
Poly forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall (v :: * -> *) a. Vector v a => (a -> Bool) -> v a -> v a
dropWhileEnd (forall a. Eq a => a -> a -> Bool
== a
0)
{-# INLINABLE toPoly #-}

toPoly' :: (Eq a, Semiring a, G.Vector v a) => v a -> Poly v a
toPoly' :: forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
v a -> Poly v a
toPoly' = forall (v :: * -> *) a. v a -> Poly v a
Poly forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall (v :: * -> *) a. Vector v a => (a -> Bool) -> v a -> v a
dropWhileEnd (forall a. Eq a => a -> a -> Bool
== forall a. Semiring a => a
zero)
{-# INLINABLE toPoly' #-}

-- | Return the leading power and coefficient of a non-zero polynomial.
--
-- >>> leading ((2 * X + 1) * (2 * X^2 - 1) :: UPoly Int)
-- Just (3,4)
-- >>> leading (0 :: UPoly Int)
-- Nothing
--
-- @since 0.3.0.0
leading :: G.Vector v a => Poly v a -> Maybe (Word, a)
leading :: forall (v :: * -> *) a. Vector v a => Poly v a -> Maybe (Word, a)
leading (Poly v a
v)
  | forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
v  = forall a. Maybe a
Nothing
  | Bool
otherwise = forall a. a -> Maybe a
Just (forall a b. (Integral a, Num b) => a -> b
fromIntegral (forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
v forall a. Num a => a -> a -> a
- Int
1), forall (v :: * -> *) a. Vector v a => v a -> a
G.last v a
v)
{-# INLINABLE leading #-}

-- | Note that 'abs' = 'id' and 'signum' = 'const' 1.
instance (Eq a, Num a, G.Vector v a) => Num (Poly v a) where
  + :: Poly v a -> Poly v a -> Poly v a
(+) = (forall a (v :: * -> *).
(Eq a, Num a, Vector v a) =>
v a -> Poly v a
toPoly forall b c a. (b -> c) -> (a -> b) -> a -> c
.) forall b c a. (b -> c) -> (a -> b) -> a -> c
. coerce :: forall a b. Coercible a b => a -> b
coerce (forall (v :: * -> *) a.
Vector v a =>
(a -> a -> a) -> v a -> v a -> v a
plusPoly @v @a forall a. Num a => a -> a -> a
(+))
  (-) = (forall a (v :: * -> *).
(Eq a, Num a, Vector v a) =>
v a -> Poly v a
toPoly forall b c a. (b -> c) -> (a -> b) -> a -> c
.) forall b c a. (b -> c) -> (a -> b) -> a -> c
. coerce :: forall a b. Coercible a b => a -> b
coerce (forall (v :: * -> *) a.
Vector v a =>
(a -> a) -> (a -> a -> a) -> v a -> v a -> v a
minusPoly @v @a forall a. Num a => a -> a
negate (-))
  * :: Poly v a -> Poly v a -> Poly v a
(*) = (forall a (v :: * -> *).
(Eq a, Num a, Vector v a) =>
v a -> Poly v a
toPoly forall b c a. (b -> c) -> (a -> b) -> a -> c
.) forall b c a. (b -> c) -> (a -> b) -> a -> c
. coerce :: forall a b. Coercible a b => a -> b
coerce (forall a. a -> a
inline (forall (v :: * -> *) a.
Vector v a =>
a
-> (a -> a -> a)
-> (a -> a -> a)
-> (a -> a -> a)
-> v a
-> v a
-> v a
karatsuba @v @a a
0 forall a. Num a => a -> a -> a
(+) (-) forall a. Num a => a -> a -> a
(*)))

  negate :: Poly v a -> Poly v a
negate (Poly v a
xs) = forall (v :: * -> *) a. v a -> Poly v a
Poly forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a b.
(Vector v a, Vector v b) =>
(a -> b) -> v a -> v b
G.map forall a. Num a => a -> a
negate v a
xs
  abs :: Poly v a -> Poly v a
abs = forall a. a -> a
id
  signum :: Poly v a -> Poly v a
signum = forall a b. a -> b -> a
const Poly v a
1
  fromInteger :: Integer -> Poly v a
fromInteger Integer
n = case forall a. Num a => Integer -> a
fromInteger Integer
n of
    a
0 -> forall (v :: * -> *) a. v a -> Poly v a
Poly forall (v :: * -> *) a. Vector v a => v a
G.empty
    a
m -> forall (v :: * -> *) a. v a -> Poly v a
Poly forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a. Vector v a => a -> v a
G.singleton a
m

  {-# INLINE (+) #-}
  {-# INLINE (-) #-}
  {-# INLINE negate #-}
  {-# INLINE fromInteger #-}
  {-# INLINE (*) #-}

-- | Note that 'times' is significantly slower than '(*)' for large polynomials,
-- because Karatsuba multiplication algorithm requires subtraction, which is not
-- provided by 'Semiring' class. Use 'timesRing' instead.
instance (Eq a, Semiring a, G.Vector v a) => Semiring (Poly v a) where
  zero :: Poly v a
zero = forall (v :: * -> *) a. v a -> Poly v a
Poly forall (v :: * -> *) a. Vector v a => v a
G.empty
  one :: Poly v a
one
    | (forall a. Semiring a => a
one :: a) forall a. Eq a => a -> a -> Bool
== forall a. Semiring a => a
zero = forall a. Semiring a => a
zero
    | Bool
otherwise = forall (v :: * -> *) a. v a -> Poly v a
Poly forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a. Vector v a => a -> v a
G.singleton forall a. Semiring a => a
one

  plus :: Poly v a -> Poly v a -> Poly v a
plus  = (forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
v a -> Poly v a
toPoly' forall b c a. (b -> c) -> (a -> b) -> a -> c
.) forall b c a. (b -> c) -> (a -> b) -> a -> c
. coerce :: forall a b. Coercible a b => a -> b
coerce (forall (v :: * -> *) a.
Vector v a =>
(a -> a -> a) -> v a -> v a -> v a
plusPoly @v @a forall a. Semiring a => a -> a -> a
plus)
  times :: Poly v a -> Poly v a -> Poly v a
times = (forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
v a -> Poly v a
toPoly' forall b c a. (b -> c) -> (a -> b) -> a -> c
.) forall b c a. (b -> c) -> (a -> b) -> a -> c
. coerce :: forall a b. Coercible a b => a -> b
coerce (forall a. a -> a
inline (forall (v :: * -> *) a.
Vector v a =>
a -> (a -> a -> a) -> (a -> a -> a) -> v a -> v a -> v a
convolution @v @a forall a. Semiring a => a
zero forall a. Semiring a => a -> a -> a
plus forall a. Semiring a => a -> a -> a
times))

  {-# INLINE zero #-}
  {-# INLINE one #-}
  {-# INLINE plus #-}
  {-# INLINE times #-}

  fromNatural :: Natural -> Poly v a
fromNatural Natural
n = if a
n' forall a. Eq a => a -> a -> Bool
== forall a. Semiring a => a
zero then forall a. Semiring a => a
zero else forall (v :: * -> *) a. v a -> Poly v a
Poly forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a. Vector v a => a -> v a
G.singleton a
n'
    where
      n' :: a
      n' :: a
n' = forall a. Semiring a => Natural -> a
fromNatural Natural
n
  {-# INLINE fromNatural #-}

instance (Eq a, Ring a, G.Vector v a) => Ring (Poly v a) where
  negate :: Poly v a -> Poly v a
negate (Poly v a
xs) = forall (v :: * -> *) a. v a -> Poly v a
Poly forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a b.
(Vector v a, Vector v b) =>
(a -> b) -> v a -> v b
G.map forall a. Ring a => a -> a
Semiring.negate v a
xs
  {-# INLINABLE negate #-}

-- | Karatsuba multiplication algorithm for polynomials over rings.
timesRing :: forall v a. (Eq a, Ring a, G.Vector v a) => Poly v a -> Poly v a -> Poly v a
timesRing :: forall (v :: * -> *) a.
(Eq a, Ring a, Vector v a) =>
Poly v a -> Poly v a -> Poly v a
timesRing = (forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
v a -> Poly v a
toPoly' forall b c a. (b -> c) -> (a -> b) -> a -> c
.) forall b c a. (b -> c) -> (a -> b) -> a -> c
. coerce :: forall a b. Coercible a b => a -> b
coerce (forall a. a -> a
inline (forall (v :: * -> *) a.
Vector v a =>
a
-> (a -> a -> a)
-> (a -> a -> a)
-> (a -> a -> a)
-> v a
-> v a
-> v a
karatsuba @v @a forall a. Semiring a => a
zero forall a. Semiring a => a -> a -> a
plus forall a. Ring a => a -> a -> a
minus forall a. Semiring a => a -> a -> a
times))
{-# INLINE timesRing #-}

dropWhileEnd
  :: G.Vector v a
  => (a -> Bool)
  -> v a
  -> v a
dropWhileEnd :: forall (v :: * -> *) a. Vector v a => (a -> Bool) -> v a -> v a
dropWhileEnd a -> Bool
p v a
xs = forall (v :: * -> *) a. Vector v a => Int -> Int -> v a -> v a
G.unsafeSlice Int
0 (Int -> Int
go (forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs)) v a
xs
  where
    go :: Int -> Int
go Int
0 = Int
0
    go Int
n = if a -> Bool
p (forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
xs (Int
n forall a. Num a => a -> a -> a
- Int
1)) then Int -> Int
go (Int
n forall a. Num a => a -> a -> a
- Int
1) else Int
n
{-# INLINE dropWhileEnd #-}

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

plusPoly
  :: G.Vector v a
  => (a -> a -> a)
  -> v a
  -> v a
  -> v a
plusPoly :: forall (v :: * -> *) a.
Vector v a =>
(a -> a -> a) -> v a -> v a -> v a
plusPoly a -> a -> a
add v a
xs v a
ys = forall a. (forall s. ST s a) -> a
runST forall a b. (a -> b) -> a -> b
$ do
  let lenXs :: Int
lenXs = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs
      lenYs :: Int
lenYs = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
ys
      lenMn :: Int
lenMn = Int
lenXs forall a. Ord a => a -> a -> a
`min` Int
lenYs
      lenMx :: Int
lenMx = Int
lenXs forall a. Ord a => a -> a -> a
`max` Int
lenYs

  Mutable v s a
zs <- forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
Int -> m (v (PrimState m) a)
MG.unsafeNew Int
lenMx
  forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenMn forall a. Num a => a -> a -> a
- Int
1] forall a b. (a -> b) -> a -> b
$ \Int
i ->
    forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
MG.unsafeWrite Mutable v s a
zs Int
i (a -> a -> a
add (forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
xs Int
i) (forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
ys Int
i))
  forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> v a -> m ()
G.unsafeCopy
    (forall (v :: * -> * -> *) a s.
MVector v a =>
Int -> Int -> v s a -> v s a
MG.unsafeSlice Int
lenMn (Int
lenMx forall a. Num a => a -> a -> a
- Int
lenMn) Mutable v s a
zs)
    (forall (v :: * -> *) a. Vector v a => Int -> Int -> v a -> v a
G.unsafeSlice  Int
lenMn (Int
lenMx forall a. Num a => a -> a -> a
- Int
lenMn) (if Int
lenXs forall a. Ord a => a -> a -> Bool
<= Int
lenYs then v a
ys else v a
xs))

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

minusPoly
  :: G.Vector v a
  => (a -> a)
  -> (a -> a -> a)
  -> v a
  -> v a
  -> v a
minusPoly :: forall (v :: * -> *) a.
Vector v a =>
(a -> a) -> (a -> a -> a) -> v a -> v a -> v a
minusPoly a -> a
neg a -> a -> a
sub v a
xs v a
ys = forall a. (forall s. ST s a) -> a
runST forall a b. (a -> b) -> a -> b
$ do
  let lenXs :: Int
lenXs = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs
      lenYs :: Int
lenYs = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
ys
      lenMn :: Int
lenMn = Int
lenXs forall a. Ord a => a -> a -> a
`min` Int
lenYs
      lenMx :: Int
lenMx = Int
lenXs forall a. Ord a => a -> a -> a
`max` Int
lenYs

  Mutable v s a
zs <- forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
Int -> m (v (PrimState m) a)
MG.unsafeNew Int
lenMx
  forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenMn forall a. Num a => a -> a -> a
- Int
1] forall a b. (a -> b) -> a -> b
$ \Int
i ->
    forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
MG.unsafeWrite Mutable v s a
zs Int
i (a -> a -> a
sub (forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
xs Int
i) (forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
ys Int
i))

  if Int
lenXs forall a. Ord a => a -> a -> Bool
< Int
lenYs
    then forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
lenXs .. Int
lenYs forall a. Num a => a -> a -> a
- Int
1] forall a b. (a -> b) -> a -> b
$ \Int
i ->
      forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
MG.unsafeWrite Mutable v s a
zs Int
i (a -> a
neg (forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
ys Int
i))
    else forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> v a -> m ()
G.unsafeCopy
      (forall (v :: * -> * -> *) a s.
MVector v a =>
Int -> Int -> v s a -> v s a
MG.unsafeSlice Int
lenYs (Int
lenXs forall a. Num a => a -> a -> a
- Int
lenYs) Mutable v s a
zs)
      (forall (v :: * -> *) a. Vector v a => Int -> Int -> v a -> v a
G.unsafeSlice  Int
lenYs (Int
lenXs forall a. Num a => a -> a -> a
- Int
lenYs) v a
xs)

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

karatsubaThreshold :: Int
karatsubaThreshold :: Int
karatsubaThreshold = Int
32

karatsuba
  :: G.Vector v a
  => a
  -> (a -> a -> a)
  -> (a -> a -> a)
  -> (a -> a -> a)
  -> v a
  -> v a
  -> v a
karatsuba :: forall (v :: * -> *) a.
Vector v a =>
a
-> (a -> a -> a)
-> (a -> a -> a)
-> (a -> a -> a)
-> v a
-> v a
-> v a
karatsuba a
zer a -> a -> a
add a -> a -> a
sub a -> a -> a
mul = v a -> v a -> v a
go
  where
    conv :: v a -> v a -> v a
conv = forall a. a -> a
inline forall (v :: * -> *) a.
Vector v a =>
a -> (a -> a -> a) -> (a -> a -> a) -> v a -> v a -> v a
convolution a
zer a -> a -> a
add a -> a -> a
mul
    go :: v a -> v a -> v a
go v a
xs v a
ys
      | Int
lenXs forall a. Ord a => a -> a -> Bool
<= Int
karatsubaThreshold Bool -> Bool -> Bool
|| Int
lenYs forall a. Ord a => a -> a -> Bool
<= Int
karatsubaThreshold
      = v a -> v a -> v a
conv v a
xs v a
ys
      | Bool
otherwise = forall a. (forall s. ST s a) -> a
runST forall a b. (a -> b) -> a -> b
$ do
        Mutable v s a
zs <- forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
Int -> m (v (PrimState m) a)
MG.unsafeNew Int
lenZs
        forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenZs forall a. Num a => a -> a -> a
- Int
1] forall a b. (a -> b) -> a -> b
$ \Int
k -> do
          let z0 :: a
z0 = if Int
k forall a. Ord a => a -> a -> Bool
< forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
zs0
                   then forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
zs0 Int
k
                   else a
zer
              z11 :: a
z11 = if Int
k forall a. Num a => a -> a -> a
- Int
m forall a. Ord a => a -> a -> Bool
>= Int
0 Bool -> Bool -> Bool
&& Int
k forall a. Num a => a -> a -> a
- Int
m forall a. Ord a => a -> a -> Bool
< forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
zs11
                   then forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
zs11 (Int
k forall a. Num a => a -> a -> a
- Int
m)
                   else a
zer
              z10 :: a
z10 = if Int
k forall a. Num a => a -> a -> a
- Int
m forall a. Ord a => a -> a -> Bool
>= Int
0 Bool -> Bool -> Bool
&& Int
k forall a. Num a => a -> a -> a
- Int
m forall a. Ord a => a -> a -> Bool
< forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
zs0
                   then forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
zs0 (Int
k forall a. Num a => a -> a -> a
- Int
m)
                   else a
zer
              z12 :: a
z12 = if Int
k forall a. Num a => a -> a -> a
- Int
m forall a. Ord a => a -> a -> Bool
>= Int
0 Bool -> Bool -> Bool
&& Int
k forall a. Num a => a -> a -> a
- Int
m forall a. Ord a => a -> a -> Bool
< forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
zs2
                   then forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
zs2 (Int
k forall a. Num a => a -> a -> a
- Int
m)
                   else a
zer
              z2 :: a
z2 = if Int
k forall a. Num a => a -> a -> a
- Int
2 forall a. Num a => a -> a -> a
* Int
m forall a. Ord a => a -> a -> Bool
>= Int
0 Bool -> Bool -> Bool
&& Int
k forall a. Num a => a -> a -> a
- Int
2 forall a. Num a => a -> a -> a
* Int
m forall a. Ord a => a -> a -> Bool
< forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
zs2
                   then forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
zs2 (Int
k forall a. Num a => a -> a -> a
- Int
2 forall a. Num a => a -> a -> a
* Int
m)
                   else a
zer
          forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
MG.unsafeWrite Mutable v s a
zs Int
k (a
z0 a -> a -> a
`add` (a
z11 a -> a -> a
`sub` (a
z10 a -> a -> a
`add` a
z12)) a -> a -> a
`add` a
z2)
        forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> m (v a)
G.unsafeFreeze Mutable v s a
zs
      where
        lenXs :: Int
lenXs = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs
        lenYs :: Int
lenYs = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
ys
        lenZs :: Int
lenZs = Int
lenXs forall a. Num a => a -> a -> a
+ Int
lenYs forall a. Num a => a -> a -> a
- Int
1

        m :: Int
m    = ((Int
lenXs forall a. Ord a => a -> a -> a
`min` Int
lenYs) forall a. Num a => a -> a -> a
+ Int
1) forall a. Bits a => a -> Int -> a
`shiftR` Int
1

        xs0 :: v a
xs0  = forall (v :: * -> *) a.
(HasCallStack, Vector v a) =>
Int -> Int -> v a -> v a
G.slice Int
0 Int
m v a
xs
        xs1 :: v a
xs1  = forall (v :: * -> *) a.
(HasCallStack, Vector v a) =>
Int -> Int -> v a -> v a
G.slice Int
m (Int
lenXs forall a. Num a => a -> a -> a
- Int
m) v a
xs
        ys0 :: v a
ys0  = forall (v :: * -> *) a.
(HasCallStack, Vector v a) =>
Int -> Int -> v a -> v a
G.slice Int
0 Int
m v a
ys
        ys1 :: v a
ys1  = forall (v :: * -> *) a.
(HasCallStack, Vector v a) =>
Int -> Int -> v a -> v a
G.slice Int
m (Int
lenYs forall a. Num a => a -> a -> a
- Int
m) v a
ys

        xs01 :: v a
xs01 = forall (v :: * -> *) a.
Vector v a =>
(a -> a -> a) -> v a -> v a -> v a
plusPoly a -> a -> a
add v a
xs0 v a
xs1
        ys01 :: v a
ys01 = forall (v :: * -> *) a.
Vector v a =>
(a -> a -> a) -> v a -> v a -> v a
plusPoly a -> a -> a
add v a
ys0 v a
ys1
        zs0 :: v a
zs0  = v a -> v a -> v a
go v a
xs0 v a
ys0
        zs2 :: v a
zs2  = v a -> v a -> v a
go v a
xs1 v a
ys1
        zs11 :: v a
zs11 = v a -> v a -> v a
go v a
xs01 v a
ys01
{-# INLINABLE karatsuba #-}

convolution
  :: G.Vector v a
  => a
  -> (a -> a -> a)
  -> (a -> a -> a)
  -> v a
  -> v a
  -> v a
convolution :: forall (v :: * -> *) a.
Vector v a =>
a -> (a -> a -> a) -> (a -> a -> a) -> v a -> v a -> v a
convolution a
zer a -> a -> a
add a -> a -> a
mul = \v a
xs v a
ys ->
  let lenXs :: Int
lenXs = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs
      lenYs :: Int
lenYs = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
ys
      lenZs :: Int
lenZs = Int
lenXs forall a. Num a => a -> a -> a
+ Int
lenYs forall a. Num a => a -> a -> a
- Int
1 in
  if Int
lenXs forall a. Eq a => a -> a -> Bool
== Int
0 Bool -> Bool -> Bool
|| Int
lenYs forall a. Eq a => a -> a -> Bool
== Int
0
  then forall (v :: * -> *) a. Vector v a => v a
G.empty
  else forall (v :: * -> *) a. Vector v a => Int -> (Int -> a) -> v a
G.generate Int
lenZs forall a b. (a -> b) -> a -> b
$ \Int
k -> forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl'
    (\a
acc Int
i -> a
acc a -> a -> a
`add` a -> a -> a
mul (forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
xs Int
i) (forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
ys (Int
k forall a. Num a => a -> a -> a
- Int
i)))
    a
zer
    [forall a. Ord a => a -> a -> a
max (Int
k forall a. Num a => a -> a -> a
- Int
lenYs forall a. Num a => a -> a -> a
+ Int
1) Int
0 .. forall a. Ord a => a -> a -> a
min Int
k (Int
lenXs forall a. Num a => a -> a -> a
- Int
1)]
{-# INLINABLE convolution #-}

-- | Create a monomial from a power and a coefficient.
--
-- @since 0.3.0.0
monomial :: (Eq a, Num a, G.Vector v a) => Word -> a -> Poly v a
monomial :: forall a (v :: * -> *).
(Eq a, Num a, Vector v a) =>
Word -> a -> Poly v a
monomial Word
_ a
0 = forall (v :: * -> *) a. v a -> Poly v a
Poly forall (v :: * -> *) a. Vector v a => v a
G.empty
monomial Word
p a
c = forall (v :: * -> *) a. v a -> Poly v a
Poly forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a. Vector v a => Int -> (Int -> a) -> v a
G.generate (forall a b. (Integral a, Num b) => a -> b
fromIntegral Word
p forall a. Num a => a -> a -> a
+ Int
1) (\Int
i -> if Int
i forall a. Eq a => a -> a -> Bool
== forall a b. (Integral a, Num b) => a -> b
fromIntegral Word
p then a
c else a
0)
{-# INLINE monomial #-}

monomial' :: (Eq a, Semiring a, G.Vector v a) => Word -> a -> Poly v a
monomial' :: forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
Word -> a -> Poly v a
monomial' Word
p a
c
  | a
c forall a. Eq a => a -> a -> Bool
== forall a. Semiring a => a
zero = forall (v :: * -> *) a. v a -> Poly v a
Poly forall (v :: * -> *) a. Vector v a => v a
G.empty
  | Bool
otherwise = forall (v :: * -> *) a. v a -> Poly v a
Poly forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a. Vector v a => Int -> (Int -> a) -> v a
G.generate (forall a b. (Integral a, Num b) => a -> b
fromIntegral Word
p forall a. Num a => a -> a -> a
+ Int
1) (\Int
i -> if Int
i forall a. Eq a => a -> a -> Bool
== forall a b. (Integral a, Num b) => a -> b
fromIntegral Word
p then a
c else forall a. Semiring a => a
zero)
{-# INLINE monomial' #-}

scaleInternal
  :: G.Vector v a
  => a
  -> (a -> a -> a)
  -> Word
  -> a
  -> v a
  -> v a
scaleInternal :: forall (v :: * -> *) a.
Vector v a =>
a -> (a -> a -> a) -> Word -> a -> v a -> v a
scaleInternal a
zer a -> a -> a
mul Word
yp a
yc v a
xs = forall a. (forall s. ST s a) -> a
runST forall a b. (a -> b) -> a -> b
$ do
  let lenXs :: Int
lenXs = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs
  Mutable v s a
zs <- forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
Int -> m (v (PrimState m) a)
MG.unsafeNew (forall a b. (Integral a, Num b) => a -> b
fromIntegral Word
yp forall a. Num a => a -> a -> a
+ Int
lenXs)
  forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. forall a b. (Integral a, Num b) => a -> b
fromIntegral Word
yp forall a. Num a => a -> a -> a
- Int
1] forall a b. (a -> b) -> a -> b
$ \Int
k ->
    forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
MG.unsafeWrite Mutable v s a
zs Int
k a
zer
  forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenXs forall a. Num a => a -> a -> a
- Int
1] forall a b. (a -> b) -> a -> b
$ \Int
k ->
    forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
MG.unsafeWrite Mutable v s a
zs (forall a b. (Integral a, Num b) => a -> b
fromIntegral Word
yp forall a. Num a => a -> a -> a
+ Int
k) (a -> a -> a
mul a
yc forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
xs Int
k)
  forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> m (v a)
G.unsafeFreeze Mutable v s a
zs
{-# INLINABLE scaleInternal #-}

-- | Multiply a polynomial by a monomial, expressed as a power and a coefficient.
--
-- >>> scale 2 3 (X^2 + 1) :: UPoly Int
-- 3 * X^4 + 0 * X^3 + 3 * X^2 + 0 * X + 0
--
-- @since 0.3.0.0
scale :: (Eq a, Num a, G.Vector v a) => Word -> a -> Poly v a -> Poly v a
scale :: forall a (v :: * -> *).
(Eq a, Num a, Vector v a) =>
Word -> a -> Poly v a -> Poly v a
scale Word
yp a
yc (Poly v a
xs) = forall a (v :: * -> *).
(Eq a, Num a, Vector v a) =>
v a -> Poly v a
toPoly forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a.
Vector v a =>
a -> (a -> a -> a) -> Word -> a -> v a -> v a
scaleInternal a
0 forall a. Num a => a -> a -> a
(*) Word
yp a
yc v a
xs

scale' :: (Eq a, Semiring a, G.Vector v a) => Word -> a -> Poly v a -> Poly v a
scale' :: forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
Word -> a -> Poly v a -> Poly v a
scale' Word
yp a
yc (Poly v a
xs) = forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
v a -> Poly v a
toPoly' forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a.
Vector v a =>
a -> (a -> a -> a) -> Word -> a -> v a -> v a
scaleInternal forall a. Semiring a => a
zero forall a. Semiring a => a -> a -> a
times Word
yp a
yc v a
xs

unscale' :: (Eq a, Euclidean a, G.Vector v a) => Word -> a -> Poly v a -> Poly v a
unscale' :: forall a (v :: * -> *).
(Eq a, Euclidean a, Vector v a) =>
Word -> a -> Poly v a -> Poly v a
unscale' Word
yp a
yc (Poly v a
xs) = forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
v a -> Poly v a
toPoly' forall a b. (a -> b) -> a -> b
$ forall a. (forall s. ST s a) -> a
runST forall a b. (a -> b) -> a -> b
$ do
  let lenZs :: Int
lenZs = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs forall a. Num a => a -> a -> a
- forall a b. (Integral a, Num b) => a -> b
fromIntegral Word
yp
  Mutable v s a
zs <- forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
Int -> m (v (PrimState m) a)
MG.unsafeNew Int
lenZs
  forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenZs forall a. Num a => a -> a -> a
- Int
1] forall a b. (a -> b) -> a -> b
$ \Int
k ->
    forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
MG.unsafeWrite Mutable v s a
zs Int
k (forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
xs (Int
k forall a. Num a => a -> a -> a
+ forall a b. (Integral a, Num b) => a -> b
fromIntegral Word
yp) forall a. Euclidean a => a -> a -> a
`quot` a
yc)
  forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> m (v a)
G.unsafeFreeze Mutable v s a
zs
{-# INLINABLE unscale' #-}

data StrictPair a b = !a :*: !b

infixr 1 :*:

fst' :: StrictPair a b -> a
fst' :: forall a b. StrictPair a b -> a
fst' (a
a :*: b
_) = a
a

-- | Evaluate the polynomial at a given point.
--
-- >>> eval (X^2 + 1 :: UPoly Int) 3
-- 10
--
-- @since 0.2.0.0
eval :: (Num a, G.Vector v a) => Poly v a -> a -> a
eval :: forall a (v :: * -> *). (Num a, Vector v a) => Poly v a -> a -> a
eval = forall (v :: * -> *) a b.
(Vector v a, Num b) =>
(a -> b -> b) -> Poly v a -> b -> b
substitute forall a. Num a => a -> a -> a
(*)
{-# INLINE eval #-}

eval' :: (Semiring a, G.Vector v a) => Poly v a -> a -> a
eval' :: forall a (v :: * -> *).
(Semiring a, Vector v a) =>
Poly v a -> a -> a
eval' = forall (v :: * -> *) a b.
(Vector v a, Semiring b) =>
(a -> b -> b) -> Poly v a -> b -> b
substitute' forall a. Semiring a => a -> a -> a
times
{-# INLINE eval' #-}

-- | Substitute another polynomial instead of 'X'.
--
-- >>> subst (X^2 + 1 :: UPoly Int) (X + 1 :: UPoly Int)
-- 1 * X^2 + 2 * X + 2
--
-- @since 0.3.3.0
subst :: (Eq a, Num a, G.Vector v a, G.Vector w a) => Poly v a -> Poly w a -> Poly w a
subst :: forall a (v :: * -> *) (w :: * -> *).
(Eq a, Num a, Vector v a, Vector w a) =>
Poly v a -> Poly w a -> Poly w a
subst = forall (v :: * -> *) a b.
(Vector v a, Num b) =>
(a -> b -> b) -> Poly v a -> b -> b
substitute (forall a (v :: * -> *).
(Eq a, Num a, Vector v a) =>
Word -> a -> Poly v a -> Poly v a
scale Word
0)
{-# INLINE subst #-}

subst' :: (Eq a, Semiring a, G.Vector v a, G.Vector w a) => Poly v a -> Poly w a -> Poly w a
subst' :: forall a (v :: * -> *) (w :: * -> *).
(Eq a, Semiring a, Vector v a, Vector w a) =>
Poly v a -> Poly w a -> Poly w a
subst' = forall (v :: * -> *) a b.
(Vector v a, Semiring b) =>
(a -> b -> b) -> Poly v a -> b -> b
substitute' (forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
Word -> a -> Poly v a -> Poly v a
scale' Word
0)
{-# INLINE subst' #-}

substitute :: (G.Vector v a, Num b) => (a -> b -> b) -> Poly v a -> b -> b
substitute :: forall (v :: * -> *) a b.
(Vector v a, Num b) =>
(a -> b -> b) -> Poly v a -> b -> b
substitute a -> b -> b
f (Poly v a
cs) b
x = forall a b. StrictPair a b -> a
fst' forall a b. (a -> b) -> a -> b
$
  forall (v :: * -> *) b a.
Vector v b =>
(a -> b -> a) -> a -> v b -> a
G.foldl' (\(b
acc :*: b
xn) a
cn -> b
acc forall a. Num a => a -> a -> a
+ a -> b -> b
f a
cn b
xn forall a b. a -> b -> StrictPair a b
:*: b
x forall a. Num a => a -> a -> a
* b
xn) (b
0 forall a b. a -> b -> StrictPair a b
:*: b
1) v a
cs
{-# INLINE substitute #-}

substitute' :: (G.Vector v a, Semiring b) => (a -> b -> b) -> Poly v a -> b -> b
substitute' :: forall (v :: * -> *) a b.
(Vector v a, Semiring b) =>
(a -> b -> b) -> Poly v a -> b -> b
substitute' a -> b -> b
f (Poly v a
cs) b
x = forall a b. StrictPair a b -> a
fst' forall a b. (a -> b) -> a -> b
$
  forall (v :: * -> *) b a.
Vector v b =>
(a -> b -> a) -> a -> v b -> a
G.foldl' (\(b
acc :*: b
xn) a
cn -> b
acc forall a. Semiring a => a -> a -> a
`plus` a -> b -> b
f a
cn b
xn forall a b. a -> b -> StrictPair a b
:*: b
x forall a. Semiring a => a -> a -> a
`times` b
xn) (forall a. Semiring a => a
zero forall a b. a -> b -> StrictPair a b
:*: forall a. Semiring a => a
one) v a
cs
{-# INLINE substitute' #-}

-- | Take the derivative of the polynomial.
--
-- >>> deriv (X^3 + 3 * X) :: UPoly Int
-- 3 * X^2 + 0 * X + 3
--
-- @since 0.2.0.0
deriv :: (Eq a, Num a, G.Vector v a) => Poly v a -> Poly v a
deriv :: forall a (v :: * -> *).
(Eq a, Num a, Vector v a) =>
Poly v a -> Poly v a
deriv (Poly v a
xs)
  | forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
xs = forall (v :: * -> *) a. v a -> Poly v a
Poly forall (v :: * -> *) a. Vector v a => v a
G.empty
  | Bool
otherwise = forall a (v :: * -> *).
(Eq a, Num a, Vector v a) =>
v a -> Poly v a
toPoly forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a b.
(Vector v a, Vector v b) =>
(Int -> a -> b) -> v a -> v b
G.imap (\Int
i a
x -> forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int
i forall a. Num a => a -> a -> a
+ Int
1) forall a. Num a => a -> a -> a
* a
x) forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a. Vector v a => v a -> v a
G.tail v a
xs
{-# INLINE deriv #-}

deriv' :: (Eq a, Semiring a, G.Vector v a) => Poly v a -> Poly v a
deriv' :: forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
Poly v a -> Poly v a
deriv' (Poly v a
xs)
  | forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
xs = forall (v :: * -> *) a. v a -> Poly v a
Poly forall (v :: * -> *) a. Vector v a => v a
G.empty
  | Bool
otherwise = forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
v a -> Poly v a
toPoly' forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a b.
(Vector v a, Vector v b) =>
(Int -> a -> b) -> v a -> v b
G.imap (\Int
i a
x -> forall a. Semiring a => Natural -> a
fromNatural (forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int
i forall a. Num a => a -> a -> a
+ Int
1)) forall a. Semiring a => a -> a -> a
`times` a
x) forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a. Vector v a => v a -> v a
G.tail v a
xs
{-# INLINE deriv' #-}

-- | Compute an indefinite integral of the polynomial,
-- setting the constant term to zero.
--
-- >>> integral (3 * X^2 + 3) :: UPoly Double
-- 1.0 * X^3 + 0.0 * X^2 + 3.0 * X + 0.0
--
-- @since 0.2.0.0
integral :: (Eq a, Fractional a, G.Vector v a) => Poly v a -> Poly v a
integral :: forall a (v :: * -> *).
(Eq a, Fractional a, Vector v a) =>
Poly v a -> Poly v a
integral (Poly v a
xs)
  | forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
xs = forall (v :: * -> *) a. v a -> Poly v a
Poly forall (v :: * -> *) a. Vector v a => v a
G.empty
  | Bool
otherwise = forall a (v :: * -> *).
(Eq a, Num a, Vector v a) =>
v a -> Poly v a
toPoly forall a b. (a -> b) -> a -> b
$ forall a. (forall s. ST s a) -> a
runST forall a b. (a -> b) -> a -> b
$ do
    Mutable v s a
zs <- forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
Int -> m (v (PrimState m) a)
MG.unsafeNew (Int
lenXs forall a. Num a => a -> a -> a
+ Int
1)
    forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
MG.unsafeWrite Mutable v s a
zs Int
0 a
0
    forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenXs forall a. Num a => a -> a -> a
- Int
1] forall a b. (a -> b) -> a -> b
$ \Int
i ->
      forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
MG.unsafeWrite Mutable v s a
zs (Int
i forall a. Num a => a -> a -> a
+ Int
1) (forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
xs Int
i forall a. Num a => a -> a -> a
* forall a. Fractional a => a -> a
recip (forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
i forall a. Num a => a -> a -> a
+ a
1))
    forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> m (v a)
G.unsafeFreeze Mutable v s a
zs
    where
      lenXs :: Int
lenXs = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs
{-# INLINABLE integral #-}

integral' :: (Eq a, Field a, G.Vector v a) => Poly v a -> Poly v a
integral' :: forall a (v :: * -> *).
(Eq a, Field a, Vector v a) =>
Poly v a -> Poly v a
integral' (Poly v a
xs)
  | forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
xs = forall (v :: * -> *) a. v a -> Poly v a
Poly forall (v :: * -> *) a. Vector v a => v a
G.empty
  | Bool
otherwise = forall a (v :: * -> *).
(Eq a, Semiring a, Vector v a) =>
v a -> Poly v a
toPoly' forall a b. (a -> b) -> a -> b
$ forall a. (forall s. ST s a) -> a
runST forall a b. (a -> b) -> a -> b
$ do
    Mutable v s a
zs <- forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
Int -> m (v (PrimState m) a)
MG.unsafeNew (Int
lenXs forall a. Num a => a -> a -> a
+ Int
1)
    forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
MG.unsafeWrite Mutable v s a
zs forall a. Semiring a => a
zero forall a. Semiring a => a
zero
    forall (t :: * -> *) (m :: * -> *) a b.
(Foldable t, Monad m) =>
t a -> (a -> m b) -> m ()
forM_ [Int
0 .. Int
lenXs forall a. Num a => a -> a -> a
- Int
1] forall a b. (a -> b) -> a -> b
$ \Int
i ->
      forall (m :: * -> *) (v :: * -> * -> *) a.
(PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
MG.unsafeWrite Mutable v s a
zs (Int
i forall a. Num a => a -> a -> a
+ Int
1) (forall (v :: * -> *) a. Vector v a => v a -> Int -> a
G.unsafeIndex v a
xs Int
i forall a. Euclidean a => a -> a -> a
`quot` forall a b. (Integral a, Ring b) => a -> b
Semiring.fromIntegral (Int
i forall a. Num a => a -> a -> a
+ Int
1))
    forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> m (v a)
G.unsafeFreeze Mutable v s a
zs
    where
      lenXs :: Int
lenXs = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs
{-# INLINABLE integral' #-}

-- | The polynomial 'X'.
--
-- > X == monomial 1 1
--
-- @since 0.2.0.0
pattern X :: (Eq a, Num a, G.Vector v a) => Poly v a
pattern $bX :: forall a (v :: * -> *). (Eq a, Num a, Vector v a) => Poly v a
$mX :: forall {r} {a} {v :: * -> *}.
(Eq a, Num a, Vector v a) =>
Poly v a -> ((# #) -> r) -> ((# #) -> r) -> r
X <- (isVar -> True)
  where X = forall a (v :: * -> *). (Eq a, Num a, Vector v a) => Poly v a
var

var :: forall a v. (Eq a, Num a, G.Vector v a) => Poly v a
var :: forall a (v :: * -> *). (Eq a, Num a, Vector v a) => Poly v a
var
  | (a
1 :: a) forall a. Eq a => a -> a -> Bool
== a
0 = forall (v :: * -> *) a. v a -> Poly v a
Poly forall (v :: * -> *) a. Vector v a => v a
G.empty
  | Bool
otherwise     = forall (v :: * -> *) a. v a -> Poly v a
Poly forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a. Vector v a => [a] -> v a
G.fromList [a
0, a
1]
{-# INLINE var #-}

isVar :: forall v a. (Eq a, Num a, G.Vector v a) => Poly v a -> Bool
isVar :: forall (v :: * -> *) a.
(Eq a, Num a, Vector v a) =>
Poly v a -> Bool
isVar (Poly v a
xs)
  | (a
1 :: a) forall a. Eq a => a -> a -> Bool
== a
0 = forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
xs
  | Bool
otherwise     = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs forall a. Eq a => a -> a -> Bool
== Int
2 Bool -> Bool -> Bool
&& v a
xs forall (v :: * -> *) a.
(HasCallStack, Vector v a) =>
v a -> Int -> a
G.! Int
0 forall a. Eq a => a -> a -> Bool
== a
0 Bool -> Bool -> Bool
&& v a
xs forall (v :: * -> *) a.
(HasCallStack, Vector v a) =>
v a -> Int -> a
G.! Int
1 forall a. Eq a => a -> a -> Bool
== a
1
{-# INLINE isVar #-}

-- | Create an identity polynomial.
pattern X' :: (Eq a, Semiring a, G.Vector v a) => Poly v a
pattern $bX' :: forall a (v :: * -> *). (Eq a, Semiring a, Vector v a) => Poly v a
$mX' :: forall {r} {a} {v :: * -> *}.
(Eq a, Semiring a, Vector v a) =>
Poly v a -> ((# #) -> r) -> ((# #) -> r) -> r
X' <- (isVar' -> True)
  where X' = forall a (v :: * -> *). (Eq a, Semiring a, Vector v a) => Poly v a
var'

var' :: forall a v. (Eq a, Semiring a, G.Vector v a) => Poly v a
var' :: forall a (v :: * -> *). (Eq a, Semiring a, Vector v a) => Poly v a
var'
  | (forall a. Semiring a => a
one :: a) forall a. Eq a => a -> a -> Bool
== forall a. Semiring a => a
zero = forall (v :: * -> *) a. v a -> Poly v a
Poly forall (v :: * -> *) a. Vector v a => v a
G.empty
  | Bool
otherwise          = forall (v :: * -> *) a. v a -> Poly v a
Poly forall a b. (a -> b) -> a -> b
$ forall (v :: * -> *) a. Vector v a => [a] -> v a
G.fromList [forall a. Semiring a => a
zero, forall a. Semiring a => a
one]
{-# INLINE var' #-}

isVar' :: forall v a. (Eq a, Semiring a, G.Vector v a) => Poly v a -> Bool
isVar' :: forall (v :: * -> *) a.
(Eq a, Semiring a, Vector v a) =>
Poly v a -> Bool
isVar' (Poly v a
xs)
  | (forall a. Semiring a => a
one :: a) forall a. Eq a => a -> a -> Bool
== forall a. Semiring a => a
zero = forall (v :: * -> *) a. Vector v a => v a -> Bool
G.null v a
xs
  | Bool
otherwise          = forall (v :: * -> *) a. Vector v a => v a -> Int
G.length v a
xs forall a. Eq a => a -> a -> Bool
== Int
2 Bool -> Bool -> Bool
&& v a
xs forall (v :: * -> *) a.
(HasCallStack, Vector v a) =>
v a -> Int -> a
G.! Int
0 forall a. Eq a => a -> a -> Bool
== forall a. Semiring a => a
zero Bool -> Bool -> Bool
&& v a
xs forall (v :: * -> *) a.
(HasCallStack, Vector v a) =>
v a -> Int -> a
G.! Int
1 forall a. Eq a => a -> a -> Bool
== forall a. Semiring a => a
one
{-# INLINE isVar' #-}