diff --git a/src/gcd.jl b/src/gcd.jl
index 6688059a..1d615355 100644
--- a/src/gcd.jl
+++ b/src/gcd.jl
@@ -79,9 +79,16 @@ end
 _coefficient_gcd(α, β) = gcd(α, β)
 _coefficient_gcd(α::AbstractFloat, β) = one(Base.promote_typeof(α, β))
 _coefficient_gcd(α, β::AbstractFloat) = one(Base.promote_typeof(α, β))
+_coefficient_gcd(α::Complex, β) = one(Base.promote_typeof(α, β))
+_coefficient_gcd(α, β::Complex) = one(Base.promote_typeof(α, β))
+_coefficient_gcd(α::Complex, β::AbstractFloat) = one(Base.promote_typeof(α, β))
+_coefficient_gcd(α::AbstractFloat, β::Complex) = one(Base.promote_typeof(α, β))
 function _coefficient_gcd(α::AbstractFloat, β::AbstractFloat)
     return one(Base.promote_typeof(α, β))
 end
+function _coefficient_gcd(α::Complex, β::Complex)
+    return one(Base.promote_typeof(α, β))
+end
 
 function Base.lcm(
     p::_APL,
@@ -306,11 +313,9 @@ function MA.operate(
     defl,
 )
     return polynomial(
-        map(terms(p)) do t
-            return term(
-                coefficient(t),
-                MA.operate(op, monomial(t), shift, defl),
-            )
+        coefficients(p),
+        map(monomials(p)) do m
+            return MA.operate(op, m, shift, defl)
         end,
     )
 end
diff --git a/test/commutative/gcd.jl b/test/commutative/gcd.jl
index c6731309..0b791639 100644
--- a/test/commutative/gcd.jl
+++ b/test/commutative/gcd.jl
@@ -291,4 +291,19 @@ end
             end
         end
     end
+    @testset "Non-concrete polynomial GCD" begin
+        Mod.@polyvar p q
+        p1 = MP.polynomial(p^2 + q^2, Number)
+        p2 = MP.polynomial(p * q, Number)
+        g = @inferred gcd(p1, p2)
+        @test isone(g)
+    end
+
+    @testset "Complex coefficient GCD" begin
+        Mod.@polyvar p q
+        p1 = MP.polynomial(p^2 + q^2, ComplexF64)
+        p2 = MP.polynomial(p * q, ComplexF64)
+        g = @inferred gcd(p1, p2)
+        @test isone(g)
+    end
 end