51Nod-1238 最小公倍数和V3(杜教筛)

$Description$
$Calculate \quad \sum_{i=1}^{n}\sum_{j=1}^{n}lcm(i,j)\quad n\leq 10^{10}$

$solution \;1$
$\sum_{i=1}^{n}\sum_{j=1}^{n}lcm(i,j)\\ =\sum_{i=1}^{n}\sum_{j=1}^{n}\frac{i\cdot j}{gcd(i,j)} \\= \sum_{d=1}^{n}d\sum_{i=1}^{\frac{n}{d}}i\sum_{i=1}^{\frac{n}{d}}j[gcd(i,j)=1] \\=\sum_{d=1}^{n}d(2\sum_{i=1}^{\frac{n}{d}}i\sum_{j=1}^{i}j\;-1) \\=\sum_{d=1}^{n}d\sum_{i=1}^{\frac{n}{d}}i^2 \cdot \varphi(i) \\ S(n)=\sum_{i=1}^{n}\sum_{d|i}d^2\cdot \varphi(d) \\ S(n)=(\frac{n\cdot (n+1)}{2})^2-\sum_{d=2}^{n}d^2\cdot S(\frac{n}{d})$
$Code$

#include<bits/stdc++.h>
using namespace std;
#define me(a,x) memset(a,x,sizeof(a))
#define IN freopen("in.txt","r",stdin);
#define STR clock_t startTime = clock();
#define END clock_t endTime = clock();cout << double(endTime - startTime) / CLOCKS_PER_SEC *1000<< "ms" << endl;
const int N=1e7+6;
const int mod=1e9+7;
const int inv2=mod+1>>1;
const int inv6=166666668;
typedef long long ll;
int prime[N],tot=0;
bool vis[N]={0};
ll f[N],phi[N];
void pre(){
  phi[1]=1;f[0]=phi[0]=0;
  for(int i=2;i<N;i++){
      if(!vis[i])prime[++tot]=i,phi[i]=i-1;
      for(int j=1;j<=tot&&i*prime[j]<N;j++){
          vis[i*prime[j]]=1;
          if(i%prime[j]==0){
              phi[i*prime[j]]=phi[i]*prime[j];
              break;
          }else phi[i*prime[j]]=phi[i]*(prime[j]-1);
      }
  }
  for(int i=1;i<N;i++){
      f[i]=(1LL*i*i%mod*phi[i]%mod+f[i-1])%mod;
  }
}
map<ll,ll>p;
ll F(ll x){
  x%=mod;
  x=x*(x+1)%mod*(2LL*x+1)%mod*inv6%mod;
  return x;
}
ll LR(ll x){
  x%=mod;
  x=x*(x+1)%mod*inv2%mod;
  return x;
}
ll cal(ll n){
  if(n<N)return f[n];
  if(p[n])return p[n];
  ll ans=LR(n);ans=ans*ans%mod;
  for(ll i=2,last;i<=n;i=last+1){
      last=n/(n/i);
      ans-=1LL*(F(last)-F(i-1)+mod)%mod*cal(n/i)%mod;
  }
  return p[n]=ans;
}
int main(){
  pre();
  ll n;scanf("%lld",&n);
  ll ans=0;
  for(ll i=1,last;i<=n;i=last+1){
      last=n/(n/i);
      ans+=1LL*(LR(last)-LR(i-1)+mod)%mod*cal(n/i)%mod;
      ans%=mod;
  }
  printf("%lld\n",ans);
}

$solution \;2$
$\sum_{i=1}^{n}\sum_{j=1}^{n}lcm(i,j)=\sum_{T=1}^{n}(\frac{\left \lfloor \frac{n}{T}\right \rfloor \cdot (\left \lfloor \frac{n}{T}\right \rfloor+1) }{2})^2\cdot \sum_{d|T}\frac{T}{d}\cdot d^2\cdot \mu(d)$ 公式推导
$S(n)=\sum_{i=1}^{n}\sum_{d|i}\frac{i}{d}\cdot d^2\cdot \mu(d)$
$f(n)=\sum_{d|n}\frac{n}{d}\cdot d^2\cdot \mu(d) \quad h(n)=n \quad g(n)=n^2 \cdot \mu(n) \quad f=h\cdot g$
$考虑消掉 \mu \; 构造 p(n)=n^2 ,\; p\cdot g=n^2\sum_{d|n}\mu(d)=[n=1]$
$f\cdot p=h\cdot p\cdot g=h\cdot I=h=n$
$\sum_{i=1}^{n}(f\cdot p)(i)=\sum_{i=1}^{n}i=\sum_{d=1}^{n}d^2 \cdot S(\frac{n}{d})$
$S(n)=\frac{n\cdot (n+1)}{2}-\sum_{d=2}^{n}d^2 \cdot S(\frac{n}{d})$
$Code$

#include<bits/stdc++.h>
using namespace std;
#define me(a,x) memset(a,x,sizeof(a))
#define IN freopen("in.txt","r",stdin);
#define STR clock_t startTime = clock();
#define END clock_t endTime = clock();cout << double(endTime - startTime) / CLOCKS_PER_SEC *1000<< "ms" << endl;
const int N=1e7+6;
const int mod=1e9+7;
const int inv2=mod+1>>1;
const int inv6=166666668;
typedef long long ll;
int prime[N],tot=0,mu[N];
bool vis[N]={0};
ll f[N];
void pre(){
  mu[1]=1;f[0]=mu[0]=0;
  for(int i=2;i<N;i++){
      if(!vis[i])prime[++tot]=i,mu[i]=-1;
      for(int j=1;j<=tot&&i*prime[j]<N;j++){
          vis[i*prime[j]]=1;
          if(i%prime[j]==0){
              mu[i*prime[j]]=0;
              break;
          }else mu[i*prime[j]]=-mu[i];
      }
  }
  for(int i=1;i<N;i++)
  for(int j=i;j<N;j+=i){
      f[j]+=1LL*j*i%mod*mu[i];
      f[j]=(f[j]+mod)%mod;
  }
  for(int i=1;i<N;i++)f[i]=(f[i]+f[i-1])%mod;
}
ll G(ll x){
  x%=mod;
  x=x*(x+1)%mod*inv2%mod;
  return x;
}
ll F(ll x){
  x%=mod;
  x=x*(x+1)%mod*(2LL*x+1)%mod*inv6%mod;
  return x;
}
map<ll,ll>p;
ll cal(ll n){
  if(n<N)return f[n];
  if(p[n])return p[n];
  ll ans=G(n);
  for(ll i=2,last;i<=n;i=last+1){
      last=n/(n/i);
      ans-=(F(last)-F(i-1)+mod)%mod*cal(n/i);
      ans=(ans+mod)%mod;
  }
  return p[n]=ans;
}
int main(){
  pre();
  ll n;scanf("%lld",&n);
  ll ans=0;
  for(ll i=1,last;i<=n;i=last+1){
      last=n/(n/i);
      ll x=G(n/i);
      x=x*x%mod;
      ans+=(cal(last)-cal(i-1)+mod)%mod*x%mod;
      ans%=mod;
  }
  printf("%lld\n",ans);
}