diff --git a/src/MutableArithmetics.jl b/src/MutableArithmetics.jl
index c92c60b..b880b7d 100644
--- a/src/MutableArithmetics.jl
+++ b/src/MutableArithmetics.jl
@@ -7,6 +7,7 @@
 module MutableArithmetics
 
 include("interface.jl")
+include("shortcuts.jl")
 
 # Test that can be used to test an implementation of the interface
 include("Test/Test.jl")
diff --git a/src/Test/int.jl b/src/Test/int.jl
index c7a8deb..89e0752 100644
--- a/src/Test/int.jl
+++ b/src/Test/int.jl
@@ -1,7 +1,6 @@
 function int_add_test(::Type{T}) where T
     @testset "add_to! / add!" begin
-        @test MA.mutability(T, MA.add_to!, T, T) isa MA.IsMutable
-        @test MA.mutability(T, MA.add!, T) isa MA.IsMutable
+        @test MA.mutability(T, +, T, T) isa MA.IsMutable
 
         t(n) = convert(T, n)
         a = t(5)
@@ -21,8 +20,7 @@ function int_add_test(::Type{T}) where T
 end
 function int_mul_test(::Type{T}) where T
     @testset "mul_to! / mul!" begin
-        @test MA.mutability(T, MA.mul_to!, T, T) isa MA.IsMutable
-        @test MA.mutability(T, MA.mul!, T) isa MA.IsMutable
+        @test MA.mutability(T, *, T, T) isa MA.IsMutable
 
         t(n) = convert(T, n)
         a = t(5)
@@ -39,12 +37,11 @@ function int_mul_test(::Type{T}) where T
         @test a == t(420)
     end
 end
-function int_muladd_test(::Type{T}) where T
-    @testset "muladd_to! / muladd! / muladd_buf_to! /muladd_buf!" begin
-        @test MA.mutability(T, MA.muladd_to!, T, T, T) isa MA.IsMutable
-        @test MA.mutability(T, MA.muladd!, T, T) isa MA.IsMutable
-        @test MA.mutability(T, MA.muladd_buf_to!, T, T, T, T) isa MA.IsMutable
-        @test MA.mutability(T, MA.muladd_buf!, T, T, T) isa MA.IsMutable
+function int_add_mul_test(::Type{T}) where T
+    @testset "add_mul_to! / add_mul! / add_mul_buf_to! / add_mul_buf!" begin
+        @test MA.mutability(T, MA.add_mul, T, T) isa MA.IsMutable
+        @test MA.mutability(T, MA.add_mul, T, T, T) isa MA.IsMutable
+        @test MA.mutability(T, MA.add_mul, T, T, T, T) isa MA.IsMutable
 
         t(n) = convert(T, n)
         a = t(5)
@@ -53,16 +50,16 @@ function int_muladd_test(::Type{T}) where T
         d = t(20)
         buf = t(24)
 
-        @test MA.muladd_to!(a, b, c, d) == t(69)
+        @test MA.add_mul_to!(a, b, c, d) == t(69)
         @test a == t(69)
         a = t(5)
-        @test MA.muladd!(b, c, d) == t(69)
+        @test MA.add_mul!(b, c, d) == t(69)
         @test b == t(69)
         b = t(9)
 
-        @test MA.muladd_buf_to!(buf, a, b, c, d) == t(69)
+        @test MA.add_mul_buf_to!(buf, a, b, c, d) == t(69)
         @test a == t(69)
-        @test MA.muladd_buf!(buf, b, c, d) == t(69)
+        @test MA.add_mul_buf!(buf, b, c, d) == t(69)
         @test b == t(69)
 
         a = t(148)
@@ -70,19 +67,19 @@ function int_muladd_test(::Type{T}) where T
         c = t(17)
         d = t(42)
         buf = t(56)
-        @test MA.muladd!(a, b, c) == t(420)
+        @test MA.add_mul!(a, b, c) == t(420)
         @test a == t(420)
         a = t(148)
-        @test MA.muladd_buf_to!(buf, d, a, b, c) == t(420)
+        @test MA.add_mul_buf_to!(buf, d, a, b, c) == t(420)
         @test d == t(420)
