text:
###############################################################################
#
#   Poly.jl : Generic polynomials over rings
#
###############################################################################

export PolynomialRing, hash, coeff, isgen, lead,
       var, truncate, mullow, reverse, shift_left, shift_right, divexact,
       pseudorem, pseudodivrem, gcd, degree, content, primpart, evaluate,
       compose, derivative, integral, resultant, discriminant, gcdx, zero, one,
       gen, length, iszero, normalise, isone, isunit, addeq!, mul!, fit!,
       setcoeff!, mulmod, powmod, invmod, lcm, divrem, mod, gcdinv, resx,
       canonical_unit, var, chebyshev_t, chebyshev_u, set_length!,
       mul_classical, mul_ks, subst, mul_karatsuba, trail,
       pow_multinomial, monomial_to_newton!, newton_to_monomial!, ismonomial,
       base_ring, parent_type, elem_type, check_parent, promote_rule,
       needs_parentheses, isnegative, show_minus_one, remove, zero!, add!,
       interpolate, sylvester_matrix

###############################################################################
#
#   Data type and parent object methods
#
###############################################################################

parent_type(::Type{Poly{T}}) where T <: RingElement = PolyRing{T}

elem_type(::Type{PolyRing{T}}) where T <: RingElement = Poly{T}

doc"""
    base_ring(R::Nemo.PolyRing)
> Return the base ring of the given polynomial ring.
"""
base_ring(R::Nemo.PolyRing{T}) where T <: RingElement = R.base_ring::parent_type(T)

doc"""
    base_ring(a::Nemo.PolyElem)
> Return the base ring of the polynomial ring of the given polynomial.
"""
base_ring(a::Nemo.PolyElem) = base_ring(parent(a))

doc"""
    parent(a::Nemo.PolyElem)
> Return the parent of the given polynomial.
"""
parent(a::Nemo.PolyElem) = a.parent

function isdomain_type(::Type{T}) where {S <: RingElement, T <: Nemo.PolyElem{S}}
   return isdomain_type(S)
end

function isexact_type(a::Type{T}) where {S <: RingElement, T <: Nemo.PolyElem{S}}
   return isexact_type(S)
end

doc"""
    var(a::Nemo.PolyRing)
> Return the internal name of the generator of the polynomial ring. Note that
> this is returned as a `Symbol` not a `String`.
"""
var(a::Nemo.PolyRing) = a.S

doc"""
    vars(a::Nemo.PolyRing)
> Return an array of the variable names for the polynomial ring. Note that
> this is returned as an array of `Symbol` not `String`.
"""
vars(a::Nemo.PolyRing) = [a.S]

function check_parent(a::Nemo.PolyElem, b::Nemo.PolyElem)
   parent(a) != parent(b) &&
                error("Incompatible polynomial rings in polynomial operation")
end

###############################################################################
#
#   Basic manipulation
#
###############################################################################

function Base.hash(a::Nemo.PolyElem, h::UInt)
   b = 0x53dd43cd511044d1%UInt
   for i in 0:length(a) - 1
      b = xor(b, xor(hash(coeff(a, i), h), h))
      b = (b << 1) | (b >> (sizeof(Int)*8 - 1))
   end
   return b
end

function setcoeff!(c::Poly{T}, n::Int, a::T) where {T <: RingElement}
   if !iszero(a) || n + 1 <= length(c)
      fit!(c, n + 1)
      c.coeffs[n + 1] = a
      c.length = max(length(c), n + 1)
      # don't normalise
   end
   return c
end

function normalise(a::Poly, n::Int)
   while n > 0 && iszero(a.coeffs[n])
      n -= 1
   end
   return n
end

length(a::Nemo.PolyElem) = a.length

doc"""
    degree(a::Nemo.PolyElem)
> Return the degree of the given polynomial. This is defined to be one less
> than the length, even for constant polynomials.
"""
degree(a::Nemo.PolyElem) = length(a) - 1

doc"""
    modulus{T <: ResElem}(a::Nemo.PolyElem{T})
> Return the modulus of the coefficients of the given polynomial.
"""
modulus(a::Nemo.PolyElem{T}) where {T <: ResElem} = modulus(base_ring(a))

coeff(a::Poly, n::Int) = n >= length(a) ? base_ring(a)(0) : a.coeffs[n + 1]

doc"""
    lead(x::Nemo.PolyElem)
> Return the leading coefficient of the given polynomial. This will be the
> nonzero coefficient of the term with highest degree unless the polynomial
> in the zero polynomial, in which case a zero coefficient is returned.
"""
lead(a::Nemo.PolyElem) = length(a) == 0 ? base_ring(a)(0) : coeff(a, length(a) - 1)

doc"""
    trail(x::Nemo.PolyElem)
> Return the trailing coefficient of the given polynomial. This will be the
> nonzero coefficient of the term with lowest degree unless the polynomial
> in the zero polynomial, in which case a zero coefficient is returned.
"""
function trail(a::Nemo.PolyElem)
   if iszero(a)
      return base_ring(a)(0)
   else
      for i = 1:length(a)
         c = coeff(a, i - 1)
         if !iszero(c)
            return c
         end
      end
      return coeff(a, length(a) - 1)
   end
end

doc"""
    zero(R::Nemo.PolyRing)
> Return the zero polynomial in the given polynomial ring.
"""
zero(R::Nemo.PolyRing) = R(0)

doc"""
    one(R::Nemo.PolyRing)
> Return the constant polynomial $1$ in the given polynomial ring.
"""
one(R::Nemo.PolyRing) = R(1)

doc"""
    gen(R::Nemo.PolyRing)
> Return the generator of the given polynomial ring.
"""
gen(R::Nemo.PolyRing) = R([zero(base_ring(R)), one(base_ring(R))])

doc"""
    iszero(a::Nemo.PolyElem)
> Return `true` if the given polynomial is zero, otherwise return `false`.
"""
iszero(a::Nemo.PolyElem) = length(a) == 0

doc"""
    isone(a::Nemo.PolyElem)
> Return `true` if the given polynomial is the constant polynomial $1$,
> otherwise return `false`.
"""
isone(a::Nemo.PolyElem) = length(a) == 1 && isone(coeff(a, 0))

doc"""
    isgen(a::Nemo.PolyElem)
> Return `true` if the given polynomial is the constant generator of its
> polynomial ring, otherwise return `false`.
"""
function isgen(a::Nemo.PolyElem)
    return length(a) == 2 && iszero(coeff(a, 0)) && isone(coeff(a, 1))
end

doc"""
    isunit(a::Nemo.PolyElem)
> Return `true` if the given polynomial is a unit in its polynomial ring,
> otherwise return `false`.
"""
isunit(a::Nemo.PolyElem) = length(a) == 1 && isunit(coeff(a, 0))

isterm(a::T) where {T <: RingElement} = true

doc"""
    isterm(a::Nemo.PolyElem)
> Return `true` if the given polynomial is has one term. This function is
> recursive, with all scalar types returning true.
"""
function isterm(a::Nemo.PolyElem)
   if !isterm(lead(a))
      return false
   end
   for i = 1:length(a) - 1
      if !iszero(coeff(a, i - 1))
         return false
      end
   end
   return true
end

ismonomial(a::T) where {T <: RingElement} = isone(a)

doc"""
    ismonomial(a::Nemo.PolyElem)
> Return `true` if the given polynomial is a monomial.
"""
function ismonomial(a::Nemo.PolyElem)
   if !ismonomial(lead(a))
      return false
   end
   for i = 1:length(a) - 1
      if !iszero(coeff(a, i - 1))
         return false
      end
   end
   return true
end

function deepcopy_internal(a::Poly{T}, dict::ObjectIdDict) where {T <: RingElement}
   coeffs = Array{T}(length(a))
   for i = 1:length(a)
      coeffs[i] = deepcopy(a.coeffs[i])
   end
   return parent(a)(coeffs)
end

###############################################################################
#
#   Canonicalisation
#
###############################################################################

canonical_unit(x::Nemo.PolyElem) = canonical_unit(lead(x))

###############################################################################
#
#   String I/O
#
###############################################################################

function show(io::IO, x::Nemo.PolyElem)
   len = length(x)
   S = var(parent(x))
   if len == 0
      print(io, base_ring(x)(0))
   else
      for i = 1:len - 1
         c = coeff(x, len - i)
         bracket = needs_parentheses(c)
         if !iszero(c)
            if i != 1 && !isnegative(c)
               print(io, "+")
            end
            if !isone(c) && (c != -1 || show_minus_one(typeof(c)))
               if bracket
                  print(io, "(")
               end
               show(io, c)
               if bracket
                  print(io, ")")
               end
               print(io, "*")
            end
            if c == -1 && !show_minus_one(typeof(c))
               print(io, "-")
            end
            print(io, string(S))
            if len - i != 1
               print(io, "^")
               print(io, len - i)
            end
         end
      end
      c = coeff(x, 0)
      bracket = needs_parentheses(c)
      if !iszero(c)
         if len != 1 && !isnegative(c)
            print(io, "+")
         end
         if bracket
            print(io, "(")
         end
         show(io, c)
         if bracket
            print(io, ")")
         end
      end
   end
end

function show(io::IO, p::Nemo.PolyRing)
   print(io, "Univariate Polynomial Ring in ")
   print(io, string(var(p)))
   print(io, " over ")
   show(io, base_ring(p))
end

needs_parentheses(x::Nemo.PolyElem) = length(x) > 1

isnegative(x::Nemo.PolyElem) = length(x) <= 1 && isnegative(coeff(x, 0))

show_minus_one(::Type{Poly{T}}) where {T <: RingElement} = show_minus_one(T)

###############################################################################
#
#   Unary operations
#
###############################################################################

doc"""
    -(a::Nemo.PolyElem)
> Return $-a$.
"""
function -(a::Nemo.PolyElem)
   len = length(a)
   z = parent(a)()
   fit!(z, len)
   for i = 1:len
      z = setcoeff!(z, i - 1, -coeff(a, i - 1))
   end
   set_length!(z, len)
   return z
end

###############################################################################
#
#   Binary operations
#
###############################################################################

doc"""
    +{T <: RingElement}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T})
> Return $a + b$.
"""
function +(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: RingElement}
   check_parent(a, b)
   lena = length(a)
   lenb = length(b)
   lenz = max(lena, lenb)
   z = parent(a)()
   fit!(z, lenz)
   i = 1
   while i <= min(lena, lenb)
      z = setcoeff!(z, i - 1, coeff(a, i - 1) + coeff(b, i - 1))
      i += 1
   end
   while i <= lena
      z = setcoeff!(z, i - 1, deepcopy(coeff(a, i - 1)))
      i += 1
   end
   while i <= lenb
      z = setcoeff!(z, i - 1, deepcopy(coeff(b, i - 1)))
      i += 1
   end
   set_length!(z, normalise(z, i - 1))
   return z
end

doc"""
    -{T <: RingElement}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T})
> Return $a - b$.
"""
function -(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: RingElement}
   check_parent(a, b)
   lena = length(a)
   lenb = length(b)
   lenz = max(lena, lenb)
   z = parent(a)()
   fit!(z, lenz)
   i = 1
   while i <= min(lena, lenb)
      z = setcoeff!(z, i - 1, coeff(a, i - 1) - coeff(b, i - 1))
      i += 1
   end
   while i <= lena
      z = setcoeff!(z, i - 1, deepcopy(coeff(a, i - 1)))
      i += 1
   end
   while i <= lenb
      z = setcoeff!(z, i - 1, -coeff(b, i - 1))
      i += 1
   end
   set_length!(z, normalise(z, i - 1))
   return z
end

function mul_karatsuba(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: RingElement}
   # we assume len(a) != 0 != lenb and parent(a) == parent(b)
   lena = length(a)
   lenb = length(b)
   m = div(max(lena, lenb) + 1, 2)
   if m < lena
      a1 = shift_right(a, m)
      a0 = truncate(a, m)
   else
      return a*truncate(b, m) + shift_left(a*shift_right(b, m), m)
   end
   if a !== b
      if m < lenb
         b1 = shift_right(b, m)
         b0 = truncate(b, m)
      else
         return b*truncate(a, m) + shift_left(b*shift_right(a, m), m)
      end
   else
      b1 = a1
      b0 = a0
   end
   z0 = a0*b0
   z2 = a1*b1
   if a !== b
      z1 = (a1 + a0)*(b1 + b0) - z2 - z0
   else
      s = a1 + a0
      z1 = s*s - z2 - z0
   end
   r = parent(a)()
   fit!(r, lena + lenb - 1)
   for i = 1:length(z0)
      r = setcoeff!(r, i - 1, coeff(z0, i - 1))
   end
   for i = length(z0) + 1:2m
      r = setcoeff!(r, i - 1, base_ring(a)())
   end
   for i = 1:length(z2)
      r = setcoeff!(r, 2m + i - 1, coeff(z2, i - 1))
   end
   for i = 1:length(z1)
      u = coeff(r, i + m - 1)
      u = addeq!(u, coeff(z1, i - 1))
      setcoeff!(r, i + m - 1, u)
   end
   return r
end

function mul_ks(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: Nemo.PolyElem}
   lena = length(a)
   lenb = length(b)
   if lena == 0 || lenb == 0
      return parent(a)()
   end
   maxa = 0
   nza = 0
   for i = 1:lena
      lenc = length(coeff(a, i - 1))
      maxa = max(lenc, maxa)
      nza += (lenc == 0 ? 0 : 1)
   end
   if a !== b
      maxb = 0
      nzb = 0
      for i = 1:lenb
         lenc = length(coeff(b, i - 1))
         maxb = max(lenc, maxb)
         nzb += (lenc == 0 ? 0 : 1)
      end
   else
      maxb = maxa
      nzb = nza
   end
   if nza*nzb < 4*max(lena, lenb)
      return mul_classical(a, b)
   end
   m = maxa + maxb - 1
   z = base_ring(base_ring(a))()
   A1 = Array{elem_type(base_ring(base_ring(a)))}(m*lena)
   for i = 1:lena
      c = coeff(a, i - 1)
      for j = 1:length(c)
         A1[(i - 1)*m + j] = coeff(c, j - 1)
      end
      for j = length(c) + 1:m
         A1[(i - 1)*m + j] = z
      end
   end
   ksa = base_ring(a)(A1)
   if a !== b
      A2 = Array{elem_type(base_ring(base_ring(a)))}(m*lenb)
      for i = 1:lenb
         c = coeff(b, i - 1)
         for j = 1:length(c)
            A2[(i - 1)*m + j] = coeff(c, j - 1)
         end
         for j = length(c) + 1:m
            A2[(i - 1)*m + j] = z
         end
      end
      ksb = base_ring(b)(A2)
   else
      ksb = ksa
   end
   p = ksa*ksb
   r = parent(a)()
   lenr = lena + lenb - 1
   fit!(r, lenr)
   for i = 1:lenr
      u = coeff(r, i - 1)
      fit!(u, m)
      for j = 1:m
         u = setcoeff!(u, j - 1, coeff(p, (i - 1)*m + j - 1))
      end
      setcoeff!(r, i - 1, u)
   end
   set_length!(r, normalise(r, lenr))
   return r
