分数代码优化

mathe

为了在复数范围而不是有限域范围内验证果树问题，我现在实现了一个“短”有理数代码，对于计算溢出，通过C++异常来处理。
结果发现这个代码的速度还是比较让人满意的，用在果树问题上，速度大概只比有限域版本慢一倍多一点。
但是其中使用的分数类代码还不是很优化，特意贴出来让大家帮忙优化一下，最后能够让速度达到有限域版本:

1. 
   
2. long long gcd(long long a, long long b)
   
3. {
   
4.         long long t;
   
5.         while(b!=0){
   
6.                 a%=b;
   
7.                 t=b;
   
8.                 b=a;
   
9.                 a=t;
   
10.         }
    
11.         return a;
    
12. }
    
13. #define MINI (-0x7fffffff)
    
14. #define MAXI 0x7fffffff
    
15. bool in_range(long long& u, long long& d)
    
16. {
    
17.         if(MINI<=u&&u<=MAXI&&MINI<=d&&d<=MAXI)
    
18.                 return true;
    
19.         long long r=gcd(u,d);
    
20.         u/=r;d/=r;
    
21.         if(MINI<=u&&u<=MAXI&&MINI<=d&&d<=MAXI)
    
22.                 return true;
    
23.         return false;
    
24. }
    
25. class MyException{
    
26.         int v;
    
27. public:
    
28.         MyException(int u):v(u){}
    
29. };
    
30. 
    
31. class number{
    
32.     int up;
    
33.     int down;
    
34. public:
    
35.     number(const number& n):up(n.up),down(n.down){}
    
36.     number(int n=0):up(n),down(1){}
    
37.     number& operator+=(const number& n);
    
38.     number& operator-=(const number& n);
    
39.     number& operator*=(const number& n);
    
40.     number& operator/=(const number& n);
    
41.     bool is_zero()const{return up==0;}
    
42.     bool is_one()const{return up==down;}
    
43.         bool is_mone()const{return up==-down;}
    
44.     bool operator==(const number& n)const{return
    
45.         n.down*(long long)up==n.up*(long long)down;}
    
46.     bool operator<(const number& n)const;
    
47.     bool operator>(const number& n)const{return n<*this;}
    
48.     bool operator<=(const number& n)const{return !(*this>n);}
    
49.     bool operator!=(const number& n)const{return !(n==*this);}
    
50.     bool operator>=(const number& n)const{return !(*this<n);}
    
51.     number operator+(const number& n)const{number m(*this);return m+=n;}
    
52.     number operator-(const number& n)const{number m(*this);return m-=n;}
    
53.     number operator*(const number& n)const{number m(*this);return m*=n;}
    
54.     number operator/(const number& n)const{number m(*this);return m/=n;}
    
55.     bool is_integer()const{return up%down==0;}
    
56.     int get_integer()const{if(is_integer())return up/down;else throw  MyException(7);}
    
57.     number& operator=(int n){up=n;down=1;return *this;}
    
58.     int get_up()const{return up;}
    
59.     int get_down()const{return down;}
    
60.         number neg()const{number x;x.up=-up;x.down=down; return x;}
    
61.     number abs()const{number x;x.up=up>=0?up:-up;x.down=down>=0?down:-down;return x;}
    
62. };
    
63. 
    
64. number& number::operator +=(const number& n)
    
65. {
    
66.         long long u=(long long)up*n.down+(long long)down*n.up;
    
67.         long long d=(long long)down*n.down;
    
68.         if(!in_range(u,d))
    
69.                 throw  MyException(1);
    
70.     up=(int)u;
    
71.     down=(int)d;
    
72.     return *this;
    
73. }
    
74. 
    
75. number& number::operator -=(const number& n)
    
76. {
    
77.         long long u=(long long)up*n.down-(long long)down*n.up;
    
78.         long long d=(long long)down*n.down;
    
79.         if(!in_range(u,d))
    
80.                 throw  MyException(2);
    
81.     up=(int)u;
    
82.     down=(int)d;
    
83.     return *this;
    
84. 
    
85. }
    
86. 
    
87. number& number::operator *=(const number& n)
    
88. {
    
89.         long long u=(long long)up*n.up;
    
90.         long long d=(long long)down*n.down;
    
91.         if(!in_range(u,d))
    
92.                 throw  MyException(3);
    
93.     up=(int)u;
    
94.     down=(int)d;
    
95.     return *this;
    
96. }
    
97. 
    
98. number& number::operator /=(const number& n)
    