-        @test MA.muladd_buf!(buf, a, b, c) == t(420)
+        @test MA.add_mul_buf!(buf, a, b, c) == t(420)
         @test a == t(420)
     end
 end
 
 function int_zero_test(::Type{T}) where T
     @testset "zero!" begin
-        @test MA.mutability(T, MA.zero!) isa MA.IsMutable
+        @test MA.mutability(T, zero, T) isa MA.IsMutable
 
         t(n) = convert(T, n)
         a = t(5)
@@ -94,7 +91,7 @@ end
 
 function int_one_test(::Type{T}) where T
     @testset "one!" begin
-        @test MA.mutability(T, MA.one!) isa MA.IsMutable
+        @test MA.mutability(T, one, T) isa MA.IsMutable
 
         t(n) = convert(T, n)
         a = t(5)
@@ -107,7 +104,7 @@ end
 const int_tests = Dict(
     "int_add" => int_add_test,
     "int_mul" => int_mul_test,
-    "int_muladd" => int_muladd_test,
+    "int_add_mul" => int_add_mul_test,
     "int_zero" => int_zero_test,
     "int_one" => int_one_test
 )
diff --git a/src/bigint.jl b/src/bigint.jl
index 22da995..b5aff91 100644
--- a/src/bigint.jl
+++ b/src/bigint.jl
@@ -1,27 +1,31 @@
-# This example contains a full implementation of the MutableArithmetics API for the BigInt datatype.
+mutability(::Type{BigInt}) = IsMutable()
 
-using MutableArithmetics
-const MA = MutableArithmetics
+# zero
+promote_operation(::typeof(zero), ::Type{BigInt}) = BigInt
+mutable_operate!(::typeof(zero), x::BigInt) = Base.GMP.MPZ.set_si!(x, 0)
 
-# zero!
-MA.mutability(::Type{BigInt}, ::typeof(MA.zero!)) = MA.IsMutable()
-MA.zero_impl!(x::BigInt) = Base.GMP.MPZ.set_si!(x, 0)
+# one
+promote_operation(::typeof(one), ::Type{BigInt}) = BigInt
+mutable_operate!(::typeof(one), x::BigInt) = Base.GMP.MPZ.set_si!(x, 1)
 
-# one!
-MA.mutability(::Type{BigInt}, ::typeof(MA.one!)) = MA.IsMutable()
-MA.one_impl!(x::BigInt) = Base.GMP.MPZ.set_si!(x, 1)
+# +
+promote_operation(::typeof(+), ::Type{BigInt}...) = BigInt
+function mutable_operate_to!(output::BigInt, ::typeof(+), a::BigInt, b::BigInt)
+    return Base.GMP.MPZ.add!(output, a, b)
+end
 
-# add_to! / add!
-MA.mutability(::Type{BigInt}, ::typeof(MA.add_to!), ::Type{BigInt}, ::Type{BigInt}) = MA.IsMutable()
-MA.add_to_impl!(x::BigInt, a::BigInt, b::BigInt) = Base.GMP.MPZ.add!(x, a, b)
-MA.add_impl!(a::BigInt, b::BigInt) = Base.GMP.MPZ.add!(a, a, b)
+# *
+promote_operation(::typeof(*), ::Type{BigInt}...) = BigInt
+function mutable_operate_to!(output::BigInt, ::typeof(*), a::BigInt, b::BigInt)
+    return Base.GMP.MPZ.mul!(output, a, b)
+end
 