end

function mul_classical(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: RingElement}
   lena = length(a)
   lenb = length(b)
   if lena == 0 || lenb == 0
      return parent(a)()
   end
   t = base_ring(a)()
   lenz = lena + lenb - 1
   d = Array{T}(lenz)
   for i = 1:lena
      d[i] = coeff(a, i - 1)*coeff(b, 0)
   end
   for i = 2:lenb
      d[lena + i - 1] = a.coeffs[lena]*coeff(b, i - 1)
   end
   for i = 1:lena - 1
      for j = 2:lenb
         t = mul!(t, coeff(a, i - 1), b.coeffs[j])
         d[i + j - 1] = addeq!(d[i + j - 1], t)
      end
   end
   z = parent(a)(d)
   set_length!(z, normalise(z, lenz))
   return z
end

doc"""
    *{T <: RingElement}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T})
> Return $a\times b$.
"""
function *(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: RingElement}
   check_parent(a, b)
   return mul_classical(a, b)
end

###############################################################################
#
#   Ad hoc binary operators
#
###############################################################################

doc"""
    *{T <: RingElem}(a::T, b::Nemo.PolyElem{T})
> Return $a\times b$.
"""
function *(a::T, b::Nemo.PolyElem{T}) where {T <: RingElem}
   len = length(b)
   z = parent(b)()
   fit!(z, len)
   for i = 1:len
      z = setcoeff!(z, i - 1, a*coeff(b, i - 1))
   end
   set_length!(z, normalise(z, len))
   return z
end

doc"""
    *(a::Union{Integer, Rational, AbstractFloat}, b::Nemo.PolyElem)
> Return $a\times b$.
"""
function *(a::Union{Integer, Rational, AbstractFloat}, b::Nemo.PolyElem)
   len = length(b)
   z = parent(b)()
   fit!(z, len)
   for i = 1:len
      z = setcoeff!(z, i - 1, a*coeff(b, i - 1))
   end
   set_length!(z, normalise(z, len))
   return z
end

doc"""
    *{T <: RingElem}(a::Nemo.PolyElem{T}, b::T)
> Return $a\times b$.
"""
*(a::Nemo.PolyElem{T}, b::T) where {T <: RingElem} = b*a

doc"""
    *(a::Nemo.PolyElem, b::Union{Integer, Rational, AbstractFloat})
> Return $a\times b$.
"""
*(a::Nemo.PolyElem, b::Union{Integer, Rational, AbstractFloat}) = b*a

###############################################################################
#
#   Powering
#
###############################################################################

function pow_multinomial(a::Nemo.PolyElem{T}, e::Int) where {T <: RingElement}
   e < 0 && throw(DomainError())
   lena = length(a)
   lenz = (lena - 1) * e + 1
   res = Array{T}(lenz)
   for k = 1:lenz
      res[k] = base_ring(a)()
   end
   d = base_ring(a)()
   first = coeff(a, 0)
   res[1] = first ^ e
   for k = 1 : lenz - 1
      u = -k
      for i = 1 : min(k, lena - 1)
         t = coeff(a, i) * res[(k - i) + 1]
         u += e + 1
         res[k + 1] = addeq!(res[k + 1], t * u)
      end
      d = addeq!(d, first)
      res[k + 1] = divexact(res[k + 1], d)
   end
   z = parent(a)(res)
   set_length!(z, normalise(z, lenz))
   return z
end

doc"""
    ^{T <: RingElement}(a::Nemo.PolyElem{T}, b::Int)
> Return $a^b$. We require $b \geq 0$.
"""
function ^(a::Nemo.PolyElem{T}, b::Int) where {T <: RingElement}
   b < 0 && throw(DomainError())
   # special case powers of x for constructing polynomials efficiently
   R = parent(a)
   if isgen(a)
      z = R()
      fit!(z, b + 1)
      z = setcoeff!(z, b, coeff(a, 1))
      for i = 1:b
         z = setcoeff!(z, i - 1, coeff(a, 0))
      end
      set_length!(z, b + 1)
      return z
   elseif length(a) == 0
      return zero(R)
   elseif length(a) == 1
      return R(coeff(a, 0)^b)
   elseif b == 0
      return one(R)
   elseif b == 1
      return deepcopy(a)
   else
      if T <: FieldElement && characteristic(base_ring(R)) == 0
         zn = 0
         while iszero(coeff(a, zn))
            zn += 1
         end
         if length(a) - zn < 8 && b > 4
             f = shift_right(a, zn)
             return shift_left(pow_multinomial(f, b), zn*b)
         end
      end
      bit = ~((~UInt(0)) >> 1)
      while (UInt(bit) & b) == 0
         bit >>= 1
      end
      z = a
      bit >>= 1
      while bit != 0
         z = z*z
         if (UInt(bit) & b) != 0
            z *= a
         end
         bit >>= 1
      end
      return z
   end
end

###############################################################################
#
#   Comparisons
#
###############################################################################

doc"""
    =={T <: RingElement}(x::Nemo.PolyElem{T}, y::Nemo.PolyElem{T})
> Return `true` if $x == y$ arithmetically, otherwise return `false`. Recall
> that power series to different precisions may still be arithmetically
> equal to the minimum of the two precisions.
"""
function ==(x::Nemo.PolyElem{T}, y::Nemo.PolyElem{T}) where {T <: RingElement}
   check_parent(x, y)
   if length(x) != length(y)
      return false
   else
      for i = 1:length(x)
         if coeff(x, i - 1) != coeff(y, i - 1)
            return false
         end
      end
   end
   return true
end

doc"""
    isequal{T <: RingElement}(x::Nemo.PolyElem{T}, y::Nemo.PolyElem{T})
> Return `true` if $x == y$ exactly, otherwise return `false`. This function is
> useful in cases where the coefficients of the polynomial are inexact, e.g.
> power series. Only if the power series are precisely the same, to the same
> precision, are they declared equal by this function.
"""
function isequal(x::Nemo.PolyElem{T}, y::Nemo.PolyElem{T}) where {T <: RingElement}
   if parent(x) != parent(y)
      return false
   end
   if length(x) != length(y)
      return false
   end
   for i = 1:length(x)
      if !isequal(coeff(x, i - 1), coeff(y, i - 1))
         return false
      end
   end
   return true
end

###############################################################################
#
#   Ad hoc comparisons
#
###############################################################################

doc"""
    =={T <: RingElem}(x::Nemo.PolyElem{T}, y::T)
> Return `true` if $x == y$.
"""
==(x::Nemo.PolyElem{T}, y::T) where T <: RingElem = ((length(x) == 0 && y == 0)
                        || (length(x) == 1 && coeff(x, 0) == y))

doc"""
    ==(x::Nemo.PolyElem, y::Union{Integer, Rational, AbstractFloat})
> Return `true` if $x == y$ arithmetically, otherwise return `false`.
"""
==(x::Nemo.PolyElem, y::Union{Integer, Rational, AbstractFloat}) = ((length(x) == 0 && base_ring(x)(y) == 0)
                        || (length(x) == 1 && coeff(x, 0) == y))

doc"""
    =={T <: RingElem}(x::T, y::Nemo.PolyElem{T})
> Return `true` if $x = y$.
"""
==(x::T, y::Nemo.PolyElem{T}) where T <: RingElem = y == x

doc"""
    ==(x::Union{Integer, Rational, AbstractFloat}, y::Nemo.PolyElem)
> Return `true` if $x == y$ arithmetically, otherwise return `false`.
"""
==(x::Union{Integer, Rational, AbstractFloat}, y::Nemo.PolyElem) = y == x

###############################################################################
#
#   Approximation
#
###############################################################################

function Base.isapprox(f::Nemo.PolyElem, g::Nemo.PolyElem; atol::Real=sqrt(eps()))
   check_parent(f, g)
   nmin = min(length(f), length(g))
   i = 1
   while i <= nmin
      if !isapprox(coeff(f, i - 1), coeff(g, i - 1); atol=atol)
         return false
      end
      i += 1
   end
   while i <= length(f)
      if !isapprox(coeff(f, i - 1), 0; atol=atol)
         return false
      end
      i += 1
   end
   while i <= length(g)
      if !isapprox(coeff(g, i - 1), 0; atol=atol)
         return false
      end
      i += 1
   end
   return true
end

function Base.isapprox(f::Nemo.PolyElem{T}, g::T; atol::Real=sqrt(eps())) where T
   return isapprox(f, parent(f)(g); atol=atol)
end

function Base.isapprox(f::T, g::Nemo.PolyElem{T}; atol::Real=sqrt(eps())) where T
   return isapprox(parent(g)(f), g; atol=atol)
end

###############################################################################
#
#   Truncation
#
###############################################################################

doc"""
    truncate(a::Nemo.PolyElem, n::Int)
> Return $a$ truncated to $n$ terms.
"""
function truncate(a::Nemo.PolyElem, n::Int)
   n < 0 && throw(DomainError())
   lena = length(a)
   if lena <= n
      return a
   end
   lenz = min(lena, n)
   z = parent(a)()
   fit!(z, lenz)
   for i = 1:lenz
      z = setcoeff!(z, i - 1, coeff(a, i - 1))
   end
   set_length!(z, normalise(z, lenz))
   return z
end

doc"""
    mullow{T <: RingElement}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}, n::Int)
> Return $a\times b$ truncated to $n$ terms.
"""
function mullow(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}, n::Int) where {T <: RingElement}
   check_parent(a, b)
   lena = length(a)
   lenb = length(b)
   if lena == 0 || lenb == 0
      return zero(parent(a))
   end
   if n < 0
      n = 0
   end
   t = base_ring(a)()
   lenz = min(lena + lenb - 1, n)
   d = Array{T}(lenz)
   for i = 1:min(lena, lenz)
      d[i] = coeff(a, i - 1)*coeff(b, 0)
   end
   if lenz > lena
      for j = 2:min(lenb, lenz - lena + 1)
          d[lena + j - 1] = coeff(a, lena - 1)*coeff(b, j - 1)
      end
   end
   for i = 1:lena - 1
      if lenz > i
         for j = 2:min(lenb, lenz - i + 1)
            t = mul!(t, coeff(a, i - 1), coeff(b, j - 1))
            d[i + j - 1] = addeq!(d[i + j - 1], t)
         end
      end
   end
   z = parent(a)(d)
   set_length!(z, normalise(z, lenz))
   return z
end

###############################################################################
#
#   Reversal
#
###############################################################################

doc"""
    reverse(x::Nemo.PolyElem, len::Int)
> Return the reverse of the polynomial $x$, thought of as a polynomial of
> the given length (the polynomial will be notionally truncated or padded with
> zeroes before the leading term if necessary to match the specified length).
> The resulting polynomial is normalised. If `len` is negative we throw a
> `DomainError()`.
"""
function reverse(x::Nemo.PolyElem, len::Int)
   len < 0 && throw(DomainError())
   r = parent(x)()
   fit!(r, len)
   for i = 1:len
      z = setcoeff!(r, i - 1, coeff(x, len - i))
   end
   set_length!(r, normalise(r, len))
   return r
end

doc"""
    reverse(x::Nemo.PolyElem)
> Return the reverse of the polynomial $x$, i.e. the leading coefficient
> of $x$ becomes the constant coefficient of the result, etc. The resulting
> polynomial is normalised.
"""
function reverse(x::Nemo.PolyElem)
   reverse(x, length(x))
end

###############################################################################
#
#   Shifting
#
###############################################################################

doc"""
    shift_left(x::Nemo.PolyElem, n::Int)
> Return the polynomial $f$ shifted left by $n$ terms, i.e. multiplied by
> $x^n$.
"""
function shift_left(f::Nemo.PolyElem, n::Int)
   n < 0 && throw(DomainError())
   if n == 0
      return f
   end
   flen = length(f)
   r = parent(f)()
   fit!(r, flen + n)
   for i = 1:n
      r = setcoeff!(r, i - 1, zero(base_ring(f)))
   end
   for i = 1:flen
      r = setcoeff!(r, i + n - 1, coeff(f, i - 1))
   end
   return r
end

doc"""
    shift_right(f::Nemo.PolyElem, n::Int)
> Return the polynomial $f$ shifted right by $n$ terms, i.e. divided by
> $x^n$.
"""
function shift_right(f::Nemo.PolyElem, n::Int)
   n < 0 && throw(DomainError())
   flen = length(f)
   if n >= flen
      return zero(parent(f))
   end
   if n == 0
      return f
   end
   r = parent(f)()
   fit!(r, flen - n)
   for i = 1:flen - n
      r = setcoeff!(r, i - 1, coeff(f, i + n - 1))
   end
   return r
end

###############################################################################
#
#   Modular arithmetic
#
###############################################################################

doc"""
    mulmod{T <: Union{Nemo.ResElem, FieldElement}}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}, d::Nemo.PolyElem{T})
> Return $a\times b \pmod{d}$.
"""
function mulmod(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}, d::Nemo.PolyElem{T}) where {T <: Union{Nemo.ResElem, FieldElement}}
   check_parent(a, b)
   check_parent(a, d)
   return mod(a*b, d)
end

doc"""
    powmod{T <: Union{Nemo.ResElem, FieldElement}}(a::Nemo.PolyElem{T}, b::Int, d::Nemo.PolyElem{T})
> Return $a^b \pmod{d}$. There are no restrictions on $b$.
"""
function powmod(a::Nemo.PolyElem{T}, b::Int, d::Nemo.PolyElem{T}) where {T <: Union{Nemo.ResElem, FieldElement}}
   check_parent(a, d)
   if length(a) == 0
      z = zero(parent(a))
   elseif length(a) == 1
      z = parent(a)(coeff(a, 0)^b)
   elseif b == 0
      z = one(parent(a))
   else
      if b < 0
         a = invmod(a, d)
         b = -b
      end
      bit = ~((~UInt(0)) >> 1)
      while (UInt(bit) & b) == 0
         bit >>= 1
      end
      z = a
      bit >>= 1
      while bit !=0
         z = mulmod(z, z, d)
         if (UInt(bit) & b) != 0
            z = mulmod(z, a, d)
         end
         bit >>= 1
      end
   end
   if length(z) >= length(d)
      z = mod(z, d)
   end
   return z
