CSU 1325 A very hard problem (莫比烏斯反演＋分塊求和優化)

_TCgogogo_發表於2015-08-23
優化


 A very hard problem

Time Limit: 3 Sec  Memory Limit: 160 MB
Submit: 266  Solved: 63
[Submit][Status][Web Board]

Description

CX老溼經常被人黑，被黑得多了，自己也就麻木了。於是經常聽到有人黑他，他都會深情地說一句：禽獸啊！

一天CX老溼突發奇想，給大家出了一個難題，並且聲稱誰能夠準確地回答出問題才能繼續黑他，否則他就要反擊了。

這個難題就是：

給出兩個數p和q，接下來q個詢問，每個詢問給出兩個數A和B，請分別求出：

一、有多少個有序數對（x，y）滿足1<=x<=A，1<=y<=B，並且gcd（x，y）為p的一個約數；

二、有多少個有序數對（x，y）滿足1<=x<=A，1<=y<=B，並且gcd（x，y）為p的一個倍數。

Input

只有一組測試資料。

第一行兩個數：p和q。（1<p<10^7 ，1<q<1000。）

接下來有q行，每行兩個數A和B。（1<A,B<10^7）

Output

輸出共q行。每行兩個數。用空格隔開。

分別表示題目描述中的兩個對應的答案。

(x,y)=(2,3)和(x,y)=(3,2)被視為兩個不同有序數對哦！

Sample Input

6 3
8 8
15 32
13 77

Sample Output

58 1
423 10
883 24

HINT

 對於64位整型請用lld，或者cin，cout。T_T

CSU_LQ

Source

CSU Monthly 2013 Oct.


題目連結：http://acm.csu.edu.cn/OnlineJudge/problem.php?id=1325


題目分析：倍數直接(a / p) * (b / p)即可，不解釋，主要是約數，開始列舉因子做，各種超時，這題有點類似HDU 4746，考慮預處理出每個因子對答案的貢獻，這裡貢獻指的是它的莫比烏斯函式的係數和，舉個的例子如果p＝6，它的因子有1，2，3，6，設F(x) = (a / x) * (b / x)，x為最大公約數，f(x)為x為最大公約數時的組數

f(1) = u(1)F(1 * 1) + u(2)F(1 * 2) + u(3)F(1 * 3) + u(4)F(1 * 4) + u(5)F(1 * 5) + u(6)F(1 * 6) + ...

f(2) = u(1)F(2 * 1) + u(2)F(2 * 2) + u(3)F(2 * 3) + ...

f(3) = u(1)F(3 * 1) + u(2)F(3 * 2) + ...

f(6) = u(1)F(6 * 1) + ...

ans = f(1) + f(2) + f(3) + f(6)，我們需要預處理出F函式前面的係數，顯然F(6)的係數為u(6) + u(3) + u(2) + u(1)，這些都是它的約數，因此預處理時通過類似篩法的nlogn複雜度即可完成對fac_msum陣列的預處理，fac_msum[i]表示，最終答案裡F(i)前面的係數

光這樣還是不夠的，依然超時，必須再加上分塊求和優化才可通過本題，總複雜度約為p + plog(p) + q*sqrt(p)，吐槽一句，單組樣例，1e7的資料，能不memset儘量不memset。。。

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <vector>
#define ll long long
using namespace std;
int const MAX = 1e7 + 5;
int p, q, pnum;
short mob[MAX];
bool noprime[MAX];
int pr[MAX], sum[MAX], fac_msum[MAX];
int fac[1000], facnum;

void Mobius()
{
    pnum = 0;
    mob[1] = 1;
    for(int i = 2; i < MAX; i++)
    {
        if(!noprime[i])
        {
            pr[pnum ++] = i;
            mob[i] = -1;
        }
        for(int j = 0; j < pnum && i * pr[j] < MAX; j++)
        {
            noprime[i * pr[j]] = true;
            if(i % pr[j] == 0)
            {
                mob[i * pr[j]] = 0;
                break;
            }
            mob[i * pr[j]] = -mob[i];
        }
    }
}

void pre()
{
    for(int i = 1; i * i <= p; i++)
    {
        if(p % i == 0)
        {   
            fac[facnum ++] = i;
            if(i * i != p)
                fac[facnum ++] = p / i;
        }
    }
    for(int i = 0; i < facnum; i++)
        for(int j = 1; j * fac[i] < MAX; j++)
            fac_msum[j * fac[i]] += mob[j];
    for(int i = 1; i < MAX; i++)
        sum[i] = sum[i - 1] + fac_msum[i];
}

ll cal(int l, int r)
{
    ll ans = 0;
    if(l > r)
        swap(l, r);
    for(int i = 1, last = 0; i <= l; i = last + 1)
    {
        last = min(l / (l / i), r / (r / i));
        ans += (ll) (l / i) * (r / i) * (sum[last] - sum[i - 1]);  
    }
    return ans;
}   

int main()
{
    Mobius();
    scanf("%d %d", &p, &q);
    pre();
    while(q --)
    {
        int a, b;
        facnum = 0;
        ll ans = 0;
        scanf("%d %d", &a, &b);
        printf("%lld %lld\n", cal(a, b), (ll) (a / p) * (b / p));
    }   
}


相關文章