Avoid Stalled SymPy Plots: How to Handle High‑Degree Complex Roots with Exact Symbols, nroots, and Domain Coloring

Learn why SymPy plots hang when mixing math with symbolic ComplexRootOf results, and fix them using exact construction, Poly.nroots, nsolve, and domain coloring.

When you try to solve and plot complex roots of a high-degree SymPy expression, a seemingly innocuous call can stall the whole process. The usual culprit is mixing numeric modules with symbolic ones and then forcing heavy symbolic roots to evaluate numerically just before plotting. Below is a compact walkthrough of what goes wrong and how to fix it, with robust alternatives for both visualization and numeric root-finding.

Problem setup

The following snippet constructs a high-degree polynomial from parameters derived with math, requests all complex roots symbolically via solveset, and then attempts to extract real and imaginary parts for plotting. This can hang before any plot appears.

from __future__ import division
from sympy import *
import numpy as np
import math
import matplotlib.pyplot as mp
u1, u2, w, phi = symbols('u1 u2 w phi')
step = 0.04 * math.pi
lag = 2.0 * math.pi / step
print(lag)
root_set = solveset(w**math.floor(lag + 0.5) + 1 - w, w, domain=Complexes)
print(root_set)
ang = np.linspace(0, 2 * np.pi, 100)
ux = np.cos(ang)
uy = np.sin(ang)
re_vals = [rv.as_real_imag()[0] for rv in root_set]
im_vals = [rv.as_real_imag()[1] for rv in root_set]
mp.figure(figsize=(6, 6))
mp.scatter(re_vals, im_vals, color='green', marker='o', label='Complex Solutions')
mp.plot(ux, uy, 'b--', alpha=0.5, label='Unit Circle Reference')
mp.xlabel('Real Axis')
mp.ylabel('Imaginary Axis')
mp.title('Solutions on Unit Circle')
mp.gca().set_aspect('equal')
mp.grid(True)
mp.legend()
mp.show()

What actually goes wrong

The core expression is w**floor(lag + 1/2) + 1 - w, where floor(lag + 1/2) evaluates to 50. Calling solveset on this returns a set of 50 roots in the complex domain. Many of them are ComplexRootOf objects, which are exact but unevaluated algebraic numbers. Converting them to floating-point for plotting requires numerical evaluation. Even evaluating a single such root can take a very long time in this case. A quick experiment illustrates the bottleneck:

rt_list = list(root_set)
# evaluating even one entry can be prohibitively slow
# rt_list[10].n()

There is also a second source of friction. Intermixing math (which produces floats) with sympy (which can preserve exactness) pushes SymPy into numeric rather than exact arithmetic early, which is rarely what you want when building symbolic expressions.

Reliable fixes

The first fix is to build the expression using SymPy objects only. That keeps exact arithmetic and symbolic operations consistent and fast to manipulate. The second fix is to avoid symbolic root extraction when you only need a visualization or numerical approximations. Below are two practical paths.

Fast visualization without computing symbolic roots

Domain coloring gives an immediate picture of zeros and poles in the complex plane without solving anything symbolically. The following uses SymPy Plotting Backend to render the magnitude pattern and overlay the unit circle.

from sympy import *
from spb import *
xv, yv, w, ang = symbols('xv yv w ang')
# exact SymPy parameters
step_q = Rational(4, 100) * pi
lag_q = 2 * pi / step_q
pow_deg = floor(lag_q + Rational(1, 2))
poly_f = w**pow_deg + 1 - w
scl = 1.05
graphics(
    domain_coloring(poly_f, (w, -scl - scl*I, scl + scl*I), coloring='k', n=500),
    line_parametric_2d(cos(ang), sin(ang), (ang, 0, 2*pi), rendering_kw={'color': 'b'}, use_cm=False, label='Unit Circle', show_in_legend=False),
    grid=False, aspect='equal', size=(10, 10)
)

In this magnitude view, zeros appear as black points. For this expression they appear distributed around the unit circle, and the plot shows all of them instantly without requiring solveset.

Numerical roots when you do need the values

If you want actual numeric approximations of all roots, use a numeric solver with well-chosen initial guesses or leverage the polynomial API. Both approaches below are fast for this case.

The first approach seeds nsolve along the unit circle and collects roots:

from sympy import *
from spb import *
w, th = symbols('w th')
step_q = Rational(4, 100) * pi
lag_q = 2 * pi / step_q
pow_deg = floor(lag_q + Rational(1, 2))
poly_f = w**pow_deg + 1 - w
root_count = 50
offset = (2 * pi / root_count) / 2
found_roots = []
for k in range(root_count):
    theta = 2 * pi / root_count * k + offset
    guess = cos(theta) + I * sin(theta)
    try:
        sol = nsolve([poly_f], [w], [complex(guess)])
        found_roots.append(sol[0, 0])
    except:
        pass
scl = 1.05
th = symbols('th')
graphics(
    domain_coloring(poly_f, (w, -scl - scl*I, scl + scl*I), coloring='k', n=500),
    line_parametric_2d(cos(th), sin(th), (th, 0, 2*pi), rendering_kw={'color': 'b'}, use_cm=False, label='Unit Circle', show_in_legend=False),
    complex_points(*found_roots, rendering_kw={'marker': '.', 'color': 'r'}),
    grid=False, aspect='equal', xlabel='Real', ylabel='Imaginary'
)

The second approach uses a polynomial representation and computes numerical roots directly:

from sympy import *
w = symbols('w')
step_q = Rational(4, 100) * pi
lag_q = 2 * pi / step_q
pow_deg = floor(lag_q + Rational(1, 2))
poly_f = w**pow_deg + 1 - w
prec_digits = int(pow_deg)
P = Poly(poly_f, w)
roots_numeric = P.nroots(n=prec_digits)
print(roots_numeric)

Why this matters

When SymPy returns ComplexRootOf objects, turning them into floats can be extremely slow for certain high-degree expressions. Trying to plot those roots by extracting real and imaginary parts forces that heavy numeric step implicitly, which looks like a hang. Building expressions with SymPy primitives and switching to numeric strategies at the right time avoids that trap entirely. A practical note from the field is that evaluation can run very fast if the quadratic term has a coefficient near, but not exactly, one and the roots stay near the unit circle.

Takeaways

Keep symbolic and numeric worlds separate: construct expressions with SymPy’s Rational, pi, and floor rather than mixing math. Don’t insist on exact root expressions if all you need is a picture or floating-point values; domain coloring and numeric root-finding are the right tools here. If you need the set of all roots, seed nsolve around the unit circle or use Poly.nroots to obtain numerical approximations quickly. This way, you retain precision where it helps and switch to numerics where it pays off.