end

doc"""
    invmod{T <: Union{Nemo.ResElem, FieldElement}}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T})
> Return $a^{-1} \pmod{d}$.
"""
function invmod(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: Union{Nemo.ResElem, FieldElement}}
   check_parent(a, b)
   g, z = gcdinv(a, b)
   if g != 1
      error("Impossible inverse in invmod")
   end
   return z
end

###############################################################################
#
#   Exact division
#
###############################################################################

doc"""
    divexact{T <: RingElement}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T})
> Return $a/b$ where the quotient is expected to be exact.
"""
function divexact(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T}) where {T <: RingElement}
   check_parent(f, g)
   iszero(g) && throw(DivideError())
   if iszero(f)
      return zero(parent(f))
   end
   lenq = length(f) - length(g) + 1
   d = Array{T}(lenq)
   for i = 1:lenq
      d[i] = zero(base_ring(f))
   end
   x = gen(parent(f))
   leng = length(g)
   while length(f) >= leng
      lenf = length(f)
      q1 = d[lenf - leng + 1] = divexact(coeff(f, lenf - 1), coeff(g, leng - 1))
      f = f - shift_left(q1*g, lenf - leng)
      if length(f) == lenf # inexact case
         set_length!(f, normalise(f, lenf - 1))
      end
   end
   q = parent(f)(d)
   set_length!(q, lenq)
   return q
end

###############################################################################
#
#   Ad hoc exact division
#
###############################################################################

doc"""
    divexact{T <: RingElem}(a::Nemo.PolyElem{T}, b::T)
> Return $a/b$ where the quotient is expected to be exact.
"""
function divexact(a::Nemo.PolyElem{T}, b::T) where {T <: RingElem}
   iszero(b) && throw(DivideError())
   z = parent(a)()
   fit!(z, length(a))
   for i = 1:length(a)
      z = setcoeff!(z, i - 1, divexact(coeff(a, i - 1), b))
   end
   set_length!(z, length(a))
   return z
end

doc"""
    divexact(a::Nemo.PolyElem, b::Union{Integer, Rational, AbstractFloat})
> Return $a/b$ where the quotient is expected to be exact.
"""
function divexact(a::Nemo.PolyElem, b::Union{Integer, Rational, AbstractFloat})
   iszero(b) && throw(DivideError())
   z = parent(a)()
   fit!(z, length(a))
   for i = 1:length(a)
      z = setcoeff!(z, i - 1, divexact(coeff(a, i - 1), b))
   end
   set_length!(z, length(a))
   return z
end

###############################################################################
#
#   Euclidean division
#
###############################################################################

doc"""
    mod{T <: Union{Nemo.ResElem, FieldElement}}(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T})
> Return $f \pmod{g}$.
"""
function mod(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T}) where {T <: Union{Nemo.ResElem, FieldElement}}
   check_parent(f, g)
   if length(g) == 0
      throw(DivideError())
   end
   if length(f) >= length(g)
      f = deepcopy(f)
      b = lead(g)
      g = inv(b)*g
      c = base_ring(f)()
      while length(f) >= length(g)
         l = -lead(f)
         for i = 1:length(g) - 1
            c = mul!(c, coeff(g, i - 1), l)
            u = coeff(f, i + length(f) - length(g) - 1)
            u = addeq!(u, c)
            f = setcoeff!(f, i + length(f) - length(g) - 1, u)
         end
         set_length!(f, normalise(f, length(f) - 1))
      end
   end
   return f
end

doc"""
    divrem{T <: Union{Nemo.ResElem, FieldElement}}(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T})
> Return a tuple $(q, r)$ such that $f = qg + r$ where $q$ is the euclidean
> quotient of $f$ by $g$.
"""
function divrem(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T}) where {T <: Union{Nemo.ResElem, FieldElement}}
   check_parent(f, g)
   if length(g) == 0
      throw(DivideError())
   end
   if length(f) < length(g)
      return zero(parent(f)), f
   end
   f = deepcopy(f)
   binv = inv(lead(g))
   g = divexact(g, lead(g))
   qlen = length(f) - length(g) + 1
   q = parent(f)()
   fit!(q, qlen)
   c = base_ring(f)()
   while length(f) >= length(g)
      q1 = lead(f)
      l = -q1
      q = setcoeff!(q, length(f) - length(g), q1*binv)
      for i = 1:length(g) - 1
         c = mul!(c, coeff(g, i - 1), l)
         u = coeff(f, i + length(f) - length(g) - 1)
         u = addeq!(u, c)
         f = setcoeff!(f, i + length(f) - length(g) - 1, u)
      end
      set_length!(f, normalise(f, length(f) - 1))
   end
   return q, f
end

doc"""
    div{T <: Union{Nemo.ResElem, FieldElement}}(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T})
> Return a tuple $q$ such that $f = qg + r$ where $q$ is the euclidean
> quotient of $f$ by $g$.
"""
function div(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T}) where {T <: Union{Nemo.ResElem, FieldElement}}
   q, r = divrem(f, g)
   return q
end

###############################################################################
#
#   Pseudodivision
#
###############################################################################

doc"""
    pseudorem{T <: RingElement}(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T})
> Return the pseudoremainder of $a$ divided by $b$. If $b = 0$ we throw a
> `DivideError()`.
"""
function pseudorem(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T}) where {T <: RingElement}
   check_parent(f, g)
   g == 0 && throw(DivideError())
   if length(f) < length(g)
      return f
   end
   k = length(f) - length(g) + 1
   b = coeff(g, length(g) - 1)
   x = gen(parent(f))
   while length(f) >= length(g)
      f = f*b - shift_left(coeff(f, length(f) - 1)*g, length(f) - length(g))
      k -= 1
   end
   return f*b^k
end

doc"""
    pseudodivrem{T <: RingElement}(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T})
> Return a tuple $(q, r)$ consisting of the pseudoquotient and pseudoremainder
> of $a$ divided by $b$. If $b = 0$ we throw a `DivideError()`.
"""
function pseudodivrem(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T}) where {T <: RingElement}
   check_parent(f, g)
   g == 0 && throw(DivideError())
   if length(f) < length(g)
      return zero(parent(f)), f
   end
   lenq = length(f) - length(g) + 1
   k = lenq
   q = parent(f)()
   fit!(q, lenq)
   b = coeff(g, length(g) - 1)
   x = gen(parent(f))
   while length(f) >= length(g)
      for i = length(f) - length(g) + 2:lenq
         q = setcoeff!(q, i - 1, coeff(q, i - 1) * b)
      end
      q = setcoeff!(q, length(f) - length(g), coeff(f, length(f) - 1))
      f = f*b - shift_left(coeff(f, length(f) - 1)*g, length(f) - length(g))
      k -= 1
   end
   while lenq > 0 && coeff(q, lenq - 1) == 0
      lenq -= 1
   end
   set_length!(q, lenq)
   s = b^k
   return q*s, f*s
end

################################################################################
#
#   Remove and valuation
#
################################################################################

#CF TODO: use squaring for fast large valuation

doc"""
    remove{T <: RingElement}(z::Nemo.PolyElem{T}, p::Nemo.PolyElem{T})
> Computes the valuation of $z$ at $p$, that is, the largest $k$ such that
> $p^k$ divides $z$. Additionally, $z/p^k$ is returned as well.
>
> See also `valuation`, which only returns the valuation.
"""
function remove(z::Nemo.PolyElem{T}, p::Nemo.PolyElem{T}) where T <: RingElement
  check_parent(z, p)
  !isexact_type(T) && error("remove requires an exact ring")
  z == 0 && error("Not yet implemented")
  flag, q = divides(z, p)
  if !flag
    return 0, z
  end
  v = 0
  qn = q
  while flag
    q = qn
    flag, qn = divides(q, p)
    v += 1
  end
  return v, q
end

doc"""
    remove{T <: Union{Nemo.ResElem, FieldElement}}(z::Nemo.PolyElem{T}, p::Nemo.PolyElem{T})
> Computes the valuation of $z$ at $p$, that is, the largest $k$ such that
> $p^k$ divides $z$. Additionally, $z/p^k$ is returned as well.
>
> See also `valuation`, which only returns the valuation.
"""
function remove(z::Nemo.PolyElem{T}, p::Nemo.PolyElem{T}) where T <: Union{Nemo.ResElem, FieldElement}
  check_parent(z, p)
  !isexact_type(T) && error("remove requires an exact ring")
  z == 0 && error("Not yet implemented")
  q, r = divrem(z, p)
  if !iszero(r)
    return 0, z
  end
  v = 0
  qn = q
  while iszero(r)
    q = qn
    qn, r = divrem(q, p)
    v += 1
  end
  return v, q
end

doc"""
    valuation{T <: RingElement}(z::Nemo.PolyElem{T}, p::Nemo.PolyElem{T})
> Computes the valuation of $z$ at $p$, that is, the largest $k$ such that
> $p^k$ divides $z$.
>
> See also `remove`, which also returns $z/p^k$.
"""
function valuation(z::Nemo.PolyElem{T}, p::Nemo.PolyElem{T}) where {T <: RingElement}
  v, _ = remove(z, p)
  return v
end

doc"""
    divides{T <: RingElement}(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T})
> Returns a pair consisting of a flag which is set to `true` if $f$ divides
> $g$ and `false` otherwise, and a polynomial $h$ such that $f = gh$ if
> such a polynomial exists. If not, the value of $h$ is undetermined.
"""
function divides(f::Nemo.PolyElem{T}, g::Nemo.PolyElem{T}) where {T <: RingElement}
   check_parent(f, g)
   !isexact_type(T) && error("divides requires an exact ring")
   if length(g) == 0
      throw(DivideError())
   end
   if length(f) == 0
      return true, parent(f)()
   end
   if length(f) < length(g)
      return false, parent(f)()
   end
   f = deepcopy(f)
   g_lead = lead(g)
   qlen = length(f) - length(g) + 1
   q = parent(f)()
   fit!(q, qlen)
   c = base_ring(f)()
   while length(f) >= length(g)
      q1 = lead(f)
      flag, d = divides(q1, g_lead)
      if !flag
         return false, parent(f)()
      end
      q = setcoeff!(q, length(f) - length(g), d)
      d = -d
      for i = 1:length(g)
         c = mul!(c, coeff(g, i - 1), d)
         u = coeff(f, i + length(f) - length(g) - 1)
         u = addeq!(u, c)
         f = setcoeff!(f, i + length(f) - length(g) - 1, u)
      end
      set_length!(f, normalise(f, length(f)))
   end
   return f == 0, q
end

doc"""
    divides{T <: RingElement}(f::Nemo.PolyElem{T}, g::T)
> Returns a pair consisting of a flag which is set to `true` if $g$ divides
> $f$ and `false` otherwise, and a polynomial $h$ such that $f = gh$ if
> such a polynomial exists. If not, the value of $h$ is undetermined.
"""
function divides(z::Nemo.PolyElem{T}, x::T) where {T <: RingElement}
   parent(x) != base_ring(z) && error("Wrong parents in divides")
   q = parent(z)()
   fit!(q, length(z))
   flag = true
   for i = 1:length(z)
      flag, c = divides(coeff(z, i - 1), x)
      if !flag
         break
      end
      q = setcoeff!(q, i - 1, c)
   end
   set_length!(q, flag ? length(z) : 0)
   return flag, q
end

###############################################################################
#
#   Content, primitive part, GCD and LCM
#
###############################################################################

function term_gcd(a::T, b::T) where {T <: RingElement}
   return gcd(a, b)
end

function term_content(a::T) where {T <: RingElement}
   return a
end

function term_gcd(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: RingElement}
   d = min(degree(a), degree(b))
   x = gen(parent(a))
   return term_gcd(coeff(a, degree(a)), coeff(b, degree(b)))*x^d
end

function term_content(a::Nemo.PolyElem{T}) where {T <: RingElement}
   for i = 1:length(a)
      c = coeff(a, i - 1)
      if !iszero(c)
         g = term_content(c)
         for j = i + 1:length(a)
            c = coeff(a, j - 1)
            if !iszero(c)
               g = term_gcd(g, term_content(c))
            end
         end
         x = gen(parent(a))
         return g*x^(i - 1)
      end
   end
   return parent(a)()
end

doc"""
    gcd{T <: RingElement}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T})
> Return a greatest common divisor of $a$ and $b$ if it exists.
"""
function gcd(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}, ignore_content::Bool = false) where {T <: RingElement}
   check_parent(a, b)
   if length(b) > length(a)
      (a, b) = (b, a)
   end
   if iszero(b)
      if a == 0
         return a
      else
         return divexact(a, canonical_unit(lead(a)))
      end
   end
   if isone(b)
      return b
   end
   if !ignore_content
      c1 = content(a)
      c2 = content(b)
      a = divexact(a, c1)
      b = divexact(b, c2)
      c = gcd(c1, c2)
   end
   lead_monomial = isterm(lead(a)) || isterm(lead(b))
   trail_monomial = isterm(trail(a)) || isterm(trail(b))
   lead_a = lead(a)
   lead_b = lead(b)
   g = one(parent(a))
   h = one(parent(a))
   while true
      d = length(a) - length(b)
      r = pseudorem(a, b)
      if r == 0
         break
      end
      if length(r) == 1
         b = one(parent(a))
         break
      end
      (a, b) = (b, divexact(r, g*h^d))
      g = lead(a)
      if d > 1
         h = divexact(g^d, h^(d - 1))
      else
         h = h^(1 - d)*g^d
      end
   end
   if !ignore_content
      if !isterm(lead(b)) && !isterm(trail(b))
         if lead_monomial # lead term monomial, so content contains rest
            d = divexact(lead(b), term_content(lead(b)))
            b = divexact(b, d)
         elseif trail_monomial # trail term is monomial, so ditto
            d = divexact(trail(b), term_content(trail(b)))
            b = divexact(b, d)
         else
            glead = gcd(lead_a, lead_b)
            if isterm(glead)
               d = divexact(lead(b), term_content(lead(b)))
               b = divexact(b, d)
            else # last ditched attempt to find easy content
               h = gcd(lead(b), glead)
               h = divexact(h, term_content(h))
               flag, q = divides(b, h)
               if flag
                  b = q
               end
            end
         end
      end
      b = divexact(b, content(b))
   end
   b = c*b
   return divexact(b, canonical_unit(lead(b)))
