ln(1+x)

ln(1+x) = x*(ln(1+x)/x)
这个近似计算公式的原理是什么？
我们知道，如果x比较小，1+x在计算机中容易丢失精度，所以采用这个近似公式
可以近似成x*u(x); u(x)= ln(1+x)/x;只要能够精确计算u(x)就可以了
在文章中指出，只要找出一个x~ 约等于x，然后能够精确计算u(x~)就可以了
可以找到一个  
x~  = 1+x-1 

第一，u(x~)来近似u(x)这个是极限相关，x趋向于零时，u(x)趋向于1，也即是说u(x)的值在x趋向于零时，与x的大小基本无关

第二，估计是计算机计算1+x时不准确，所以需要找一个计算机能够准确计算尺来的x~来近似取代x
第三，x~=1+x-1得出来的值，再1+x~时可逆，也就是说1+x~的值是精缺，不会由于x~较小而导致误差ln(1+x)当 x 很小时直接计算会丢失精度，但是没关系。我们变通一下，我们直接计算出来的ln(1+x) 存在 y 变量中。
那么： y 约等于 ln(1+x)    
y / x 的结果与 ln(1+x) / x 的结果相差可能会很大。
  
我们可以根据 y ，凑个 x' 使得 y = ln(1+x') 这个式子很精确的成立。那么
y/x' 与 ln(1+x')/x' 相等。
  
而我们又知道 当 x 和 x' 都很小且相差不大时，
ln(1+x')/x' 与 ln(1+x) / x 的结果相差很小。

=====================================================

 
  double gsl_log1p (const double x)  
  {  
  volatile double y, z;  
  y = 1 + x;  
   z = y - 1;  
   return log(y) - (z - x) / y ;  /* cancels errors with IEEE arithmetic */  
}   

当x比较小时，1 + x = 1 + x1 + x2, 其中x2是 1 + x 时精度丢失部分。比如说，在MATLAB里，
  
>> x=pi/10^12
  
x =
  
     3.141592653589793e-012
  
>> y=1+x
  
y =
  
     1.000000000003142e+000
  
>> z=y-1
  
z =
  
     3.141487070479343e-012
  
y相当于 1 + x1, z相当于x1, x - z 相当于x2. 在 x2 = 0 处展开：
  
ln(a+x) = ln(a) + x/a - x^2/(2a^2) + O(x^3)
  
即:
  
ln(1+x1+x2) ≈ ln(1+x1) + x2/(1+x1) = ln(y) + (x-z)/y
  
回到MATLAB里那个例子:
  
>> log(y)
  
ans =
  
     3.141487070474409e-012
  
>> log(y)+(x-z)/y
  
ans =
  
     3.141592653584858e-012
  
而 http://www.wolframalpha.com/input/?i=ln%281%2Bpi%2F10^12%29 给出的结果大概是
  
3.14159265358485843626211E-12