montgomery() 函數計算 t1^t2 % p 的值

 1 int montgomery(int t1, int t2, int p)
 2 {
 3     int ans = 1;
 4     while (t2) {
 5         if (t2 & 1)
 6             ans = ans * t1 % p;
 7         t2 >>= 1;
 8         t1 = t1 * t1 % p;
 9     }
10     return ans;
11 }
12 
13 int C(int dn, int up, int
 
 p)
14 {
15     if (dn == up)
16         return 1;
17     if (dn < up)
18         return 0;
19     if (up > up - dn)
20         up = up - dn;
21     int i, t1 = 1, t2 = 1;
22     for (i = 1; i <= up; ++i) {
23         t1 = t1 * i % p;
24         t2 = t2 * (dn - i + 1) % p;
25     }
26     return
 
 t2 * montgomery(t1, p - 2, p) % p;
27 }
28 
29 int lucas(int dn, int up, int p)
30 {
31     if (up == 0)
32         return 1;
33         return C(dn % p, up % p, p) * lucas(dn / p, up / p, p) % p;
34 }

樸素盧卡斯偽代碼