2个给定数字之间的分数密度
Per*_*dor 5 c++ algorithm math fractions
我正在尝试对一个简单的Fraction类做一些分析,我想要一些数据来比较那个类型doubles.
问题
对,知道我正在寻找一些好的方法来获得2个数字之间的分数密度.馏分是基本上2点的整数(例如pair< long, long>),以及之间的密度s和t在该范围内可表示数的量.并且它需要是在O(1)或非常快的情况下完成的精确或非常好的近似.
为了使它更简单,假设我想要s和t之间的所有数字(不是分数)a/b,其中0 <= s <= a/b <t <= M,0 <= a,b < = M(b> 0,a和b是整数)
例
如果我的分数是一个只计数为6(M = 6)的数据类型,并且我希望密度在0和1之间,那么答案将是12.这些数字是:
0, 1/6, 1/5, 1/4, 1/3, 2/5, 1/2, 3/5, 2/3, 3/4, 4/5, 5/6.
我已经想到了什么
一种非常天真的方法是循环所有可能的分数,并计算那些不能简化的分数.就像是:
long fractionsIn(double s, double t){
long density = 0;
long M = LONG_MAX;
for(int d = 1; d < floor(M/t); d++){
for(int n = ceil(d*s); n < M; n++){
if( gcd(n,d) == 1 )
density++;
}
}
return density;
}
但是gcd()很慢,所以它不起作用.我也尝试做一些数学但我无法做任何好事.
解
感谢@ m69的答案,我为此制作了以下代码Fraction = pair<Long,Long>:
//this should give the density of fractions between first and last, or less.
double fractionsIn(unsigned long long first, unsigned long long last){
double pi = 3.141592653589793238462643383279502884;
double max = LONG_MAX; //i can't use LONG_MAX directly
double zeroToOne = max/pi * max/pi * 3; // = approx. amount of numbers in Farey's secuence of order LONG_MAX.
double res = 0;
if(first == 0){
res = zeroToOne;
first++;
}
for(double i = first; i < last; i++){
res += zeroToOne/(i * i+1);
if(i == i+1)
i = nextafter(i+1, last); //if this happens, i might not count some fractions, but i have no other choice
}
return floor(res);
}
主要变化是nextafter,这对于大数字很重要(1e17)
结果
正如我在开始时解释,我是想比较Fractions有double.这是结果Fraction = pair<Long,Long>(这里我得到了双打密度):
Density between 0,1: | 1,2 | 1e6,1e6+1 | 1e14,1e14+1 | 1e15-1,1e15 | 1e17-10,1e17 | 1e19-10000,1e19 | 1e19-1000,1e19
Doubles: 4607182418800017408 | 4503599627370496 | 8589934592 | 64 | 8 | 1 | 5 | 0
Fraction: 2.58584e+37 | 1.29292e+37 | 2.58584e+25 | 2.58584e+09 | 2.58584e+07 | 2585 | 1 | 0
密度介于0和1之间
如果表示分数的整数在0~M范围内,那么值0(包括)和1(不包括)之间的分数密度为:
M: 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31
0~(1): 1 2 4 6 10 12 18 22 28 32 42 46 58 64 72 80 96 102 120 128 140 150 172 180 200 212 230 242 270 278 308 ...
这是OEIS上的序列A002088.如果向下滚动到公式部分,您将找到有关如何近似的信息,例如:
Φ(Ñ)=(3÷π 2)× Ñ 2 + O [ Ñ ×(LN Ñ)2/3 ×(LN LN Ñ)4/3 ]
(不幸的是,没有更详细地介绍O [x]部分中涉及的常数.请参阅下面关于近似质量的讨论.)
跨范围分布
从0到1的间隔包含可以用最多M的数字表示的唯一分数总数的一半; 例如,这是M = 15时的分布(即4位整数):
0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16
72 36 12 6 4 2 2 2 1 1 1 1 1 1 1 1
总共144个独特的分数.如果查看M的不同值的序列,您将看到此序列中的步骤收敛:
0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16
1: 1 1
2: 2 1 1
3: 4 2 1 1
4: 6 3 1 1 1
5: 10 5 2 1 1 1
6: 12 6 2 1 1 1 1
7: 18 9 3 2 1 1 1 1
8: 22 11 4 2 1 1 1 1 1
9: 28 14 5 2 2 1 1 1 1 1
10: 32 16 5 3 2 1 1 1 1 1 1
11: 42 21 7 4 2 2 1 1 1 1 1 1
12: 46 23 8 4 2 2 1 1 1 1 1 1 1
13: 58 29 10 5 3 2 2 1 1 1 1 1 1 1
14: 64 32 11 5 4 2 2 1 1 1 1 1 1 1 1
15: 72 36 12 6 4 2 2 2 1 1 1 1 1 1 1 1
密度不仅在分数总数的0到1的一半之间,而且1到2之间的密度是四分之一,2到3之间的密度接近第十二,依此类推.
随着M值的增加,0-1,1-2,2-3范围内的分数分布收敛于:
1/2, 1/4, 1/12, 1/24, 1/40, 1/60, 1/84, 1/112, 1/144, 1/180, 1/220, 1/264 ...
这个序列可以通过1/2开始计算,然后:
0-1: 1/2 x 1/1 = 1/2
1-2: 1/2 x 1/2 = 1/4
2-3: 1/4 x 1/3 = 1/12
3-4: 1/12 x 2/4 = 1/24
4-5: 1/24 x 3/5 = 1/40
5-6: 1/40 x 4/6 = 1/60
6-7: 1/60 x 5/7 = 1/84
7-8: 1/84 x 6/8 = 1/112
8-9: 1/112 x 7/9 = 1/144 ...
您当然可以直接计算这些值中的任何一个,而无需中间的步骤:
0-1: 1/2
6-7: 1/2 x 1/6 x 1/7 = 1/84
(另请注意,分布序列的后半部分由1组成;这些都是除以1的整数.)
在给定间隔内逼近密度
使用OEIS页面上提供的公式,您可以计算或近似区间0-1中的密度,并乘以2这是可以表示为分数的唯一值的总数.
给定两个值s和t,然后可以计算和求和间隔s~s + 1,s + 1~s + 2,... t-1~t中的密度,或者使用插值来获得更快但是近似值不太精确.
例
假设我们使用的是10位整数,能够表示0到1023之间的值.使用从OEIS页面链接的这个表,我们发现0~1之间的密度是318452,并且总分数是636904 .
如果我们想在间隔s~t = 100~105中找到密度:
100~101: 1/2 x 1/100 x 1/101 = 1/20200 ; 636904/20200 = 31.53
101~102: 1/2 x 1/101 x 1/102 = 1/20604 ; 636904/20604 = 30.91
102~103: 1/2 x 1/102 x 1/103 = 1/21012 ; 636904/21012 = 30.31
103~104: 1/2 x 1/103 x 1/104 = 1/21424 ; 636904/21424 = 29.73
104~105: 1/2 x 1/104 x 1/105 = 1/21840 ; 636904/21840 = 29.16
舍入这些值得出总和:
32 + 31 + 30 + 30 + 29 = 152
蛮力算法给出了这样的结果:
32 + 32 + 30 + 28 + 28 = 150
因此,对于这个低M值和小间隔只有5个值,我们的偏差为1.33%.如果我们在第一个和最后一个值之间使用了线性插值:
100~101: 31.53
104~105: 29.16
average: 30.345
total: 151.725 -> 152
我们已达到相同的价值.对于较大的区间,所有密度的总和可能会更接近实际值,因为舍入误差将相互抵消,但线性插值的结果可能会变得不那么准确.对于更大的M值,计算的密度应该与实际值收敛.
Φ(n)的逼近质量
使用这个简化的公式:
Φ(Ñ)=(3÷π 2)× Ñ 2
结果几乎总是小于实际值,但是对于n≥182,它们在1%以内,对于n≥1880,它们在0.1%之内,对于n≥19494,它们在0.01%之内.我建议硬编码较低范围(第一个这里可以找到50,000个值,然后使用从近似值足够好的点开始的简化公式.
这是一个简单的代码示例,其中Φ(n)的前182个值是硬编码的.分布序列的近似似乎增加了与Φ(n)的近似值相似的误差,因此应该可以得到合适的近似值.代码简单地遍历区间s~t中的每个整数并对分数求和.为了加速代码并仍然获得良好的结果,您应该计算区间中几个点的分数,然后使用某种非线性插值.
function fractions01(M) {
var phi = [0,1,2,4,6,10,12,18,22,28,32,42,46,58,64,72,80,96,102,120,128,140,150,172,180,200,212,230,242,270,278,308,
324,344,360,384,396,432,450,474,490,530,542,584,604,628,650,696,712,754,774,806,830,882,900,940,964,1000,
1028,1086,1102,1162,1192,1228,1260,1308,1328,1394,1426,1470,1494,1564,1588,1660,1696,1736,1772,1832,1856,
1934,1966,2020,2060,2142,2166,2230,2272,2328,2368,2456,2480,2552,2596,2656,2702,2774,2806,2902,2944,3004,
3044,3144,3176,3278,3326,3374,3426,3532,3568,3676,3716,3788,3836,3948,3984,4072,4128,4200,4258,4354,4386,
4496,4556,4636,4696,4796,4832,4958,5022,5106,5154,5284,5324,5432,5498,5570,5634,5770,5814,5952,6000,6092,
6162,6282,6330,6442,6514,6598,6670,6818,6858,7008,7080,7176,7236,7356,7404,7560,7638,7742,7806,7938,7992,
8154,8234,8314,8396,8562,8610,8766,8830,8938,9022,9194,9250,9370,9450,9566,9654,9832,9880,10060];
if (M < 182) return phi[M];
return Math.round(M * M * 0.30396355092701331433 + M / 4); // experimental; see below
}
function fractions(M, s, t) {
var half = fractions01(M);
var frac = (s == 0) ? half : 0;
for (var i = (s == 0) ? 1 : s; i < t && i <= M; i++) {
if (2 * i < M) {
var f = Math.round(half / (i * (i + 1)));
frac += (f < 2) ? 2 : f;
}
else ++frac;
}
return frac;
}
var M = 1023, s = 100, t = 105;
document.write(fractions(M, s, t));将Φ(n)的近似值与50,000个第一值的列表进行比较表明,添加M÷4是公式第二部分的可行替代品; 我没有对n的较大值进行测试,因此请谨慎使用.
蓝色:简化公式.红色:改进了简化公式.
分布近似的质量
将M = 1023的结果与蛮力算法的结果进行比较,误差实际上很小,从不超过-7或+6,而在205~206之间,它们被限制在-1~ + 1.然而,范围的大部分(57~1024)每个整数具有少于100个分数,并且在区间171~1024中,每个整数仅有10个分数或更少.这意味着-1或+1的小错误和舍入误差会对结果产生很大影响,例如:
interval: 241 ~ 250
fractions/integer: 6
approximation: 5
total: 50 (instead of 60)
为了改善每个整数的分数很少的区间的结果,我建议将上述方法与该范围的最后部分的单独方法相结合:
范围的最后部分的替代方法
如已经提到的,并且在代码示例中实现,范围的后半部分M÷2~M,每个整数具有1个分数.另外,间隔M÷3~M÷2有2; 区间M÷4~M÷3有4.这当然是Φ(n)序列:
M/2 ~ M : 1
M/3 ~ M/2: 2
M/4 ~ M/3: 4
M/5 ~ M/4: 6
M/6 ~ M/5: 10
M/7 ~ M/6: 12
M/8 ~ M/7: 18
M/9 ~ M/8: 22
M/10 ~ M/9: 28
M/11 ~ M/10: 32
M/12 ~ M/11: 42
M/13 ~ M/12: 46
M/14 ~ M/13: 58
M/15 ~ M/14: 64
M/16 ~ M/15: 72
M/17 ~ M/16: 80
M/18 ~ M/17: 96
M/19 ~ M/18: 102 ...
在这些间隔之间,一个整数可以具有不同数量的分数,具体取决于M的确切值,例如:
interval fractions
202 ~ 203 10
203 ~ 204 10
204 ~ 205 9
205 ~ 206 6
206 ~ 207 6
区间204~205位于区间之间的边缘,因为M÷5 = 204.6; 它有6 + 3 = 9个分数,因为M模5是3.如果M是1022或1024而不是1023,它将有8或10个分数.(这个例子很简单,因为5是素数;见下文.)
同样,我建议使用Φ(n)的硬编码值来计算范围最后部分的分数.如果使用上面列出的前17个值,则这将覆盖范围的每个整数少于100个分数的部分,这样可以减少舍入误差低于1%的影响.前56个值将给出0.1%,前182个值为0.01%.
与Φ(n)的值一起,您可以对每个模数值的边缘间隔的分数进行硬编码,例如:
modulo: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
M/ 2 1 2
M/ 3 2 3 4
M/ 4 4 5 5 6
M/ 5 6 7 8 9 10
M/ 6 10 11 11 11 11 12
M/ 7 12 13 14 15 16 17 18
M/ 8 18 19 19 20 20 21 21 22
M/ 9 22 23 24 24 25 26 26 27 28
M/10 28 29 29 30 30 30 30 31 31 32
M/11 32 33 34 35 36 37 38 39 40 41 42
M/12 42 43 43 43 43 44 44 45 45 45 45 46
M/13 46 47 48 49 50 51 52 53 54 55 56 57 58
M/14 58 59 59 60 60 61 61 61 61 62 62 63 63 64
M/15 64 65 66 66 67 67 67 68 69 69 69 70 70 71 72
M/16 72 73 73 74 74 75 75 76 76 77 77 78 78 79 79 80
M/17 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96
M/18 96 97 97 97 97 98 98 99 99 99 99 100 100 101 101 101 101 102
| 归档时间： |
|
| 查看次数： |
303 次 |
| 最近记录： |