给定 \(n\) 和长度为 \(n\) 的正整数序列 \(a\)。
定义一次操作为，选择一个数 \(a_i\) 和质数 \(p\)，将 \(a_i\) 赋值为 \(a_i\times p\) 或 \(\frac{a_i}{p}\)，必须满足赋值后 \(a_i\) 为正整数。
定义区间 \(l,r\) 的代价 \(val(l,r)\) 为，使 \(a_l,a_{l+1},\cdots,a_r\) 全都相等的最小操作次数。
请求出 \(\sum_{l=1}^n \sum_{r=l}^n val(l,r)\)，对 \(10^9+7\) 取模。
对于 \(100\%\) 的数据，\(1\le n\le 2\times 10^5\)，\(1\le a_i\le 10^6\)。
口胡题解：
先转化一下式子：
假设 \(b_i\) 为 \(a_i\) 的 \(p\) 次方,这有两个结论：
和
那么其中 \(f\) 计算方法为：
我们设 \(sum_r=\sum^{r}_{i=1}([b_i\le k]-[b_i> k])\) 考虑在 \(p,k\) 确定的情况下，也就是 \(sum\) 确定的情况下快速计算 \(\sum^{n}_{l=1}\sum^{n}_{r=l}f(l,r,p,k)\).
考虑先对 \(sum\) 排序，前面的 \((r-l+1)\) 当成区间长度，这个求和是好做的，我们记为 \(len\)，此时我们可以把式子变成：
换个方式，就是：
这样计算复杂度就降低了，考虑现在复杂度：发现是 \(O(Pn\log n\log V)\)，其中 \(P\) 是质数个数，算出来不超过 \(8\times 10^4\)，考虑优化。
接下来是真口胡了，我不知道具体怎么做。
看看怎么优化这个东西，因为每个数都由 \(O(\log V)\) 个质数构成，所以发现 \(sum_i\) 的非 \(-1\) 的变化次数是 \(O(\frac{n\log V}{P})\) 的。也就是说我们先给每个数 \(-i\)，那么就只有 \(O(\frac{n\log V}{P})\) 个值了，那么直接找出所有的变化点，相当于把 \(sum_i\) 加权，接着做这个式子，现在复杂度就是 \(O(n\log^2V\log(\log\frac{n\log V}{P}))\)，能过了。
官方题解：
每个质数分别计算，问题变成每次将一个数加 1 或减 1，使所有数都相等的最小操作次数。
可以发现变为中位数次数最小，考虑计算。
可能要自行脑补一档部分分，只有 \(01\) 怎么做？
那就是区间答案就是 \(min\) 0,1 的个数。
可以从前到后扫描线，将 0 看做 \(-1\)，利用前缀和来计算。
用线段树可以做到 \(O(nlogn)\)，但考虑每次前缀和的变化量为 1，可以开桶统计每种前缀和的个数，做到 \(O(n)\) 。
对于一般情况，枚举分界 \(x\) ，将 \(\le x\) 的数看做 \(0\)，其余看做 \(1\)，然后将所有 \(x\) 的答案加起来，可以发现这样是对的。
于是 \(a_i=2^k\) 就做完了，但对于一般情况，不能每次都 \(O(n)\) 计算。
发现非 \(0\) 数的总和是 \(nlogV\) 的，考虑如过以某个点为端点的区间中位数一定为零，那么可以直接计算，称这些位置是无效的。
设非 \(0\) 的数有 \(m\) 个，那么有效位置不超过 \(3m\) 个。
考虑如何求出有效的位置，枚举非零数，然后向左向右扩展，需要记录前后缀的某个最值。
std 实现的比较随意，时间复杂度没有仔细分析，但不超过 \(O(nlog^2V)\)。
#include<bits/stdc++.h>
using namespace std;
const int N=1e6+5,P=1e9+7,i2=P+1>>1;
int n,ans,c[N],d[N<<1];
vector<int>p[N];
bool vis[N],ban[N];
inline int md(int x){return x>=P?x-P:x;
}
int calc(vector<pair<int,int> >&a,int lim){int res=0;for(int i=0,j,s;j=i,i<a.size();i=j+1){c[s=1]=(a[i].second>=lim?1:-1);while(j+1<a.size()&&a[j+1].first==a[j].first+1)++j,c[++s]=(a[j].second>=lim?1:-1);int cl=0,s1=0,s2=0,c1=0,c2=1;d[n]=1;for(int k=1;k<=s;++k){c[k]+=c[k-1],cl=md(cl+k);if(c[k]<c[k-1]){s1=(1ll*(P-c[k])*d[c[k]+n]+s1)%P;c1-=d[c[k]+n];s2=(1ll*(c[k]+P)*d[c[k]+n]+s2)%P;c2+=d[c[k]+n];}else{s2=(1ll*(P-c[k-1])*d[c[k-1]+n]+s2)%P;c2-=d[c[k-1]+n];s1=(1ll*(c[k-1]+P)*d[c[k-1]+n]+s1)%P;c1+=d[c[k-1]+n];}res=(1ll*(P-c[k])*c1+s1+1ll*(c[k]+P)*c2+P-s2+res+cl)%P;s2=md(s2+md(c[k]+P)),++d[c[k]+n],++c2;}for(int i=0;i<=s;++i)d[c[i]+n]=0;}res=1ll*res*i2%P;return res;
}
int solve(vector<pair<int,int> >&a){sort(a.begin(),a.end());int sz=a.size(),mx=1;for(int i=0,o=1e9;i<sz;++i){o=min(o,2*i-a[i].first);vis[a[i].first]=true,mx=max(mx,a[i].second);for(int j=a[i].first+1,ed=(i==sz-1)?n:a[i+1].first-1;j<=ed;++j){if(2*i-j+1<o)break;vis[j]=true;a.push_back(make_pair(j,0));}}for(int i=sz-1,o=-1e9;~i;--i){o=max(o,2*i-a[i].first);for(int j=a[i].first-1,ed=i?a[i-1].first+1:1;j>=ed;--j){if(2*i-j-1>o||vis[j])break;vis[j]=true;a.push_back(make_pair(j,0));}}sort(a.begin(),a.end());int res=0;for(int i=1;i<=mx;++i)res=md(res+calc(a,i));for(int i=0,j;j=i,i<a.size();i=j+1){while(j+1<a.size()&&a[j+1].first==a[j].first+1)++j;for(int k=i,w;k<=j;++k)if(a[k].second){w=(1ll*(a[i].first-1)*(n-a[k].first+1)+1ll*(n-a[j].first)*a[k].first-1ll*(a[i].first-1)*(n-a[j].first))%P;res=(1ll*w*a[k].second+res)%P;}}for(auto v:a)vis[v.first]=false;return res;
}
int main(){freopen("unnamed.in","r",stdin);freopen("unnamed.out","w",stdout);scanf("%d",&n);for(int i=1,x;i<=n;++i){scanf("%d",&x);p[x].push_back(i);}for(int i=2;i<=1000000;++i)if(!ban[i]){long long x=i;vector<pair<int,int> >tmp;for(int j=i;j<=1000000;j+=i)ban[j]=true;for(int j=1;x<=1000000;x*=i,++j)for(int k=x;k<=1000000;k+=x)if((k/x)%i)for(auto v:p[k])tmp.push_back(make_pair(v,j));if(!tmp.empty())ans=md(ans+solve(tmp));}printf("%d\n",ans);return 0;
}