polynomial: introduce a helper function to divide polynomial by ring element.

[ginac.git] / ginac / polynomial / upoly.hpp

diff --git a/ginac/polynomial/upoly.hpp b/ginac/polynomial/upoly.hpp
index 3564ae02bc3f80a5cf27c8ccad408c075e86abbf..2e9092f1d47a24ad1888a6af5a881f353f8682c4 100644 (file)
--- a/ginac/polynomial/upoly.hpp
+++ b/ginac/polynomial/upoly.hpp
@@ -8,6 +8,7 @@
 #include <stdexcept>
 #include <iterator>
 #include <limits>
+#include "ring_traits.hpp"
 #include "debug.hpp"
 #include "compiler.h"
 
@@ -113,6 +114,36 @@ operator*=(T& p, const typename T::value_type& c)
        return p;
 }
 
+/// Divide the polynomial @a p by the ring element @a c, put the result
+/// into @a r. If @a p is not divisible by @a c @a r is undefined.
+template<typename T> bool divide(T& r, const T& p, const typename T::value_type& c)
+{
+       if (p.empty()) {
+               r.clear();
+               return true;
+       }
+       if (r.size() != p.size())
+               r.resize(p.size());
+       bool divisible = true;
+       for (std::size_t i = p.size(); i-- != 0; ) {
+               divisible = div(r[i], p[i], c);
+               if (!divisible)
+                       break;
+       }
+       return divisible;
+}
+
+template<typename T> bool divide(T& p, const typename T::value_type& c)
+{
+       if (p.empty())
+               return true;
+       T r(p.size());
+       bool divisible = divide(r, p, c);
+       if (divisible)
+               swap(p, r);
+       return divisible;
+}
+
 // Convert Z[x] -> Z/p[x]
 
 static void