米缇

首先有：

\[\sigma_k(xy)=\sum_{a|x} \sum_{b|y} [\gcd(a,b)=1] (a\frac{y}{b})^k \\ \]

大概是枚举最简分数保证不重

\[\sum_{i=1}^n \sum_{j=1}^n\sigma_k(ij) \\ \sum_{i=1}^n \sum_{j=1}^n \sum_{x|i} \sum_{y|j} [\gcd(x,y)=1] (x\frac{j}{y})^k \\ \sum_{x=1}^n \sum_{y=1}^n [\gcd(x,y)=1](\frac{x}{y})^k \sum_{i=1}^{n/x}\sum_{j=1}^{n/y} (jy)^k \\ \sum_{x=1}^n \sum_{y=1}^n [\gcd(x,y)=1]x^k \lfloor \frac{n}{x}\rfloor\sum_{j=1}^{n/y} j^k \\ \sum_{d=1}^n \mu(d) \sum_{x=1}^{n/d}\sum_{y=1}^{n/d} (xd)^k \lfloor \frac{n}{xd} \rfloor \sum_{j=1}^{n/yd }j^k \\ \sum_{d=1}^n \mu(d)d^k \sum_{x=1}^{n/d} x^k \lfloor \frac{n}{xd} \rfloor \sum_{y=1}^{n/d}\sum_{j=1}^{n/yd }j^k \\ \sum_{d=1}^n \mu(d)d^k \left(\sum_{x=1}^{n/d} x^k \lfloor \frac{n}{xd} \rfloor \right)^2 \]

令 \(f(n)=\mu(n)n^k,g(n)=\sum_{i=1}^n i^k \lfloor \frac{n}{d} \rfloor\) ,

\[\sum_{d=1}^n f(d) g^2(\lfloor \frac{n}{d} \rfloor) \]

---

\(f(n)\): 使用杜教筛

\[f: \mu *id_k \\ g: id_k\\ h: n^k \epsilon(n)\\ \]

---

\(g(n)\)：

\[g(n)=\sum_{i=1}^n i^k \lfloor \frac{n}{d} \rfloor \]

---

\[\delta g = \sum_{i=1}^{n-1} [i|n] i^k + n^k=\sigma_k(n) \\ \Rightarrow g(n)=\sum_{i=1}^n\sigma_k(i) \]

线性筛一部分，剩下的用原始式子整除分块

---

自然数幂和使用第二类斯特林数

\[\sum_{i=1}^ni^k \\ \sum_{i=0}^n i^k - [k=0] \\ \sum_{i=0}^n \sum_{j=0}^k \begin{Bmatrix}k\\j\end{Bmatrix} \binom{i}{j}j! - [k=0] \\ \sum_{j=0}^k \begin{Bmatrix}k\\j\end{Bmatrix} j! \sum_{i=0}^{n} \binom{i}{j} - [k=0] \\ \sum_{j=0}^k \begin{Bmatrix}k\\j\end{Bmatrix} j! \sum_{i=0}^{n-j} \binom{j+i}{j} - [k=0] \\ \sum_{j=0}^k \begin{Bmatrix}k\\j\end{Bmatrix} j! \sum_{i=0}^{n-j} \binom{j+i}{i} - [k=0] \\ \sum_{j=0}^k \begin{Bmatrix}k\\j\end{Bmatrix} j! \binom{n+1}{n-j} - [k=0] \\ \sum_{j=0}^k \begin{Bmatrix}k\\j\end{Bmatrix} j! \binom{n+1}{j+1} - [k=0] \\ \]

所有计算部分都要使用记忆化

#include <cstdio>
#include <unordered_map>
using namespace std;
#define ll long long

const int MAXN = 1e7 , MAXK = 7000 , Mod = 1e9 + 7;
inline int Add( int x , int y ) { x += y; return x >= Mod ? x - Mod : x; }
inline int Sub( int x , int y ) { x -= y; return x < 0 ? x + Mod : x; }
inline int Mul( int x , int y ) { return 1ll * x * y % Mod; }
inline int Qkpow( int x , int po ) { int p = 1; for( ; po ; po >>= 1 , x = Mul( x , x ) ) if( po & 1 ) p = Mul( p , x ); return p; }
inline int Inv( int x ) { return Qkpow( x , Mod - 2 ); }

