基于查找表的Arcsine CORDIC实现问题：输入-1时输出偏差过大

问题描述 #

以下代码尝试用CORDIC算法实现反正弦功能，但输入为-1时，预期输出应为-90度（即-π/2弧度），实际输出仅为-56.44度左右，偏差明显：

import math
import numpy as np
import matplotlib.pyplot as plt

# Define constants
K = 0.6072529350088813  # K value for CORDIC
NUM_ITERATIONS = 10     # Number of iterations for CORDIC
bitw = 15

# Function to compute arcsine using CORDIC
def arcsin_cordic(angle_ratio):
    x = K                                   # Initial value for x
    y = 0                         # Initial value for y
    z = 0                                   # Initial value for z (angle accumulator)
    inpLUK = np.zeros(bitw)
    for j in range(bitw):
        tmp = math.atan(2 ** -j)
        inpLUK[j] = tmp

    # CORDIC iterations
    for i in range(NUM_ITERATIONS):
        if y < angle_ratio:
            d = 1
        else:
            d = -1

        x_new = x - d * y * 2**(-i)
        y = y + d * x * 2**(-i)
        z = z + d * inpLUK[i]
        x = x_new

    return y

angle_array = np.arange(-1, 1, 0.15)
for angle_ratio in angle_array:
    arcsin_value_cordic = arcsin_cordic(angle_ratio)

错误原因分析 #

• 返回值完全错误：函数返回的是y，但y的作用是逼近输入的angle_ratio（即sin(z)），而真正的反正弦结果是累积角度z。
• 初始值设置错误：CORDIC的总缩放因子为K，若初始x=K，迭代后y的最大绝对值仅能达到K≈0.607，永远无法逼近输入的-1，导致角度累积不足。
• 迭代次数不足：仅10次迭代的精度不足以让结果收敛到正确值。

修正后的代码 #

import math
import numpy as np
import matplotlib.pyplot as plt

# Define constants
K = 0.6072529350088813  # CORDIC总缩放因子
NUM_ITERATIONS = 15     # 增加迭代次数提升精度

# Function to compute arcsine using CORDIC
def arcsin_cordic(angle_ratio):
    # 预计算各次迭代的旋转角度（弧度）
    inpLUK = [math.atan(2 ** -j) for j in range(NUM_ITERATIONS)]
    
    # 初始值：x设为1/K，抵消CORDIC的缩放效应，迭代后x≈cos(z), y≈sin(z)
    x = 1.0 / K
    y = 0.0
    z = 0.0  # 角度累积器，最终存储arcsin(angle_ratio)的弧度值
    
    for i in range(NUM_ITERATIONS):
        # 判断旋转方向：让y向目标angle_ratio逼近
        d = 1 if y < angle_ratio else -1
        
        # CORDIC迭代公式
        x_new = x - d * y * (2 ** (-i))
        y_new = y + d * x * (2 ** (-i))
        z_new = z + d * inpLUK[i]
        
        x, y, z = x_new, y_new, z_new
    
    # 返回弧度值，如需角度则乘以180/math.pi
    return z

# 测试输入-1
angle_ratio = -1
arcsin_rad = arcsin_cordic(angle_ratio)
arcsin_deg = arcsin_rad * 180 / math.pi
print(f"输入{angle_ratio}，输出弧度：{arcsin_rad:.4f}，输出角度：{arcsin_deg:.2f}°")

# 批量测试
angle_array = np.arange(-1, 1, 0.15)
for angle_ratio in angle_array:
    arcsin_rad = arcsin_cordic(angle_ratio)
    arcsin_deg = arcsin_rad * 180 / math.pi
    # 对比math库的结果
    math_deg = math.asin(angle_ratio) * 180 / math.pi
    print(f"输入{angle_ratio:.2f}，CORDIC结果：{arcsin_deg:.2f}°，Math库结果：{math_deg:.2f}°")

修正效果 #

• 输入-1时，输出角度接近-90°，与预期一致。
• 所有输入的结果与math.asin的偏差在可接受范围内，精度随迭代次数增加而提升。

内容的提问来源于stack exchange，提问作者Pournami Prasad