Advertisement
cosenza987

ss

Nov 5th, 2023
769
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 18.89 KB | None | 0 0
  1. const int BASE_DIGITS = 9;
  2. const int BASE = 1000000000;
  3.  
  4. struct BigInt {
  5.     int sign;
  6.     vector<int> a;
  7.  
  8.     // -------------------- Constructors --------------------
  9.     // Default constructor.
  10.     BigInt() : sign(1) {}
  11.  
  12.     // Constructor from long long.
  13.     BigInt(long long v) { *this = v; }
  14.     BigInt& operator=(long long v) {
  15.         sign = 1;
  16.         if (v < 0) {
  17.             sign = -1;
  18.             v = -v;
  19.         }
  20.         a.clear();
  21.         for (; v > 0; v = v / BASE) a.push_back(v % BASE);
  22.         return *this;
  23.     }
  24.  
  25.     // Initialize from string.
  26.     BigInt(const string& s) { read(s); }
  27.  
  28.     // -------------------- Input / Output --------------------
  29.     void read(const string& s) {
  30.         sign = 1;
  31.         a.clear();
  32.         int pos = 0;
  33.         while (pos < (int)s.size() && (s[pos] == '-' || s[pos] == '+')) {
  34.             if (s[pos] == '-') sign = -sign;
  35.             ++pos;
  36.         }
  37.         for (int i = s.size() - 1; i >= pos; i -= BASE_DIGITS) {
  38.             int x = 0;
  39.             for (int j = max(pos, i - BASE_DIGITS + 1); j <= i; j++)
  40.                 x = x * 10 + s[j] - '0';
  41.             a.push_back(x);
  42.         }
  43.         trim();
  44.     }
  45.     friend istream& operator>>(istream& stream, BigInt& v) {
  46.         string s;
  47.         stream >> s;
  48.         v.read(s);
  49.         return stream;
  50.     }
  51.  
  52.     friend ostream& operator<<(ostream& stream, const BigInt& v) {
  53.         if (v.sign == -1 && !v.isZero()) stream << '-';
  54.         stream << (v.a.empty() ? 0 : v.a.back());
  55.         for (int i = (int)v.a.size() - 2; i >= 0; --i)
  56.             stream << setw(BASE_DIGITS) << setfill('0') << v.a[i];
  57.         return stream;
  58.     }
  59.  
  60.     // -------------------- Comparison --------------------
  61.     bool operator<(const BigInt& v) const {
  62.         if (sign != v.sign) return sign < v.sign;
  63.         if (a.size() != v.a.size())
  64.             return a.size() * sign < v.a.size() * v.sign;
  65.         for (int i = ((int)a.size()) - 1; i >= 0; i--)
  66.             if (a[i] != v.a[i]) return a[i] * sign < v.a[i] * sign;
  67.         return false;
  68.     }
  69.  
  70.     bool operator>(const BigInt& v) const { return v < *this; }
  71.     bool operator<=(const BigInt& v) const { return !(v < *this); }
  72.     bool operator>=(const BigInt& v) const { return !(*this < v); }
  73.     bool operator==(const BigInt& v) const {
  74.         return !(*this < v) && !(v < *this);
  75.     }
  76.     bool operator!=(const BigInt& v) const { return *this < v || v < *this; }
  77.  
  78.     // Returns:
  79.     // 0 if |x| == |y|
  80.     // -1 if |x| < |y|
  81.     // 1 if |x| > |y|
  82.     friend int __compare_abs(const BigInt& x, const BigInt& y) {
  83.         if (x.a.size() != y.a.size()) {
  84.             return x.a.size() < y.a.size() ? -1 : 1;
  85.         }
  86.  
  87.         for (int i = ((int)x.a.size()) - 1; i >= 0; --i) {
  88.             if (x.a[i] != y.a[i]) {
  89.                 return x.a[i] < y.a[i] ? -1 : 1;
  90.             }
  91.         }
  92.         return 0;
  93.     }
  94.  
  95.     // -------------------- Unary operator - and operators +-
  96.     // --------------------
  97.     BigInt operator-() const {
  98.         BigInt res = *this;
  99.         if (isZero()) return res;
  100.  
  101.         res.sign = -sign;
  102.         return res;
  103.     }
  104.  
  105.     // Note: sign ignored.
  106.     void __internal_add(const BigInt& v) {
  107.         if (a.size() < v.a.size()) {
  108.             a.resize(v.a.size(), 0);
  109.         }
  110.         for (int i = 0, carry = 0; i < (int)max(a.size(), v.a.size()) || carry;
  111.             ++i) {
  112.             if (i == (int)a.size()) a.push_back(0);
  113.  
  114.             a[i] += carry + (i < (int)v.a.size() ? v.a[i] : 0);
  115.             carry = a[i] >= BASE;
  116.             if (carry) a[i] -= BASE;
  117.         }
  118.     }
  119.  
  120.     // Note: sign ignored.
  121.     void __internal_sub(const BigInt& v) {
  122.         for (int i = 0, carry = 0; i < (int)v.a.size() || carry; ++i) {
  123.             a[i] -= carry + (i < (int)v.a.size() ? v.a[i] : 0);
  124.             carry = a[i] < 0;
  125.             if (carry) a[i] += BASE;
  126.         }
  127.         this->trim();
  128.     }
  129.  
  130.     BigInt operator+=(const BigInt& v) {
  131.         if (sign == v.sign) {
  132.             __internal_add(v);
  133.         }
  134.         else {
  135.             if (__compare_abs(*this, v) >= 0) {
  136.                 __internal_sub(v);
  137.             }
  138.             else {
  139.                 BigInt vv = v;
  140.                 swap(*this, vv);
  141.                 __internal_sub(vv);
  142.             }
  143.         }
  144.         return *this;
  145.     }
  146.  
  147.     BigInt operator-=(const BigInt& v) {
  148.         if (sign == v.sign) {
  149.             if (__compare_abs(*this, v) >= 0) {
  150.                 __internal_sub(v);
  151.             }
  152.             else {
  153.                 BigInt vv = v;
  154.                 swap(*this, vv);
  155.                 __internal_sub(vv);
  156.                 this->sign = -this->sign;
  157.             }
  158.         }
  159.         else {
  160.             __internal_add(v);
  161.         }
  162.         return *this;
  163.     }
  164.  
  165.     // Optimize operators + and - according to
  166.     // https://stackoverflow.com/questions/13166079/move-semantics-and-pass-by-rvalue-reference-in-overloaded-arithmetic
  167.     template <typename L, typename R>
  168.     typename std::enable_if<std::is_convertible<L, BigInt>::value&&
  169.         std::is_convertible<R, BigInt>::value&&
  170.         std::is_lvalue_reference<R&&>::value,
  171.         BigInt>::type friend
  172.         operator+(L&& l, R&& r) {
  173.         BigInt result(std::forward<L>(l));
  174.         result += r;
  175.         return result;
  176.     }
  177.     template <typename L, typename R>
  178.     typename std::enable_if<std::is_convertible<L, BigInt>::value&&
  179.         std::is_convertible<R, BigInt>::value&&
  180.         std::is_rvalue_reference<R&&>::value,
  181.         BigInt>::type friend
  182.         operator+(L&& l, R&& r) {
  183.         BigInt result(std::move(r));
  184.         result += l;
  185.         return result;
  186.     }
  187.  
  188.     template <typename L, typename R>
  189.     typename std::enable_if<std::is_convertible<L, BigInt>::value&&
  190.         std::is_convertible<R, BigInt>::value,
  191.         BigInt>::type friend
  192.         operator-(L&& l, R&& r) {
  193.         BigInt result(std::forward<L>(l));
  194.         result -= r;
  195.         return result;
  196.     }
  197.  
  198.     // -------------------- Operators * / % --------------------
  199.     friend pair<BigInt, BigInt> divmod(const BigInt& a1, const BigInt& b1) {
  200.         assert(b1 > 0);  // divmod not well-defined for b < 0.
  201.  
  202.         long long norm = BASE / (b1.a.back() + 1);
  203.         BigInt a = a1.abs() * norm;
  204.         BigInt b = b1.abs() * norm;
  205.         BigInt q = 0, r = 0;
  206.         q.a.resize(a.a.size());
  207.  
  208.         for (int i = a.a.size() - 1; i >= 0; i--) {
  209.             r *= BASE;
  210.             r += a.a[i];
  211.             long long s1 = r.a.size() <= b.a.size() ? 0 : r.a[b.a.size()];
  212.             long long s2 =
  213.                 r.a.size() <= b.a.size() - 1 ? 0 : r.a[b.a.size() - 1];
  214.             long long d = ((long long)BASE * s1 + s2) / b.a.back();
  215.             r -= b * d;
  216.             while (r < 0) {
  217.                 r += b, --d;
  218.             }
  219.             q.a[i] = d;
  220.         }
  221.  
  222.         q.sign = a1.sign * b1.sign;
  223.         r.sign = a1.sign;
  224.         q.trim();
  225.         r.trim();
  226.         auto res = make_pair(q, r / norm);
  227.         if (res.second < 0) res.second += b1;
  228.         return res;
  229.     }
  230.     BigInt operator/(const BigInt& v) const { return divmod(*this, v).first; }
  231.  
  232.     BigInt operator%(const BigInt& v) const { return divmod(*this, v).second; }
  233.  
  234.     void operator/=(int v) {
  235.         assert(v > 0);  // operator / not well-defined for v <= 0.
  236.         if (llabs(v) >= BASE) {
  237.             *this /= BigInt(v);
  238.             return;
  239.         }
  240.         if (v < 0) sign = -sign, v = -v;
  241.         for (int i = (int)a.size() - 1, rem = 0; i >= 0; --i) {
  242.             long long cur = a[i] + rem * (long long)BASE;
  243.             a[i] = (int)(cur / v);
  244.             rem = (int)(cur % v);
  245.         }
  246.         trim();
  247.     }
  248.  
  249.     BigInt operator/(int v) const {
  250.         assert(v > 0);  // operator / not well-defined for v <= 0.
  251.  
  252.         if (llabs(v) >= BASE) {
  253.             return *this / BigInt(v);
  254.         }
  255.         BigInt res = *this;
  256.         res /= v;
  257.         return res;
  258.     }
  259.     void operator/=(const BigInt& v) { *this = *this / v; }
  260.  
  261.     long long operator%(long long v) const {
  262.         assert(v > 0);  // operator / not well-defined for v <= 0.
  263.         assert(v < BASE);
  264.         int m = 0;
  265.         for (int i = a.size() - 1; i >= 0; --i)
  266.             m = (a[i] + m * (long long)BASE) % v;
  267.         return m * sign;
  268.     }
  269.  
  270.     void operator*=(int v) {
  271.         if (llabs(v) >= BASE) {
  272.             *this *= BigInt(v);
  273.             return;
  274.         }
  275.         if (v < 0) sign = -sign, v = -v;
  276.         for (int i = 0, carry = 0; i < (int)a.size() || carry; ++i) {
  277.             if (i == (int)a.size()) a.push_back(0);
  278.             long long cur = a[i] * (long long)v + carry;
  279.             carry = (int)(cur / BASE);
  280.             a[i] = (int)(cur % BASE);
  281.             // asm("divl %%ecx" : "=a"(carry), "=d"(a[i]) : "A"(cur),
  282.             // "c"(base));
  283.             /*
  284.              int val;
  285.              __asm {
  286.              lea esi, cur
  287.              mov eax, [esi]
  288.              mov edx, [esi+4]
  289.              mov ecx, base
  290.              div ecx
  291.              mov carry, eax
  292.              mov val, edx;
  293.              }
  294.              a[i] = val;
  295.              */
  296.         }
  297.         trim();
  298.     }
  299.  
  300.     BigInt operator*(int v) const {
  301.         if (llabs(v) >= BASE) {
  302.             return *this * BigInt(v);
  303.         }
  304.         BigInt res = *this;
  305.         res *= v;
  306.         return res;
  307.     }
  308.  
  309.     // Convert BASE 10^old --> 10^new.
  310.     static vector<int> convert_base(const vector<int>& a, int old_digits,
  311.         int new_digits) {
  312.         vector<long long> p(max(old_digits, new_digits) + 1);
  313.         p[0] = 1;
  314.         for (int i = 1; i < (int)p.size(); i++) p[i] = p[i - 1] * 10;
  315.         vector<int> res;
  316.         long long cur = 0;
  317.         int cur_digits = 0;
  318.         for (int i = 0; i < (int)a.size(); i++) {
  319.             cur += a[i] * p[cur_digits];
  320.             cur_digits += old_digits;
  321.             while (cur_digits >= new_digits) {
  322.                 res.push_back((long long)(cur % p[new_digits]));
  323.                 cur /= p[new_digits];
  324.                 cur_digits -= new_digits;
  325.             }
  326.         }
  327.         res.push_back((int)cur);
  328.         while (!res.empty() && !res.back()) res.pop_back();
  329.         return res;
  330.     }
  331.  
  332.     void fft(vector<complex<double> >& a, bool invert) const {
  333.         int n = (int)a.size();
  334.  
  335.         for (int i = 1, j = 0; i < n; ++i) {
  336.             int bit = n >> 1;
  337.             for (; j >= bit; bit >>= 1) j -= bit;
  338.             j += bit;
  339.             if (i < j) swap(a[i], a[j]);
  340.         }
  341.  
  342.         for (int len = 2; len <= n; len <<= 1) {
  343.             double ang = 2 * 3.14159265358979323846 / len * (invert ? -1 : 1);
  344.             complex<double> wlen(cos(ang), sin(ang));
  345.             for (int i = 0; i < n; i += len) {
  346.                 complex<double> w(1);
  347.                 for (int j = 0; j < len / 2; ++j) {
  348.                     complex<double> u = a[i + j];
  349.                     complex<double> v = a[i + j + len / 2] * w;
  350.                     a[i + j] = u + v;
  351.                     a[i + j + len / 2] = u - v;
  352.                     w *= wlen;
  353.                 }
  354.             }
  355.         }
  356.         if (invert)
  357.             for (int i = 0; i < n; ++i) a[i] /= n;
  358.     }
  359.  
  360.     void multiply_fft(const vector<int>& a, const vector<int>& b,
  361.         vector<int>& res) const {
  362.         vector<complex<double> > fa(a.begin(), a.end());
  363.         vector<complex<double> > fb(b.begin(), b.end());
  364.         int n = 1;
  365.         while (n < (int)max(a.size(), b.size())) n <<= 1;
  366.         n <<= 1;
  367.         fa.resize(n);
  368.         fb.resize(n);
  369.  
  370.         fft(fa, false);
  371.         fft(fb, false);
  372.         for (int i = 0; i < n; ++i) fa[i] *= fb[i];
  373.         fft(fa, true);
  374.  
  375.         res.resize(n);
  376.         long long carry = 0;
  377.         for (int i = 0; i < n; ++i) {
  378.             long long t = (long long)(fa[i].real() + 0.5) + carry;
  379.             carry = t / 1000;
  380.             res[i] = t % 1000;
  381.         }
  382.     }
  383.  
  384.     BigInt mul_simple(const BigInt& v) const {
  385.         BigInt res;
  386.         res.sign = sign * v.sign;
  387.         res.a.resize(a.size() + v.a.size());
  388.         for (int i = 0; i < (int)a.size(); ++i)
  389.             if (a[i])
  390.                 for (int j = 0, carry = 0; j < (int)v.a.size() || carry; ++j) {
  391.                     long long cur =
  392.                         res.a[i + j] +
  393.                         (long long)a[i] * (j < (int)v.a.size() ? v.a[j] : 0) +
  394.                         carry;
  395.                     carry = (int)(cur / BASE);
  396.                     res.a[i + j] = (int)(cur % BASE);
  397.                 }
  398.         res.trim();
  399.         return res;
  400.     }
  401.  
  402.     typedef vector<long long> vll;
  403.  
  404.     static vll karatsubaMultiply(const vll& a, const vll& b) {
  405.         int n = a.size();
  406.         vll res(n + n);
  407.         if (n <= 32) {
  408.             for (int i = 0; i < n; i++)
  409.                 for (int j = 0; j < n; j++) res[i + j] += a[i] * b[j];
  410.             return res;
  411.         }
  412.  
  413.         int k = n >> 1;
  414.         vll a1(a.begin(), a.begin() + k);
  415.         vll a2(a.begin() + k, a.end());
  416.         vll b1(b.begin(), b.begin() + k);
  417.         vll b2(b.begin() + k, b.end());
  418.  
  419.         vll a1b1 = karatsubaMultiply(a1, b1);
  420.         vll a2b2 = karatsubaMultiply(a2, b2);
  421.  
  422.         for (int i = 0; i < k; i++) a2[i] += a1[i];
  423.         for (int i = 0; i < k; i++) b2[i] += b1[i];
  424.  
  425.         vll r = karatsubaMultiply(a2, b2);
  426.         for (int i = 0; i < (int)a1b1.size(); i++) r[i] -= a1b1[i];
  427.         for (int i = 0; i < (int)a2b2.size(); i++) r[i] -= a2b2[i];
  428.  
  429.         for (int i = 0; i < (int)r.size(); i++) res[i + k] += r[i];
  430.         for (int i = 0; i < (int)a1b1.size(); i++) res[i] += a1b1[i];
  431.         for (int i = 0; i < (int)a2b2.size(); i++) res[i + n] += a2b2[i];
  432.         return res;
  433.     }
  434.  
  435.     BigInt mul_karatsuba(const BigInt& v) const {
  436.         vector<int> a6 = convert_base(this->a, BASE_DIGITS, 6);
  437.         vector<int> b6 = convert_base(v.a, BASE_DIGITS, 6);
  438.         vll a(a6.begin(), a6.end());
  439.         vll b(b6.begin(), b6.end());
  440.         while (a.size() < b.size()) a.push_back(0);
  441.         while (b.size() < a.size()) b.push_back(0);
  442.         while (a.size() & (a.size() - 1)) a.push_back(0), b.push_back(0);
  443.         vll c = karatsubaMultiply(a, b);
  444.         BigInt res;
  445.         res.sign = sign * v.sign;
  446.         long long carry = 0;
  447.         for (int i = 0; i < (int)c.size(); i++) {
  448.             long long cur = c[i] + carry;
  449.             res.a.push_back((int)(cur % 1000000));
  450.             carry = cur / 1000000;
  451.         }
  452.         res.a = convert_base(res.a, 6, BASE_DIGITS);
  453.         res.trim();
  454.         return res;
  455.     }
  456.  
  457.     void operator*=(const BigInt& v) { *this = *this * v; }
  458.     BigInt operator*(const BigInt& v) const {
  459.         if (a.size() * v.a.size() <= 1000111) return mul_simple(v);
  460.         if (a.size() > 500111 || v.a.size() > 500111) return mul_fft(v);
  461.         return mul_karatsuba(v);
  462.     }
  463.  
  464.     BigInt mul_fft(const BigInt& v) const {
  465.         BigInt res;
  466.         res.sign = sign * v.sign;
  467.         multiply_fft(convert_base(a, BASE_DIGITS, 3),
  468.             convert_base(v.a, BASE_DIGITS, 3), res.a);
  469.         res.a = convert_base(res.a, 3, BASE_DIGITS);
  470.         res.trim();
  471.         return res;
  472.     }
  473.  
  474.     // -------------------- Misc --------------------
  475.     BigInt abs() const {
  476.         BigInt res = *this;
  477.         res.sign *= res.sign;
  478.         return res;
  479.     }
  480.     void trim() {
  481.         while (!a.empty() && !a.back()) a.pop_back();
  482.         if (a.empty()) sign = 1;
  483.     }
  484.  
  485.     bool isZero() const { return a.empty() || (a.size() == 1 && !a[0]); }
  486.  
  487.     friend BigInt gcd(const BigInt& a, const BigInt& b) {
  488.         return b.isZero() ? a : gcd(b, a % b);
  489.     }
  490.     friend BigInt lcm(const BigInt& a, const BigInt& b) {
  491.         return a / gcd(a, b) * b;
  492.     }
  493.  
  494.     friend BigInt sqrt(const BigInt& a1) {
  495.         BigInt a = a1;
  496.         while (a.a.empty() || a.a.size() % 2 == 1) a.a.push_back(0);
  497.  
  498.         int n = a.a.size();
  499.  
  500.         int firstDigit = (int)sqrt((double)a.a[n - 1] * BASE + a.a[n - 2]);
  501.         int norm = BASE / (firstDigit + 1);
  502.         a *= norm;
  503.         a *= norm;
  504.         while (a.a.empty() || a.a.size() % 2 == 1) a.a.push_back(0);
  505.  
  506.         BigInt r = (long long)a.a[n - 1] * BASE + a.a[n - 2];
  507.         firstDigit = (int)sqrt((double)a.a[n - 1] * BASE + a.a[n - 2]);
  508.         int q = firstDigit;
  509.         BigInt res;
  510.  
  511.         for (int j = n / 2 - 1; j >= 0; j--) {
  512.             for (;; --q) {
  513.                 BigInt r1 =
  514.                     (r - (res * 2 * BigInt(BASE) + q) * q) * BigInt(BASE) *
  515.                     BigInt(BASE) +
  516.                     (j > 0 ? (long long)a.a[2 * j - 1] * BASE + a.a[2 * j - 2]
  517.                         : 0);
  518.                 if (r1 >= 0) {
  519.                     r = r1;
  520.                     break;
  521.                 }
  522.             }
  523.             res *= BASE;
  524.             res += q;
  525.  
  526.             if (j > 0) {
  527.                 int d1 =
  528.                     res.a.size() + 2 < r.a.size() ? r.a[res.a.size() + 2] : 0;
  529.                 int d2 =
  530.                     res.a.size() + 1 < r.a.size() ? r.a[res.a.size() + 1] : 0;
  531.                 int d3 = res.a.size() < r.a.size() ? r.a[res.a.size()] : 0;
  532.                 q = ((long long)d1 * BASE * BASE + (long long)d2 * BASE + d3) /
  533.                     (firstDigit * 2);
  534.             }
  535.         }
  536.  
  537.         res.trim();
  538.         return res / norm;
  539.     }
  540. };
  541.  
  542. BigInt fib_slow(int n) {
  543.     BigInt a(0);
  544.     BigInt b(1);
  545.     for (int i = 2; i <= n; ++i) {
  546.         BigInt t = a;
  547.         a = b;
  548.         b += t;
  549.     }
  550.     return b;
  551. }
  552.  
  553. struct fraction {
  554.     static int precision;
  555.     BigInt x, y;
  556.     fraction(BigInt x_ = 0, BigInt y_ = 1) : x(x_), y(y_) {}
  557.     fraction& operator +=(fraction a) {
  558.         fraction ans(this->x * a.y + this->y * a.x, this->y * a.y);
  559.         BigInt g = gcd(ans.x, ans.y);
  560.         ans.x /= g;
  561.         ans.y /= g;
  562.         return *this = ans;
  563.     }
  564.     friend ostream& operator << (ostream& os, fraction const& a) {
  565.         BigInt f = a.x / a.y;
  566.         vector<BigInt> decimal;
  567.         BigInt remainder = a.x % a.y;
  568.         for(int i = 0; i < precision; i++) {
  569.             remainder *= 10;
  570.             if(i == precision - 1) {
  571.                 if(((remainder % a.y) * 10) > a.y * 5) {
  572.                     decimal.push_back((remainder + a.y - 1) / a.y);
  573.                 } else {
  574.                     decimal.push_back(remainder / a.y);
  575.                 }
  576.             } else {
  577.                 decimal.push_back(remainder / a.y);
  578.             }
  579.             remainder = remainder % a.y;
  580.         }
  581.         int carry = 0;
  582.         for(int i = precision - 1; i >= 0; i--) {
  583.             decimal[i] += carry;
  584.             carry = 0;
  585.             if(decimal[i] >= 10) {
  586.                 carry = 1;
  587.                 decimal[i] = decimal[i] % 10;
  588.             }
  589.         }
  590.         if(carry) {
  591.             f += 1;
  592.         }
  593.         cout << f << ".";
  594.         for(auto e : decimal) {
  595.             cout << e;
  596.         }
  597.         return os;
  598.     }
  599. };
  600.  
  601. int fraction::precision = 6;
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement