fork download
copy 
  1. /// @title: Count non negative integer tuple (a, b, c) satisfied (sum <= S) and (product <= T)
  2. /// @testing: for large S <= 1e18, T <= 1e12
  3. /// > in O(T^(5/9)) time complexity = 830ms codeforces = 310ms ideone
  4. /// > in O(const) memory complexity
  5. /// > in O(3700Mb) code memory in 200 lines, 4444 characters
  6. ///
  7. /// @date:    23/08/2021
  8. /// @link:    https://i...content-available-to-author-only...e.com/yMi8tu
  9. /// @author:  SPyofgame
  10. /// @lisence: free lisence
  11.  
  12. #pragma GCC target ("avx2")
  13. #pragma GCC optimization ("O3")
  14. #pragma GCC optimization ("unroll-loops")
  15.  
  16. #include <iostream>
  17. #include <cmath>
  18.  
  19. typedef long long ll;
  20. const int MOD = 1e9 + 7;
  21.  
  22. ///====*====*====*====*====*====*====*====*====*====*====*====*====*====*====*====*====*====*====*====
  23.  
  24. /// Utility function
  25. ll min(ll a, ll b) { return a < b ? a : b; }
  26. ll max(ll a, ll b) { return a > b ? a : b; }
  27. ll square(ll x) { return x * x; }
  28.  
  29. /// Sum of (1 + 2 + ... + x)
  30. ll natural_sum(ll x)
  31. {
  32.     return (x <= 0) ? 0 : x * (x + 1) / 2;
  33. }
  34.  
  35. /// Sum of (l + (l + 1) + ... + r)
  36. ll natural_sum(ll l, ll r)
  37. {
  38.     return (l > r) ? 0 : natural_sum(r) - natural_sum(l - 1);
  39. }
  40.  
  41. /// Sigma(p = y -> x) floor(n / p) in O(cbrt(n))
  42. /// @link: https://a...content-available-to-author-only...v.org/pdf/1206.3369.pdf
  43. /// @author: Richard Sladkey
  44. ll fastsumdiv(ll y, ll x, ll n)
  45. {
  46.     if (y > x) return 0;
  47.  
  48.     ll S = 0;
  49.     ll B = n / (x + 1);
  50.     ll E = n % (x + 1);
  51.     ll D = n / x - B;
  52.     ll G = B - x * D;
  53.     ll d = 0;
  54.  
  55.     for (; x >= y; --x)
  56.     {
  57.         E += G;
  58.         if (E >= x)
  59.         {
  60.             D += 1, G -= x, E -= x;
  61.             if (E >= x)
  62.             {
  63.                 D += 1, G -= x, E -= x;
  64.                 if (E >= x) break; /// not likely to happen more
  65.             }
  66.         }
  67.         else if (E < 0)
  68.         {
  69.             D -= 1, G += x, E += x;
  70.         }
  71.  
  72.         G += 2 * D, B += D, S += B;
  73.     }
  74.  
  75.     E = n % (x + 1);
  76.     D = n / x - B;
  77.     G = B - x * D;
  78.     for (; x >= y; --x)
  79.     {
  80.         E += G;
  81.         d = E / x;
  82.         D += d;
  83.         E -= x * d;
  84.         G += 2 * D - x * d, B += D, S += B;
  85.     }
  86.  
  87.     for (; x >= y; --x)
  88.     {
  89.         S += n / x;
  90.     }
  91.  
  92.     return S;
  93. }
  94.  
  95. /// you can modify to -> Bignum / Modulo / Overflow
  96. void add(ll &res, ll val)
  97. {
  98.     val %= MOD;
  99.     res += val;
  100.     if (res >= MOD) res -= MOD;
  101. }
  102.  
  103. /// you can modify to -> Bignum / Modulo / Overflow
  104. void sub(ll &res, ll val)
  105. {
  106.     val %= MOD;
  107.     res -= val;
  108.     if (res < 0) res += MOD;
  109. }
  110.  
  111. /// Count (a, b, c) satisfied
  112. /// { min(a, b, c) >= 0
  113. /// { a + b + c <= S
  114. /// { a * b * c <= T
  115. /// > O(T^(5/9)) complexity
  116. ///
  117. /// @proof: https://m...content-available-to-author-only...e.com/questions/4230187/faster-algorithm-for-counting-non-negative-tripplea-b-c-satisfied-a-b-c
  118. /// @author: SPyofgame
  119. ll solve_ABC(ll S, ll T)
  120. {
  121.     ll cbT = cbrt(T); /// Not very accuracy
  122.     while (cbT * cbT * cbT < T) ++cbT;
  123.     while (cbT * cbT * cbT > T) --cbT;
  124.  
  125.  
  126.  
  127. /// [1] 0 <= a < b < c && a <= cbrt(T) -> cnt1 * 6
  128.     ll cnt1 = 0;
  129.     add(cnt1, (S / 2) * S - 2 * natural_sum(S / 2)); /// a = 0
  130.  
  131.     ll k = 1;
  132.     for (ll a = 1, upa = min(S, cbT); a <= upa; ++a) /// a > 0 -> count(b < c)
  133.     {
  134.         ll SSS = S - a;
  135.         ll TTT = T / a;
  136.         ll KKK = min(SSS / 2, min((ll)sqrt(TTT), TTT / 2 + 1));
  137.  
  138.         /// Binary search for max k satisfy (S - a - k) <= (T / a / k)
  139.         ll k = a;
  140.         for (ll l = a + 1, r = KKK; l <= r; )
  141.         {
  142.             ll m = (l + r) >> 1;
  143.             ll v = TTT / m;
  144.             if (SSS - m <= v)
  145.             {
  146.                 k = m;
  147.                 l = m + 1;
  148.             }
  149.             else
  150.             {
  151.                 r = m - 1;
  152.             }
  153.         }
  154.  
  155.         sub(cnt1, natural_sum(a + 1, KKK));
  156.         add(cnt1, SSS * (k - a) - natural_sum(a + 1, k));
  157.         add(cnt1, fastsumdiv(k + 1, KKK, TTT));
  158.     }
  159.  
  160.  
  161.  
  162. /// [2] 0 <= a < b = c && a <= cbrt(T) -> cnt2 * 3
  163.     ll cnt2 = 0;
  164.     add(cnt2, S / 2); /// a = 0
  165.     for (ll a = 1, upa = min(S, cbT); a <= upa; ++a) /// a > 0 -> count(b = c)
  166.     {
  167.         add(cnt2, max(0LL, min((S - a) / 2, ll(floor(sqrt(T / a)))) - a));
  168.     }
  169.  
  170.  
  171.  
  172. /// [3] 0 <= a = b < c && a <= cbrt(T) -> cnt3 * 3
  173.     ll cnt3 = 0;
  174.     add(cnt3, S); /// a = 0
  175.     for (ll a = 1, upa = min(S / 2, cbT); a <= upa; ++a) /// a > 0 -> count(a < c)
  176.     {
  177.         add(cnt3, max(0LL, min(S - a - a, T / a / a) - a));
  178.     }
  179.  
  180.  
  181.  
  182. /// [4] 0 <= a = b = c && a <= cbrt(T) -> cnt4 * 1
  183.     ll cnt4 = 0;
  184.     add(cnt4, min(S / 3, cbT) + 1); /// a = b = c >= 0
  185.  
  186.  
  187.  
  188. /// Final result: total counting
  189.     ll res = 0;
  190.     add(res, cnt1 * 6 + cnt2 * 3 + cnt3 * 3 + cnt4 * 1);
  191.     return res;
  192. }
  193.  
  194. /// Main function
  195. int main()
  196. {
  197. 	/// Input any number 0 <= S, T <= 1e18
  198.     std::cout << solve_ABC(1000000000000000000LL, 1000000000000LL);
  199.     return 0;
  200. }
  201.  
Success #stdin #stdout 0.31s 5640KB
comments (?)
stdin
copy 
Standard input is empty
stdout
copy 
533160243
https://ideone.com/yMi8tu
Count non negative integer tuple (a, b, c) satisfied (sum <= S) and (product <= T)
language:
C++14 (gcc 8.3)
created:
2 years ago
visibility:
public

Share or Embed source code

Discover > Sphere Engine API

The brand new service which powers Ideone!

Discover > IDE Widget

Widget for compiling and running the source code in a web browser!