如何求解无解析解的复方程根？附Python复函数求解场景

问题背景与现有实现 #

我定义了由以下函数构成的复函数，并通过连分式离散化构建了待求解方程：

基础函数定义 #

from sympy import *

def alpha(n, rho, l):
    return n**2 + (2*rho + 2)*n + 2*rho + 1

def beta(n, rho, l):
    return -(2*n**2 + (8*rho + 2)*n + 8*rho**2 + 4*rho + l*(l+1) -3)

def gamma(n, rho, l):
    return n**2 + 4*rho*n + 4*rho**2 - 4

连分式离散化 #

def ContinuedFrac(n, rho, l):
    cf = 0
    for i in reversed(range(1,n+1)):
        cf = cf + beta(i, rho, l)/alpha(i, rho, l)
        cf = (-gamma(i, rho, l)/alpha(i, rho, l)) / cf
    return cf

待求解方程构建 #

omega = symbols("ω")
n = 5
rho_guess = -I*omega
l=2

f = ContinuedFrac(n, rho_guess, l) + beta(0, rho_guess, l)/alpha(0, I*omega, l)

该方程无解析解，仅存在无穷多数值复根。目前我用lambdfy结合mpmath的findroot只能根据初始值求解单个根：

import mpmath as mp

omega = symbols("ω")

def findomega(n, rho_guess, l):
    return lambdify(omega, ContinuedFrac(n, rho_guess, l) + beta(0, rho_guess, l)/alpha(0, rho_guess , l))

f = findomega(15, -I*omega, 3)
mp.findroot(f, 1.16-0.56j)

# 输出: mpc(real='1.1652715765041771', imag='-0.56260091739349236')

我需要直接获取多个omega的复根列表，而非单个根。

批量获取复根的解决方案 #

由于方程无解析解，只能通过网格搜索+数值求解+去重的方式批量获取复根，具体实现如下：

完整代码示例 #

import mpmath as mp
from sympy import *

# 基础函数复用
def alpha(n, rho, l):
    return n**2 + (2*rho + 2)*n + 2*rho + 1

def beta(n, rho, l):
    return -(2*n**2 + (8*rho + 2)*n + 8*rho**2 + 4*rho + l*(l+1) -3)

def gamma(n, rho, l):
    return n**2 + 4*rho*n + 4*rho**2 - 4

def ContinuedFrac(n, rho, l):
    cf = 0
    for i in reversed(range(1,n+1)):
        cf = cf + beta(i, rho, l)/alpha(i, rho, l)
        cf = (-gamma(i, rho, l)/alpha(i, rho, l)) / cf
    return cf

# 生成可数值计算的函数
omega = symbols("ω")
def get_target_func(n, l):
    rho_guess = -I*omega
    expr = ContinuedFrac(n, rho_guess, l) + beta(0, rho_guess, l)/alpha(0, rho_guess, l)
    return lambdify(omega, expr, modules="mpmath")

# 参数设置
n_val = 15
l_val = 3
target_func = get_target_func(n_val, l_val)

# 定义复平面搜索范围（可根据实际根分布调整）
real_range = (0, 2)    # 实部范围
imag_range = (-1, 0)   # 虚部范围
grid_step = 0.2        # 网格步长，步长越小搜索越精细
tolerance = 1e-6       # 根的去重阈值

found_roots = []

# 遍历网格点作为初始值求解
for real_part in mp.arange(real_range[0], real_range[1], grid_step):
    for imag_part in mp.arange(imag_range[0], imag_range[1], grid_step):
        init_guess = mp.mpc(real_part, imag_part)
        try:
            # 求解根
            root = mp.findroot(target_func, init_guess)
            # 检查是否为已发现的近似根
            is_unique = True
            for existing_root in found_roots:
                if mp.norm(root - existing_root) < tolerance:
                    is_unique = False
                    break
            if is_unique:
                found_roots.append(root)
        except (mp.NoConvergence, ValueError):
            # 跳过不收敛的初始值
            continue

# 输出结果
print("找到的复根列表：")
for idx, root in enumerate(found_roots, 1):
    print(f"根{idx}: 实部={root.real:.6f}, 虚部={root.imag:.6f}")

关键说明 #

1. 网格搜索：在复平面的目标区域生成均匀网格点，每个点作为findroot的初始猜测值，覆盖可能存在根的区域
2. 根去重：不同初始值可能收敛到同一个根，通过计算根之间的欧氏距离判断是否为同一根，避免重复记录
3. 异常处理：部分初始值可能无法收敛（比如函数奇点附近），用try-except捕获异常并跳过
4. 参数调整：
   • 可根据已知根的分布缩小real_range和imag_range，提高搜索效率
   • 减小grid_step可发现更多根，但会增加计算时间
   • tolerance需匹配findroot的求解精度，避免误判重复根
5. 并行加速（可选）：如果搜索范围大，可使用multiprocessing模块并行处理每个网格点，大幅缩短计算时间

内容的提问来源于stack exchange，提问作者Projetos Programação