diff --git a/ChangeLog b/ChangeLog
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,5 +1,8 @@
 ChangeLog
 
+v0.2.2.0 2014/03/04
+	(2014/03/04) PS1 - added high-precision quadratic equation solver.
+
 v0.2.1.1 2014/02/05
 	(2014/02/04) PS1 - improved efficiency of square root.
 	(2014/02/04) PS1 - added profiling options.
diff --git a/Data/Number/FixedPrec.hs b/Data/Number/FixedPrec.hs
--- a/Data/Number/FixedPrec.hs
+++ b/Data/Number/FixedPrec.hs
@@ -24,6 +24,7 @@
   
   -- * Other operations
   fractional,
+  solve_quadratic,
   log_double
   ) where
 
@@ -88,6 +89,28 @@
       a = 2^(b `div` 2)
       b = hibit (fromIntegral n)
 
+-- | Find the ceiling of the larger solution of a quadratic
+-- equation. Specifically, given the polynomial /p/(/x/) = /x/² + /bx/
+-- + /c/, where /b/ and /c/ are integers, find the smallest integer
+-- /x/ ≥ -/b/\/2 satisfying /p/(/x/) ≥ 0, if /b/^2 - 4/c/ ≥ 0.
+-- 
+-- This is done using integer arithmetic, so there are no rounding
+-- errors. It generalizes 'intsqrt'.
+intquad :: Integer -> Integer -> Maybe Integer
+intquad b c
+  | b^2 - 4*c < 0 = Nothing
+  | x1^2 + b*x1 + c >= 0 = Just x1
+  | otherwise = Just (iterate x0)
+  where
+    iterate x
+      | px <= 0 && px + 2*x+1+b > 0 = x+1
+      | otherwise = iterate ((x^2 - c) `div` (2*x + b))
+      where
+        px = x^2 + b*x + c
+    x1 = -(b `div` 2)
+    x0 = x1 + 2^(h `div` 2)
+    h = hibit (b^2 - 4*c)
+
 -- ----------------------------------------------------------------------
 -- ** Other general-purpose functions
             
@@ -100,7 +123,7 @@
 -- In other words, let /n/ = ⌊log[sub /b/] /x/⌋ and 
 -- /r/ = /x/ /b/[sup −/n/]. This can be more efficient than 'floor'
 -- ('logBase' /b/ /x/) depending on the type; moreover, it also works
--- for exact types such as 'Rational' and 'EReal'.
+-- for exact types such as 'Rational' and 'QRootTwo'.
 floorlog :: (Fractional b, Ord b) => b -> b -> (Integer, b)
 floorlog b x 
     | x <= 0            = error "floorlog: argument not positive"
@@ -274,6 +297,22 @@
 fractional a@(F x) = F (x `mod` one) where
   p = getprec a
   one = (decshiftL p 1)
+
+-- | Solve the quadratic equation /x/^2 + /bx/ + /c/ = 0 with maximal
+-- possible precision, using a numerically stable method. Return the
+-- pair (/x1/, /x2/) of solutions with /x1/ <= /x2/, or 'Nothing' if no
+-- solution exists.
+-- 
+-- This is far more precise, and probably more efficient, than naively
+-- using the quadratic formula.
+solve_quadratic :: (Precision e) => FixedPrec e -> FixedPrec e -> Maybe (FixedPrec e, FixedPrec e)
+solve_quadratic b c = do
+  let p = getprec b + 3
+      b' = floor (b * 10^p)
+      c' = floor (c * 100^p)
+  x2' <- intquad b' c'
+  let x1' = -b' - x2'
+  return (fromInteger x1' / 10^p, fromInteger x2' / 10^p)
 
 -- ----------------------------------------------------------------------
 -- ** Power series
diff --git a/fixedprec.cabal b/fixedprec.cabal
--- a/fixedprec.cabal
+++ b/fixedprec.cabal
@@ -1,6 +1,6 @@
 name:           fixedprec
 -- Don't forget to update the ChangeLog
-version:        0.2.1.1
+version:        0.2.2.0
 license:        BSD3
 cabal-version:  >= 1.8
 build-type:	Simple