浮点的数字分布
作者：concreteHAM
什么是浮点数，不用我多说，这里我们要讨论的是规格化的任意进制浮点数的前导数字的概率分布。
在《计算机程序设计艺术》第二卷中做了非常深入的讨论，这里我从中精炼出要点。
例如：
⒉345 E 67
这是一个十进制规格化浮点数，前导数字就是2。
就只有一个“随机”的浮点数而言，讨论其分布式没有意义的，我们要讨论的是充分多个“随机”数进行的一系列运算后产生的浮点结果的前导数字分布。
假设现在有一巨大的浮点数集，依此对数集中每个浮点数都乘以2，其中有一个十进制浮点数F，它的前导数字是1，那么它底数可能的值范围就是1.000…～1.999…，乘以一个数字2，那么它的底数就变成2.000…～3.999…，很明显乘以2以前前导数字是1的浮点个数与现在前导数字是2、3的浮点个数相同。以此我们接下来分析。
对于一个b进制的浮点数，它的前导数字x范围就是0 < x < b，设f(x）是上述数集的前导数字的概率密度函数（注：是密度函数），那么它在前导数字u和v之间（0<u<v<b）的概率就是：
∫[u,v]f(x)dx ⑴
由前面所述的，对于一个小增量Δx，f(x）必须满足这样一个公式：
f⑴Δx = x*f(x）Δx ⑵
很明显：
f(x) = f⑴/x ⑶
两边在[1,b]之间进行积分，等号左边必定为1，右边等于f⑴ln(b）：
1 = f⑴ln(b) ⑷
得：f⑴= 1/ln(b) 带入⑶中：
f(x) = 1/(x*ln(b))
那么利用⑴式得：
∫[u,v]1/(x*ln(b))dx
=ln(v/u) / ln(b)⑸
这就是求前导数字的概率分布函数。
例如b = 10进制时，前导数字为1的概率就是：
= ln（（1+1）/1） / ln⑽
≈ 0.301
前导数字为9的概率就是：
= ln（（9+1）/9） / ln⑽
≈0.0458
以下是一个测试程序（Mathematica软件）：
T[n_,b_]:=Block[{res={},ran,i,a},
For[i=1,i<b,i++;
res=Append[res,0]
];
For[i=0,i<n,i++;
ran=Random[]*Random[]*Random[]; 充分打乱模拟实际运算中的浮点数
ran=Log[b,ran];
a=Floor[b^(ran-Floor[ran])]; 取出前导数字
res[[a]]++ 对前导数字个数统计
];
Return[res]
]
执行T[100000,10]，以10进制测试100000个浮点数，得到一个分布：
{30149,18821,13317,9674,7688,6256,5306,4655,4134}
和理论值相当接近。
关于如何取出前导数字如下：
设原浮点数为a*10^e >= 0
其中a为底数范围（-10,-1]∪[1,10），e为指数，此时我们分离出底数部分和指数部分。
我们先对数求以10为底的对数从而分离出指数部分：
lg(a*10^e)=lg(a)+lg（10^e)=e+lg(a），因为a∈[1,10），所以lg(a）∈[0,1），而lg（10^e)=e是整数，所以lg(a*10^e）∈[e,e+1），因此我们可以通过向负无穷方向取整，也就是说取小于等于lg(a*10^e）的最大整数，lg(a）的值就等于lg(a*10^e）的值减去它取整后的数就可以了，10^(lg(a））也就是底数a，a像负无穷取整就是前导数字，如何得到e就不用多说了很简单。