-# mul_to! / mul!
-MA.mutability(::Type{BigInt}, ::typeof(MA.mul_to!), ::Type{BigInt}, ::Type{BigInt}) = MA.IsMutable()
-MA.mul_to_impl!(x::BigInt, a::BigInt, b::BigInt) = Base.GMP.MPZ.mul!(x, a, b)
-MA.mul_impl!(a::BigInt, b::BigInt) = Base.GMP.MPZ.mul!(a, a, b)
-
-# muladd_to! / muladd! / muladd_buf!
-MA.muladd_to_impl!(dest::BigInt, b::BigInt, c::BigInt, d::BigInt) = Base.GMP.MPZ.add!(dest, Base.GMP.MPZ.mul!(BigInt(), c, d), b)
-MA.muladd_impl!(a::BigInt, b::BigInt, c::BigInt) = Base.GMP.MPZ.add!(a, a, Base.GMP.MPZ.mul!(BigInt(), b, c))
-MA.muladd_buf_impl!(buf::BigInt, a::BigInt, b::BigInt, c::BigInt) = Base.GMP.MPZ.add!(a, Base.GMP.MPZ.mul!(buf, b, c))
+# add_mul
+function mutable_operate_to!(output::BigInt, ::typeof(add_mul), args::BigInt...)
+    return mutable_buffered_operate_to!(BigInt(), output, add_mul, args...)
+end
+function mutable_buffered_operate_to!(buffer::BigInt, output::BigInt, ::typeof(add_mul),
+                                      a::BigInt, args::BigInt...)
+    mutable_operate_to!(buffer, *, args...)
+    return mutable_operate_to!(output, +, a, buffer)
+end
diff --git a/src/interface.jl b/src/interface.jl
index 50c929d..d51a7fe 100644
--- a/src/interface.jl
+++ b/src/interface.jl
@@ -1,93 +1,15 @@
-# These functions suffixed by `_impl!` fail if the first argument cannot be modified to be the result
-
 # Example of mutable types that can implement this API: BigInt, Array, JuMP.AffExpr, MultivariatePolynomials.AbstractPolynomial
-# `..._impl!` functions are similar to `JuMP.add_to_expression`
-# `...!` functions are similar to `JuMP.destructive_add` and `MOI.Utilities.operate!`
-
-"""
-    add_to_impl!(a, b, c)
-
-Write the result of the sum of `b` and `c` to `a`.
-"""
-function add_to_impl! end
-
-"""
-    add_impl!(a, b)
-
-Write the result of the sum of `a` and `b` to `a`.
-"""
-function add_impl! end
-# Fallback
-add_impl!(a, b) = add_to_impl!(a, a, b)
-
-"""
-    mul_to_impl!(a, b, c)
-
-Write the result of the product of `b` and `c` to `a`.
-"""
-function mul_to_impl! end
-
-"""
-    mul_impl!(a, b)
-
-Write the result of the product of `a` and `b` to `a`.
-"""
-function mul_impl! end
-# Fallback
-mul_impl!(a, b) = mul_to_impl!(a, a, b)
-
-"""
-    muladd_to_impl!(a, b, c, d)
-
-Write the result of `muladd(c, d, b)` to `a`.
-"""
-function muladd_to_impl! end
-
-"""
-    muladd_buf_to_impl!(buf, a, b, c, d)
-
-Write the result of `muladd(c, d, b)` to `a`, possibly modifying `buf`.
-"""
-function muladd_buf_to_impl! end
-# Fallback
-function muladd_buf_to_impl!(buf, a, b, c, d)
-    mul_to_impl!(buf, c, d)
-    return add_to_impl!(a, b, buf)
-end
-
-"""
-    muladd_impl!(a, b, c)
-
-Write the result of `muladd(b, c, a)` to `a`.
-"""
-function muladd_impl! end
-# Fallback
-muladd_impl!(a, b, c) = muladd_to_impl!(a, a, b, c)
-
-"""
-    muladd_buf_impl!(buf, a, b, c)
+# `mutable_operate!(add_mul, ...)` is similar to `JuMP.add_to_expression(...)`
+# `operate!(add_mul, ...)` is similar to `JuMP.destructive_add(...)`
+# `operate!` is similar to `MOI.Utilities.operate!`
 
-Write the result of `muladd(b, c, a)` to `a`, possibly modifying `buf`.
 """
-function muladd_buf_impl! end
-# No fallback as it depends on the type of arguments whether we should
-# redirect to `muladd_buf_to_impl!` or `muladd_impl!`
+    promote_operation(op::Function, ArgsTypes::Type...)
 
-
-"""
-    zero_impl!(a)
-
-Write the result of `zero(a)` to `a`.
-"""
-function zero_impl! end
-
-"""
-    one_impl!(a)
-
-Write the result of `one(a)` to `a`.
+Returns the type returned to the call `operate(op, args...)` where the types of
+the arguments `args` are `ArgsTypes`.
 """
