diff --git a/src/types.jl b/src/types.jl
index f94477e6..ec2b7f11 100644
--- a/src/types.jl
+++ b/src/types.jl
@@ -164,7 +164,8 @@ end
 """Gives the full harmonic ansatz used to construct `eom`."""
 function _show_ansatz(eom::HarmonicEquation)
     output = ""
-    for nat_var in get_variables(eom.natural_equation)
+    vars = unique(getfield.(eom.variables, :natural_variable))
+    for nat_var in vars
         # the Hopf variable (limit cycle frequency) does not contribute a term
         harm_vars = filter(
             x -> isequal(nat_var, x.natural_variable) && x.type !== "Hopf", eom.variables
diff --git a/test/krylov.jl b/test/krylov.jl
index 406907c6..4f1cae94 100644
--- a/test/krylov.jl
+++ b/test/krylov.jl
@@ -19,6 +19,11 @@ add_harmonic!(diff_eom, x, ω) # x will rotate at ω
 harmonic_eq1 = get_krylov_equations(diff_eom; order=1)
 harmonic_eq2 = get_krylov_equations(diff_eom; order=2)
 
+@testset "show method" begin
+    print_variable = HB._show_ansatz(harmonic_eq1)
+    @test !(occursin("xˍt(t)", print_variable))
+end
+
 fixed = (ω0 => 1.0, γ => 0.005, α => 1.0, η => 0, F => 0.0, ψ => 0.0, θ => 0.0)
 varied = (ω => range(0.99, 1.01, 5), λ => range(1e-6, 0.05, 5))