UVA 10375 选择与除法 （唯一分解性定理）

The binomial coefficient C(m,n) is defined as
C(m,n) =
m! (m−n)! n! Given four natural numbers p, q, r, and s, compute the the result of dividing C(p,q) by C(r,s).
Input
Input consists of a sequence of lines. Each line contains four non-negative integer numbers giving values for p, q, r, and s, respectively, separated by a single space. All the numbers will be smaller than 10,000 with p ≥ q and r ≥ s.
Output
For each line of input, print a single line containing a real number with 5 digits of precision in the fraction, giving the number as described above. You may assume the result is not greater than 100,000,000.
Sample Input
10 5 14 9 93 45 84 59 145 95 143 92 995 487 996 488 2000 1000 1999 999 9998 4999 9996 4998
Sample Output
0.12587 505606.46055 1.28223 0.48996 2.00000 3.99960

#include<iostream>
#include<algorithm>
#include<string>
#include<cmath>//int dx[4]={0,0,-1,1};int dy[4]={-1,1,0,0};
#include<set>//int gcd(int a,int b){return b?gcd(b,a%b):a;}
#include<queue>
#include<cmath>
#include<stack>
#include<string.h>
#include<stdlib.h>
#include<cstdio>
#define mod 1e9+7
#define ll long long
#define maxn 10005
#define MAX 500005
#define ms memset
using namespace std;
#pragma comment(linker, "/STACK:1024000000,1024000000") ///在c++中是防止暴栈用的
vector<int> prime;
void get_prime()
{
    int check[maxn];
    memset(check,0,sizeof(check));
    for(int i=2;i<maxn;i++)
    {
        if(check[i])  continue;
        prime.push_back(i);
       /// cout<<i<<endl;
        for(int j=i*i;j<maxn;j+=i)
            check[j]=1;
    }
}
/*
题目大意：给定四个数求式子。

整数范围是在10000这样，
所以用素数筛过比较好，
即答案用素数的乘积表示，
相除和相乘都用素数表的加减代替，
这样思路就清晰了很多。
紫书上的例题。
*/
int e[maxn];
void add_integer(int n,int d)
{
    for(int i=0;i<prime.size();i++)
    {
        while(n%prime[i]==0)
        {
            n/=prime[i];
            e[i]+=d;
        }
        if(n==1) break;
    }
}

void add_factorial(int n,int d)
{
    //memset(e,0,sizeof(e));
    for(int i=1;i<=n;i++)
        add_integer(i,d);
}

int p,q,r,s;
int main()
{
    get_prime();
    while(cin>>p>>q>>r>>s)
    {
        memset(e,0,sizeof(e));
        add_factorial(p,1);
        add_factorial(q,-1);
        add_factorial(p-q,-1);
        add_factorial(r,-1);
        add_factorial(s,1);
        add_factorial(r-s,1);
        double ans=1;
        for(int i=0;i<prime.size();i++)
            ans*=pow(prime[i],e[i]);
        printf("%.5f\n",ans);
    }
    return 0;
}