-function one_impl! end
-
+function promote_operation end
 
 # Define Traits
 abstract type MutableTrait end
@@ -97,150 +19,129 @@ struct NotMutable <: MutableTrait end
 """
     mutability(T::Type, ::typeof(op), args::Type...)::MutableTrait
 
-Return `IsMutable` to indicate that `op(a::T, ::args[1], ...)` returns `a`.
-That is, the result of the operation is stored in `a` and then `a` is returned.
-Equivalently, returns whether `op_impl` is supported.
+Return `IsMutable` to indicate an object of type `T` can be modified to be
+equal to `op(args...)`.
 """
-mutability(::Type, op, args::Type...) = NotMutable()
+function mutability(T::Type, op, args::Type...)
+    if mutability(T) isa IsMutable && promote_operation(op, args...) == T
+        return IsMutable()
+    else
+        return NotMutable()
+    end
+end
 mutability(x, op, args...) = mutability(typeof(x), op, typeof.(args)...)
+mutability(::Type) = NotMutable()
 
-"""
-    add_to!(a, b, c)
-
-Return the sum of `b` and `c`, possibly modifying `a`.
-"""
-function add_to! end
-function add_to!(a, b, c)
-    add_to!(a, b, c, mutability(a, add_to!, b, c))
+function mutable_operate_to_fallback(::NotMutable, output, op::Function, args...)
+    throw(ArgumentError("Cannot call `mutable_operate_to!($output, $op, $(args...))` as `$output` cannot be modifed to equal the result of the operation. Use `operate!` or `operate_to!` instead which returns the value of the result (possibly modifying the first argument) to write generic code that also works when the type cannot be modified."))
 end
-# generic fallbacks
-add_to!(a, b, c, ::NotMutable) = b + c
-add_to!(a, b, c, ::IsMutable) = add_to_impl!(a, b, c)
 
+function mutable_operate_to_fallback(::IsMutable, op::Function, args...)
+    error("`mutable_operate_to!($op, $(args...))` is not implemented yet.")
+end
 
 """
-    add!(a, b)
+    mutable_operate_to!(output, op::Function, args...)
 
-Return the sum of `a` and `b`, possibly modifying `a`.
+Modify the value of `output` to be equal to the value of `op(args...)`. Can
+only be called if `mutability(output, op, args...)` returns `true`.
 """
-function add! end
-add!(a, b) = add!(a, b, mutability(a, add!, b))
-# generic fallbacks
-add!(a, b, ::NotMutable) = a + b
-add!(a, b, ::IsMutable) = add_impl!(a, b)
-mutability(T::Type, ::typeof(add!), U::Type) = mutability(T, add_to!, T, U)
+function mutable_operate_to!(output, op::Function, args...)
+    mutable_operate_fallback(mutability(output, op, args...), output, op, args...)
+end
 
 """
-    mul_to!(a, b, c)
+    mutable_operate!(op::Function, args...)
 
-Return the product of `b` and `c`, possibly modifying `a`.
+Modify the value of `args[1]` to be equal to the value of `op(args...)`. Can
+only be called if `mutability(args[1], op, args...)` returns `true`.
 """
-function mul_to! end
-function mul_to!(a, b, c)
-    mul_to!(a, b, c, mutability(a, mul_to!, b, c))
+function mutable_operate!(op::Function, args...)
+    mutable_operate_to!(args[1], op, args...)
 end
-# generic fallbacks
-mul_to!(a, b, c, ::NotMutable) = b * c
-mul_to!(a, b, c, ::IsMutable) = mul_to_impl!(a, b, c)
 
 """
-    mul!(a, b)
+    mutable_buffered_operate_to!(buffer, output, op::Function, args...)
 
-Return the product of `a` and `b`, possibly modifying `a`.
+Modify the value of `output` to be equal to the value of `op(args...)`,
+possibly modifying `buffer`. Can only be called if
+`mutability(output, op, args...)` returns `true`.
 """
