題目定義的函式 \(a\) 更常見的記法為 \(\Omega\)(OEIS A008472,容易說明其為加性函式),後文亦使用 \(\Omega\)。
為了好看,後文用 \(a_i\) 代指題目中的 \(m_i\)。
題意即轉化為:
給定整數 \(n, a_1, k\) 和長為 \(n\) 的整數序列 \(w\),設
\[a_i = \max\limits_{1 \leqslant j < i} \left\{ \frac {a_j} {\Omega(\operatorname{lcm}(w_i, w_j)) + \Omega(\gcd(w_i, w_j))} \right\} + k \qquad 2 \leqslant i \leqslant n \]求
\[\sum\limits_{i = 1}^n a_i \]\(1 \leqslant n \leqslant 182376, 2 \leqslant w_i \leqslant 182376\)。
對 \(2 \leqslant i \leqslant n\),由(因為 \(xy\) 和 \(\operatorname{lcm}(x, y)\) 的質因子相同)
\[\Omega(xy) = \Omega(\operatorname{lcm}(x, y))
\]
得
\[a_i = \max\limits_{1 \leqslant j < i} \left\{ \frac {a_j} {\Omega(w_i w_j) + \Omega(\gcd(w_i, w_j))} \right\} + k
\]
由(容斥質因子可得)
\[\Omega(xy) = \Omega(x) + \Omega(y) - \Omega(\gcd(x, y))
\]
得
\[\begin{align*}
a_i &= \max\limits_{1 \leqslant j < i} \left\{ \frac {a_j} {\Omega(w_i) + \Omega(w_j) - \Omega(\gcd(w_i, w_j)) + \Omega(\gcd(w_i, w_j))} \right\} + k \\
&= \max\limits_{1 \leqslant j < i} \left\{ \frac {a_j} {\Omega(w_i) + \Omega(w_j)} \right\} + k \\
\end{align*}
\]
發現 \(\Omega(w_i)\) 是定值,但位置在分母,透過取倒數讓它移動到分子,即
\[\begin{align*}
a_i &= \max\limits_{1 \leqslant j < i} \left\{ \frac {a_j} {\Omega(w_i) + \Omega(w_j)} \right\} + k \\
&= \frac 1 {\min\limits_{1 \leqslant j < i} \left\{ \frac {\Omega(w_i) + \Omega(w_j)} {a_j} \right\}} + k \\
\end{align*}
\]
考慮求
\[\min\limits_{1 \leqslant j < i} \left\{ \frac {\Omega(w_i) + \Omega(w_j)} {a_j} \right\}
\]
發現這可以看作關於 \(\Omega(w_i)\) 的一次函式,斜率為 \(\frac 1 {a_j}\),截距為 \(\frac {\Omega(w_j)} {a_j}\)。
因為 \(m_i\) 一定是正數,所以可能成為答案的一次函式的影像構成一個上凸包,李超線段樹或者平衡樹維護凸包即可。
時間複雜度 \(O(n \log w)\),\(w\) 為值域。
基於李超線段樹的參考實現如下,注意本題對精度要求較高。
#include <iomanip>
#include <iostream>
using namespace std;
typedef long double ld;
const int lim = 500 * 365;
const ld inf = 1e18;
const ld eps = 1e-15;
bool vis[lim + 5];
int pr[lim + 5], tail;
int Omega[lim + 5];
int n, k;
int w[lim];
ld a[lim];
struct line {
ld k, b;
inline ld get(int x) const { return k * x + b; }
} d[800005];
static inline void build(int s, int t, int p) {
d[p].b = inf;
if (s == t)
return;
int mid = (s + t) >> 1;
build(s, mid, p << 1);
build(mid + 1, t, p << 1 | 1);
}
static inline void insert(line nd, int s, int t, int p) {
int mid = (s + t) >> 1;
line l1 = nd;
line &l2 = d[p];
if (l1.get(mid) - l2.get(mid) < -eps)
swap(l1, l2);
if (l1.get(s) - l2.get(s) < -eps)
insert(l1, s, mid, p << 1);
if (l1.get(t) - l2.get(t) < -eps)
insert(l1, mid + 1, t, p << 1 | 1);
}
static inline ld query(int x, int s, int t, int p) {
ld ret = d[p].get(x);
if (s == t)
return ret;
int mid = (s + t) >> 1;
if (x <= mid)
return min(ret, query(x, s, mid, p << 1));
else
return min(ret, query(x, mid + 1, t, p << 1 | 1));
}
static inline line make_line(int x) { return {1. / a[x], Omega[w[x]] / a[x]}; }
static inline void solve() {
cin >> n >> a[1] >> k;
for (int i = 1; i <= n; ++i)
cin >> w[i];
build(0, 2e5, 1);
insert(make_line(1), 0, 2e5, 1);
for (int i = 2; i <= n; ++i) {
ld inv = query(Omega[w[i]], 0, 2e5, 1);
a[i] = 1. / inv + k;
insert(make_line(i), 0, 2e5, 1);
}
ld ans = 0;
for (int i = 1; i <= n; ++i)
ans += a[i];
cout << fixed << setprecision(12) << ans << endl;
}
signed main() {
#ifndef ONLINE_JUDGE
freopen("P11014.in", "r", stdin);
#endif
ios::sync_with_stdio(false);
cin.tie(0);
cout.tie(0);
Omega[1] = 0;
for (int i = 2; i <= lim; ++i) {
if (!vis[i]) {
pr[++tail] = i;
Omega[i] = i;
}
for (int j = 1; j <= tail && i * pr[j] <= lim; ++j) {
vis[i * pr[j]] = 1;
if (i % pr[j] == 0) {
Omega[i * pr[j]] = Omega[i];
break;
}
Omega[i * pr[j]] = Omega[i] + pr[j];
}
}
solve();
return 0;
}