add

Author: KumaCS

docs/fps/sampling-points-shift.md

@@ -0,0 +1,22 @@
+以下の問題を $O((N+M)\log(N+M))$ 時間で解く．
+
+> 次数 $N$ 未満の多項式 $f(x)$ について $f(0),f(1),\dots,f(N-1)$ が与えられる．
+> $f(c),f(c+1),\dots,f(c+M-1)$ を列挙せよ．
+
+## アルゴリズム
+
+下降階乗冪を用いて $f(x)=\sum_{i=0}^{N-1}a_ix^{\underline{i}}$ とおくと $0\leq i\lt N$ に対し
+
+$$f(k)=k!\sum_{i=0}^{k}\frac{a_i}{(k-i)!}$$
+
+が成り立つため，$a_i$ は
+
+$$a_i=[x^i]e^{-x}\sum_{k=0}^{N-1}\frac{f(k)}{k!}x^k$$
+
+として計算できる．
+
+また下降階乗冪について $(a+b)^{\underline{n}}=\sum_{i=0}^{n}\binom{n}{i}a^{\underline{i}}b^{\underline{n-i}}$ が成り立つため，$f(x+c)=\sum_{i=0}^{N-1}b_ix^{\underline{i}}$ とおくと
+
+$$b_i=\frac{1}{i!}\sum_{j}(i+j)!a_{i+j}\cdot\frac{c^{\underline{j}}}{j!}$$
+
+が成り立つ．あとは $a_i$ の計算の逆をすればよい．

fps/relaxed.hpp

@@ -200,6 +200,6 @@ class RelaxedSqrt {
 };
 
 /**
- * @brief Relaxed 畳み込み
+ * @brief Relaxed 
  * @docs docs/fps/relaxed.md
  */

fps/sampling-points-shift.hpp

@@ -0,0 +1,32 @@
+#pragma once
+#include "modint/factorial.hpp"
+#include "fps/formal-power-series.hpp"
+
+// f(0),f(1),...,f(n-1) -> f(c),...,f(c+m-1)
+template <class mint>
+vector<mint> SamplingPointsShift(const vector<mint>& f, mint c, int m) {
+  using fps = FormalPowerSeries<mint>;
+  using fact = Factorial<mint>;
+  int n = f.size();
+  fact::reserve(m);
+  fps f1(n), ei(n);
+  for (int i = 0; i < n; i++) f1[i] = f[i] * fact::fact_inv(i);
+  for (int i = 0; i < n; i++) ei[i] = fact::fact_inv(i) * (i % 2 ? -1 : 1);
+  f1 *= ei;
+  for (int i = n; i < f1.size(); i++) f1[i] = 0;
+  for (int i = 0; i < n; i++) f1[i] *= fact::fact(i);
+  fps g(n, 1);
+  for (int i = 1; i < n; i++) g[i] = g[i - 1] * (c + 1 - i) * fact::inv(i);
+  g = g.middle_product(f1);
+  for (int i = 0; i < n; i++) g[i] *= fact::fact_inv(i);
+  fps e(m);
+  for (int i = 0; i < m; i++) e[i] = fact::fact_inv(i);
+  g *= e;
+  g.resize(m);
+  for (int i = 0; i < m; i++) g[i] *= fact::fact(i);
+  return g;
+}
+/**
+ * @brief 評価点シフト
+ * @docs docs/fps/sampling-points-shift.md
+ */

math/stern-brocot-tree.hpp

@@ -24,7 +24,7 @@ struct SternBrocotTreeNode {
   }
   SternBrocotTreeNode(pair<T, T> p) : SternBrocotTreeNode(p.first, p.second) {}
   SternBrocotTreeNode(const vector<T> seq_) {
-    for (auto &v : seq_) {
+    for (auto& v : seq_) {
       assert(v != 0);
       if (v > 0)
         go_right(v);
@@ -53,10 +53,10 @@ struct SternBrocotTreeNode {
   }
   T depth() const {
     T d = 0;
-    for (auto &v : seq) d += abs(v);
+    for (auto& v : seq) d += abs(v);
     return d;
   }
-  static Node lca(const Node &x, const Node &y) {
+  static Node lca(const Node& x, const Node& y) {
     Node res;
     for (int i = 0; i < min(x.seq.size(), y.seq.size()); i++) {
       T d1 = x.seq[i], d2 = y.seq[i];
@@ -110,7 +110,7 @@ struct SternBrocotTreeNode {
       }
     }
   }
-};  // namespace SternBrocotTree
+};
 
 /**
  * @brief Stern-Brocot Tree

verify/fps/LC_division_of_polynomials.test.cpp

@@ -0,0 +1,21 @@
+#define PROBLEM "https://judge.yosupo.jp/problem/division_of_polynomials"
+
+#include "template/template.hpp"
+#include "modint/modint.hpp"
+using mint = ModInt<998244353>;
+#include "fps/fps-ntt-friendly.hpp"
+using fps = FormalPowerSeries<mint>;
+
+int main() {
+  int n, m;
+  in(n, m);
+  fps f(n), g(m);
+  in(f, g);
+  fps q = f / g;
+  q.shrink();
+  fps r = f - g * q;
+  r.shrink();
+  out(q.size(), r.size());
+  out(q);
+  out(r);
+}

verify/fps/LC_shift_of_sampling_points_of_polynomial.test.cpp

@@ -0,0 +1,18 @@
+#define PROBLEM "https://judge.yosupo.jp/problem/shift_of_sampling_points_of_polynomial"
+
+#include "template/template.hpp"
+#include "modint/modint.hpp"
+using mint = ModInt<998244353>;
+#include "fps/fps-ntt-friendly.hpp"
+using fps = FormalPowerSeries<mint>;
+#include "fps/sampling-points-shift.hpp"
+
+int main() {
+  int n, m;
+  mint c;
+  in(n, m, c);
+  vector<mint> f(n);
+  in(f);
+  auto g = SamplingPointsShift(f, c, m);
+  out(g);
+}

verify/fps/LC_sqrt_of_formal_power_series.relaxed.test.cpp

@@ -13,10 +13,14 @@ int main() {
   in(n);
   fps a(n);
   in(a);
-  if (ModSqrt(a[0].val(), mint::get_mod()) == -1) {
+  if (a[0] == 0) {
+    a = FpsSqrt(a);
+    if (a.empty())
+      out(-1);
+    else
+      out(a);
+  } else if (ModSqrt(a[0].val(), mint::get_mod()) == -1) {
     out(-1);
-  } else if (a[0] == 0) {
-    out(FpsSqrt(a));
   } else {
     fps b(n);
     RelaxedSqrt<mint> sqrt;