2019牛客暑期多校训练营（第一场）——B.Integration——（推公式）

最新推荐文章于 2020-07-13 16:44:07 发布

通信男神杨丽斌

最新推荐文章于 2020-07-13 16:44:07 发布

阅读量238

点赞数 1

分类专栏：数论

本文链接：https://blog.csdn.net/Q755100802/article/details/96447702

版权

数论专栏收录该内容

15 篇文章 0 订阅

订阅专栏

链接：https://ac.nowcoder.com/acm/contest/881/B
来源：牛客网

Integration

时间限制：C/C++ 2秒，其他语言4秒
空间限制：C/C++ 524288K，其他语言1048576K
64bit IO Format: %lld

题目描述

Bobo knows that ∫∞011+x2 dx=π2.∫0∞11+x2 dx=π2.

Given n distinct positive integers a1,a2,…,ana1,a2,…,an, find the value of
1π∫∞01∏ni=1(a2i+x2) dx.1π∫0∞1∏i=1n(ai2+x2) dx.

It can be proved that the value is a rational number pqpq.
Print the result as (p⋅q−1)mod(109+7)(p⋅q−1)mod(109+7).

输入描述:

The input consists of several test cases and is terminated by end-of-file.

The first line of each test case contains an integer n.
The second line contains n integers a1,a2,…,ana1,a2,…,an.

* 1≤n≤1031≤n≤103
* 1≤ai≤1091≤ai≤109
* {a1,a2,…,an}{a1,a2,…,an} are distinct.
* The sum of n2n2 does not exceed 107107.

输出描述:

For each test case, print an integer which denotes the result.

示例1

输入

复制

输出

复制

500000004
250000002
83333334

这题血亏，早早推出来公式，脑袋一抽每次都取逆元再乘，被卡的明明白白。。。

其实就是将乘式化成和式，就可以用积分的线性公式了

考虑两个式子相乘：

$\frac{1}{x^{2}+ai^{2}}\cdot \frac{1}{x^{2}+aj^{2}}$

很明显可以裂项

$(\frac{1}{aj^{2}-ai^{2}})\cdot (\frac{1}{x^{2}+ai^{2}}-\frac{1}{x^{2}+aj^{2}})$

而题目都告诉我们了，

$\int \frac{1}{x^{2}+ai^{2}}dx=\frac{1}{ai}\int\frac{1}{1+(\frac{x}{ai})^{2}}d\frac{x}{ai}=\frac{\pi }{2\cdot ai}$

所以再乘上系数就是答案。对于两个的系数是这样，那么对于三个的呢？

那就是

**$\frac{1}{x^{2}+ai^{2}}\frac{1}{x^{2}+ak^{2}}+\frac{1}{x^{2}+aj^{2}}\frac{1}{x^{2}+ak^{2}}$**

那还是两个式子相乘的裂，再乘上原来的系数就好啦

你可以再推一下，反正最后，对于每个乘式i，裂项后它的系数就是

$\frac{1}{ai}\cdot \prod_{j=1,j!=i}^{n}\frac{1}{aj^{2}-ai ^{2}}$

再乘以1/2的逆元，和起来就是答案

（不要每次取逆元！乘完再取！）

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <iostream>
#include <queue>
#include <cmath>
#define int long long
using namespace std;

const int mod=1e9+7;
const int maxn=1e3+10;


//扩展GCD
int ex_gcd(int a,int b,int &x,int &y){
   if(b==0){x = 1;y = 0;return a;}
   int g = ex_gcd(b,a%b,x,y);
   int temp = x;
   x = y;
   y = temp - a/b*y;
   return g;
}

//扩展欧几里得的另一种写法
void ex_gcd(int a,int b,int &d,int &x,int &y){
    //求解ax+by=gcd(a,b)的一组解
    if(!b){
        d=a,x=1,y=0;
    }
    else{
        ex_gcd(b,a%b,d,y,x);
        y-=x*(a/b);
    }
}
//逆元
int inv(int a,int mod){
   int X,Y;
   int g = ex_gcd(a,mod,X,Y);
   if(g!=1)return -1;
   return (X%mod + mod)%mod;
}
int n;
int a[maxn];

int get_num(int x){
    int ans=1;
    for(int i=1;i<=n;i++){
        if(i==x)continue;
        int temp=(a[i]*a[i]%mod-a[x]*a[x]%mod+mod)%mod;
        ans=ans*temp%mod;
    }
    return inv(ans,mod);
}

int inva[maxn];

signed main(){
    while(scanf("%lld",&n)!=EOF){
        for(int i=1;i<=n;i++){
            scanf("%lld",&a[i]);
            inva[i]=inv(a[i],mod);
        }
        //sort(a+1,a+n+1);
        int ans=0;
        int inv2=inv(2,mod);
        for(int i=1;i<=n;i++){
            int now=get_num(i);
            now=now*inva[i]%mod;
            ans=(ans+now*inv2)%mod;
        }
        printf("%lld\n",ans);
    }
    return 0;
}