Nhận xét: Vì ~4~ số ~a_1,a_2,b_1,b_2~ đôi một khác nhau nên số bộ số ~a_1,a_2,b_1,b_2~ sao cho ~a_1 \times a_2 > b_1 \times b_2~ bằng số bộ số ~a_1,a_2,b_1,b_2~ sao cho ~a_1 \times a_2 < b_1 \times b_2~.

Vì vậy, ta sử dụng bao hàm loại trừ bằng cách đếm số bộ số ~a_1,a_2,b_1,b_2~ sao cho ~a_1 \times a_ 2 = b_1 \times b_2~, sử dụng tổng số bộ ~4~ có thể, trừ đi và chia cho ~2~ sẽ là đáp án.

Ta có: ~a_1 \times a_2 = b_1 \times b_2~ ~\leftrightarrow~ ~\frac{a_1}{b_1} = \frac{b_2}{a_2}~

Đặt ~\frac{a_1}{b_1} = \frac{b_2}{a_2} = \frac{x}{y}~, với ~\frac{x}{y}~ là một phân số tối giản.

Bài toán: đếm số bộ ~a_1,a_2,b_1,b_2~ sao cho ~\frac{a_1}{b_1} = \frac{b_2}{a_2} = \frac{x}{y}~ và ~4~ số phân biệt nhau.

Ta thấy: ~\frac{x}{y}~ là phân số tối giản, khi và chỉ khi ~x~ và ~y~ nguyên tố cùng nhau. Có ~2~ trường hợp:

• ~x < y~, thì số lượng số ~x~ sao cho ~\frac{x}{y}~ là phân số tối giản là ~phi(y)~. Với ~y~ xác định, ~y > x~ nên số phân số ~\frac{a}{b}~ sao cho ~a, b \leq N~ và ~\frac{a}{b} = \frac{x}{y}~ luôn xác định, vì vậy có thể đếm số bộ số trong ~O(1)~.

• ~x > y~, thì số lượng số ~y~ sao cho ~\frac{x}{y}~ là phân số tối giản là ~phi(x)~. Với ~x~ xác định thì cũng sẽ tính được tương tự như trường hợp trên.

Nhưng như vậy thì sẽ có một trường hợp bị sai: Nếu ~\frac{a_1}{b_1} = \frac{b_2}{a_2}~ ~(1)~ thì sẽ có thể xảy ra trường hợp ~b_1 = b_2~, không thoả mãn điều kiện đề bài. Ta sẽ cần loại trừ số bộ ~4~ như vậy.

Sau đó, thế ~b_2 = b_1~ vào ~(1)~, ta có: ~\frac{a_1}{b_1} = \frac{b_1}{a_2}~ ~\leftrightarrow~ ~a_1 \times a_2 = {b_1}^2~

Lúc này ta sẽ có một bài toán như sau: Đếm số cặp số ~(a_1,a_2)~ sao cho tích ~a_1 \times a_2~ là số chính phương. Đây là bài khá cổ điển và có thể giải bằng sàng nguyên tố.

Code mẫu

#include <bits/stdc++.h>

using namespace std;

#define fo(i, a, b) for(int i = a; i <= b; i++)
#define _fo(i, a, b) for(int i = a; i >= b; i--)
#define foa(i, a) for (auto &i : a)
#define sz(a) ((int) a.size())
#define all(a) begin(a), end(a)
#define fi first
#define se second
#define pb(x) push_back(x)
#define mk(x, y) make_pair(x, y)

typedef long long ll;
typedef pair<ll, ll> pll;
typedef vector<ll> vl;
typedef pair<int, int> pii;
typedef vector<int> vi;

const int LOG = 22; 
const int MAX = 1e7+5;
const int MOD = 1e9+7;
const int SQRT = 400;
const ll INF = 2e9;
const ll lon = 1e18;

const ll inv2 = 500000004;

ll add(ll a, ll b) { return (a+b) % MOD; }
ll sub(ll a, ll b) { return (a-b+MOD) % MOD; }
ll mul(ll a, ll b) { return (a % MOD)*(b % MOD) % MOD; }
ll pick(ll k) { return mul(mul(k, k-1), inv2); }

int n;
int sieve[MAX], phi[MAX], mask[MAX], sqr[MAX];

void precompute() {
    for(ll p = 2; p <= n; p++) {
        if(sieve[p] == 0) {
            sieve[p] = p; 
            for(ll i = p*p; i <= n; i += p) {
                if(sieve[i] == 0) sieve[i] = p;
            }
        }
    }

    phi[1] = 1;
    for(ll i = 2; i <= n; i++) {
        ll prev = i/sieve[i];
        if(prev % sieve[i] == 0) phi[i] = phi[prev]*sieve[i];
        else phi[i] = phi[prev]*(sieve[i]-1);
    }

    mask[1] = 1;
    for(ll i = 2; i <= n; i++) {
        ll temp = sieve[i]*sieve[i];
        if(i % temp == 0) mask[i] = mask[i/temp];
        else mask[i] = mask[i/sieve[i]]*sieve[i];
    }
}

int solve() {
    int equal = 0;
    fo(i, 2, n) equal = add(equal, mul(phi[i], pick(n/i)));
    fo(i, 1, n) {
        equal = sub(equal, sqr[mask[i]]);
        sqr[mask[i]]++; 
    } 
    return mul(sub(mul(pick(n), pick(n-2)), equal), inv2);
}

signed main() {
    ios_base::sync_with_stdio(0); cin.tie(0); cout.tie(0);

    cin >> n;
    if(n < 4) cout << 0;
    else {
        precompute();               
        cout << solve();
    }
}