题意：

有这样一个有关最大公约数的函数: 函数 f(x, y):

{
     c=0
     当 y>0:
     {
          c +=1
          t = x % y
          x = y
          y = t
      }
      返回 c * x * x
}

给出三个正整数n,m,p，你需要计算:

$$ \sum_{i=1}^{n} \sum_{j=1}^{m} \left \lfloor \frac{i*j}{f(i,j))} \right \rfloor $$

n <= 666,666,666, m <= 666, p <= 666,666,666。

思路：

打表找规律。

但是找规律也要按照基本法

观察到m比较小，对于固定的j,容易看出f(i,j)和f(i+j*k,k)是等价的。

比赛的时候没做出来，因为纠结取整的问题...

解决办法竟然是....通过循环节观察orz

转载一篇靠谱的题解：

一开始，我自己假设先不考虑c。那么就变成了ΣΣi/gcd*j/gcd=Σj/gcd*Σi/gcd，如此一来，由于m比较小，我就可以枚举j，然后对应求出j所有的因子作为gcd，gcd确定之后再根据容斥来统计i/gcd的和。具体统计方法和15年沈阳regional的frog那题类似，用n的因子来进行暴力的容斥。但是很显然这样子很难把c的影响带进来，而且这里的c还要向下取整，更加的麻烦。

于是打表找规律，首先很容易知道f(i,j)=f(i+kj,j)。根据这个，在不考虑向下取整的情况下，对于同一个j，我们就可以列出一个等差数列，其中首项是ij/f(i,j)，公差为jj/f(i,j)。但是这里要考虑这个向下取整。我们设i为模9为7的数j为9，可以打出如下i*j/f(i,j)的表，括号内为相邻值的差：

可以发现，每c组i*j/f(i,j)是一个循环节，也就是说可以看作c个等差数列，然后对于每一个等差数列，它的首项我们可以暴力算出，而公差也很容易求出。利用等差数列求和公式可以很快速的计算出结果。复杂度的话，我们需要枚举j和在j剩余系下的i，然后还有c个等差数列，复杂度为O(m^2logN)在接受范围之内。具体见代码：

以及。。膜真的很浪费生命啊？

一开始取了5个%，TLE,4个%,50%概率AC,3个%就比较稳得过了...

一个%大概300ms orz

/* ***********************************************
Author :111qqz
Created Time :2017年11月02日 星期四 09时38分56秒
File Name :5970.cpp
************************************************ */

#include <bits/stdc++.h>
#define PB push_back
#define fst first
#define sec second
#define lson l,m,rt<<1
#define rson m+1,r,rt<<1|1
#define ms(a,x) memset(a,x,sizeof(a))
typedef long long LL;
#define pi pair < int ,int >
#define MP make_pair

using namespace std;
const double eps = 1E-8;
const int dx4[4]={1,0,0,-1};
const int dy4[4]={0,-1,1,0};
const int inf = 0x3f3f3f3f;
LL n,m,p;
LL f(LL x,LL y,LL &c)
{
    LL t;
    while (y>0)
    {
    c++;
    t = x % y;
    x = y;
    y = t;
    }
    return c * x * x;
}
int main()
{
    #ifndef  ONLINE_JUDGE 
    freopen("./in.txt","r",stdin);
  #endif
    int T;
    cin>>T;
    while (T--)
    {
        scanf("%lld%lld%lld",&n,&m,&p);
        LL ans = 0 ;
        for ( int j = 1 ; j <= m ; j++)
        {
        for ( int i = 1 ; i <= j && i <= n ; i++)
        {
            LL c=0,fun = f(i,j,c);
            for ( int k = 0 ; k < c;  k++)
            {
            if (i+k*j>n) break;
            LL a1 = (i+k*j)*j/fun ;
            LL d = c*j*j/fun;
            LL num = (n-i-k*j)/(c*j) + 1;
            ans = (ans + num * a1 % p + (num-1)*(num)/2%p*d )%p;
            }
        }
        }
        printf("%lld\n",ans);
    }                                                                         
  #ifndef ONLINE_JUDGE  
  fclose(stdin);
  #endif
    return 0;
}