题目链接

https://lydsy.com/JudgeOnline/problem.php?id=2154

https://lydsy.com/JudgeOnline/problem.php?id=2693

题解

i=1nj=1mlcm(i,j)=i=1nj=1mijgcd(i,j)=T=1min(n,m)S(nT,mT)TdTμ(d)d \begin{aligned} & \sum_{i=1}^n \sum_{j=1}^m \mathrm{lcm}(i,j)\\ = & \sum_{i=1}^n \sum_{j=1}^m \frac{ij}{\gcd(i,j)}\\ = & \sum_{T=1}^{\min(n,m)}S(\lfloor\frac{n}{T}\rfloor,\lfloor\frac{m}{T}\rfloor)T\sum_{d|T}\mu(d)d \end{aligned}

其中
S(x,y)=i=1xj=1yij=x(x+1)2y(y+1)2 \begin{aligned} S(x,y)= & \sum_{i=1}^x \sum_{j=1}^y ij\\ = & \frac{x(x+1)}{2}\frac{y(y+1)}{2} \end{aligned}
考虑
g(T)=TdTμ(d)d g(T)= T\sum_{d|T}\mu(d) d
满足以下条件
g(T)={1T=1T(1T)TPg(pT)={g(p)g(T)pP,pTpg(T)pP,pT g(T)= \begin{cases} 1 & T=1\\ T(1-T) & T\in \mathbb P \end{cases}\\ g(pT)= \begin{cases} g(p)g(T) & p\in \mathbb P,p\nmid T\\ pg(T) & p\in \mathbb P,p\mid T \end{cases}
因此可以线筛求出gg,答案可以整除分块。

代码

BZOJ 2154

#include <cstdio>
#include <algorithm>

int read()
{
  int x=0,f=1;
  char ch=getchar();
  while((ch<'0')||(ch>'9'))
    {
      if(ch=='-')
        {
          f=-f;
        }
      ch=getchar();
    }
  while((ch>='0')&&(ch<='9'))
    {
      x=x*10+ch-'0';
      ch=getchar();
    }
  return x*f;
}

const int maxn=10000000;
const int mod=20101009;

int p[maxn+10],prime[maxn+10],cnt,f[maxn+10],sum[maxn+10];

int getprime()
{
  p[1]=1;
  f[1]=1;
  for(int i=2; i<=maxn; ++i)
    {
      if(!p[i])
        {
          prime[++cnt]=i;
          f[i]=1-i;
        }
      for(int j=1; (j<=cnt)&&(i*prime[j]<=maxn); ++j)
        {
          p[i*prime[j]]=1;
          if(i%prime[j]==0)
            {
              f[i*prime[j]]=f[i];
              break;
            }
          f[i*prime[j]]=1ll*f[i]*(mod+1-prime[j])%mod;
        }
    }
  for(int i=1; i<=maxn; ++i)
    {
      sum[i]=(sum[i-1]+1ll*f[i]*i)%mod;
    }
  return 0;
}

int n,m,ans;

int S(int x)
{
  return (1ll*x*(x+1)/2)%mod;
}

int main()
{
  getprime();
  n=read();
  m=read();
  for(int l=1,r; l<=std::min(n,m); l=r+1)
    {
      r=std::min(n/(n/l),m/(m/l));
      ans=(ans+1ll*(sum[r]-sum[l-1]+mod)*S(n/l)%mod*S(m/l))%mod;
    }
  printf("%d\n",ans);
  return 0;
}

BZOJ 2693

#include <cstdio>
#include <algorithm>

int read()
{
  int x=0,f=1;
  char ch=getchar();
  while((ch<'0')||(ch>'9'))
    {
      if(ch=='-')
        {
          f=-f;
        }
      ch=getchar();
    }
  while((ch>='0')&&(ch<='9'))
    {
      x=x*10+ch-'0';
      ch=getchar();
    }
  return x*f;
}

const int maxn=10000000;
const int mod=100000009;

int p[maxn+10],prime[maxn+10],cnt,f[maxn+10],sum[maxn+10];

int getprime()
{
  p[1]=1;
  f[1]=1;
  for(int i=2; i<=maxn; ++i)
    {
      if(!p[i])
        {
          prime[++cnt]=i;
          f[i]=1-i;
        }
      for(int j=1; (j<=cnt)&&(i*prime[j]<=maxn); ++j)
        {
          p[i*prime[j]]=1;
          if(i%prime[j]==0)
            {
              f[i*prime[j]]=f[i];
              break;
            }
          f[i*prime[j]]=1ll*f[i]*(mod+1-prime[j])%mod;
        }
    }
  for(int i=1; i<=maxn; ++i)
    {
      sum[i]=(sum[i-1]+1ll*f[i]*i)%mod;
    }
  return 0;
}

int T,n,m,ans;

int S(int x)
{
  return (1ll*x*(x+1)/2)%mod;
}

int main()
{
  getprime();
  T=read();
  while(T--)
    {
      n=read();
      m=read();
      ans=0;
      for(int l=1,r; l<=std::min(n,m); l=r+1)
        {
          r=std::min(n/(n/l),m/(m/l));
          ans=(ans+1ll*(sum[r]-sum[l-1]+mod)*S(n/l)%mod*S(m/l))%mod;
        }
      printf("%d\n",ans);
    }
  return 0;
}