fork download
  1. #include <cstdio>
  2. #include <vector>
  3.  
  4. using int64 = long long;
  5.  
  6. namespace solve {
  7. const int N = 1000000 + 10;
  8. std::vector<int> divs[N];
  9. std::vector<std::pair<int, int>> edges[N];
  10. int p[N], m;
  11. int64 f[N], mu[N];
  12. void init() {
  13. for (int i = 1; i < N; ++i) {
  14. for (int j = i; j < N; j += i) {
  15. divs[j].push_back(i);
  16. }
  17. }
  18. mu[1] = 1;
  19. m = 0;
  20. for (int i = 2; i < N; ++i) {
  21. if(!p[i]) p[m++] = i, mu[i] = -1;
  22. for (int j = 0; i * p[j] < N && j < m; ++j) {
  23. p[i * p[j]] = 1;
  24. if (i % p[j]) mu[i * p[j]] = -mu[i];
  25. else {
  26. mu[i * p[j]] = 0;
  27. break;
  28. }
  29. }
  30. }
  31. }
  32. void run(int n) {
  33. for (int i = 1; i <= n; ++i) {
  34. edges[i].clear();
  35. f[i] = n / i;
  36. }
  37. int64 ret = 0, A = 0, B = 0;
  38. for (int i = 1; i <= n; ++i) {
  39. ret += mu[i] * f[i] * f[i] * f[i];
  40. }
  41. std::vector<int> from, to, wei;
  42. m = 0;
  43. std::vector<int> deg(n + 1);
  44. for (int i = 1; i <= n; ++i) {
  45. if (!mu[i]) continue;
  46. for (int j = i; j <= n; j += i) {
  47. if (mu[j] && mu[j / i]) {
  48. for (auto &&k: divs[j / i]) {
  49. int x = i * k, y = j / k;
  50. if (x >= y) continue;
  51. deg[x]++;
  52. deg[y]++;
  53. B += (mu[x] * f[y] + mu[y] * f[x]) * f[j] * f[j];
  54. from.push_back(x);
  55. to.push_back(y);
  56. wei.push_back(f[j]);
  57. ++m;
  58. }
  59. }
  60. }
  61. }
  62. for (int i = 0; i < m; ++i) {
  63. int u = from[i], v = to[i], w = wei[i];
  64. if (deg[u] < deg[v]) edges[u].emplace_back(v, w);
  65. else edges[v].emplace_back(u, w);
  66. }
  67. std::vector<int> mark(n + 1, -1);
  68. for (int a = 1; a <= n; ++a) {
  69. for (auto &&e: edges[a]) mark[e.first] = e.second;
  70. for (auto &&e: edges[a]) {
  71. int b = e.first, w1 = e.second;
  72. for (auto &&ee: edges[b]) {
  73. int c = ee.first, w2 = ee.second;
  74. if (mark[c] != -1) {
  75. A += mu[a] * mu[b] * mu[c] * w1 * w2 * mark[c];
  76. }
  77. }
  78. }
  79. for (auto &&e: edges[a]) mark[e.first] = -1;
  80. }
  81. ret += A * 6 + B * 3;
  82. printf("%lld\n", ret);
  83. }
  84. }
  85.  
  86. int main() {
  87. solve::init();
  88. int n;
  89. scanf("%d", &n);
  90. solve::run(n);
  91. return 0;
  92. }
  93.  
Success #stdin #stdout 6.31s 400528KB
stdin
1000000
stdout
286747058490057124