weyl_diff() [2/3]

Nterm * WeylAlgebra::weyl_diff ( const ring_elem c, const_exponents expf, const int * derivatives, const Nterm * g ) const

protected

Definition at line 339 of file weylalg.cpp.

{
  // This isn't really differentiation, but it is close.
  // It is the inner loop of the multiplication routine for the Weyl algebra.
  // Returns: sum of d*[n,derivative]*c*n*m/(derivative,derivatives) e_i, for
  // each
  // term d*n*e_i of v which satisfies: x part of n is >= derivatives,
  // and where the multiplication and division of monomials is in the
  // commutative
  // monoid.
 
  Nterm head;
  head.next = nullptr;
  Nterm *result = &head;
 
  int i;
  exponents_t exp = newarray_atomic(int, _nderivatives);
  exponents_t deriv_exp = newarray_atomic(int, nvars_);
  exponents_t result_exp = newarray_atomic(int, nvars_);
  for (i = 0; i < nvars_; i++) deriv_exp[i] = 0;
  if (_homogeneous_weyl_algebra)
    {
      int sum = 0;
      for (i = 0; i < _nderivatives; i++)
        {
          sum += 2 * derivatives[i];
          deriv_exp[_derivative[i]] = derivatives[i];
          deriv_exp[_commutative[i]] = derivatives[i];
        }
      deriv_exp[_homog_var] = -sum;
    }
  else
    for (i = 0; i < _nderivatives; i++)
      {
        deriv_exp[_derivative[i]] = derivatives[i];
        deriv_exp[_commutative[i]] = derivatives[i];
      }
 
  for (Nterm& t : ring_elem(g))
    {
      // This first part checks whether the x-part of t->monom is divisible by
      // 'derivatives'.  If so, true is returned, and the resulting monomial is
      // set.
      M_->to_expvector(t.monom, result_exp);
      extractCommutativePart(result_exp, exp);
      if (divides(derivatives, exp))
        {
          ring_elem a = diff_coefficients(c, derivatives, exp);
          if (K_->is_zero(a))
            {
              K_->remove(a);
              continue;
            }
          ring_elem b = K_->mult(a, t.coeff);
          K_->remove(a);
          if (K_->is_zero(b))
            {
              K_->remove(b);
              continue;
            }
          // Now compute the new monomial:
          Nterm *tm = new_term();
          tm->coeff = b;
          for (int i2 = 0; i2 < nvars_; i2++)
            result_exp[i2] += expf[i2] - deriv_exp[i2];
          M_->from_expvector(result_exp, tm->monom);
 
          // Append to the result
          result->next = tm;
          result = tm;
        }
    }
  freemem(exp);
  freemem(result_exp);
  result->next = nullptr;
  return head.next;
}

References _commutative, _derivative, _homog_var, _homogeneous_weyl_algebra, _nderivatives, Nterm::coeff, diff_coefficients(), divides(), extractCommutativePart(), freemem(), PolynomialRing::K_, PolynomialRing::M_, Nterm::monom, PolyRing::new_term(), newarray_atomic, Nterm::next, PolynomialRing::nvars_, and result().

Referenced by mult_by_term().