Resultant patch

[Ralf Stephan]

Please include. This works with any polynomial.

Best,
ralf

diff -u --recursive GiNaC-1.2.1-orig/doc/tutorial/ginac.texi GiNaC-1.2.1/doc/tutorial/ginac.texi
--- GiNaC-1.2.1-orig/doc/tutorial/ginac.texi	Mon Mar 15 18:47:17 2004
+++ GiNaC-1.2.1/doc/tutorial/ginac.texi	Sun Jun 20 15:45:19 2004
@@ -4371,7 +4399,7 @@
 original polynomial.
 
 
- at subsection GCD and LCM
+ at subsection GCD, LCM, resultant
 @cindex GCD
 @cindex LCM
 @cindex @code{gcd()}
@@ -4408,7 +4436,37 @@
 @}
 @end example
 
+ at cindex resultant
+ at cindex @code{resultant()}
+
+The resultant of two expressions only makes sense with polynomials.
+It is always computed with respect to a specific symbol within the
+expressions. The function has the interface
+
+ at example
+ex resultant(const ex & a, const ex & b, const symbol & s);
+ at end example
 
+Resultants are commutative. The following example computes the resultant 
+of two expressions with respect to both @code{x} and @code{y}:
+
+ at example
+#include <ginac/ginac.h>
+using namespace GiNaC;
+
+int main()
+@{
+    symbol x("x"), y("y");
+
+    ex e1 = x+pow(y,2), e2 = 2*pow(x,3)-1; // x+y^2, 2*x^3-1
+    ex r;
+    
+    r = resultant (e1, e2, x); 
+    // -> 1+2*y^6
+    r = resultant (e1, e2, y); 
+    // -> 1-4*x^3+4*x^6
+@}
+ at end example
 @subsection Square-free decomposition
 @cindex square-free decomposition
 @cindex factorization
diff -u --recursive GiNaC-1.2.1-orig/ginac/normal.cpp GiNaC-1.2.1/ginac/normal.cpp
--- GiNaC-1.2.1-orig/ginac/normal.cpp	Thu Jan  8 16:06:50 2004
+++ GiNaC-1.2.1/ginac/normal.cpp	Sun Jun 20 11:58:11 2004
@@ -2374,5 +2374,33 @@
 		return e;
 }
 
+/** Resultant of two expressions e1,e2 with respect to symbol s.
+ *  Method: Compute determinant of Sylvester matrix of e1,e2,s.  */
+ex resultant (const ex &e1, const ex &e2, const symbol &s)
+{
+  const int h1 = e1.degree (s);
+  const int l1 = e1.ldegree (s);
+  const int h2 = e2.degree (s);
+  const int l2 = e2.ldegree (s);
+
+  const int msize = h1 + h2;
+  matrix m (msize, msize);
+
+  for (int l = h1; l >= l1; --l)
+  {
+    const ex e = e1.coeff (s, l);
+    for (int k = 0; k < h2; ++k)
+      m (k, k+h1-l) = e;
+  }
+  for (int l = h2; l >= l2; --l)
+  {
+    const ex e = e2.coeff (s, l);
+    for (int k = 0; k < h1; ++k)
+      m (k+h2, k+h2-l) = e;
+  }
+
+  return m.determinant();
+}
+
 
 } // namespace GiNaC
diff -u --recursive GiNaC-1.2.1-orig/ginac/normal.h GiNaC-1.2.1/ginac/normal.h
--- GiNaC-1.2.1-orig/ginac/normal.h	Thu Jan  8 16:06:50 2004
+++ GiNaC-1.2.1/ginac/normal.h	Sun Jun 20 11:57:22 2004
@@ -66,6 +66,9 @@
 // Collect common factors in sums.
 extern ex collect_common_factors(const ex & e);
 
+// Resultant of two polynomials e1,e2 with respect to symbol s
+extern ex resultant (const ex& e1, const ex& e2, const symbol& s);
+
 } // namespace GiNaC
 
 #endif // ndef __GINAC_NORMAL_H__