在Ada中的二次方程

Question

我刚刚开始学习Ada并尝试编写一些代码。不过缺点是，它的语法和函数与C++有很大区别，因此我必须强行将各种东西塞进去才能让这个程序正常运行。

我的问题是，是否有更好的方法来执行我在这里所做的计算？

   IF(B < 0.0) THEN
      B := ABS(B);
      X1 := (B / 2.0) + Sqrt( (B / 2.0) ** 2.0 + ABS(C));
      X2 := (B / 2.0) - Sqrt( (B / 2.0) ** 2.0 + ABS(C));
   ELSE
      X1 := -(B / 2.0) + Sqrt( (B / 2.0) ** 2.0 - C);
      X2 := -(B / 2.0) - Sqrt( (B / 2.0) ** 2.0 - C);
   END IF;

我遇到了一些关于负数的问题，所以我使用了IF语句并使用ABS（）将其变成正数。但奇怪的是，这对于另一种情况完美地起作用，这很奇怪...

- starcorn

---

关于前两行代码 - 如果你已经知道B是负数，最好避免使用abs()函数。可以使用B:=-B来代替。即使编译器很聪明并且可以内联处理这些内容。 - DarenW

Answer 1

解决二次方程并不像大多数人想象的那样简单。

解决 a x^2 + b x + c = 0 的标准公式为：

delta = b^2 - 4 a c
x1 = (-b + sqrt(delta)) / (2 a)   (*)
x2 = (-b - sqrt(delta)) / (2 a)

但是当 4 a c << b^2 时，计算 x1 需要减去接近的数字，这会导致精度丢失。因此，您可以改用以下方法：

delta as above
x1 = 2 c / (-b - sqrt(delta))     (**)
x2 = 2 c / (-b + sqrt(delta))

这个方法可以得到更好的x1，但是x2有和之前x1一样的问题。

因此，正确计算根的方式是

q = -0.5 (b + sign(b) sqrt(delta))

使用x1 = q / a和x2 = c / q非常高效。如果您想处理delta为负数或系数为复数的情况，则必须使用复数算术（这也很棘手）。

编辑：使用Ada代码：

DELTA := B * B - 4.0 * A * C;

IF(B > 0.0) THEN
    Q := -0.5 * (B + SQRT(DELTA));
ELSE
    Q := -0.5 * (B - SQRT(DELTA));
END IF;

X1 := Q / A;
X2 := C / Q;

- Alexandre C.

---

很好的解释！这是一个很好的例子，说明标准教科书公式可能不是编写代码的最佳指南。 - DarenW

风格：保留字通常为小写，标识符为大写，'if'语句中不需要括号，例如if B > 0.0 then和Delta := etc。 - debater

问题：我们在这里使用浮点数吗？如果是这样，那么接近值之间的减法是否真的会有精度问题？ - debater

@辩手：是的。1.0000001 - 1 = 0.0000001，你已经失去了8个有效数字。如果你还不相信，请尝试使用上述代码和高中公式，其中a=1，b=10000000，c=1。 - Alexandre C.

Answer 2

给定ax² + bx + c = 0，二次方程公式可以给出解为x = (-b +/- sqrt(b²-4ac) ) / 2a。判别式d = b²-4ac为正数时有实根，为负数时有非零虚部的根（即非实复数），当根是二重根时，判别式为0。

因此，这个Ada代码将是：

D := B ** 2.0 - 4.0 * A * C;
IF D >= 0.0 THEN
  X1 := (-B + Sqrt(D)) / (2.0 * A);
  X2 := (-B - Sqrt(D)) / (2.0 * A);
ELSE
  -- Deal with the fact that the result is a non-real complex number.
END IF;

注意：我对Ada有点生疏，但这应该接近正确的语法。

- andand

---

关闭；所有数字都需要是实数（第一行为4.0，第二行为0.0，第三行和第四行为2.0）。此外，在条件表达式周围不必加括号；if D <= 0.0 then。 - Simon Wright

Answer 3

二次公式为x = ( -b +/- sqrt ( b ** 2 - 4*a*c ) ) / ( 2 * a )

我猜测a的值为1。

所以，x = -( b/2 ) +/- sqrt ( ( ( b ** 2 ) / 4 ) - c )

计算d = ( b ** 2 ) * 0.25 - c然后检查它的符号。

如果d的符号为负，则有复数根；按您希望处理它们。

如果b恰好为负，用+ abs( c )替换- c会得到垃圾结果。

通常乘以0.5或0.25比除以2.0或4.0更好。

- Pete Kirkham

Answer 4

虽然我不懂Ada，但我看到以下可以优化的地方：

1. 在IF指令的第一个分支中，您已经知道B是负数。因此，您可以使用B := -B而不是B := ABS(B)。或者更好的方法是：在第一个分支中使用-B代替B。
2. 您正在四次使用子表达式B/2.0。将其分配给辅助变量B_2（或者如果您不想再使用另一个变量，则将其再次分配给B）可能更有效（并且更清晰）。
   还有，sqrt被计算了两次。将其分配给辅助变量可以节省运行时间（并使读者明确地知道完全相同的子表达式被使用了两次）。
3. 然后，使用B_2*B_2可能比使用**2.0更快；如果Ada中有专用的平方函数，则更好。

- Curd

Answer 5

对我来说，这个问题更多是与数值算法而非Ada语言有关。 对于数值计算，通常必须（如果不总是）参考参考/学术论文。

这种问题总是让我想起这个： https://en.wikipedia.org/wiki/Fast_inverse_square_root 只有在“做数学运算”或找到一些解决你问题的论文时，才能找到以下技巧。

float Q_rsqrt( float number )
{
    long i;
    float x2, y;
    const float threehalfs = 1.5F;

    x2 = number * 0.5F;
    y  = number;
    i  = * ( long * ) &y;                       // evil floating point bit level hacking
    i  = 0x5f3759df - ( i >> 1 );               // what the...? 
    y  = * ( float * ) &i;
    y  = y * ( threehalfs - ( x2 * y * y ) );   // 1st iteration
//  y  = y * ( threehalfs - ( x2 * y * y ) );   // 2nd iteration, this can be removed

    return y;
}

附注：正如维基百科文章所指出的那样，这种实现方式对于大多数平台来说可能已经过时了。

- LoneWanderer