关于阶乘问题
本文主要探討數(shù)學(xué)中的階乘問題。比如求階乘最后非零位的位數(shù)，計算階乘的值，等等。
題目：http://acm.hdu.edu.cn/showproblem.php?pid=1042
題目：http://www.51nod.com/onlineJudge/questionCode.html#!problemId=1059
分析：題意比較簡單，都是給定一個數(shù)，求階乘的值。第二題的數(shù)比較大。
???? 首先如果直接相乘肯定是不行的，稍微高級一點的做法就是采用萬進位，這樣就能過第一道題了。如果要暴 ???? 力，可以將階乘寫成冪的形式（先進行素數(shù)篩選）。用Java寫了一個，能過這兩個題，只不過效率太低。
代碼： import java.math.BigInteger; import java.util.Arrays; import java.util.Scanner;public class Main {public final static int N = 100005;public static int[] p = new int[N];public static boolean[] prime = new boolean[N];public static int cnt;@SuppressWarnings("resource")public static void main(String[] args){isprime();Scanner cin = new Scanner(System.in);int n = cin.nextInt();BigInteger ans = BigInteger.ONE;for(int i = 0; i < cnt; i++){if(p[i] > n) break;ans = ans.multiply(BigInteger.valueOf(p[i]).pow(getFactor(p[i], n)));}StringBuffer str = new StringBuffer("");char[] ch = ans.toString().toCharArray();for(int i = 0; i < ch.length; i++){if(i % 1000 == 0 && i > 0){System.out.println(str.toString());str.delete(0, str.length());}str.append(ch[i]);}if(str.length() > 0){System.out.println(str.toString());}}public static void isprime(){cnt = 0;Arrays.fill(prime, true);for(int i = 2; i < N; i++){if(prime[i]){p[cnt++] = i;for(int j = i + i; j < N; j += i){prime[j] = false;}}}}public static int getFactor(int p, int n){int ans = 0;while(n != 0){ans += n / p;n /= p;}return ans;} }
這樣的效率確實太低下，必須想更好的方法。有位大神的代碼效率比較高，使用了FFT優(yōu)化，我對其進行了修
改并加上了一些注釋。如下
代碼：
#include <iostream> #include <algorithm> #include <string.h> #include <stdio.h> #include <math.h>using namespace std;//FFT模板開始 const double PI = acos(-1.0);struct Virt {double r, i;Virt(double r = 0.0, double i = 0.0){this->r = r;this->i = i;}Virt operator + (const Virt &x){return Virt(r + x.r, i + x.i);}Virt operator - (const Virt &x){return Virt(r - x.r, i - x.i);}Virt operator * (const Virt &x){return Virt(r * x.r - i * x.i, i * x.r + r * x.i);} };//雷德算法--倒位序 void Rader(Virt F[], int len) {int j = len >> 1;for(int i = 1; i < len - 1; i++){if(i < j) swap(F[i], F[j]);int k = len >> 1;while(j >= k){j -= k;k >>= 1;}if(j < k) j += k;} }//FFT實現(xiàn) void FFT(Virt F[], int len, int on) {Rader(F, len);for(int h = 2; h <= len; h <<= 1) //分治后計算長度為h的DFT{Virt wn( cos(-on * 2 * PI / h), sin(-on * 2 * PI / h)); //單位復(fù)根e^(2*PI/m)用歐拉公式展開for(int j = 0; j < len; j += h){Virt w(1, 0); //旋轉(zhuǎn)因子for(int k = j; k < j + h / 2; k++){Virt u = F[k];Virt t = w * F[k + h / 2];F[k] = u + t; //蝴蝶合并操作F[k + h / 2] = u - t;w = w * wn; //更新旋轉(zhuǎn)因子}}}if(on == -1)for(int i = 0; i < len; i++)F[i].r /= len; }//求卷積 void Conv(Virt a[], Virt b[], int n) {FFT(a, n, 1);FFT(b, n, 1);for(int i = 0; i < n; i++)a[i] = a[i] * b[i];FFT(a, n, -1); } //FFT模板結(jié)束/***********************///代碼正文 typedef long long LL;const int N = 1000005; const int maxn = 500000; const int BITS = 1000; const int LEN = 3;Virt x1[N], x2[N];int num[maxn]; int temp[maxn]; int tmp[20][maxn];char buf[maxn]; int ans[maxn]; char rel[maxn * LEN]; LL bt[6]= {1, 10, 100, 1000, 10000, 100000};void Solve(int x[], int &L1, int y[], int L2) {int L = 1;while(L < 2 * L1 || L < 2 * L2) L <<= 1;for(int i = 0; i < L; i++){if(i < L1) x1[i].r = x[i];else x1[i].r = 0;if(i < L2) x2[i].r = y[i];else x2[i].r = 0;x1[i].i = x2[i].i = 0;}Conv(x1, x2, L);LL over = 0;L1 = 0;for(int i = 0; i < L; i++){over += x1[i].r + 0.5;x[i] = over % BITS;if(over) L1 = i;over /= BITS;}if(over) x[++L1] = over;++L1; }//將buf字符串轉(zhuǎn)化為int數(shù)組 void CharToInt(char *buf, int x[], int &len) {int bLen = strlen(buf);//求出余項的值int val = 0;for(int i = 0; i < bLen % LEN; i++)val = val * 10 + buf[i] - '0';len = 0;if(val) x[len++] = val;for(int i = bLen % LEN; i < bLen; i += LEN){val = 0;for (int j = 0; j < LEN; j++)val = val * 10 + buf[i + j] - '0';x[len++] = val;}//反轉(zhuǎn)字符串for(int i = 0; i < len / 2; i++)swap(x[i], x[len - 1 - i]); }void PutStr(int x[], int n, char *buf, int &len, int BITS_LEN) {int j = 0;sprintf(buf, "%d", x[n - 1]);while(buf[j]) j++;for(int i = n - 2; i >= 0; i--){for(int k = 0; k < BITS_LEN; k++){buf[j + k] = x[i] / bt[BITS_LEN - 1 - k] + '0';x[i] %= bt[BITS_LEN - 1 - k];}j += BITS_LEN;}buf[j] = 0;len = j; }//每1000位一千位地輸出 void Print(char *buf, int len) {for(int i = 0; i < len; i += 1000){char ch = buf[i + 1000];buf[i + 1000] = 0;printf("%s\n", buf + i);buf[i + 1000] = ch;} }void ConvertNum(int *num5, int len5, int *num, int &len) {PutStr(num5, len5, buf, len, 5);CharToInt(buf, num, len); }void Run(int l, int r, int *num, int &m) {memset(temp, 0, (r - l + 1) * sizeof(int));memset(num, 0, (r - l + 1) * sizeof(int));temp[0] = 1;m = 0;for(int i = l; i <= r; i++){LL top = 0;for(int j = 0; j <= m; j++){top += (LL)(temp[j]) * i;temp[j] = top % 100000;top /= 100000;}if(top)temp[++m] = top;}ConvertNum(temp, m + 1, num, m); }//分治法 void Divide(int dept, int ans[], int &n, int l, int r) {if (r - l <= 200){Run(l, r, ans, n);return ;}int mid = (l + r) >> 1;int rlen;Divide(dept + 1, ans, n, l, mid);Divide(dept + 1, tmp[dept] + mid * 2 + 2, rlen, mid + 1, r);Solve(ans, n, tmp[dept] + mid * 2 + 2, rlen); }int main() {int n, m;while (scanf("%d", &n) != EOF){Divide(1, ans, m, 1, n);PutStr(ans, m, rel, n, LEN);Print(rel, n);}return 0; }
總結(jié)