-function mul! end
-mul!(a, b) = mul!(a, b, mutability(a, mul!, b))
-# generic fallbacks
-mul!(a, b, ::NotMutable) = a * b
-mul!(a, b, ::IsMutable) = mul_impl!(a, b)
-mutability(T::Type, ::typeof(mul!), U::Type) = mutability(T, mul_to!, T, U)
+function mutable_buffered_operate_to! end
 
 """
-    muladd_to!(a, b, c, d)
+    mutable_buffered_operate!(buffer, op::Function, args...)
 
-Return `muladd(c, d, b)`, possibly modifying `a`.
+Modify the value of `args[1]` to be equal to the value of `op(args...)`,
+possibly modifying `buffer`. Can only be called if
+`mutability(args[1], op, args...)` returns `true`.
 """
-function muladd_to! end
-function muladd_to!(a, b, c, d)
-    muladd_to!(a, b, c, d, mutability(a, muladd_to!, b, c, d))
-end
-# generic fallbacks
-muladd_to!(a, b, c, d, ::NotMutable) = muladd(c, d, b)
-muladd_to!(a, b, c, d, ::IsMutable) = muladd_to_impl!(a, b, c, d)
-function mutability(S::Type, ::typeof(muladd_to!), T::Type, U::Type, V::Type)
-    return mutability(S, add_to!, T, typeof(zero(U) * zero(V)))
+function mutable_buffered_operate!(buffer, op::Function, args...)
+    mutable_buffered_operate_to!(buffer, args[1], op, args...)
 end
 
 """
-    muladd!(a, b, c)
+    operate_to!(output, op::Function, args...)
 
-Return `muladd(b, c, a)`, possibly modifying `a`.
+Returns the value of `op(args...)`, possibly modifying `output`.
 """
-function muladd! end
-function muladd!(a, b, c)
-    muladd!(a, b, c, mutability(a, muladd!, b, c))
-end
-# generic fallbacks
-muladd!(a, b, c, ::NotMutable) = muladd(b, c, a)
-muladd!(a, b, c, ::IsMutable) = muladd_impl!(a, b, c)
-function mutability(S::Type, ::typeof(muladd!), T::Type, U::Type)
-    return mutability(S, add!, typeof(zero(T) * zero(U)))
+function operate_to!(output, op::Function, args...)
+    return operate_to_fallback!(mutability(output, op, args...), output, op, args...)
 end
 
+function operate_to_fallback!(::NotMutable, output, op::Function, args...)
+    return op(args...)
+end
+function operate_to_fallback!(::IsMutable, output, op::Function, args...)
+    return mutable_operate_to!(output, op, args...)
+end
 
 """
-    muladd_buf_to!(buf, a, b, c, d)
+    operate!(op::Function, args...)
 
-Return `muladd(c, d, b)`, possibly modifying `a` and `buf`.
+Returns the value of `op(args...)`, possibly modifying `args[1]`.
 """
-function muladd_buf_to! end
-function muladd_buf_to!(buf, a, b, c, d)
-    muladd_buf_to!(buf, a, b, c, d, mutability(a, muladd_to!, b, c, d))
+function operate!(op::Function, args...)
+    return operate_fallback!(mutability(args[1], op, args...), op, args...)
 end
-# generic fallbacks
-muladd_buf_to!(buf, a, b, c, d, ::NotMutable) = muladd(c, d, b)
-muladd_buf_to!(buf, a, b, c, d, ::IsMutable) = muladd_buf_to_impl!(buf, a, b, c, d)
-function mutability(S::Type, ::typeof(muladd_buf_to!), T::Type, U::Type, V::Type, W::Type)
-    if mutability(S, add_to!, U, typeof(zero(V) * zero(W))) isa IsMutable && mutability(T, add_to!, U, typeof(zero(V) * zero(W))) isa IsMutable
-        return IsMutable()
-    end
-    return NotMutable()
+
+function operate_fallback!(::NotMutable, op::Function, args...)
+    return op(args...)
+end
+function operate_fallback!(::IsMutable, op::Function, args...)
+    return mutable_operate!(op, args...)
 end
 
 """
-    muladd_buf!(buf, a, b, c)
+    buffered_operate_to!(buffer, output, op::Function, args...)
 
-Return `muladd(b, c, a)`, possibly modifying `a` and `buf`.
+Returns the value of `op(args...)`, possibly modifying `buffer` and `output`.
 """