end

function gcd(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: Union{Nemo.ResElem, FieldElement}}
   check_parent(a, b)
   if length(a) > length(b)
      (a, b) = (b, a)
   end
   if iszero(b)
      if iszero(a)
         return(a)
      else
         d = lead(a)
         return divexact(a, d)
      end
   end
   g = gcd(content(a), content(b))
   a = divexact(a, g)
   b = divexact(b, g)
   while !iszero(a)
      (a, b) = (mod(b, a), a)
   end
   b = g*b
   d = lead(b)
   return divexact(b, d)
end

doc"""
    lcm{T <: RingElement}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T})
> Return a least common multiple of $a$ and $b$ if it exists.
"""
function lcm(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: RingElement}
   check_parent(a, b)
   return a*divexact(b, gcd(a, b))
end

doc"""
    content(a::Nemo.PolyElem)
> Return the content of $a$, i.e. the greatest common divisor of its
> coefficients.
"""
function content(a::Nemo.PolyElem)
   z = coeff(a, 0)
   for i = 2:length(a)
      z = gcd(z, coeff(a, i - 1))
   end
   if z == 0
      return z
   else
      return divexact(z, canonical_unit(z))
   end
end

doc"""
    primpart(a::Nemo.PolyElem)
> Return the primitive part of $a$, i.e. the polynomial divided by its content.
"""
function primpart(a::Nemo.PolyElem)
   d = content(a)
   if d == 0
      return 0
   else
      return divexact(a, d)
   end
end

###############################################################################
#
#   Evaluation/composition
#
###############################################################################

doc"""
    evaluate{T <: RingElement}(a::Nemo.PolyElem{T}, b::T)
> Evaluate the polynomial $a$ at the value $b$ and return the result.
"""
function evaluate(a::Nemo.PolyElem, b::T) where {T <: RingElement}
   i = length(a)
   R = base_ring(a)
   if i == 0
       return zero(R)
   end
   if i > 25
      return subst(a, b)
   end
   z = R(coeff(a, i - 1))
   while i > 1
      i -= 1
      z = z*b + R(coeff(a, i - 1))
   end
   return z
end

doc"""
    compose(a::Nemo.PolyElem, b::Nemo.PolyElem)
> Compose the polynomial $a$ with the polynomial $b$ and return the result,
> i.e. return $a\circ b$.
"""
function compose(a::Nemo.PolyElem, b::Nemo.PolyElem)
   i = length(a)
   R = parent(a)
   if i == 0
       return zero(R)
   end
   if i*length(b) > 25
      return subst(a, b)
   end
   z = R(coeff(a, i - 1))
   while i > 1
      i -= 1
      z = z*b + coeff(a, i - 1)
   end
   return z
end

###############################################################################
#
#   Derivative
#
###############################################################################

doc"""
    derivative(a::Nemo.PolyElem)
> Return the derivative of the polynomial $a$.
"""
function derivative(a::Nemo.PolyElem)
   if a == 0
      return zero(parent(a))
   end
   len = length(a)
   z = parent(a)()
   fit!(z, len - 1)
   for i = 1:len - 1
      z = setcoeff!(z, i - 1, i*coeff(a, i))
   end
   set_length!(z, normalise(z, len - 1))
   return z
end

###############################################################################
#
#   Integral
#
###############################################################################

doc"""
    integral{T <: Union{Nemo.ResElem, FieldElement}}(x::Nemo.PolyElem{T})
> Return the integral of the polynomial $a$.
"""
function integral(x::Nemo.PolyElem{T}) where {T <: Union{Nemo.ResElem, FieldElement}}
   len = length(x)
   p = parent(x)()
   fit!(p, len + 1)
   p = setcoeff!(p, 0, zero(base_ring(x)))
   for i = 1:len
      p = setcoeff!(p, i, divexact(coeff(x, i - 1), base_ring(x)(i)))
   end
   len += 1
   while len > 0 && coeff(p, len - 1) == 0 # FIXME: cannot use normalise here
      len -= 1
   end
   set_length!(p, len)
   return p
end

###############################################################################
#
#   Resultant
#
###############################################################################

# Dichotomous Lazard, computes Se0 from Sd0 and Sd1. See the paper,
# "Optimizations of the subresultant algorithm" by Lionel Ducos, J. Pure and
# Appl. Algebra 2000.
function subresultant_lazard(Sd0::Nemo.PolyElem{T}, Sd1::Nemo.PolyElem{T}) where T <: RingElement
   n = length(Sd0) - length(Sd1) - 1
   if n == 0
      return Sd1
   end
   x = lead(Sd1)
   y = lead(Sd0)
   a = 1 << (ndigits(n, 2) - 1) # floor(log_2(a))
   c = x
   n = n - a
   while a != 1
      a = a >> 1
      c = divexact(c*c, y)
      if n >= a
         c = divexact(c*x, y)
         n -= a
      end
   end
   return divexact(c*Sd1, y)
end

# Ducos optimised calculation of Se1. See the paper, "Optimizations of the
# subresultant algorithm" by Lionel Ducos, J. Pure and Appl. Algebra 2000.
function subresultant_ducos(A::Nemo.PolyElem{T}, Sd1::Nemo.PolyElem{T}, Se0::Nemo.PolyElem{T}, sd::T) where T <: RingElement
   d1 = length(A)
   e1 = length(Sd1)
   cd1 = lead(Sd1)
   se = lead(Se0)
   D = parent(A)()
   fit!(D, d1 - 1)
   for j = 0:e1 - 2
      setcoeff!(D, j, se*coeff(A, j))
   end
   set_length!(D, normalise(D, e1 - 1))
   Hj = parent(A)()
   fit!(Hj, e1)
   setcoeff!(Hj, e1 - 1, se)
   Hj -= Se0
   D += coeff(A, e1 - 1)*Hj
   for j = e1:d1 - 2
      Hj = shift_left(Hj, 1)
      Hj -= divexact(coeff(Hj, e1 - 1)*Sd1, cd1)
      D += coeff(A, j)*Hj
   end
   D = divexact(D, lead(A))
   Hj = shift_left(Hj, 1)
   r = divexact((Hj + D)*cd1 - coeff(Hj, e1 - 1)*Sd1, sd)
   return iseven(d1 - e1) ? -r : r
end

doc"""
    resultant{T <: RingElement}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T})
> Return the resultant of the $a$ and $b$.
"""
# See the paper, "Optimizations of the subresultant algorithm" by Lionel
# Ducos, J. Pure and Appl. Algebra 2000.
function resultant_ducos(p::Nemo.PolyElem{T}, q::Nemo.PolyElem{T}) where {T <: RingElement}
   check_parent(p, q)
   if length(p) == 0 || length(q) == 0
      return zero(base_ring(p))
   end
   sgn = 1
   if length(p) < length(q)
      p, q = q, p
      if iseven(length(p)) && iseven(length(q))
         sgn = -sgn
      end
   end
   lp = length(p)
   lq = length(q)
   if lq == 1
      return coeff(q, 0)^(lp - 1)
   end
   c1 = content(p)
   c2 = content(q)
   p = divexact(p, c1)
   q = divexact(q, c2)
   sd = lead(q)^(lp - lq)
   Sd0 = parent(p)()
   A = q
   B = pseudorem(p, -A)
   while true
      d1 = length(A)
      e1 = length(B)
      if e1 == 0
         return zero(base_ring(p))
      end
      Sd1 = B
      delta = d1 - e1
      if delta > 1
         if length(Sd0) == 0
            C = divexact(lead(B)^(delta - 1)*B, sd^(delta - 1))
         else
            C = subresultant_lazard(Sd0, Sd1)
         end
      else
         C = B
      end
      if e1 == 1
         return coeff(C, 0)*c1^(lq - 1)*c2^(lp - 1)*sgn
      end
      B = subresultant_ducos(A, Sd1, C, sd)
      Sd0 = C
      Sd1 = B
      A = C
      sd = lead(A)
   end
end

# details can be found in, "Optimizations of the subresultant algorithm" by
# Lionel Ducos, J. Pure and Appl. Algebra 2000. Note, the resultant is
# the constant coefficient of S_0 (aka S_00 in other sources)
function resultant_subresultant(p::Nemo.PolyElem{T}, q::Nemo.PolyElem{T}) where {T <: RingElement}
   check_parent(p, q)
   if length(p) == 0 || length(q) == 0
      return zero(base_ring(p))
   end
   sgn = 1
   if length(p) < length(q)
      p, q = q, p
      if iseven(length(p)) && iseven(length(q))
         sgn = -sgn
      end
   end
   lp = length(p)
   lq = length(q)
   if lq == 1
      return coeff(q, 0)^(lp - 1)
   end
   s = lead(q)^(lp - lq)
   S = parent(p)()
   A = q
   B = pseudorem(p, -q)
   while true
      d1 = length(A)
      e1 = length(B)
      if e1 == 0
         return zero(base_ring(p))
      end
      S = B
      delta = d1 - e1
      if delta > 1
         C = divexact(lead(B)^(delta - 1)*B, s^(delta - 1))
         S = C
      else
         C = B
      end
      if e1 == 1
         return coeff(S, 0)*sgn
      end
      B = divexact(pseudorem(A, -B), s^delta*lead(A))
      A = C
      s = lead(A)
   end
end

function resultant_lehmer(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: Union{Nemo.ResElem, FieldElement}}
   crossover = 40
   R = base_ring(a)
   check_parent(a, b)
   if length(a) == 0 || length(b) == 0
      return zero(base_ring(a))
   end
   sgn = 1
   if length(a) < length(b)
      a, b = b, a
      if iseven(length(a)) && iseven(length(b))
         sgn = -sgn
      end
   end
   lA = lenA = length(a)
   lB = lenB = length(b)
   if lenB == 1
      return coeff(b, 0)^(lA - 1)
   end
   c1 = content(a)
   c2 = content(b)
   A = divexact(a, c1)
   B = divexact(b, c2)
   s = R(1)
   while lenB > crossover/2 + 1
      shift = max(lenA - crossover, 0)
      a = shift_right(A, shift)
      b = shift_right(B, shift)
      u1, v1 = R(1), R(0)
      u2, v2 = R(0), R(1)
      lena = lenA - shift
      lenb = lenB - shift
      if lenb > crossover/2 + 1
         A = truncate(A, shift)
         B = truncate(B, shift)
         while lenb > crossover/2 + 1
            if iseven(lena + shift) && iseven(lenb + shift)
               sgn = -sgn
            end
            (q, b), a = divrem(a, b), b
            u1, u2 = u2, u1 - q*u2
            v1, v2 = v2, v1 - q*v2
            s *= lead(a)^(lena - length(b))
            lena = lenb
            lenb = length(b)
         end
         A, B = u1*A + v1*B + shift_left(a, shift), u2*A + v2*B + shift_left(b, shift)
      else
         if iseven(lenA) && iseven(lenB)
               sgn = -sgn
         end
         B, A = mod(A, B), B
         s *= lead(A)^(lenA - length(B))
      end
      lenA = length(A)
      lenB = length(B)
      if lenB == 0
         return zero(base_ring(a)), parent(A)(1), parent(A)(1)
      end
   end
   while lenB > 1
      if iseven(lenA) && iseven(lenB)
         sgn = -sgn
      end
      B, A = mod(A, B), B
      s *= lead(A)^(lenA - length(B))
      lenA = lenB
      lenB = length(B)
      if lenB == 0
         return zero(base_ring(A))
      end
   end
   s *= lead(B)^(lenA - 1)
   return c1^(lB - 1)*c2^(lA - 1)*s*sgn
end

function sylvester_matrix(p::Nemo.PolyElem{T}, q::Nemo.PolyElem{T}) where T <: RingElement
   check_parent(p, q)
   R = base_ring(p)
   if length(p) == 0 || length(q) == 0
      return zero_matrix(R, 0, 0)
   end
   m = degree(p)
   n = degree(q)
   M = zero_matrix(R, m + n, m + n)
   for i = 1:n
      for j = m:-1:0
         M[i, m - j + i] = coeff(p, j)
      end
   end
   for i = 1:m
      for j = n:-1:0
         M[i + n, n - j + i] = coeff(q, j)
       end
   end
   return M
end

function resultant_sylvester(p::Nemo.PolyElem{T}, q::Nemo.PolyElem{T}) where T <: RingElement
   check_parent(p, q)
   R = base_ring(p)
   if length(p) == 0 || length(q) == 0
      return R(0)
   end
   return det_df(sylvester_matrix(p, q))
end

function resultant(p::Nemo.PolyElem{T}, q::Nemo.PolyElem{T}) where {T <: RingElement}
  R = parent(p)
  if !isexact_type(T)
     return resultant_sylvester(p, q)
  end
  try
     return resultant_ducos(p, q)
  catch
     return resultant_sylvester(p, q)
  end
end

