给定三角形边长平方，如何精确计算内角？数值精度问题求解

已知三角形的边长平方为a, b, c > 0，需要精确计算其内角。现有两种参考方法：Kahan基于未开平方边长的常用方法，以及依赖边向量的主流方法。目前尝试用余弦定理实现计算：

float angle = acos( (b+a-c) / (2.*sqrt(b*a)) );

但部分非退化输入会返回0.0，存在精度问题。

测试案例说明 #

输入由边向量的平方长度生成，测试案例如下：

// 三角形顶点为(0,0) (1,1) (1+e,1)
float e = 1e-7;
float a = ((1+e)-1)*((1+e)-1);
float b = (1+e)*(1+e) + 1*1;
float c = 1*1 + 1*1;

该三角形的真实角度约为：

4.9999997529193436e-08 2.3561944901923448 0.7853981133974508

不同方法的计算结果对比 #

• 直接对a,b,c开平方后应用Kahan方法，得到结果：

  0            3.14159  0
  

• 应用余弦定理计算，得到结果：

  0            2.35619  0.785398
  

  这个结果相对更优，但原本大于0的角度被计算为0。
• 若将a,b,c按非2的幂次缩放，余弦定理的计算结果偏差极大：

  0            2.02856  1.11303
  

待解决疑问 #

1. 是否存在基于边长平方直接计算内角的更精准数值方法？
2. 是不是float精度不足以满足这类计算需求？
3. 会不会是float类型存储的边长平方违反了三角不等式？

测试环境为macOS下执行clang++ main.cpp -g && ./a.out，对应实现代码为测试用C++代码。

内容的提问来源于stack exchange，提问作者Alec Jacobson