<style> .reveal .slides { text-align: left; font-size:30px; } </style> # 快速冪 && 質數篩 ---- ## 快速冪 * 快速計算 $x^y$ ---- 原本我們在計算 $x^y$ 的時候呢，有一個最直觀且簡短的寫法，就是跑一個時間複雜度為 $O(y)$ 的迴圈 ```cpp= int ans=1; for(int i=1;i<=y;i++){ ans*=x; } cout<<ans<<"\n"; ``` 但是如果今天 $y$ 太大的時候，就會需要優化這個程式，需要讓他更快速的解決這個問題，於是就衍生出了快速冪的這個算法。 ---- 先把 y 換成二進位 以 $2^{13}$ 為例 $13$ $= 1101$ $= 1\cdot 2^3 + 1\cdot 2^2 + 0\cdot 2^1 + 1\cdot 2^0$ $= 8 + 4 + 1$ $\to 2^{13} = 2^8 * 2^4 * 2^1$ ---- 因此，我們只要求出 $2^8, 2^4, 2^1$ 就好 也就是只需要求出每個目標的二的冪次在做相加，也就是 $\log N$ 次 ---- ## 實作 一樣以 $2^{13}$ 為例 $1101$ 就從最右邊的 bit 開始 判斷最右邊的 bit 是不是 $1$ 如果是 $1$ 則把答案乘上去 接著每次把 $x$ 平方 接著把整個右移一個，變成 $110$ 因此一樣變成判斷最右邊的 bit 即可 $再右移\to 11$ $再右移\to 1$ 一直做到 y 右移到 0 為止 ---- ## 過程 | x | y | 原本答案 | 計算後答案 | | -------- | -------- | -------- | -------- | | $2^1$ | 1101 | $2^0$ | $2^1$ | | $2^2$ | 110 | $2^1$ | $2^1$ | | $2^4$ | 11 | $2^1$ | $2^5$ | | $2^8$ | 1 | $2^5$ | $2^{13}$ | ---- ## 程式碼 ```cpp= long long mypow(long long x,long long y){ long long ans = 1; while(y){ if( y&1 ) ans = ans * x; x = x * x; //每次把自己平方 y >>= 1; //每次右移一格 } return ans; } ``` ---- ## 溢位處理 如果數字太大的話，很有可能會超過 long long 的範圍 因此許多題目會要我們把結果 mod 一個數字 ---- ## 加入 mod 後的程式碼 ```cpp= const ll MOD = 1e9 + 7; long long mypow(long long x,long long y){ long long ans = 1; while(y){ if( y&1 ) ans = (ans * x) % MOD; x = (x * x) % MOD; //每次把自己平方 y >>= 1; //每次右移一格 } return ans; } ``` ---- ## mod與乘法的交換律 為什麼我們在實作的過程中可以邊 mod 邊乘? 在 $(A \times B) \% m$ 的情況下 $A = (C_1 \times m) + R_1$, $B = (C_2 \times m) + R_2$ $(A \times B) \% m = ((C_1 \times m) \times (C_2 \times m)) +((C_1 \times m) \times R_2) +((C_2 \times m) \times R_1)$ $+ (R_1 \times R_2) \% m$ 我們可以發現上方的一行皆可以被 $m$ 整除，因此我們得到這樣的結論 $(A \times B) \% m = ((A \% m) \times (B \% m)) \% m$ ---- 換成 2 進位，只需要計算 $\log y$ 次 當要計算 $x^y$ ，$y$ 到 $10^{18}$ 的時候 $\log 10^{18} = 60$ 只需做 60 次就可以算出 $x 的 10^{18}$ 次方 --- ## 質數篩法 ---- 對於範圍內的數字使用快速的方法判斷它是否為質數。 ---- ## 方法一 窮舉 詢問數字 $N$ 是否為質數? 窮舉每個小於 $N$ 的數字，判斷該數字是否整除 $N$ 若整除代表枚舉的該數字是 $N$ 的因數 以 $7$ 為例子，我們要跑一個迴圈從 $2$ 跑到 $6$ 並且用每一個數字去除除看。 最差情況 : 判斷的數字是質數，因此要跑 $2\sim N-1$ 總共 $N-2$ 個數字 (除了 $1$ 跟 $N$ ) ---- 簡單的快速判斷，暴力跑過每個數字 ```cpp int n; cin>>n; for(int i=2;i<n;i++){ if(n%i==0){ //若整除則此數字絕對不是質數 cout<<"NO\n"; break; //跳出迴圈 } if(i==n-1){ //若到n-1(也就是迴圈的最後一次) 都沒有跑出去代表他是質數 cout<<"YES\n"; break; } } ``` ---- ## 針對方法一的優化 可以再檢查更少次?降低該方法的時間複雜度? 由於因數性質的關係，所以對於一個數字的檢查，其實只需要從 $2$ 跑到 $\sqrt{N}$ 即可。 所以只需要檢查區間 $[2,\sqrt{N}]$ 的數字就好 那我們的寫法會是從 $2$ 一路跑到 $\sqrt N$ ，而其中由於小數點常常造成不可預期的錯誤，因此我們的慣用寫法會把根號寫在另外一邊，改成使用平方的。 ---- 使用暴力判斷判斷此數字 $n$ 是否為質數 ```cpp void isPrime(int n){ bool flag=1; for(int i=2;i*i<=n;i++){ //檢查範圍只需要到根號n if(n%i==0){ //若整除則此數字絕對不是質數 flag=0; break; } } if(flag) cout<<"YES\n"; else cout<<"NO\n"; } ``` ---- ## 小技巧 在迴圈內請用 $i\times i \le n$ 而非 $i < \sqrt{n}+1$ $\sqrt{n}$ 容易有浮點數誤差 ---- ## 時間複雜度分析 針對一個數字詢問是否為質數，需要跑 $2$ 直到 $\sqrt n$ 因此時間複雜度為 $O(\sqrt n)$ ---- 倘若今天有多筆詢問呢? 當針對範圍內的數字有多筆詢問時，用一樣的方法容易造成超時。 如果使用的是剛剛的方式，根據我們的時間複雜度分析，會變成 $Q \cdot \sqrt{N}$ 因此當數字一大時就很容易 Time Limit Exceed ---- ## Sieve of Eratosthenes ## 質數篩法 ---- ## Sieve of Eratosthenes 也就是建立 "質數表" 想法: 最小的質數開始，所有質數的倍數都一定不是質數 因此我們就先用一個陣列，儲存每個數字是不是質數 一開始先把所有大於1的數字當成質數 再從2依序把所有質數的倍數設成非質數 ---- ## Sieve of Eratosthenes  ---- ## 演算法步驟 從 $2$ 開始跑到 $N$ 判斷當前數字是不是質數 (是否被更改過) 如果是質數(尚未被更改過) 把所有質數的倍數設成非質數(更改) ---- ## 核心程式碼 ```cpp= //n為會跑到的最大值 bool isprime[1000005]; //紀錄每個數字是否是質數 vector<int> prime; // 儲存範圍內所有的質數 isprime[1]=1; // 1 代表此數字不是質數，否則為 0 //先將所有大於1數字設為質數(0) for(int i=2;i<=n;i++){ if(!isprime[i]){ //如果為質數 prime.push_back(i); for(long long j=2;i*j<=n;j++){ //記得容易會爆int的話要設long long //所有質數的倍數設成非質數 isprime[i*j] = 1; } } } ``` ---- ## 時間複雜度 主程式 ```cpp= vector<int> prime; void sieve(){ for(int i=2;i<=n;i++){ if(!isprime[i]){ prime.push_back(i); //把質數記錄進質數表裡面 for(long long j=2;i*j<=n;j++){ //主要迴圈 判定每個非質數 isprime[i*j] = 1; } } } } ``` 次數為 $\frac{n}{2} + \frac{n}{3} + ... + \frac{n}{n}$ 為調和級數 $= n \log n$ 時間複雜度為 $O(n \log n)$ ---- ## 針對質數篩的優化 在內迴圈裡的倍數可以直接從 $i$ 倍開始 因為小於 $i$ 的倍數之前都跑過了 ```cpp= vector<int> prime; void sieve(){ for(int i=2;i<=n;i++){ if(!isprime[i]){ prime.push_back(i); for(long long j=i; i*j<=n; j++){ //這裡的j可以直接從i開始 isprime[i*j] = 1; } } } } ``` 複雜度為 $O(N\log\log N)$ ---- ## 質數篩的應用 對於這段程式 ```cpp= vector<int> prime; bool isprime[N]; void sieve(){ for(int i=2;i<=n;i++){ if(!isprime[i]){ prime.push_back(i); for(long long j=i; i*j<=n; j++){ isprime[i*j] = 1; } } } } ``` 若我們把第六行的替代的元素換成 $i$ 是什麼意思呢? 也就是說 isprime$[i]$ 若等於 $0$ 代表他是質數，若不等於 $0$ 則代表其不是質數且 isprime$[i]$ 為 $i$ 的因數。 ---- 找到質數就往上乘覆蓋 $factor$ 的倍數 ```cpp= void sieve(){ for(int i=2;i<=n;i++){ if(!factor[i]){ prime.push_back(i); for(long long j=i; i*j<=n; j+=i){ factor[j] = i; //紀錄當前數字i遇到的因數 } } } } ``` ---- 如果我們建好質數表，每次只要查質數表 就可以直接判斷是不是質數 yeah :bread: ---- ## 質因數分解 那我們有質數表之後，除了判斷一個數字是否為質數外，還可以用來質因數分解。 首先可以想到的是，利用質數表去進行枚舉，那我們只需要好好維護他的邊界即可，裡面就當整除時代表此數字為因數，把此數字記錄下來並除掉。 ---- ## 程式碼實現 ```cpp void factorize(int x){ vector<int> factors; //質因數分解的答案存在這 //for 迴圈判斷當前遍歷到的質數還在範圍裡 for(int i = 0; i < primeSize && prime[i] * prime[i] <= x; i++){ while(x % prime[i] == 0) //若當前枚舉到的質數是此數的因數 factors.push_back(prime[i]), x/=prime[i]; //紀錄答案並更新 x } if(x > 1) factors.push_back(x); //最後記得把自己加進去，若他本身是很大的質數則有可能沒算到 for(int i=0;i<factors.size();i++){ //遍歷並輸出答案 cout<<factors[i]<<" \n"[i==factors.size()-1]; } } ``` ---- 另外還有一種很直觀的方法，既然表內 $factor[i]$ 的值都是 $i$ 的因數，那我們在質因數分解 $n$ 的時候，只需要不斷除以 $factor[n]$ ，並且把中途的數字記錄下來即可。 ---- ## 程式碼實現 ```cpp void factorize(int x){ //需要質因數分解的數字 vector<int> factors; while(x != 1){ factors.push_back(factor[x]); //記錄所有因數 x/=factor[x]; //把他除以自己的因數 } for(int i=0;i<factors.size();i++){ //最後factors裡面即為此數字質因數分解的結果 cout<<factors[i]<<" \n"[i==factors.size()-1]; } } ``` ---- ## 兩種方法分析 * 範圍問題 對於第一種方法來說，他會遍歷到 $prime[i]*prime[i]$，因此值域上可以處理到 $N*N$ ($N$為最大會遍歷到的範圍) 而對於第二種方法，由於會戳到 $factor[x]$ ， 因此需要讓 $x$ 的大小小於 $factor$ 的大小 (質數表大小) * 時間複雜度分析 對於第一種方法，不難看出他的時間複雜度就是 $O(\sqrt x)$ (迴圈遍歷) 對於第二種方法，最糟糕的情況是所有質因數都是 $2$ ，因此 $x$ 會不斷的除 $2$ 直到結束，因此時間複雜度為 $O(\log x)$ 針對以上兩種分析就可以看出這兩種方法的優劣以及該如何選擇。 ---- 之後會上到更快搜尋質數的方法 miller-rabin 有興趣再自己去看 --- ## 補充 因數的數量最多只有 $n^{\frac{1}{3}}$ 個 --- ## 練習時間&&本週題單 https://vjudge.net/contest/585204