author Ralf Stephan Thu, 31 Dec 2015 15:36:46 +0000 (16:36 +0100) committer Richard Kreckel Thu, 31 Dec 2015 15:36:46 +0000 (16:36 +0100)
Coeff is too general as long as we have only ints as exponents.
This patch fixes a potentially cubic slowdown.

index cc756ac43e9c70b0cb3bfec01af080aa8a63ed90..8f267e8de4d49c0062d5d1ed910aea1c61f53cff 100644 (file)
@@ -823,11 +823,11 @@ ex pseries::mul_series(const pseries &other) const

// Series multiplication
epvector new_seq;
-       int a_max = degree(var);
-       int b_max = other.degree(var);
-       int a_min = ldegree(var);
-       int b_min = other.ldegree(var);
-       int cdeg_min = a_min + b_min;
+       const int a_max = degree(var);
+       const int b_max = other.degree(var);
+       const int a_min = ldegree(var);
+       const int b_min = other.ldegree(var);
+       const int cdeg_min = a_min + b_min;
int cdeg_max = a_max + b_max;

int higher_order_a = std::numeric_limits<int>::max();
@@ -836,18 +836,30 @@ ex pseries::mul_series(const pseries &other) const
higher_order_a = a_max + b_min;
if (is_order_function(other.coeff(var, b_max)))
higher_order_b = b_max + a_min;
-       int higher_order_c = std::min(higher_order_a, higher_order_b);
+       const int higher_order_c = std::min(higher_order_a, higher_order_b);
if (cdeg_max >= higher_order_c)
cdeg_max = higher_order_c - 1;
-
+
+       std::map<int, ex> rest_map_a, rest_map_b;
+       for (const auto& it : seq)
+               rest_map_a[ex_to<numeric>(it.coeff).to_int()] = it.rest;
+
+       if (other.var.is_equal(var))
+               for (const auto& it : other.seq)
+                       rest_map_b[ex_to<numeric>(it.coeff).to_int()] = it.rest;
+
for (int cdeg=cdeg_min; cdeg<=cdeg_max; ++cdeg) {
ex co = _ex0;
// c(i)=a(0)b(i)+...+a(i)b(0)
for (int i=a_min; cdeg-i>=b_min; ++i) {
-                       ex a_coeff = coeff(var, i);
-                       ex b_coeff = other.coeff(var, cdeg-i);
-                       if (!is_order_function(a_coeff) && !is_order_function(b_coeff))
-                               co += a_coeff * b_coeff;
+                       const auto& ita = rest_map_a.find(i);
+                       if (ita == rest_map_a.end())
+                               continue;
+                       const auto& itb = rest_map_b.find(cdeg-i);
+                       if (itb == rest_map_b.end())
+                               continue;
+                       if (!is_order_function(ita->second) && !is_order_function(itb->second))
+                               co += ita->second * itb->second;
}
if (!co.is_zero())
new_seq.push_back(expair(co, numeric(cdeg)));