Skip to content

Commit

Permalink
quotientsByGCD
Browse files Browse the repository at this point in the history
  • Loading branch information
stla committed Apr 26, 2024
1 parent f86bfdc commit f03b483
Show file tree
Hide file tree
Showing 2 changed files with 27 additions and 10 deletions.
22 changes: 12 additions & 10 deletions scripts/script6.hs
Original file line number Diff line number Diff line change
@@ -1,11 +1,13 @@
import Math.Algebra.Hspray
import Data.Ratio
x = lone 1 :: Spray Rational
y = lone 2 :: Spray Rational
z = lone 3 :: Spray Rational
poly = ((2%3) *^ (x^**^3 ^*^ y ^*^ z) ^+^ x^**^2) ^*^ ((7%4) *^ (x ^*^ y ^*^ z))
p1 = x^**^2 ^*^ y^**^2
p2 = x ^+^ (y ^*^ z^**^2)
r = bbDivision poly [p1, p2]
-- prettySpray show "x" poly
-- "(7 % 4) * x^(3, 1, 1) + (7 % 6) * x^(4, 2, 2)"
x = qlone 1
y = qlone 2

sprayA = (x^**^3 ^-^ 4 *^ x^**^2 ^+^ 5 *^ x) ^*^ y^**^2
sprayB = (x^**^2 ^-^ 6 *^ x) ^*^ y

(qA, qB) = quotientsByGCD sprayA sprayB

g = gcdSpray sprayA sprayB
((qA', _), (qB', _)) = (sprayDivision sprayA g, sprayDivision sprayB g)


15 changes: 15 additions & 0 deletions src/Math/Algebra/Hspray.hs
Original file line number Diff line number Diff line change
Expand Up @@ -217,6 +217,7 @@ module Math.Algebra.Hspray
, isPolynomialOf
, bombieriSpray
, collinearSprays
, quotientsByGCD
) where
import qualified Algebra.Additive as AlgAdd
import qualified Algebra.Differential as AlgDiff
Expand Down Expand Up @@ -1746,6 +1747,20 @@ collinearSprays spray1 spray2 = r *^ spray2 == spray1

-- division stuff -------------------------------------------------------------

-- | quotients of two univariate sprays by their gcd
quotientsByGCD ::
(Eq a, AlgField.C a) => Spray a -> Spray a -> (Spray a, Spray a)
quotientsByGCD sprayA sprayB =
go (sprayA, sprayB) (unitSpray, zeroSpray) (zeroSpray, unitSpray)
where
go (oldr, r) (olds, s) (oldt, t)
| isZeroSpray r = (t, s)
| otherwise =
go (r, remainder) (s, olds ^-^ quo ^*^ s) (t, oldt ^-^ quo ^*^ t)
where
(quo, remainder) = sprayDivision oldr r


-- | index of the maximum of a list
-- maxWithIndex :: Ord a => [a] -> (Int, a)
-- maxWithIndex = maximumBy (comparing snd) . zip [0 .. ]
Expand Down

0 comments on commit f03b483

Please sign in to comment.