-function muladd_buf! end
-function muladd_buf!(buf, a, b, c)
-    muladd_buf!(buf, a, b, c, mutability(a, muladd!, b, c))
-end
-# generic fallbacks
-muladd_buf!(buf, a, b, c, ::NotMutable) = muladd(b, c, a)
-muladd_buf!(buf, a, b, c, ::IsMutable) = muladd_buf_impl!(buf, a, b, c)
-function mutability(S::Type, ::typeof(muladd_buf!), T::Type, U::Type, V::Type)
-    return mutability(S, add_to!, T, typeof(zero(U) * zero(V)))
+function buffered_operate_to!(buffer, output, op::Function, args...)
+    return buffered_operate_to_fallback!(mutability(output, op, args...),
+                                         buffer, output, op, args...)
 end
 
-"""
-    zero!(a)
+function buffered_operate_to_fallback!(::NotMutable, buffer, output, op::Function, args...)
+    return op(args...)
+end
+function buffered_operate_to_fallback!(::IsMutable, buffer, output, op::Function, args...)
+    return mutable_buffered_operate_to!(buffer, output, op, args...)
+end
 
-Return `zero(a)`, possibly modifying `a`.
 """
-function zero! end
-zero!(x) = zero!(x, mutability(x, zero!))
-# fallback
-zero!(x, ::NotMutable) = zero(x)
-zero!(x, ::IsMutable) = zero_impl!(x)
+    buffered_operate!(buffer, op::Function, args...)
 
+Returns the value of `op(args...)`, possibly modifying `buffer`.
 """
-    one!(a)
+function buffered_operate!(buffer, op::Function, args...)
+    return buffered_operate_fallback!(mutability(args[1], op, args...),
+                                      buffer, op, args...)
+end
 
-Return `one(a)`, possibly modifying `a`.
-"""
-function one! end
-one!(x) = one!(x, mutability(x, one!))
-# fallback
-one!(x, ::NotMutable) = one(x)
-one!(x, ::IsMutable) = one_impl!(x)
+function buffered_operate_fallback!(::NotMutable, buffer, op::Function, args...)
+    return op(args...)
+end
+function buffered_operate_fallback!(::IsMutable, buffer, op::Function, args...)
+    return mutable_buffered_operate!(buffer, op, args...)
+end
diff --git a/src/linear_algebra.jl b/src/linear_algebra.jl
index 2d568c6..17945f8 100644
--- a/src/linear_algebra.jl
+++ b/src/linear_algebra.jl
@@ -1,11 +1,11 @@
 import LinearAlgebra
 
-function mutability(::Type{<:Vector}, ::typeof(mul_to!),
-                    ::Type{<:AbstractVecOrMat}, ::Type{<:AbstractVector})
-    return IsMutable() # Assume the element type of the first vector is of correct type which is the case if it is called from `mul`
+mutability(::Type{<:Vector}) = IsMutable()
+function promote_operation(::typeof(*), ::Type{<:AbstractMatrix{S}}, ::Type{<:AbstractVector{T}}) where {S, T}
+    return Vector{Base.promote_op(LinearAlgebra.matprod, S, T)}
 end
-function mul_to_impl!(C::Vector, A::AbstractVecOrMat, B::AbstractVector)
-    if mutability(eltype(C), muladd!, eltype(A), eltype(B)) isa NotMutable
+function mutable_operate_to!(C::Vector, ::typeof(*), A::AbstractMatrix, B::AbstractVector)
+    if mutability(eltype(C), add_mul, eltype(C), eltype(A), eltype(B)) isa NotMutable
         return LinearAlgebra.mul!(C, A, B)
     end
     # If `mutability(S, muladd!, T, U)` is `NotMutable`, we might as well redirect to `LinearAlgebra.mul!(C, A, B)`
@@ -33,19 +33,19 @@ function mul_to_impl!(C::Vector, A::AbstractVecOrMat, B::AbstractVector)
         b = B[k]
         for i = 1:mA
             # `C[i] = muladd_buf!(mul_buffer, C[i], A[aoffs + i], b)`