int prnum , prime[ MAXN + 5 ] , lowp[ MAXN + 5 ] , loc[ MAXN + 5 ] , sk[ MAXN + 5 ] , pk[ MAXN + 5 ] , mu[ MAXN + 5 ];
int sf[ MAXN + 5 ] , g[ MAXN + 5 ] , spk[ MAXN + 5 ];
int k , iv[ MAXN + 5 ] , fac[ MAXN + 5 ] , S2[ MAXK + 5 ][ MAXK + 5 ];
bool vis[ MAXN + 5 ];
void sieve( ) {
	fac[ 0 ] = 1;
	for( int i = 1 ; i <= MAXN ; i ++ ) fac[ i ] = Mul( fac[ i - 1 ] , i );
	iv[ 1 ] = 1;
	for( int i = 2 ; i <= MAXN ; i ++ ) iv[ i ] = Mul( iv[ Mod % i ] , Mod - Mod / i );
	
	S2[ 0 ][ 0 ] = 1;
	for( int i = 1 ; i <= MAXK ; i ++ )
		for( int j = 1 ; j <= i ; j ++ )
			S2[ i ][ j ] = Add( S2[ i - 1 ][ j - 1 ] , Mul( S2[ i - 1 ][ j ] , j ) );
	
	sk[ 1 ] = 1; mu[ 1 ] = 1; pk[ 1 ] = 1;
	for( int i = 2 ; i <= MAXN ; i ++ ) {
		if( !vis[ i ] ) {
			prime[ ++ prnum ] = i; lowp[ i ] = i; loc[ i ] = 1;
			sk[ i ] = Add( 1 , pk[ i ] = Qkpow( i , k ) ); mu[ i ] = Mod - 1;
		}
		for( int j = 1 ; j <= prnum && 1ll * i * prime[ j ] <= MAXN ; j ++ ) {
			vis[ i * prime[ j ] ] = 1; pk[ i * prime[ j ] ] = Mul( pk[ i ] , pk[ prime[ j ] ] );
			if( i % prime[ j ] == 0 ) {
				lowp[ i * prime[ j ] ] = lowp[ i ] * prime[ j ]; loc[ i * prime[ j ] ] = loc[ i ] + 1;
				if( lowp[ i * prime[ j ] ] == i * prime[ j ] ) {
					if( k != 0 ) {
						int p = Qkpow( prime[ j ] , k );
						sk[ i * prime[ j ] ] = Mul( Qkpow( p , loc[ i * prime[ j ] ] + 1 ) - 1 , Inv( p - 1 ) );
					}
					else sk[ i * prime[ j ] ] = loc[ i * prime[ j ] ] + 1;
				}
				else sk[ i * prime[ j ] ] = Mul( sk[ i * prime[ j ] / lowp[ i * prime[ j ] ] ] , sk[ lowp[ i * prime[ j ] ] ] );
				break;
			}
			lowp[ i * prime[ j ] ] = prime[ j ]; loc[ i * prime[ j ] ] = 1;
			sk[ i * prime[ j ] ] = Mul( sk[ i ] , sk[ prime[ j ] ] );
			mu[ i * prime[ j ] ] = Mod - mu[ i ];
		}
	}
	for( int i = 1 ; i <= MAXN ; i ++ ) {
		sf[ i ] = Add( sf[ i - 1 ] , Mul( mu[ i ] , pk[ i ] ) );
		g[ i ] = Add( g[ i - 1 ] , sk[ i ] ); 
		spk[ i ] = Add( spk[ i - 1 ] , pk[ i ] );
	}
}

unordered_map< ll , int > mSpk;
int Spk( ll n ) {
	if( n <= MAXN ) return spk[ n ]; if( mSpk[ n ] ) return mSpk[ n ];
	if( k == 0 ) return n % Mod;
	int Ans = 0 , c = 1;
	for( int i = 0 ; i <= k ; i ++ ) {
		c = Mul( c , Mul( ( n - i + 1 ) % Mod , iv[ i + 1 ] ) );
		Ans = Add( Ans , Mul( Mul( S2[ k ][ i ] , fac[ i ] ) , c ) );
	}
	return mSpk[ n ] = Ans;
}

unordered_map< ll , int > mSF;
int SF( ll n ) {
	if( n <= MAXN ) return sf[ n ]; if( mSF[ n ] ) return mSF[ n ];
	int Ans = 1;
	for( ll l = 2 , r ; l <= n ; l = r + 1 ) {
		r = n / ( n / l );
		Ans = Sub( Ans , Mul( Sub( Spk( r ) , Spk( l - 1 ) ) , SF( n / l ) ) );
	}
	return mSF[ n ] = Ans;
}

unordered_map< ll , int > mG;
int G( ll n ) {
	if( n <= MAXN ) return g[ n ]; if( mG[ n ] ) return mG[ n ];
	int Ans = 0;
	for( ll l = 1 , r ; l <= n ; l = r + 1 ) {
		r = n / ( n / l );
		Ans = Add( Ans , Mul( Sub( Spk( r ) , Spk( l - 1 ) ) , n / l % Mod ) );
	}
	return mG[ n ] = Ans;
}

ll n;
int main( ) {
	scanf("%lld %d",&n,&k); sieve();
	
	int Ans = 0;
	for( ll l = 1 , r ; l <= n ; l = r + 1 ) {
		r = n / ( n / l );
		Ans = Add( Ans , Mul( Sub( SF( r ) , SF( l - 1 ) ) , Mul( G( n / l ) , G( n / l ) ) ) ); 
	}
	printf("%d\n", Ans );
	return 0;
}