99. {
    
100.         long long u=(long long)up*n.down;
     
101.         long long d=(long long)down*n.up;
     
102.         if(!in_range(u,d))
     
103.                 throw  MyException(4);
     
104.         if(d==0)throw  MyException(5);
     
105.         if(d<0){
     
106.                 d=-d;u=-u;
     
107.         }
     
108.     up=(int)u;
     
109.     down=(int)d;
     
110.     return *this;
     
111. }
     
112. 
     
113. bool number::operator <(const number &n)const
     
114. {
     
115.     return (long long)up*n.down<n.up*(long long)down;
     
116. }
     

复制代码

wayne

想问一下mathe，你的那个 gcd函数里面的while如果改成这样，效率会提高吗

1. while(b!=0){
   
2.                  t=b;
   
3.                  b=a%b;
   
4.                  a=t;
   
5.          }
   

复制代码

mathe

这种变化应该没有区别。现在的编译器对于局部优化可以做得很不错的。

KeyTo9_Fans

所做的操作都是很必要的，一点优化的余地都没有了，怎么可能提速2倍以上？

如果硬要做点改动的话，可以在in_range函数里把这3行代码加在最前面：

1. 
   
2. if((u&15)==0&&(d&15)==0){u>>=4;d>>=4}
   
3. if((u&3)==0&&(d&3)==0){u>>=2;d>>=2}
   
4. if((u&1)==0&&(d&1)==0){u>>=1;d>>=1}
   

复制代码

添加理由：

1.有64%的可约分数含有公因子2，有52%的可约分数只需把因子2约去就不可约了（对于要进来检查的分数来说，这个比例可能更大）。
2.约去因子2只要用移位操作就可以完成。
3.移位操作花费的代价很小。
4.gcd操作花费的代价很大，应尽量少用。

上述改动实际上是试图通过盲目地增加一些移位操作来减少gcd的调用次数。

花费大量的移位操作的代价来换取一次gcd的调用是否值得？

如果这个交易是亏本的，那么速度反而会下降。

即使是盈利的，速度改善的程度也是很有限的，能提速10%就很不错了，不可能提速2倍以上。

mathe

呵呵，最终提高多少不重要，有多少是多少吧。
不过这里有很多汇编高手，也许他们可以帮不少忙的。
到星期一我可以提供一个完整的版本和测试数据。

gxqcn

分数的四则运算看似简单，实际上有许多技巧而变得复杂，尤其是如何通分、结合率的使用等。
举个简单的例子，计算：$1/2+1/6+1/12+1/20+1/30+1/42+1/56+1/72+1/90$

4# KeyTo9_Fans 建议 in_range 中采用移位运算，
但其变量都是有符号的，恐怕不大适合。
是不是可以采用新的数据结构，分别定义：符号、分子和分母（均为无符号型）

求最大公约数最好采用 Josef Stein 算法，以避免除法指令。

KeyTo9_Fans

回6#：

测了几个数据，似乎没有任何问题。

1. 
   
2. #include<stdio.h>
   
3. 
   
4. int i;
   
5. __int64 _i;
   
6. 
   
7. int main()
   
8. {
   
9.         //test1
   
10.         i=-128;
    
11.         printf("%d\n",i&15);                //output: 0
    
12.         printf("%d\n",i>>4);                //output: -8
    
13.         //test2
    
14.         i=-37;
    
15.         printf("%d\n",i&1);                        //output: 1
    
16.         printf("%d\n",i&3);                        //output: 3
    
17.         //test3
    
18.         i=-2147483648;
    
19.         printf("%d\n",i&65535);                //output: 0
    
20.         printf("%d\n",i>>16);                //output: -32768
    
21.         printf("%d\n",i>>31);                //output: -1
    
22.         //test4
    
23.         _i=-2222222222222222;
    
24.         printf("%I64d\n",_i&1);                //output: 0
    
25.         printf("%I64d\n",_i>>1);        //output: -1111111111111111
    
26.         printf("%I64d\n",_i&3);                //output: 2
    
27.         return 0;
    
28. }
    

复制代码

其中，-128的32个比特位为

11111111 11111111 11111111 10000000

后面连续出现7个0，和+128的性质相同：

10000000

所以&1、&3、&15的运算结果和+128是完全一样的。

在默认情况下，>>运算是算术右移，符号会保留。

所以((-128)>>4)的比特位为：

11111111 11111111 11111111 11111000

即 ((-128)>>4) = -8

与除以16的效果完全相同。

最特殊负数-2147483648的&运算与>>运算也没有问题。

所以那3行代码应该是没问题的。

gxqcn

楼上说的是对的。
我搞大数计算多了，对无符号类型变量有偏爱。

无心人

我觉得long long拖了后腿
用汇编吧

mathe

我也觉得long long是个问题。
不过用汇编也会有汇编的问题。
这两天我女儿生病了，现在我暂时没有时间提供测试版本了