-            muladd_buf_impl!(mul_buffer, C[i], A[aoffs + i], b)
+            mutable_buffered_operate!(mul_buffer, add_mul, C[i], A[aoffs + i], b)
         end
     end
     end # @inbounds
     return C
 end
 
-function mul(A::AbstractVecOrMat{T}, B::AbstractVector{S}) where {T, S}
-    TS = Base.promote_op(LinearAlgebra.matprod, T, S)
-    C = similar(B, TS, axes(A,1))
+function mul(A::AbstractMatrix{S}, B::AbstractVector{T}) where {T, S}
+    U = Base.promote_op(LinearAlgebra.matprod, S, T)
+    C = similar(B, U, axes(A, 1))
     # C now contains only undefined values, we need to fill this with actual zeros
     for i in eachindex(C)
-        @inbounds C[i] = zero(TS)
+        @inbounds C[i] = zero(U)
     end
     return mul_to!(C, A, B)
 end
diff --git a/src/shortcuts.jl b/src/shortcuts.jl
new file mode 100644
index 0000000..e3bc0ba
--- /dev/null
+++ b/src/shortcuts.jl
@@ -0,0 +1,85 @@
+"""
+    add_to!(a, b, c)
+
+Return the sum of `b` and `c`, possibly modifying `a`.
+"""
+add_to!(output, args...) = operate_to!(output, +, args...)
+
+"""
+    add!(a, b, ...)
+
+Return the sum of `a`, `b`, ..., possibly modifying `a`.
+"""
+add!(args...) = operate!(+, args...)
+
+"""
+    mul_to!(a, b, c)
+
+Return the product of `b` and `c`, possibly modifying `a`.
+"""
+mul_to!(output, args...) = operate_to!(output, *, args...)
+
+"""
+    mul!(a, b, ...)
+
+Return the product of `a`, `b`, ..., possibly modifying `a`.
+"""
+mul!(args...) = operate!(*, args...)
+
+"""
+    add_mul(a, args...)
+
+Return `a + *(args...)`. Note that `add_mul(a, b, c) = muladd(b, c, a)`.
+"""
+function add_mul end
+add_mul(a, b, c) = muladd(b, c, a)
+
+function promote_operation(::typeof(add_mul), T::Type, args::Type...)
+    return promote_operation(+, T, promote_operation(*, args...))
+end
+
+"""
+    add_mul_to!(output, args...)
+
+Return `add_mul(args...)`, possibly modifying `output`.
+"""
+add_mul_to!(output, args...) = operate_to!(output, add_mul, args...)
+
+"""
+    add_mul!(args...)
+
+Return `add_mul(args...)`, possibly modifying `args[1]`.
+"""
+add_mul!(args...) = operate!(add_mul, args...)
+
+"""
+    add_mul_buf_to!(buffer, output, args...)
+
+Return `add_mul(args...)`, possibly modifying `output` and `buffer`.
+"""
+function add_mul_buf_to!(buffer, output, args...)
+    buffered_operate_to!(buffer, output, add_mul, args...)
+end
+
+"""
+    add_mul_buf!(buffer, args...)
+
+Return `add_mul(args...)`, possibly modifying `args[1]` and `buffer`.
+"""
+function add_mul_buf!(buffer, args...)
+    buffered_operate!(buffer, add_mul, args...)
+end
+
+"""
+    zero!(a)
+
+Return `zero(a)`, possibly modifying `a`.
+"""
+zero!(a) = operate!(zero, a)
+
+"""
+    one!(a)
+
+Return `one(a)`, possibly modifying `a`.
+"""
+one!(a) = operate!(one, a)
diff --git a/test/int.jl b/test/int.jl
index 2384e88..1062c89 100644
--- a/test/int.jl
+++ b/test/int.jl
@@ -24,35 +24,35 @@ end
     @test MA.mul!(a, b) == 420
 end
 
