(擴充套件)歐幾里得演算法、素性測試、埃式篩法、區間篩法、快速冪運算
阿新 • • 發佈:2019-01-03
來自挑戰程式設計競賽2.6 數學問題的解題竅門
1.歐幾里得演算法
求解最大公約數,時間複雜度在O(log max(a,b))以內,可以看出,輾轉相除法是非常高效的
int gcd(int a,int b)
{
return (b==0)?a:gcd(b,a%b);
}
2.擴充套件歐幾里得演算法
求解方程a*x+b*y=gcd(a,b),a、b、x、y均為整數,時間複雜度和輾轉相除法是相同的,函式返回gcd(a,b)。
int gcd(int a,int b) { return (b==0)?a:gcd(b,a%b); } int extgcd(int a,int b,int& x,int& y) { int d=a; if(b!=0){ d=extgcd(b,a%b,y,x); y-=(a/b)*x; } else{ x=1; y=0; } return d; }
3.素數測試
//素性測試O(√n) bool is_prime(int n) { for(int i=2;i*i<=n;i++){ if(n%i==0) return false; } return n!=1;//n等於1是例外 } //約數列舉O(√n) vector<int> divisor(int n) { vector<int> res; for(int i=2;i*i<=n;i++){ if(n%i==0){ res.push_back(i); if(i!=n/i) res.push_back(n/i); } } return res; } //整數分解O(√n)素數因子 map<int,int> prime_factor(int n) { map<int,int> res; for(int i=2;i*i<=n;i++){ while(n%i==0){ ++res[i]; n/=i; } } if(n!=1) res[n]=1; return res; }
其中map第一個int是n的素數因子,對應的第二個int是這個因子的個數
2.埃式篩法
給定一個正整數n,請問n以內有多少個素數。
首先將2到n的所有整數記下來,其中最小的整數2是素數,將表中2的倍數劃去。表中最小的整數是3,它不能被更小的數整除,所以是素數,將表中3的倍數劃去,以此類推,表中剩餘的最小數字是m,則m就是素數。這樣反覆操作即可列舉得到n以內的所有素數,時間複雜度為O(nloglogn),可以看作是線性時間。
int prime[maxn]; bool is_prime[maxn];//is_prime[i]是true表示i是素數 //返回n以內素數的個數 int sieve(int n) { int p=0; for(int i=0;i<=n;i++) is_prime[i]=true; is_prime[0]=is_prime[1]=false; for(int i=2;i<=n;i++){ if(is_prime[i]){ prime[p++]=i; for(int j=2*i;j<=n;j+=i) is_prime[j]=false; } } return p; }
3.區間篩法
給定整數a、b,求區間[a,b)內有多少個素數
b以內的合數的最小質因數一定不超過√b。如果有√b以內的素數表,就可以把埃式篩法用到區間[a,b)上了。也就是說,先做好[2,√b)和[a,b)的表,然後從[2,√b)表中篩選素數的同時,也將其倍數從[a,b)的表中劃去,最後剩下的就是[a,b)內的素數了。
typedef long long ll;
bool is_prime[1000000+5];
bool is_prime_small[1000000+5];
//對區間[a,b)內的整數執行篩法。is_prime[i-a]=true→i是素數
void segment_sieve(ll a,ll b)
{
for(int i=0;(ll)i*i<b;i++) is_prime_small[i]=true;
for(int i=0;i<b-a;i++) is_prime[i]=true;
for(int i=2;(ll)i*i<b;i++){
if(is_prime_small[i]){
for(int j=2*i;(ll)j*j<b;j+=i) is_prime_small[i]=false;//篩[2,√b)
for(ll j=max(2ll,(a+i-1)/i)*i;j<b;j+=i) is_prime[j-a]=false;//篩[a,b)
//2LL是2的長整數形式
//((a+i-1)/i)*i是滿足>=a&&%i==0的離a最近的數
//也可以寫成(a%i==0)?a:(a/i+1)*i
}
}
}
4.快速冪運算
typedef long long ll;
ll mod_pow(ll x,ll n,ll mod)
{
ll res=1;
while(n>0){
if(n&1) res=res*x%mod;
x=x*x%mod;
n>>=1;
}
return res;
}