( ext{Description})
( ext{Solution})
现在真的干啥啥不行，吃饭第一名。
方法壹
其实就是 7.3.3.用法叁，可以得到：
[sum_{T=1}^n ext{sum}left(frac{n}{T}
ight) imes T^2sum_{d|T}muleft(frac{T}{d}
ight) imes d
]
如果是普通莫比乌斯反演题，到这里就可以预处理 (+) 分块做了。
但这题不行。
观察 (sum_{d|T}muleft(frac{T}{d} ight) imes d) 的结构，其实就是 (mu* ext{id}=varphi)。
那我们只需要预处理 (i^2 imes varphi(i)) 的前缀和，显然这是个积性函数，试试杜教筛。
令 (f(i)=i^2 imes varphi(i))，为了卷积时抵消 (i^2) 项，令 (g(i)=i^2)。
则：
[h(n)=(f*g)(n)=sum_{d|n}f(d)gleft(frac{n}{d}
ight)
]
[=sum_{d|n}d^2 imes varphi(d) imes left(frac{n}{d}
ight)^2
]
[=n^2sum_{d|n}varphi(d)1left(frac{n}{d}
ight)
]
由于 (varphi*1= ext{id})，则：
[=n^2 imes n=n^3
]
带入杜教筛的套路式，做 ( ext{sum}) 的分块套上杜教筛套路式的分块即可。
方法贰
当然可以
[sum_{i=1}^nsum_{j=1}^ni imes j imes ext{id}(gcd(i,j))
]
[=sum_{i=1}^nsum_{j=1}^ni imes j imes sum_{d|i,d|j}varphi(d)
]
由于 (gcd) 是满足同时为 (i,j) 的因数的最大的数，所以 (d) 同时为 (i,j) 的因数就一定是 (gcd) 的因数，就可以巧妙地去掉 (gcd)。
有：
[=sum_{d=1}^n varphi(d)sum_{d|i}sum_{d|j}i imes j
]
[=sum_{d=1}^nvarphi(d) imes d^2left(sum_{i=1}^{n/d}i
ight)^2
]
然后和 方法壹 一样分块再杜教筛即可。
( ext{Code})
#include <cstdio>
#define rep(i,_l,_r) for(register signed i=(_l),_end=(_r);i<=_end;++i)
#define fep(i,_l,_r) for(register signed i=(_l),_end=(_r);i>=_end;--i)
#define erep(i,u) for(signed i=head[u],v=to[i];i;i=nxt[i],v=to[i])
#define efep(i,u) for(signed i=Head[u],v=to[i];i;i=nxt[i],v=to[i])
#define print(x,y) write(x),putchar(y)
template <class T> inline T read(const T sample) {
T x=0; int f=1; char s;
while((s=getchar())>'9'||s<'0') if(s=='-') f=-1;
while(s>='0'&&s<='9') x=(x<<1)+(x<<3)+(s^48),s=getchar();
return x*f;
}
template <class T> inline void write(const T x) {
if(x<0) return (void) (putchar('-'),write(-x));
if(x>9) write(x/10);
putchar(x%10^48);
}
template <class T> inline T Max(const T x,const T y) {if(x>y) return x; return y;}
template <class T> inline T Min(const T x,const T y) {if(x<y) return x; return y;}
template <class T> inline T fab(const T x) {return x>0?x:-x;}
template <class T> inline T gcd(const T x,const T y) {return y?gcd(y,x%y):x;}
template <class T> inline T lcm(const T x,const T y) {return x/gcd(x,y)*y;}
template <class T> inline T Swap(T &x,T &y) {x^=y^=x^=y;}
#include <map>
using namespace std;
typedef long long ll;
const int maxn=4e6;
int mod,inv6,p[maxn>>2],pc,phi[maxn+5];
bool is[maxn+5];
ll n;
map <ll,int> mp;
int qkpow(int x,int y) {
int r=1;
while(y) {
if(y&1) r=1ll*r*x%mod;
x=1ll*x*x%mod; y>>=1;
}
return r;
}
int calc2(ll x) {
x=x%mod;
return x*(x+1)%mod*(2*x+1)%mod*inv6%mod;
}
int calc3(ll x) {
x=x%mod;
return x*(x+1)/2%mod*(x*(x+1)/2%mod)%mod;
}
void init() {
phi[1]=1;
rep(i,2,maxn) {
if(!is[i]) p[++pc]=i,phi[i]=i-1;
rep(j,1,pc) {
if(i*p[j]>maxn) break;
is[i*p[j]]=1;
if(i%p[j]==0) {phi[i*p[j]]=phi[i]*p[j]; break;}
phi[i*p[j]]=phi[i]*(p[j]-1);
}
}
rep(i,1,maxn) phi[i]=(0ll+phi[i-1]+1ll*i*i%mod*phi[i]%mod)%mod;
}
int S(ll x) {
if(x<=maxn) return phi[x]%mod;
if(mp.count(x)) return mp[x];
int ret=0; ll r;
ret=calc3(x);
for(ll l=2;l<=x;l=r+1) {
r=(x/(x/l));
ret=(0ll+ret-1ll*(0ll+calc2(r)-calc2(l-1)+mod)%mod*S(x/l)%mod+mod)%mod;
}
mp.insert(make_pair(x,ret));
return ret;
}
int Work() {
ll r; int ret=0;
for(ll l=1;l<=n;l=r+1) {
r=(n/(n/l));
ret=(ret+1ll*calc3(n/l)*(0ll+S(r)-S(l-1)+mod)%mod)%mod;
}
return ret;
}
int main() {
mod=read(9),n=read(9ll); inv6=qkpow(6,mod-2);
init();
print(Work(),'
');
return 0;
}