-@testset "muladd_to! / muladd! / muladd_buf_to! /muladd_buf!" begin
-    @test MA.mutability(Int, MA.muladd_to!, Int, Int, Int) isa MA.NotMutable
-    @test MA.mutability(Int, MA.muladd!, Int, Int) isa MA.NotMutable
-    @test MA.mutability(Int, MA.muladd_buf_to!, Int, Int, Int, Int) isa MA.NotMutable
-    @test MA.mutability(Int, MA.muladd_buf!, Int, Int, Int) isa MA.NotMutable
+@testset "add_mul_to! / add_mul! / add_mul_buf_to! /add_mul_buf!" begin
+    @test MA.mutability(Int, MA.add_mul_to!, Int, Int, Int) isa MA.NotMutable
+    @test MA.mutability(Int, MA.add_mul!, Int, Int) isa MA.NotMutable
+    @test MA.mutability(Int, MA.add_mul_buf_to!, Int, Int, Int, Int) isa MA.NotMutable
+    @test MA.mutability(Int, MA.add_mul_buf!, Int, Int, Int) isa MA.NotMutable
     a = 5
     b = 9
     c = 3
     d = 20
     buf = 24
-    @test MA.muladd_to!(a, b, c, d) == 69
-    @test MA.muladd!(b, c, d) == 69
-    @test MA.muladd_buf_to!(buf, a, b, c, d) == 69
-    @test MA.muladd_buf!(buf, b, c, d) == 69
+    @test MA.add_mul_to!(a, b, c, d) == 69
+    @test MA.add_mul!(b, c, d) == 69
+    @test MA.add_mul_buf_to!(buf, a, b, c, d) == 69
+    @test MA.add_mul_buf!(buf, b, c, d) == 69
     a = 148
     b = 16
     c = 17
-    @test MA.muladd!(a, b, c) == 420
+    @test MA.add_mul!(a, b, c) == 420
     a = 148
     b = 16
     c = 17
     buf = 56
-    @test MA.muladd_buf!(buf, a, b, c) == 420
+    @test MA.add_mul_buf!(buf, a, b, c) == 420
     a = 148
     b = 16
     c = 17
     d = 42
     buf = 56
-    @test MA.muladd_buf_to!(buf, d, a, b, c) == 420
+    @test MA.add_mul_buf_to!(buf, d, a, b, c) == 420
 end
 
 @testset "zero!" begin
diff --git a/test/matmul.jl b/test/matmul.jl
index 4a9adb1..2e5f5ee 100644
--- a/test/matmul.jl
+++ b/test/matmul.jl
@@ -7,19 +7,11 @@
 		@test MA.mul(A, x) == BigInt[3; 3; 3]
 		@test MA.mul_to!(y, A, x) == BigInt[3; 3; 3] && y == BigInt[3; 3; 3]
 
-		MA.mutability(::Type{BigInt}, ::typeof(MA.zero!)) = MA.IsMutable()
-		MA.zero_impl!(x::BigInt) = Base.GMP.MPZ.set_si!(x, 0)
-		MA.mutability(::Type{BigInt}, ::typeof(MA.mul_to!), ::Type{BigInt}, ::Type{BigInt}) = MA.IsMutable()
-		MA.mul_to_impl!(x::BigInt, a::BigInt, b::BigInt) = Base.GMP.MPZ.mul!(x, a, b)
-		MA.mutability(::Type{BigInt}, ::typeof(MA.add_to!), ::Type{BigInt}, ::Type{BigInt}) = MA.IsMutable()
-		MA.add_to_impl!(x::BigInt, a::BigInt, b::BigInt) = Base.GMP.MPZ.add!(x, a, b)
-		MA.muladd_buf_impl!(buf::BigInt, a::BigInt, b::BigInt, c::BigInt) = Base.GMP.MPZ.add!(a, Base.GMP.MPZ.mul!(buf, b, c))
-
 		A = BigInt[1 1 1; 1 1 1; 1 1 1]
 		x = BigInt[1; 1; 1]
 		y = BigInt[0; 0; 0]
 
-		@test MA.mutability(y, MA.mul_to!, A, x) isa MA.IsMutable
+		@test MA.mutability(y, *, A, x) isa MA.IsMutable
 		@test MA.mul(A, x) == BigInt[3; 3; 3]
 		@test MA.mul_to!(y, A, x) == BigInt[3; 3; 3] && y == BigInt[3; 3; 3]
 		@test_throws DimensionMismatch MA.mul(BigInt[1 1; 1 1], BigInt[])