function resultant_euclidean(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where T <: Union{Nemo.ResElem, FieldElement}
   check_parent(a, b)
   if length(a) == 0 || length(b) == 0
      return zero(base_ring(a))
   end
   sgn = 1
   if length(a) < length(b)
      a, b = b, a
      if iseven(length(a)) && iseven(length(b))
         sgn = -sgn
      end
   end
   la = lena = length(a)
   lb = lenb = length(b)
   if lenb == 1
      return coeff(b, 0)^(la - 1)
   end
   c1 = content(a)
   c2 = content(b)
   A = divexact(a, c1)
   B = divexact(b, c2)
   s = base_ring(A)(1)
   la = lena = length(A)
   lb = lenb = length(B)
   while lenb > 1
      if iseven(lena) && iseven(lenb)
         sgn = -sgn
      end
      B, A = mod(A, B), B
      s *= lead(A)^(lena - length(B))
      lena = lenb
      lenb = length(B)
      if lenb == 0
         return zero(base_ring(a))
      end
   end
   s *= lead(B)^(lena - 1)
   return c1^(lb - 1)*c2^(la - 1)*s*sgn
end

function resultant(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: Union{Nemo.ResElem, FieldElement}}
   try
      return resultant_euclidean(a, b)
   catch
      return resultant_sylvester(a, b)
   end
end

###############################################################################
#
#   Discriminant
#
###############################################################################

doc"""
    discriminant(a::Nemo.PolyElem)
> Return the discrimnant of the given polynomial.
"""
function discriminant(a::Nemo.PolyElem)
   d = derivative(a)
   z = resultant(a, d)
   if length(a) - length(d) == 1
      z = divexact(z, lead(a))
   else
      z = z*lead(a)^(length(a) - length(d) - 2)
   end
   mod4 = (length(a) + 3) % 4 # degree mod 4
   return mod4 == 2 || mod4 == 3 ? -z : z
end

###############################################################################
#
#   RESX
#
###############################################################################

doc"""
    resx{T <: RingElement}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T})
> Return a tuple $(r, s, t)$ such that $r$ is the resultant of $a$ and $b$ and
> such that $r = a\times s + b\times t$.
"""
function resx(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: RingElement}
   check_parent(a, b)
   sgn = 1
   swap = false
   if length(a) < length(b)
      a, b = b, a
      swap = true
      if iseven(length(a)) && iseven(length(b))
         sgn = -sgn
      end
   end
   lena = length(a)
   lenb = length(b)
   if lenb == 0
      return zero(base_ring(a)), zero(parent(a)), zero(parent(a))
   end
   (lena <= 1 && lenb <= 1) && error("Constant polynomials in resx")
   c1 = content(a)
   c2 = content(b)
   A = divexact(a, c1)
   B = divexact(b, c2)
   g = one(base_ring(a))
   h = one(base_ring(a))
   u1, u2 = one(parent(a)), zero(parent(a))
   v1, v2 = zero(parent(a)), one(parent(a))
   while lenb > 1
      d = lena - lenb
      if iseven(lena) && iseven(lenb)
         sgn = -sgn
      end
      (Q, B), A = pseudodivrem(A, B), B
      lena = lenb
      lenb = length(B)
      if lenb == 0
         return zero(base_ring(a)), zero(parent(a)), zero(parent(a))
      end
      s = h^d
      B = divexact(B, g*s)
      t = lead(A)^(d + 1)
      u2, u1 = divexact(u1*t - Q*u2, g*s), u2
      v2, v1 = divexact(v1*t - Q*v2, g*s), v2
      g = lead(A)
      h = divexact(h*g^d, s)
   end
   s = divexact(h*lead(B)^(lena - 1), h^(lena - 1))
   res = c1^(length(b) - 1)*c2^(length(a) - 1)*s*sgn
   if length(b) > 1
      u2 *= c1^(length(b) - 2)*c2^(length(a) - 1)*sgn
   else
      u2 *= c2^(length(a) - 1)*sgn
      u2 = divexact(u2, c1)
   end
   if length(a) > 1
      v2 *= c1^(length(b) - 1)*c2^(length(a) - 2)*sgn
   else
      v2 *= c1^(length(b) - 1)*sgn
      v2 = divexact(v2, c2)
   end
   if lena != 2
      if lena > 1
         d1 = lead(B)^(lena - 2)
         d2 = h^(lena - 2)
         u2 = divexact(u2*d1, d2)
         v2 = divexact(v2*d1, d2)
      else
         u2 = divexact(u2*h, lead(B))
         v2 = divexact(v2*h, lead(B))
      end
   end
   if swap
      u2, v2 = v2, u2
   end
   return res, u2, v2
end

###############################################################################
#
#   GCDX
#
###############################################################################

doc"""
    gcdx{T <: Union{Nemo.ResElem, FieldElement}}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T})
> Return a tuple $(g, s, t)$ such that $g$ is the greatest common divisor of
> $a$ and $b$ and such that $r = a\times s + b\times t$.
"""
function gcdx(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: Union{Nemo.ResElem, FieldElement}}
   check_parent(a, b)
   !isexact_type(T) && error("gcdx requires exact Bezout domain")
   if length(a) == 0
      if length(b) == 0
         return zero(parent(a)), zero(parent(a)), zero(parent(a))
      else
         d = lead(b)
         return divexact(b, d), zero(parent(a)), divexact(one(parent(a)), d)
      end
   end
   if length(b) == 0
      d = lead(a)
      return divexact(a, d), divexact(one(parent(a)), d), zero(parent(a))
   end
   swap = false
   if length(a) < length(b)
      a, b = b, a
      swap = true
   end
   c1 = content(a)
   c2 = content(b)
   A = divexact(a, c1)
   B = divexact(b, c2)
   u1, u2 = parent(a)(inv(c1)), zero(parent(a))
   v1, v2 = zero(parent(a)), parent(a)(inv(c2))
   while length(B) > 0
      (Q, B), A = divrem(A, B), B
      u2, u1 = u1 - Q*u2, u2
      v2, v1 = v1 - Q*v2, v2
   end
   if swap
      u1, v1 = v1, u1
   end
   d = gcd(c1, c2)
   A, u1, v1 = d*A, d*u1, d*v1
   d = lead(A)
   return divexact(A, d), divexact(u1, d), divexact(v1, d)
end

doc"""
    gcdinv{T <: Union{Nemo.ResElem, FieldElement}}(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T})
> Return a tuple $(g, s)$ such that $g$ is the greatest common divisor of $a$
> and $b$ and such that $s = a^{-1} \pmod{b}$. This function is useful for
> inverting modulo a polynomial and checking that it really was invertible.
"""
function gcdinv(a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: Union{Nemo.ResElem, FieldElement}}
   check_parent(a, b)
   R = base_ring(a)
   if length(a) == 0
      if length(b) == 0
         return zero(parent(a)), zero(parent(a))
      else
         d = lead(b)
         return divexact(b, d), zero(parent(a))
      end
   end
   if length(b) == 0
      d = lead(a)
      return divexact(a, d), inv(d)
   end
   if length(a) < length(b)
      a, b = b, a
      u1, u2 = zero(parent(a)), one(parent(a))
   else
      u1, u2 = one(parent(a)), zero(parent(a))
   end
   lena = length(a)
   lenb = length(b)
   c1 = content(a)
   c2 = content(b)
   A = divexact(a, c1)
   B = divexact(b, c2)
   u1 *= inv(c1)
   u2 *= inv(c2)
   while lenb > 0
      d = lena - lenb
      (Q, B), A = divrem(A, B), B
      lena = lenb
      lenb = length(B)
      u2, u1 = u1 - Q*u2, u2
   end
   d = gcd(c1, c2)
   A, u1 = d*A, d*u1
   d = lead(A)
   return divexact(A, d), divexact(u1, d)
end

###############################################################################
#
#   Newton representation
#
###############################################################################

doc"""
    monomial_to_newton!{T <: RingElement}(P::Array{T, 1}, roots::Array{T, 1})
> Converts a polynomial $p$, given as an array of coefficients, in-place
> from its coefficients given in the standard monomial basis to the Newton
> basis for the roots $r_0, r_1, \ldots, r_{n-2}$. In other words, this
> determines output coefficients $c_i$ such that
> $$c_0 + c_1(x-r_0) + c_2(x-r_0)(x-r_1) + \ldots + c_{n-1}(x-r_0)(x-r_1)\cdots(x-r_{n-2})$$
> is equal to the input polynomial.
"""
function monomial_to_newton!(P::Array{T, 1}, roots::Array{T, 1}) where {T <: RingElement}
   n = length(roots)
   if n > 0
      R = parent(roots[1])
      t = R()
      for i = 1:n - 1
         for j = n - 1:-1:i
            t = mul!(t, P[j + 1], roots[i])
            P[j] = addeq!(P[j], t)
         end
      end
   end
   return
end

doc"""
    newton_to_monomial!{T <: RingElement}(P::Array{T, 1}, roots::Array{T, 1})
> Converts a polynomial $p$, given as an array of coefficients, in-place
> from its coefficients given in the Newton basis for the roots
> $r_0, r_1, \ldots, r_{n-2}$ to the standard monomial basis. In other words,
> this evaluates
> $$c_0 + c_1(x-r_0) + c_2(x-r_0)(x-r_1) + \ldots + c_{n-1}(x-r_0)(x-r_1)\cdots(x-r_{n-2})$$
> where $c_i$ are the input coefficients given by $p$.
"""
function newton_to_monomial!(P::Array{T, 1}, roots::Array{T, 1}) where {T <: RingElement}
   n = length(roots)
   if n > 0
      R = parent(roots[1])
      t = R()
      for i = n - 1:-1:1
         d = -roots[i]
         for j = i:n - 1
            t = mul!(t, P[j + 1], d)
            P[j] = addeq!(P[j], t)
         end
      end
   end
   return
end

###############################################################################
#
#   Interpolation
#
###############################################################################

doc"""
    interpolate{T <: RingElement}(S::Nemo.PolyRing, x::Array{T, 1}, y::Array{T, 1})
> Given two arrays of values $xs$ and $ys$ of the same length $n$, find
> the polynomial $f$ in the polynomial ring $R$ of length at most $n$ such that
> $f$ has the value $ys$ at the points $xs$. The values in the arrays $xs$ and
> $ys$ must belong to the base ring of the polynomial ring $R$. If no such
> polynomial exists, an exception is raised.
"""
function interpolate(S::Nemo.PolyRing, x::Array{T, 1}, y::Array{T, 1}) where {T <: RingElement}
   length(x) != length(y) && error("Array lengths don't match in interpolate")
   !isdomain_type(T) && error("Generic interpolation requires a domain type")
   n = length(x)
   if n == 0
      return S()
   elseif n == 1
      return S(y[1])
   end
   R = base_ring(S)
   parent(y[1]) != R && error("Polynomial ring does not match inputs")
   P = Array{T}(n)
   for i = 1:n
      P[i] = deepcopy(y[i])
   end
   for i = 2:n
      t = P[i - 1]
      for j = i:n
         p = P[j] - t
         q = x[j] - x[j - i + 1]
         t = P[j]
         P[j] = divexact(p, q) # division is exact over domain (Lipson, 1971)
      end
   end
   newton_to_monomial!(P, x)
   r = S(P)
   set_length!(r, normalise(r, n))
   return r
end

function interpolate(S::Nemo.PolyRing, x::Array{T, 1}, y::Array{T, 1}) where {T <: ResElem}
   length(x) != length(y) && error("Array lengths don't match in interpolate")
   n = length(x)
   if n == 0
      return S()
   elseif n == 1
      return S(y[1])
   end
   R = base_ring(S)
   parent(y[1]) != R && error("Polynomial ring does not match inputs")
   P = Array{T}(n)
   for i = 1:n
      P[i] = deepcopy(y[i])
   end
   for i = 2:n
      t = P[i - 1]
      for j = i:n
         p = P[j] - t
         q = x[j] - x[j - i + 1]
         t = P[j]
         P[j] = p*inv(q) # must have invertible q for now
      end
   end
   newton_to_monomial!(P, x)
   r = S(P)
   set_length!(r, normalise(r, n))
   return r
end

###############################################################################
#
#   Special functions
#
###############################################################################

function chebyshev_t_pair(n::Int, x::Nemo.PolyElem)
   if n == 0
      return one(parent(x)), x
   elseif n == 1
      return x, one(parent(x))
   elseif n < 0
      a, b = chebyshev_t_pair(1 - n, x)
      return b, a
   elseif iseven(n)
      a, b = chebyshev_t_pair(n >> 1, x)
      return 2*(a*a) - 1, 2*(a*b) - x
   else
      a, b = chebyshev_t_pair((n >> 1) + 1, x)
      return 2*(a*b) - x, 2*(b*b) - 1
   end
end

doc"""
    chebyshev_t(n::Int, x::Nemo.PolyElem)
> Return the Chebyshev polynomial of the first kind $T_n(x)$, defined by
> $T_n(x) = \cos(n \cos^{-1}(x))$.
"""
function chebyshev_t(n::Int, x::Nemo.PolyElem)
   if n == 0
      return one(parent(x))
   elseif n == 1
      return x
   elseif n < 0
      return chebyshev_t(-n, x)
   elseif iseven(n)
      a = chebyshev_t(n >> 1, x)
      return 2*(a*a) - 1
   else
      a, b = chebyshev_t_pair((n >> 1) + 1, x)
      return 2*(a*b) - x
   end
end

function chebyshev_u_pair(n::Int, x::Nemo.PolyElem)
   if n == 0
      return one(parent(x)), zero(parent(x))
   elseif n == 1
      return 2*x, one(parent(x))
   elseif n == -1
      return zero(parent(x)), -one(parent(x))
   elseif n < -1
      a, b = chebyshev_u_pair(-1 - n, x)
      return -b, -a
   elseif iseven(n)
      a, b = chebyshev_u_pair(n >> 1, x)
      return (a+b)*(a - b), 2*b*(a-x*b)
   else
      a, b = chebyshev_u_pair(n >> 1, x)
      return 2*a*(x*a - b), (a + b)*(a - b)
   end
end

doc"""
    chebyshev_u(n::Int, x::Nemo.PolyElem)
> Return the Chebyshev polynomial of the first kind $U_n(x)$, defined by
> $(n+1) U_n(x) = T'_{n+1}(x)$.
"""
function chebyshev_u(n::Int, x::Nemo.PolyElem)
   if n == 0
      return one(parent(x))
   elseif n == 1
      return 2*x
   elseif n == -1
      return zero(parent(x))
   elseif n < -1
      return -chebyshev_u(-2 - n, x)
   elseif iseven(n)
      a, b = chebyshev_u_pair(n >> 1, x)
      return (a + b)*(a - b)
   else
      a, b = chebyshev_u_pair(n >> 1, x)
      return 2*a*(x*a - b)
   end
end

###############################################################################
#
#   Unsafe functions
#
###############################################################################

function set_length!(c::Poly{T}, n::Int) where T <: RingElement
   if n < c.length
      for i = n + 1:c.length
         c.coeffs[i] = zero!(c.coeffs[i])
      end
   end
   c.length = n
end

function fit!(c::Poly{T}, n::Int) where {T <: RingElement}
   if length(c.coeffs) < n
      t = c.coeffs
      c.coeffs = Array{T}(n)
      for i = 1:length(c)
         c.coeffs[i] = t[i]
      end
      for i = length(c) + 1:n
         c.coeffs[i] = zero(base_ring(c))
      end
   end
   return nothing
end

function zero!(c::Poly{T}) where {T <: RingElement}
   set_length!(c, 0)
   return c
end

function mul!(c::Nemo.PolyElem{T}, a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: RingElement}
   lena = length(a)
   lenb = length(b)

   if lena == 0 || lenb == 0
      set_length!(c, 0)
   else
      if a === c
         a = deepcopy(a)
      end
      if b === c
         b = deepcopy(b)
      end

      t = base_ring(a)()

      lenc = lena + lenb - 1
      fit!(c, lenc)

      for i = 1:lena
         c.coeffs[i] = mul!(c.coeffs[i], coeff(a, i - 1), coeff(b, 0))
      end

      for i = 2:lenb
         c.coeffs[lena + i - 1] = mul!(c.coeffs[lena + i - 1], coeff(a, lena - 1), coeff(b, i - 1))
      end

      for i = 1:lena - 1
         for j = 2:lenb
            t = mul!(t, coeff(a, i - 1), coeff(b, j - 1))
            c.coeffs[i + j - 1] = addeq!(c.coeffs[i + j - 1], t)
         end
      end

      set_length!(c, normalise(c, lenc))
   end
   return c
end

function addeq!(c::Nemo.PolyElem{T}, a::Nemo.PolyElem{T}) where {T <: RingElement}
   lenc = length(c)
   lena = length(a)
   len = max(lenc, lena)
   fit!(c, len)
   for i = 1:lena
      c.coeffs[i] = addeq!(c.coeffs[i], coeff(a, i - 1))
   end
   set_length!(c, normalise(c, len))
   return c
end

function add!(c::Nemo.PolyElem{T}, a::Nemo.PolyElem{T}, b::Nemo.PolyElem{T}) where {T <: RingElement}
   lena = length(a)
   lenb = length(b)
   len = max(lena, lenb)
   fit!(c, len)
   i = 1
   while i <= min(lena, lenb)
      c.coeffs[i] = add!(c.coeffs[i], coeff(a, i - 1), coeff(b, i - 1))
      i += 1
   end
   while i <= lena
      c = setcoeff!(c, i - 1, coeff(a, i - 1))
      i += 1
   end
   while i <= lenb
      c = setcoeff!(c, i - 1, coeff(b, i - 1))
      i += 1
   end
   set_length!(c, normalise(c, len))
   return c
end

function addmul!(z::Nemo.PolyElem{T}, x::Nemo.PolyElem{T}, y::Nemo.PolyElem{T}, c::Nemo.PolyElem{T}) where {T <: RingElement}
   c = mul!(c, x, y)
   z = addeq!(z, c)
   return z
end

###############################################################################
#
#   Random elements
#
###############################################################################

function rand(S::Nemo.PolyRing, deg_range::UnitRange{Int}, v...)
   R = base_ring(S)
   f = S()
   x = gen(S)
   for i = 0:rand(deg_range)
      f += rand(R, v...)*x^i
   end
   return f
end

###############################################################################
#
#   Promotion rules
#
###############################################################################

promote_rule(::Type{Poly{T}}, ::Type{Poly{T}}) where T <: RingElement = Poly{T}

function promote_rule(::Type{Poly{T}}, ::Type{U}) where {T <: RingElement, U <: RingElement}
   promote_rule(T, U) == T ? Poly{T} : Union{}
end

###############################################################################
#
#   Polynomial substitution
#
###############################################################################

doc"""
    subst{T <: RingElement}(f::Nemo.PolyElem{T}, a::Any)
> Evaluate the polynomial $f$ at $a$. Note that $a$ can be anything, whether
> a ring element or not.
"""
function subst(f::Nemo.PolyElem{T}, a::Any) where {T <: RingElement}
   S = parent(a)
   n = degree(f)
   if n < 0
      return S()
   elseif n == 0
      return coeff(f, 0)*S(1)
   elseif n == 1
      return coeff(f, 0)*S(1) + coeff(f, 1)*a
   end
   d1 = isqrt(n)
   d = div(n, d1)
   A = powers(a, d)
   s = coeff(f, d1*d)*A[1]
   for j = 1:min(n - d1*d, d - 1)
      c = coeff(f, d1*d + j)
      if !iszero(c)
         s += c*A[j + 1]
      end
   end
   for i = 1:d1
      s *= A[d + 1]
      s += coeff(f, (d1 - i)*d)*A[1]
      for j = 1:min(n - (d1 - i)*d, d - 1)
         c = coeff(f, (d1 - i)*d + j)
         if !iszero(c)
            s += c*A[j + 1]
         end
      end
   end
   return s
end

###############################################################################
#
#   Parent object call overload
#
###############################################################################

function (a::PolyRing{T})(b::RingElement) where {T <: RingElement}
   return a(base_ring(a)(b))
end

function (a::PolyRing{T})() where {T <: RingElement}
   z = Poly{T}()
   z.parent = a
   return z
end

function (a::PolyRing{T})(b::Union{Integer, Rational, AbstractFloat}) where {T <: RingElement}
   z = Poly{T}(base_ring(a)(b))
   z.parent = a
   return z
end

function (a::PolyRing{T})(b::T) where {T <: RingElement}
   parent(b) != base_ring(a) && error("Unable to coerce to polynomial")
   z = Poly{T}(b)
   z.parent = a
   return z
end

function (a::PolyRing{T})(b::Nemo.PolyElem{T}) where {T <: RingElement}
   parent(b) != a && error("Unable to coerce polynomial")
   return b
end

function (a::PolyRing{T})(b::Array{T, 1}) where T <: RingElement
   R = base_ring(a)
   for i = 1:length(b)
      b[i] = R(b[i])
   end
   z = Poly{T}(b)
   z.parent = a
   return z
end

function (a::PolyRing{T})(b::Array{S, 1}) where {S <: RingElement, T <: RingElement}
   R = base_ring(a)
   len = length(b)
   entries = Array{T}(len)
   for i = 1:length(b)
      entries[i] = R(b[i])
   end
   z = Poly{T}(entries)
   z.parent = a
   return z
end

# Functions to remove ambiguities on julia 0.7
function (a::PolyRing{T})(b::T) where {T <: Rational}
   parent(b) != base_ring(a) && error("Unable to coerce to polynomial")
   z = Poly{T}(b)
   z.parent = a
   return z
end

function (a::PolyRing{T})(b::T) where {T <: AbstractFloat}
   parent(b) != base_ring(a) && error("Unable to coerce to polynomial")
   z = Poly{T}(b)
   z.parent = a
   return z
end

function (a::PolyRing{T})(b::T) where {T <: Integer}
   parent(b) != base_ring(a) && error("Unable to coerce to polynomial")
   z = Poly{T}(b)
   z.parent = a
   return z
end

###############################################################################
#
#   PolynomialRing constructor
#
###############################################################################

doc"""
    PolynomialRing(R::Nemo.Ring, s::AbstractString; cached::Bool = true)
> Given a base ring `R` and string `s` specifying how the generator (variable)
> should be printed, return a tuple `S, x` representing the new polynomial
> ring $S = R[x]$ and the generator $x$ of the ring. By default the parent
> object `S` will depend only on `R` and `x` and will be cached. Setting the
> optional argument `cached` to `false` will prevent the parent object `S` from
> being cached.
"""
function PolynomialRing(R::Nemo.Ring, s::AbstractString; cached::Bool = true)
   S = Symbol(s)
   T = elem_type(R)
   parent_obj = PolyRing{T}(R, S, cached)

   return parent_obj, parent_obj([R(0), R(1)])
end


text:
function getTetraMeshNumbering(Tetras)

# tetrahedron connectivity
t = Tetras

nc = size(t,1)
# global edge vertex numbers
edges = zeros(Int64,6*nc,2)

cnt = 1
for i=1:nc
	edges[cnt,1] = t[i,1]
	edges[cnt,2] = t[i,2]; cnt += 1

	edges[cnt,1] = t[i,1]
	edges[cnt,2] = t[i,3]; cnt += 1

	edges[cnt,1] = t[i,1]
	edges[cnt,2] = t[i,4]; cnt += 1

	edges[cnt,1] = t[i,2]
	edges[cnt,2] = t[i,3]; cnt += 1

	edges[cnt,1] = t[i,2]
	edges[cnt,2] = t[i,4]; cnt += 1

	edges[cnt,1] = t[i,3]
	edges[cnt,2] = t[i,4]; cnt += 1
end

# collate duplicate edges
edges = unique(sortrows(edges), 1)

# face degrees of freedom
faces = zeros(Int64,4*nc,3)
cnt = 1
for i=1:nc
	faces[cnt,1] = t[i,1]
	faces[cnt,2] = t[i,2]
	faces[cnt,3] = t[i,3]; cnt += 1

	faces[cnt,1] = t[i,1]
	faces[cnt,2] = t[i,2]
	faces[cnt,3] = t[i,4]; cnt += 1

	faces[cnt,1] = t[i,1]
	faces[cnt,2] = t[i,3]
	faces[cnt,3] = t[i,4]; cnt += 1

	faces[cnt,1] = t[i,2]
	faces[cnt,2] = t[i,3]
	faces[cnt,3] = t[i,4]; cnt += 1
end

# global face vertex numbers
faces = unique(sortrows(faces), 1)

return edges, faces

end

ext: jl

max_stars_repo_name: JuliaInv/JTetra

max_issues_repo_name: JuliaInv/JTetra

max_forks_repo_name: JuliaInv/JTetra


text:
using RadialBasisFunctionModels
using Test

test_files = ["README.jl", "test_machine_wrapper.jl", "test_radial_funcs.jl", "test_rbf_models.jl"]

@testset "RadialBasisFunctionModels.jl" begin
    for fn in test_files
        include(fn)
    end
end


text:
# https://www.codewars.com/kata/515f51d438015969f7000013

function pyramid(n)
  if n <= 0
    []
  else
    [[1 for _ in (1:a)] for a in (1:n)]
  end
end

max_stars_repo_name: hiljusti/codewars-solutions

max_issues_repo_name: hiljusti/codewars-solutions

max_forks_repo_name: hiljusti/codewars-solutions


text:
using Revise
using Gadfly
using ANOVADDPTest
using DataFrames
using Distributions
using Random

function generate_plot_normal(N, G, distinct_groups)
    rng = MersenneTwister(1)
    m = NormalDDP(rng, N, G)
    x = rand(rng, 1:G, N)
    y = randn(rng, N)
    for i in 1:N
        x[i] in distinct_groups && (y[i] += 1)
    end
    data = NormalData(x, y)
    pred = NormalData(expandgrid(1:G,  range(-2.0, stop = 2.0, length = 50))...)
    ch = train(rng, m, data, pred)
    df = DataFrame(x = pred.x, y = pred.y, f = mean(ch.f))
    return plot(df, x = :y, y = :f, color = :x, Geom.line, Scale.color_discrete())
end

generate_plot_normal(1000, 3, [])
generate_plot_normal(1000, 3, [2])
generate_plot_normal(1000, 3, [3])
generate_plot_normal(1000, 3, [2, 3])


text:
__precompile__()

module DiffBase

using StaticArrays
using Compat

include("results.jl")
include("testfuncs.jl")
include("rules.jl")

end # module


text:
import Test
import PredictMD


text:
using TestRunner
using FactCheck
sampleTests = quote
  using FactCheck

  facts("First facts group") do
    @fact 1 --> 1 "First group first test"
    @fact 2 --> 1 "First group second failing test"
  end

  facts("Second facts group") do
    @fact 1 --> 1 "Second group first test"
    @fact 2 --> 2 #Nameless fact
  end

  facts() do #Nameless group
    @fact 1 --> 1 #Nameless fact
    @fact 2 --> 1 #Failing nameless fact in nameless group
  end
  facts() do
  end
end

if VERSION >= v"0.5-"
  walk(f, leaf) = f(leaf)
  function walk(f, ex::Expr)
      ex = copy(ex)
      ex.args = map((x)->walk(f,x), ex.args)
      ex
  end

  remove_file_name(node::LineNumberNode) = LineNumberNode(:none, node.line)
  remove_file_name(exp) = exp
  sampleTests = walk(remove_file_name, sampleTests)
end

using TestRunner: TestStructureNode, RootNode, FactsCollectionNode, FactNode, ContextNode

include("testComparisons.jl")

sampleTestsFilePath = "FileWithSampleTests.jl"
sampleExceptionFailingTestFilePath = "FileWithExceptionFailingTest.jl"

facts("Test handling tests") do
  context("File loading") do
      @fact TestRunner._get_file_content(sampleTestsFilePath) --> sampleTests
  end

  context("Structure parsing") do
   @fact TestRunner._get_tests_structure(sampleTests) --> Vector{TestStructureNode}([
                                                        FactsCollectionNode(6,"First facts group",
                                                        [FactNode(7, "First group first test"),FactNode(8, "First group second failing test")]),
                                                        FactsCollectionNode(11,"Second facts group",
                                                        [FactNode(12,"Second group first test"),FactNode(13, "")]),
                                                        FactsCollectionNode(16,"",
                                                        [FactNode(17,""),FactNode(18,"")]),
                                                        FactsCollectionNode(20,"",[])
                                                      ])

  @fact get_tests_structure(sampleTestsFilePath) --> RootNode(Vector{TestStructureNode}([
                                                       FactsCollectionNode(6,"First facts group",
                                                       [FactNode(7, "First group first test"),FactNode(8, "First group second failing test")]),
                                                       FactsCollectionNode(11,"Second facts group",
                                                       [FactNode(12,"Second group first test"),FactNode(13, "")]),
                                                       FactsCollectionNode(16,"",
                                                       [FactNode(17,""),FactNode(18,"")]),
                                                       FactsCollectionNode(20,"",[])
                                                     ])) "File structure parsing"
 end

  context("Tests running") do
    @fact run_all_tests(sampleTestsFilePath) --> RootNode( Vector{TestStructureNode}([
                                                         FactsCollectionNode(6, "First facts group",
                                                          [
                                                           FactNode(7, "First group first test", test_success, "    Success :: (line:7) :: Tests running :: First group first test :: fact was true\n      Expression: 1 --> 1\n        Expected: 1\n        Occurred: 1"),
                                                           FactNode(8, "First group second failing test", test_failure, "    Failure :: (line:8) :: Tests running :: First group second failing test :: fact was false\n      Expression: 2 --> 1\n        Expected: 1\n        Occurred: 2")
                                                          ]),
                                                         FactsCollectionNode(11,"Second facts group",
                                                          [
                                                           FactNode(12, "Second group first test", test_success, "    Success :: (line:12) :: Tests running :: Second group first test :: fact was true\n      Expression: 1 --> 1\n        Expected: 1\n        Occurred: 1"),
                                                           FactNode(13, "", test_success, "    Success :: (line:13) :: Tests running :: fact was true\n      Expression: 2 --> 2\n        Expected: 2\n        Occurred: 2")
                                                          ]),
                                                         FactsCollectionNode(16, "",
                                                          [
                                                           FactNode(17, "", test_success, "    Success :: (line:17) :: Tests running :: fact was true\n      Expression: 1 --> 1\n        Expected: 1\n        Occurred: 1"),
                                                           FactNode(18, "", test_failure, "    Failure :: (line:18) :: Tests running :: fact was false\n      Expression: 2 --> 1\n        Expected: 1\n        Occurred: 2")
                                                          ]),
                                                         FactsCollectionNode(20,"",[])
                                                       ])) "Tests running"
  end

  context("Exception handling") do
    message = "TestName\nDoge\n in f at $(abspath(sampleExceptionFailingTestFilePath)):2"
    @pending run_all_tests(sampleExceptionFailingTestFilePath) --> Vector{TestStructureNode}([FactsCollectionNode(3,"",[FactNode(4,"TestName", false, message)])])
  end

  context("Children acquiring") do
    @fact children(FactNode(42, "Test")) --> Vector{TestStructureNode}()
    @fact children(FactsCollectionNode(17,"", [FactNode(18,""),FactNode(19,"")])) --> Vector{TestStructureNode}([FactNode(18,""),FactNode(19,"")])
  end
end

facts("JSON parsing tests") do
  include("TestRunnerJSON.tests.jl")
end

FactCheck.exitstatus()


text:
### A Pluto.jl notebook ###
# v0.15.1

using Markdown
using InteractiveUtils

# ╔═╡ 90facaa4-158c-11ec-3a0f-2b30580dc6b7
begin
	using CSV, DataFrames, PlotlyJS, WebIO
end

# ╔═╡ 40de9078-38b7-4238-b7a4-fcb1a440d753
begin
	using PlotlyBase, HypertextLiteral

	# Custom abstract-trace array to string function
	TraceData(data::Array{<:AbstractTrace}) = join(["[",join(string.(data),","),"]"])

	function Base.show(io::IO, mimetype::MIME"text/html", p::PlotlyBase.Plot)
		show(io,mimetype,@htl("""
		<script src='https://cdn.plot.ly/plotly-latest.min.js'></script>
		<div style="height: auto">
			<script>
			let parentnode = currentScript.parentElement

		  Plotly.plot(parentnode,$(HypertextLiteral.JavaScript(TraceData(p.data))),$(HypertextLiteral.JavaScript(string(p.layout))),{responsive: true});
			</script>
		</div>
	"""))
	end
end

# ╔═╡ 360b9e44-dff2-40c2-ab68-818f50447d76
begin
	d1 = CSV.read("ft2632_hanning_620.txt", DataFrame, delim="\t");
	d2 = CSV.read("ft2767_hanning_620.txt", DataFrame, delim="\t");
	d3 = CSV.read("ft2768_hanning_620.txt", DataFrame, delim="\t");
end

# ╔═╡ 8659283f-0618-4af9-b78a-948cf8755a00
begin
	p1 = plot(d1, x=:freq2632, y=:ft2632, Layout(title="Benzene", xaxis_title="Frequency / MHz", yaxis_title="Intensity", paper_bgcolor="#fff"))
	p2 = plot(d2, x=:freq2767, y=:ft2767, Layout(title="Benzene + O2", xaxis_title="Frequency / MHz", yaxis_title="Intensity", paper_bgcolor="#fff"))
	p3 = plot(d3, x=:freq2768, y=:ft2768, Layout(title="Benzene + N2", xaxis_title="Frequency / MHz", yaxis_title="Intensity", paper_bgcolor="#fff"))
	combined = [p1; p2; p3]
	relayout!(combined, showlegend=false, sharex=true)
end

# ╔═╡ c76c0762-da6b-4f3d-9d68-5bfdcf0be3fc
combined

# ╔═╡ 6d525819-4724-42ee-8934-d1d66e7f9eb6
PlotlyJS.savefig(combined, "discharge-plots.html")

# ╔═╡ c1ada4ac-6b51-4490-b5c0-64d724ae88a3
begin
	joint = make_subplots(rows=3, cols=1, shared_xaxes=true)
	add_trace!(joint, scatter(d1, x=:freq2632, y=:ft2632, name="Benzene"), row=1, col=1)
	add_trace!(joint, scatter(d2, x=:freq2767, y=:ft2767, name="Benzene + O2"), row=2, col=1)
	add_trace!(joint, scatter(d3, x=:freq2768, y=:ft2768, name="Benzene + N2"), row=3, col=1)
	relayout!(joint, xaxis_title="Frequency / MHz")
end

# ╔═╡ ad5b05e7-b300-4de5-8abe-52a44ab01c70
PlotlyJS.savefig(joint, "durr.html")

# ╔═╡ bef30ffd-d141-4520-9b34-5b93a83e0d1a

# ╔═╡ 00000000-0000-0000-0000-000000000001
PLUTO_PROJECT_TOML_CONTENTS = """
[deps]
CSV = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b"
DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
HypertextLiteral = "ac1192a8-f4b3-4bfe-ba22-af5b92cd3ab2"
PlotlyBase = "a03496cd-edff-5a9b-9e67-9cda94a718b5"
PlotlyJS = "f0f68f2c-4968-5e81-91da-67840de0976a"
WebIO = "0f1e0344-ec1d-5b48-a673-e5cf874b6c29"

[compat]
CSV = "~0.9.1"
DataFrames = "~1.2.2"
HypertextLiteral = "~0.9.0"
PlotlyBase = "~0.8.16"
PlotlyJS = "~0.18.7"
WebIO = "~0.8.15"
"""

# ╔═╡ 00000000-0000-0000-0000-000000000002
PLUTO_MANIFEST_TOML_CONTENTS = """
# This file is machine-generated - editing it directly is not advised

[[ArgTools]]
uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f"

[[Artifacts]]
uuid = "56f22d72-fd6d-98f1-02f0-08ddc0907c33"

[[AssetRegistry]]
deps = ["Distributed", "JSON", "Pidfile", "SHA", "Test"]
git-tree-sha1 = "b25e88db7944f98789130d7b503276bc34bc098e"
uuid = "bf4720bc-e11a-5d0c-854e-bdca1663c893"
version = "0.1.0"

[[Base64]]
uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f"

[[BinDeps]]
deps = ["Libdl", "Pkg", "SHA", "URIParser", "Unicode"]
git-tree-sha1 = "1289b57e8cf019aede076edab0587eb9644175bd"
uuid = "9e28174c-4ba2-5203-b857-d8d62c4213ee"
version = "1.0.2"

[[Blink]]
deps = ["Base64", "BinDeps", "Distributed", "JSExpr", "JSON", "Lazy", "Logging", "MacroTools", "Mustache", "Mux", "Reexport", "Sockets", "WebIO", "WebSockets"]
git-tree-sha1 = "08d0b679fd7caa49e2bca9214b131289e19808c0"
uuid = "ad839575-38b3-5650-b840-f874b8c74a25"
version = "0.12.5"

[[CSV]]
deps = ["CodecZlib", "Dates", "FilePathsBase", "Mmap", "Parsers", "PooledArrays", "SentinelArrays", "Tables", "Unicode", "WeakRefStrings"]
git-tree-sha1 = "c907e91e253751f5840135f4c9deb1308273338d"
uuid = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b"
version = "0.9.1"

[[CodecZlib]]
deps = ["TranscodingStreams", "Zlib_jll"]
git-tree-sha1 = "ded953804d019afa9a3f98981d99b33e3db7b6da"
uuid = "944b1d66-785c-5afd-91f1-9de20f533193"
version = "0.7.0"

[[ColorSchemes]]
deps = ["ColorTypes", "Colors", "FixedPointNumbers", "Random"]
git-tree-sha1 = "9995eb3977fbf67b86d0a0a0508e83017ded03f2"
uuid = "35d6a980-a343-548e-a6ea-1d62b119f2f4"
version = "3.14.0"

[[ColorTypes]]
deps = ["FixedPointNumbers", "Random"]
git-tree-sha1 = "024fe24d83e4a5bf5fc80501a314ce0d1aa35597"
uuid = "3da002f7-5984-5a60-b8a6-cbb66c0b333f"
version = "0.11.0"

[[Colors]]
deps = ["ColorTypes", "FixedPointNumbers", "Reexport"]
git-tree-sha1 = "417b0ed7b8b838aa6ca0a87aadf1bb9eb111ce40"
uuid = "5ae59095-9a9b-59fe-a467-6f913c188581"
version = "0.12.8"

[[Compat]]
deps = ["Base64", "Dates", "DelimitedFiles", "Distributed", "InteractiveUtils", "LibGit2", "Libdl", "LinearAlgebra", "Markdown", "Mmap", "Pkg", "Printf", "REPL", "Random", "SHA", "Serialization", "SharedArrays", "Sockets", "SparseArrays", "Statistics", "Test", "UUIDs", "Unicode"]
git-tree-sha1 = "6071cb87be6a444ac75fdbf51b8e7273808ce62f"
uuid = "34da2185-b29b-5c13-b0c7-acf172513d20"
version = "3.35.0"

[[Crayons]]
git-tree-sha1 = "3f71217b538d7aaee0b69ab47d9b7724ca8afa0d"
uuid = "a8cc5b0e-0ffa-5ad4-8c14-923d3ee1735f"
version = "4.0.4"

[[DataAPI]]
git-tree-sha1 = "bec2532f8adb82005476c141ec23e921fc20971b"
uuid = "9a962f9c-6df0-11e9-0e5d-c546b8b5ee8a"
version = "1.8.0"

[[DataFrames]]
deps = ["Compat", "DataAPI", "Future", "InvertedIndices", "IteratorInterfaceExtensions", "LinearAlgebra", "Markdown", "Missings", "PooledArrays", "PrettyTables", "Printf", "REPL", "Reexport", "SortingAlgorithms", "Statistics", "TableTraits", "Tables", "Unicode"]
git-tree-sha1 = "d785f42445b63fc86caa08bb9a9351008be9b765"
uuid = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
version = "1.2.2"

[[DataStructures]]
deps = ["Compat", "InteractiveUtils", "OrderedCollections"]
git-tree-sha1 = "7d9d316f04214f7efdbb6398d545446e246eff02"
uuid = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"
version = "0.18.10"

[[DataValueInterfaces]]
git-tree-sha1 = "bfc1187b79289637fa0ef6d4436ebdfe6905cbd6"
uuid = "e2d170a0-9d28-54be-80f0-106bbe20a464"
version = "1.0.0"

[[Dates]]
deps = ["Printf"]
uuid = "ade2ca70-3891-5945-98fb-dc099432e06a"

[[DelimitedFiles]]
deps = ["Mmap"]
uuid = "8bb1440f-4735-579b-a4ab-409b98df4dab"

[[Distributed]]
deps = ["Random", "Serialization", "Sockets"]
uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b"

[[DocStringExtensions]]
deps = ["LibGit2"]
git-tree-sha1 = "a32185f5428d3986f47c2ab78b1f216d5e6cc96f"
uuid = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
version = "0.8.5"

[[Downloads]]
deps = ["ArgTools", "LibCURL", "NetworkOptions"]
uuid = "f43a241f-c20a-4ad4-852c-f6b1247861c6"

[[FilePathsBase]]
deps = ["Dates", "Mmap", "Printf", "Test", "UUIDs"]
git-tree-sha1 = "0f5e8d0cb91a6386ba47bd1527b240bd5725fbae"
uuid = "48062228-2e41-5def-b9a4-89aafe57970f"
version = "0.9.10"

[[FileWatching]]
uuid = "7b1f6079-737a-58dc-b8bc-7a2ca5c1b5ee"

[[FixedPointNumbers]]
deps = ["Statistics"]
git-tree-sha1 = "335bfdceacc84c5cdf16aadc768aa5ddfc5383cc"
uuid = "53c48c17-4a7d-5ca2-90c5-79b7896eea93"
version = "0.8.4"

[[Formatting]]
deps = ["Printf"]
git-tree-sha1 = "8339d61043228fdd3eb658d86c926cb282ae72a8"
uuid = "59287772-0a20-5a39-b81b-1366585eb4c0"
version = "0.4.2"

[[FunctionalCollections]]
deps = ["Test"]
git-tree-sha1 = "04cb9cfaa6ba5311973994fe3496ddec19b6292a"
uuid = "de31a74c-ac4f-5751-b3fd-e18cd04993ca"
version = "0.5.0"

[[Future]]
deps = ["Random"]
uuid = "9fa8497b-333b-5362-9e8d-4d0656e87820"

[[HTTP]]
deps = ["Base64", "Dates", "IniFile", "Logging", "MbedTLS", "NetworkOptions", "Sockets", "URIs"]
git-tree-sha1 = "60ed5f1643927479f845b0135bb369b031b541fa"
uuid = "cd3eb016-35fb-5094-929b-558a96fad6f3"
version = "0.9.14"

[[Hiccup]]
deps = ["MacroTools", "Test"]
git-tree-sha1 = "6187bb2d5fcbb2007c39e7ac53308b0d371124bd"
uuid = "9fb69e20-1954-56bb-a84f-559cc56a8ff7"
version = "0.2.2"

[[HypertextLiteral]]
git-tree-sha1 = "72053798e1be56026b81d4e2682dbe58922e5ec9"
uuid = "ac1192a8-f4b3-4bfe-ba22-af5b92cd3ab2"
version = "0.9.0"

[[IniFile]]
deps = ["Test"]
git-tree-sha1 = "098e4d2c533924c921f9f9847274f2ad89e018b8"
uuid = "83e8ac13-25f8-5344-8a64-a9f2b223428f"
version = "0.5.0"

[[InteractiveUtils]]
deps = ["Markdown"]
uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240"

[[InvertedIndices]]
git-tree-sha1 = "bee5f1ef5bf65df56bdd2e40447590b272a5471f"
uuid = "41ab1584-1d38-5bbf-9106-f11c6c58b48f"
version = "1.1.0"

[[IteratorInterfaceExtensions]]
git-tree-sha1 = "a3f24677c21f5bbe9d2a714f95dcd58337fb2856"
uuid = "82899510-4779-5014-852e-03e436cf321d"
version = "1.0.0"

[[JLLWrappers]]
deps = ["Preferences"]
git-tree-sha1 = "642a199af8b68253517b80bd3bfd17eb4e84df6e"
uuid = "692b3bcd-3c85-4b1f-b108-f13ce0eb3210"
version = "1.3.0"

[[JSExpr]]
deps = ["JSON", "MacroTools", "Observables", "WebIO"]
git-tree-sha1 = "bd6c034156b1e7295450a219c4340e32e50b08b1"
uuid = "97c1335a-c9c5-57fe-bc5d-ec35cebe8660"
version = "0.5.3"

[[JSON]]
deps = ["Dates", "Mmap", "Parsers", "Unicode"]
git-tree-sha1 = "8076680b162ada2a031f707ac7b4953e30667a37"
uuid = "682c06a0-de6a-54ab-a142-c8b1cf79cde6"
version = "0.21.2"

[[Kaleido_jll]]
deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"]
git-tree-sha1 = "2ef87eeaa28713cb010f9fb0be288b6c1a4ecd53"
uuid = "f7e6163d-2fa5-5f23-b69c-1db539e41963"
version = "0.1.0+0"

[[LaTeXStrings]]
git-tree-sha1 = "c7f1c695e06c01b95a67f0cd1d34994f3e7db104"
uuid = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f"
version = "1.2.1"

[[Lazy]]
deps = ["MacroTools"]
git-tree-sha1 = "1370f8202dac30758f3c345f9909b97f53d87d3f"
uuid = "50d2b5c4-7a5e-59d5-8109-a42b560f39c0"
version = "0.15.1"

[[LibCURL]]
deps = ["LibCURL_jll", "MozillaCACerts_jll"]
uuid = "b27032c2-a3e7-50c8-80cd-2d36dbcbfd21"

[[LibCURL_jll]]
deps = ["Artifacts", "LibSSH2_jll", "Libdl", "MbedTLS_jll", "Zlib_jll", "nghttp2_jll"]
uuid = "deac9b47-8bc7-5906-a0fe-35ac56dc84c0"

[[LibGit2]]
deps = ["Base64", "NetworkOptions", "Printf", "SHA"]
uuid = "76f85450-5226-5b5a-8eaa-529ad045b433"

[[LibSSH2_jll]]
deps = ["Artifacts", "Libdl", "MbedTLS_jll"]
uuid = "29816b5a-b9ab-546f-933c-edad1886dfa8"

[[Libdl]]
uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb"

[[LinearAlgebra]]
deps = ["Libdl"]
uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"

[[Logging]]
uuid = "56ddb016-857b-54e1-b83d-db4d58db5568"

[[MacroTools]]
deps = ["Markdown", "Random"]
git-tree-sha1 = "0fb723cd8c45858c22169b2e42269e53271a6df7"
uuid = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09"
version = "0.5.7"

[[Markdown]]
deps = ["Base64"]
uuid = "d6f4376e-aef5-505a-96c1-9c027394607a"

[[MbedTLS]]
deps = ["Dates", "MbedTLS_jll", "Random", "Sockets"]
git-tree-sha1 = "1c38e51c3d08ef2278062ebceade0e46cefc96fe"
uuid = "739be429-bea8-5141-9913-cc70e7f3736d"
version = "1.0.3"

[[MbedTLS_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "c8ffd9c3-330d-5841-b78e-0817d7145fa1"

[[Missings]]
deps = ["DataAPI"]
git-tree-sha1 = "2ca267b08821e86c5ef4376cffed98a46c2cb205"
uuid = "e1d29d7a-bbdc-5cf2-9ac0-f12de2c33e28"
version = "1.0.1"

[[Mmap]]
uuid = "a63ad114-7e13-5084-954f-fe012c677804"

[[MozillaCACerts_jll]]
uuid = "14a3606d-f60d-562e-9121-12d972cd8159"

[[Mustache]]
deps = ["Printf", "Tables"]
git-tree-sha1 = "36995ef0d532fe08119d70b2365b7b03d4e00f48"
uuid = "ffc61752-8dc7-55ee-8c37-f3e9cdd09e70"
version = "1.0.10"

[[Mux]]
deps = ["AssetRegistry", "Base64", "HTTP", "Hiccup", "Pkg", "Sockets", "WebSockets"]
git-tree-sha1 = "82dfb2cead9895e10ee1b0ca37a01088456c4364"
uuid = "a975b10e-0019-58db-a62f-e48ff68538c9"
version = "0.7.6"

[[NetworkOptions]]
uuid = "ca575930-c2e3-43a9-ace4-1e988b2c1908"

[[Observables]]
git-tree-sha1 = "3469ef96607a6b9a1e89e54e6f23401073ed3126"
uuid = "510215fc-4207-5dde-b226-833fc4488ee2"
version = "0.3.3"

[[OrderedCollections]]
git-tree-sha1 = "85f8e6578bf1f9ee0d11e7bb1b1456435479d47c"
uuid = "bac558e1-5e72-5ebc-8fee-abe8a469f55d"
version = "1.4.1"

[[Parameters]]
deps = ["OrderedCollections", "UnPack"]
git-tree-sha1 = "2276ac65f1e236e0a6ea70baff3f62ad4c625345"
uuid = "d96e819e-fc66-5662-9728-84c9c7592b0a"
version = "0.12.2"

[[Parsers]]
deps = ["Dates"]
git-tree-sha1 = "438d35d2d95ae2c5e8780b330592b6de8494e779"
uuid = "69de0a69-1ddd-5017-9359-2bf0b02dc9f0"
version = "2.0.3"

[[Pidfile]]
deps = ["FileWatching", "Test"]
git-tree-sha1 = "1be8660b2064893cd2dae4bd004b589278e4440d"
uuid = "fa939f87-e72e-5be4-a000-7fc836dbe307"
version = "1.2.0"

[[Pkg]]
deps = ["Artifacts", "Dates", "Downloads", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"]
uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"

[[PlotlyBase]]
deps = ["ColorSchemes", "Dates", "DelimitedFiles", "DocStringExtensions", "JSON", "LaTeXStrings", "Logging", "Parameters", "Pkg", "REPL", "Requires", "Statistics", "UUIDs"]
git-tree-sha1 = "7eb4ec38e1c4e00fea999256e9eb11ee7ede0c69"
uuid = "a03496cd-edff-5a9b-9e67-9cda94a718b5"
version = "0.8.16"

[[PlotlyJS]]
deps = ["Base64", "Blink", "DelimitedFiles", "JSExpr", "JSON", "Kaleido_jll", "Markdown", "Pkg", "PlotlyBase", "REPL", "Reexport", "Requires", "WebIO"]
git-tree-sha1 = "ec6bc7270269be2d565b272116ca74ca2f8bd9ab"
uuid = "f0f68f2c-4968-5e81-91da-67840de0976a"
version = "0.18.7"

[[PooledArrays]]
deps = ["DataAPI", "Future"]
git-tree-sha1 = "a193d6ad9c45ada72c14b731a318bedd3c2f00cf"
uuid = "2dfb63ee-cc39-5dd5-95bd-886bf059d720"
version = "1.3.0"

[[Preferences]]
deps = ["TOML"]
git-tree-sha1 = "00cfd92944ca9c760982747e9a1d0d5d86ab1e5a"
uuid = "21216c6a-2e73-6563-6e65-726566657250"
version = "1.2.2"

[[PrettyTables]]
deps = ["Crayons", "Formatting", "Markdown", "Reexport", "Tables"]
git-tree-sha1 = "0d1245a357cc61c8cd61934c07447aa569ff22e6"
uuid = "08abe8d2-0d0c-5749-adfa-8a2ac140af0d"
version = "1.1.0"

[[Printf]]
deps = ["Unicode"]
uuid = "de0858da-6303-5e67-8744-51eddeeeb8d7"

[[REPL]]
deps = ["InteractiveUtils", "Markdown", "Sockets", "Unicode"]
uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb"

[[Random]]
deps = ["Serialization"]
uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"

[[Reexport]]
git-tree-sha1 = "45e428421666073eab6f2da5c9d310d99bb12f9b"
uuid = "189a3867-3050-52da-a836-e630ba90ab69"
version = "1.2.2"

[[Requires]]
deps = ["UUIDs"]
git-tree-sha1 = "4036a3bd08ac7e968e27c203d45f5fff15020621"
uuid = "ae029012-a4dd-5104-9daa-d747884805df"
version = "1.1.3"

[[SHA]]
uuid = "ea8e919c-243c-51af-8825-aaa63cd721ce"

[[SentinelArrays]]
deps = ["Dates", "Random"]
git-tree-sha1 = "54f37736d8934a12a200edea2f9206b03bdf3159"
uuid = "91c51154-3ec4-41a3-a24f-3f23e20d615c"
version = "1.3.7"

[[Serialization]]
uuid = "9e88b42a-f829-5b0c-bbe9-9e923198166b"

[[SharedArrays]]
deps = ["Distributed", "Mmap", "Random", "Serialization"]
uuid = "1a1011a3-84de-559e-8e89-a11a2f7dc383"

[[Sockets]]
uuid = "6462fe0b-24de-5631-8697-dd941f90decc"

[[SortingAlgorithms]]
deps = ["DataStructures"]
git-tree-sha1 = "b3363d7460f7d098ca0912c69b082f75625d7508"
uuid = "a2af1166-a08f-5f64-846c-94a0d3cef48c"
version = "1.0.1"

[[SparseArrays]]
deps = ["LinearAlgebra", "Random"]
uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"

[[Statistics]]
deps = ["LinearAlgebra", "SparseArrays"]
uuid = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"

[[TOML]]
deps = ["Dates"]
uuid = "fa267f1f-6049-4f14-aa54-33bafae1ed76"

[[TableTraits]]
deps = ["IteratorInterfaceExtensions"]
git-tree-sha1 = "c06b2f539df1c6efa794486abfb6ed2022561a39"
uuid = "3783bdb8-4a98-5b6b-af9a-565f29a5fe9c"
version = "1.0.1"

[[Tables]]
deps = ["DataAPI", "DataValueInterfaces", "IteratorInterfaceExtensions", "LinearAlgebra", "TableTraits", "Test"]
git-tree-sha1 = "368d04a820fe069f9080ff1b432147a6203c3c89"
uuid = "bd369af6-aec1-5ad0-b16a-f7cc5008161c"
version = "1.5.1"

[[Tar]]
deps = ["ArgTools", "SHA"]
uuid = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e"

[[Test]]
deps = ["InteractiveUtils", "Logging", "Random", "Serialization"]
uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

[[TranscodingStreams]]
deps = ["Random", "Test"]
git-tree-sha1 = "216b95ea110b5972db65aa90f88d8d89dcb8851c"
uuid = "3bb67fe8-82b1-5028-8e26-92a6c54297fa"
version = "0.9.6"

[[URIParser]]
deps = ["Unicode"]
git-tree-sha1 = "53a9f49546b8d2dd2e688d216421d050c9a31d0d"
uuid = "30578b45-9adc-5946-b283-645ec420af67"
version = "0.4.1"

[[URIs]]
git-tree-sha1 = "97bbe755a53fe859669cd907f2d96aee8d2c1355"
uuid = "5c2747f8-b7ea-4ff2-ba2e-563bfd36b1d4"
version = "1.3.0"

[[UUIDs]]
deps = ["Random", "SHA"]
uuid = "cf7118a7-6976-5b1a-9a39-7adc72f591a4"

[[UnPack]]
git-tree-sha1 = "387c1f73762231e86e0c9c5443ce3b4a0a9a0c2b"
uuid = "3a884ed6-31ef-47d7-9d2a-63182c4928ed"
version = "1.0.2"

[[Unicode]]
uuid = "4ec0a83e-493e-50e2-b9ac-8f72acf5a8f5"

[[WeakRefStrings]]
deps = ["DataAPI", "Parsers"]
git-tree-sha1 = "4a4cfb1ae5f26202db4f0320ac9344b3372136b0"
uuid = "ea10d353-3f73-51f8-a26c-33c1cb351aa5"
version = "1.3.0"

[[WebIO]]
deps = ["AssetRegistry", "Base64", "Distributed", "FunctionalCollections", "JSON", "Logging", "Observables", "Pkg", "Random", "Requires", "Sockets", "UUIDs", "WebSockets", "Widgets"]
git-tree-sha1 = "adc25e225bc334c7df6eec3b39496edfc451cc38"
uuid = "0f1e0344-ec1d-5b48-a673-e5cf874b6c29"
version = "0.8.15"

[[WebSockets]]
deps = ["Base64", "Dates", "HTTP", "Logging", "Sockets"]
git-tree-sha1 = "f91a602e25fe6b89afc93cf02a4ae18ee9384ce3"
uuid = "104b5d7c-a370-577a-8038-80a2059c5097"
version = "1.5.9"

[[Widgets]]
deps = ["Colors", "Dates", "Observables", "OrderedCollections"]
git-tree-sha1 = "eae2fbbc34a79ffd57fb4c972b08ce50b8f6a00d"
uuid = "cc8bc4a8-27d6-5769-a93b-9d913e69aa62"
version = "0.6.3"

[[Zlib_jll]]
deps = ["Libdl"]
uuid = "83775a58-1f1d-513f-b197-d71354ab007a"

[[nghttp2_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d"

[[p7zip_jll]]
deps = ["Artifacts", "Libdl"]
uuid = "3f19e933-33d8-53b3-aaab-bd5110c3b7a0"
"""

# ╔═╡ Cell order:
# ╠═90facaa4-158c-11ec-3a0f-2b30580dc6b7
# ╠═40de9078-38b7-4238-b7a4-fcb1a440d753
# ╠═360b9e44-dff2-40c2-ab68-818f50447d76
# ╠═8659283f-0618-4af9-b78a-948cf8755a00
# ╠═c76c0762-da6b-4f3d-9d68-5bfdcf0be3fc
# ╠═6d525819-4724-42ee-8934-d1d66e7f9eb6
# ╠═c1ada4ac-6b51-4490-b5c0-64d724ae88a3
# ╠═ad5b05e7-b300-4de5-8abe-52a44ab01c70
# ╠═bef30ffd-d141-4520-9b34-5b93a83e0d1a
# ╟─00000000-0000-0000-0000-000000000001
# ╟─00000000-0000-0000-0000-000000000002


text:
function SingleLayerModelMappingFunction()
    return SingleLayerModelMappingFunction(())
end

function SingleLayerModelMappingFunction(arg0::jdouble)
    return SingleLayerModelMappingFunction((jdouble,), arg0)
end

function mapping_factor(obj::SingleLayerModelMappingFunction, arg0::RealFieldElement)
    return jcall(obj, "mappingFactor", RealFieldElement, (RealFieldElement,), arg0)
end

function mapping_factor(obj::SingleLayerModelMappingFunction, arg0::jdouble)
    return jcall(obj, "mappingFactor", jdouble, (jdouble,), arg0)
end