在知乎上看到一个有趣的问题: 如何证明这个数列无界? 在此记录一下简单的做法,然后把这篇博客记录以后遇到的一些特殊数列。

证明:当 $x$ 是无理数时,$f_n(x) = \sum_{i=1}^n (-1)^{\lfloor ix \rfloor}$ 无界

$\lfloor x \rfloor$ 表示 $x$ 的整数部分, $0 \leq \lbrace x \rbrace \doteq x- \lfloor x \rfloor<1$ 表示 $x$ 的小数部分

为了更顺畅的介绍这个数列,我们先来点前戏,不然插入的时候就不那么顺滑了 0.0 不要笑,我没有开车 0.0

连分数

它和 Farey 序列,无理数的有理逼近,密切相关,还是由于求 Pell 方程的快速方法

我们将一个整数序列 $a_0,a_1, \cdots, a_n$ 构成的分数叫做连分数,记作

其中 $a_0$ 是整数, $a_i (i>0)$ 是正整数。

我们定义:

$p_0 = a_0,p_1 = a_1a_0 +1,p_n = a_np_{n-1} + p_{n-2} (n>1)$

$q_0 = 1,q_1 = a_1,q_n=a_nq_{n-1}+q_{n-2} (n>1)$

那么数学归纳法易证:

并且还有

  • $[a_0,a_1,\cdots, a_n] = [a_0,a_1,\cdots, a_{n-1} + \frac{1}{a_n}]$
  • $p_{n+1}q_n - q_{n+1}p_n = (-1)^n$
  • 从而 $(p_n,p_{n+1}) = (q_n,q_{n+1}) = (p_n,q_n) = 1$
  • $x_{n+2} - x_n = \frac{(-1)^n a_{n+2}}{q_n q_{n+2}}$
  • 从而 $x_0 < x_2< x_4 < \cdots < x_5 < x_3 < x_1$

对任意一个数 $x$ 构造连分数

  1. 初始 $i = 0$
  2. $a_i = \lfloor x \rfloor$ 取 $x$ 的整数部分。
  3. $x - a_i == 0$ 结束
  4. 取 $x = \frac{1}{x-a_i}$
  5. ++$i$ 回到步骤 1

从构造中,我们看出:连分数长度有限当且仅当 $x$ 是有理数,此时,最后一项得到的连分数就是 $x$。

当 $x$ 是无理数时, $x_{2n} < x < x_{2n+1}$ ,又因为

且 $x_{2n}$ 单调递增,$x_{2n+1}$ 单调递减。所以 $\lim _{n \to \infty} x_n = x$。并且

补充定义

显然,$x = \frac{p_n’}{q_n’} = \frac{p_n}{q_n} + \frac{(-1)^n}{q_n q_{n+1}’}$

唯一分解

对于任意给定 $x = [a_0,a_1,\cdots,a_n,\cdots]$,对任意 $N$,有如下唯一分解:

其中 $0 \leq c_i \leq a_{i+1}, c_m(N,x) >0$ 且对任意 $0 \leq j \leq m(N,x)$,$\sum_{i=0} ^ j c_i q_i < q_{j+1}$ 。

$ f_n(x) = \sum_{i=1}^n (-1)^{\lfloor ix \rfloor} $ 何时有界

由于 $f_n(x+2) = f_n(x) $, 所以我们仅考虑 $x \in [0,2) $ , 在不引起混淆时,简记成 $f_n $

当 $x$ 是有理数时,记 $ x = \frac{p}{q},\quad (p,q)=1$

则 $f_{2q+n} = f_{2q} + f_n$, 从而 $a_n$ 有界当且仅当 $f_{2q} \neq 0 $

  • 当 $q$ 为奇数, $f_{2q} = 0$,即 $f_n$ 有界

  • 当 $p$ 为偶数时,此时 $q$ 为奇数,且 $(\frac{p}{2},q)=1$,所以当 $i$ 跑遍 $[1, q]$ 时,$i\frac{p}{2} \mod q$ 跑遍 $[0,q-1]$,即

    从而 $f_{2q} = 2f_q = 2$ 推出 $f_n$ 无界。

对任意无理数 $x$,我们证明 $g_n(x) = f_n(2x)$ 无界,从而 $f_n(x)$ 无界

目前没理解论文中的做法,太繁琐了

显然,我们只需考虑 $0< x <1$

我们考虑 $x$ 的连分数 $[0,a_1,a_2,\cdots,a_n,\cdots]$,由于 $(p_n,q_n)=(q_n,q_{n+1})=1$ (目前不知道有什么简单的操作)

这个问题可以推广到更一般的情形:单栏阅读 双栏阅读

里面定理一挺有意思的,虽然长但是很清晰。

另一个相关问题:https://www.zhihu.com/question/392014769

是否存在 $N$ 使得, $|\sin n|>\frac{1}{n}$ 对所有 $n>N$ 成立。(dna049寨森 Lambda-CDM 建议修改)

存在无穷多个正整数 $x_n$ 使得 $|\sin x_n| < \frac{2 \pi}{x_n}$

引理: 对任意无理数 $a$ , 存在无穷多个 $x_n$,使得 $\min(\{a x_n\},\{-a x_n\}) < \frac{2}{x_n}$

我们考虑 $a$ 的连分数 $[a_0,a_1,\cdots,a_n,\cdots]$, 我们知道 $|a - \frac{p_n}{q_n} | < \frac{1}{q_n ^2}$,由于 $(p_n,q_n) = 1$, 所以存在 $u_n,v_n$ 使得 $p_n u_n + q_n v_n = 1$, 此时 $p_n(u_n+tq_n) + q_n(v_n - tp_n) = 1$,所以我们不妨假设 $q_n < u_n < 2q_n$。所以

我们仅考虑 $n$ 为奇数数的情况,此时 $-\frac{1}{q_n^2} < a - \frac{p_n}{q_n} < 0$,所以 $-\frac{1}{q_n} <(a-\frac{p_n}{q_n})u_n + \frac{1}{q_n} <\frac{1}{q_n}$。而 $x_n = u_{2n+1} > q_{2n+1}$ 即为所求。

实际上 $n$ 为偶数时,也可以取做,只是此时 $u_n$ 要换成 $3q_n - u_n$

取引理中 $a = \frac{1}{\pi}$,由于 $|\sin(x)|$ 是周期为 $\pi$ 的偶函数。所以 $|\sin(x_n)| = |\sin(\{ax_n\}\pi)| = |\sin(\{-ax_n\}\pi)| < \frac{2\pi}{x_n}$

寨森 Lambda-CDM 启发,可以不用引理直接证明:

存在无穷多个正整数 $m$,使得 $|\sin m| < \frac{\pi}{m}$

Proof:考虑 $\frac{1}{\pi}$ 的连分数 $[a_0,a_1,\cdots,a_n,\cdots]$,我们有 $|\frac{1}{\pi} - \frac{p_n}{q_n} | < \frac{1}{q_n ^2}$,所以

这个问题还关联这 Irrationality Measure,菲尔兹奖级别的工作!

Irrationality Measure

对于给定是实数 $x$, 定义

显然考虑连分数,我们知道 当 $x$ 是有理数时,$u(x) = 1$,无理数时,$\mu(x) \geq 2$。

Roth 证明,当 $x$ 是代数数时 $\mu(x) = 2$,因此获得 Field 奖。

$\mu(L) = \infty$,其中 $L = \sum_{n=1} ^{\infty} 10^{-n!}$ 为刘维尔(Joseph Liouville)数

我们下面证明:当 $u>2$ 时,$A = \{ x \mid 0< |x -\frac{p}{q}| < \frac{1}{q^u}\}$,则 $A$ 的外(lebesgue)测度 $m^{\star}(A) = 0$,从而(lebesgue)测度 $m(A)=0$

由于 当$1\leq p \leq q$时, $(\frac{p}{q}-\frac{1}{q^u},\frac{p}{q} + \frac{1}{q^u})$ 只有可数个,记作$I_1, I_2,\cdots I_n,\cdots$ ,按照定义,$A \cap (0,1) \subseteq \overline\lim I_n$, 另一方面 $\sum_{i=1}^{\infty} I_n = \sum_{q=1}^{\infty} \frac{q+1}{q^u} < +\infty$,从而 $m^{\star}(A \cup (0,1)) = 0$,同理可证了$m^{\star}(A \cup (n,n+1))=0, n \in \mathbb{Z}$,从而

从而$m^{\star}(A) = m(A) = 0$。其中 $m^{\star}(I) = \inf \{u|u=\sum_{k=1}^{\infty} |I_k|,\; \cup_{k=1}^ {\infty} I_k \supset m ,I_k \text{ 是开矩形 } \}$

sagemath 数值测试(Jupyter Notebook 上运行)

1
2
3
4
5
6
7
8
9
x = continued_fraction(1/pi).convergents()
for i in range(1,32):
n = x[i].denominator()
print(N(abs(sin(n))-pi/n),"\t= ",abs(sin(n))-pi/n ) # 后来精度就崩了

x = continued_fraction(1/pi).convergents()
for i in range(1,32):
n = x[i].denominator()
print(i,N(abs(sin(n))*n),"\t= ",abs(sin(n))*n) # 后来精度就崩了

$1 = \sum_{s \in S} \frac{1}{s^2-1}$, 其中 $S = \{a^m | a > 1, m>1, a,m \in \mathbb{N} \}$

方法:On a Series of Goldbach and Euler

其中 $B = \{n \in \mathbb{N} \mid n = k^m , k>1 \}$,$A= \{n \in \mathbb{N} \mid n >1 \} - B$

$f(n,m) = \sum_{i=0}^{i \leq n-im} C_{n-im}^i$

洛谷题目,具体题号就不便说了。

数据范围 $1 \leq n \leq 10^{18} , m \leq 100$ 。我们注意考察 $f(n,m)$ 的意义!假设我们要上楼梯,每次只能上 1 步或者 $m+1$ 步,那么 $f(n,m)$ 就是方案数,因此,我们显然有

当然了我们不考虑意义直接用 $f(n,m)-f(n-1,m)$ 也可以得到这个公式。

无法过题的 $O(m \log m \log n)$ 代码,原因是没法找一个合适的基底

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
#include<bits/stdc++.h>
using namespace std;
using LL=long long;
using BI=__int128;
const LL mod=1e9+7;
const int M =10004;
BI a[4*M],r[4*M],ans[2*M];
template<typename T,typename U>
T powmod(T x,U n,T p){
T r(1);
while(n){
if(n&1) r=r*x%p;
n>>=1; x=x*x%p;
}
return r;
}
void bitreverse(BI *x,int len){
for(int i=0,j=0;i!=len;++i){
if(i>j) swap(x[i],x[j]);
for(int l=len>>1;(j^=l)<l;l>>=1);
}
}
const BI FM = BI(29)<<57|1,gg=3;
// the mod must NTT-friendly or (len+1)*mod^2 < FM
void ntt(BI *x,int len,bool isInverse=false){
g = powmod(gg,(FM-1)/len,FM);
if(isInverse){
g=powmod(g,len-1,FM);
BI invlen = powmod(BI(len),FM-2,FM);
for(int i=0;i!=len;++i){
x[i]=x[i]*invlen%FM;
}
}
bitreverse(x,len);
for(int half=1,step=2;half!=len;half<<=1,step<<=1){
BI wn = powmod(g,len/step,FM);
for(int i=0;i<len;i+=step){
BI w(1);
for(int j = i;j<i+half;++j){
BI t=(w*x[j+half])%FM;
x[j+half]=(FM-t+x[j])%FM;
x[j]=(x[j]+t)%FM;
w = w*wn%FM;
}
}
}
}
void stand(BI *a, int n){
for(int i=2*n;i>n;--i){
a[i-1]=(a[i-1]+a[i])%mod;
a[i-n-1]=(a[i-n-1]+a[i])%mod;
a[i] = 0;
}
}
void square(BI *a, int n){
int len = 1<<(32-__builtin_clz(2*n+1));
ntt(a,len);
for(int i=0;i<len;++i){
a[i]=a[i]*a[i]%FM;
}
ntt(a,len,1);
stand(a,n);
}
void mul(BI *a,BI *b,int na,int nb){
int len = 1<<(32-__builtin_clz((na+nb)+1));
ntt(a,len);ntt(b,len);
for(int i=0;i<len;++i){
a[i] = a[i]*b[i]%FM;
}
ntt(b,len,1);ntt(a,len,1);
stand(a,na);
}
void initC(int m){
for(int i=0;i<=m;++i){
ans[i] = 1;
}
for(int i=m+1;i<=2*m;++i){
ans[i] = (ans[i-1]+ans[i-m-1])%mod;
}
memset(a,0,sizeof(a));
memset(r,0,sizeof(r));
r[0]=a[1]=1;
}
LL getans(LL n,int m){
initC(m);
while(n){
if(n&1) mul(r,a,m,m);
n>>=1; square(a,m);
}
LL ret = 0;
for(int i=0;i<=m;++i){
ret+=(r[i]*ans[i+m])%mod;
}
return ret%mod;
}

int main() {
//freopen("in","r",stdin);
//freopen("m","w",stdout);
std::ios::sync_with_stdio(false);std::cin.tie(nullptr);
LL n;
int m;
while(cin>>n>>m){
if(m==0){
cout<<powmod(2LL,n,mod)<<endl;
}else{
if(n<=m){
cout<<1<<endl;
}else{
cout<<getans(n-m,m)<<endl;
}
}
}
return 0;
}

求 $m$ 阶线性递推关系式第 $n$ 项

参考:2013/03/02 线性递推关系和矩阵乘法

上面文献的意义:矩阵加速已经成为过去式了,多项式加速时代的到来,NTT 加速多项式的时代到来。

设 $a_n = c_{m-1} a_{n-1} + c_{m-2} a_{n-2} + c_0 a_{n-m}$,给定了初值 $a_1,a_2,\cdots,a_m$ 的情况下,求 $a_n$

我们记

并且显然

我们就可以用矩阵乘法在 $O(m^3 \log n)$ 解决这个问题了。

但线性递推关系式可以不用矩阵优化到 $O(m^2 \log n)$ :

显然 $A$ 的特征多项式 $f(x) = x^m - c_1 x^{m-1}- \cdots -c_m$,所以 $A^n = c_{m-1} A^{n-1} + \cdots c_0 I$,也就是说 $A^n$ 可以由 $A^{m-1},A^{m-2},\cdots,I$ 线性表出,而 $A^{n_1+n_2} = A^{n1} A^{n_2}$ 然后运算就是卷积运算,相当于多项式乘法!也就是说我们可以在 $O(m^2 \log n)$ 时间复杂度,求出 $A^n = b_0 I + b_1 A + \cdots b_{m-1} A^{m-1}$, 注意到

所以我们只需预处理出 $a_1,\cdots,a_{2m-1}$ 这 $2m-1$ 个数即可。代码见我博客的模板

这也提供了一般 $m$ 阶矩阵的 $n$ 次方的一个 $O(m^3+m^2\log n)$ 的算法!!!卧槽!我也太帅了吧。

最后如果递推关系中仅有常数个 $c_i$ 不为 0,此时还能用 NTT(数论快速变换) 利用 多项式带模除法: Division with remainder 的 $O(m \log m)$ 算法,上述算法可以优化到 $O(m\log m \log n + m^2)$,(暂时不知道如何去掉 $m^2$)但是写法太复杂。搞定!参考 这里,不要涉及矩阵

注意到一次乘法之后,会变成 $I,A,\cdots A^{2m-1}$ 的线性组合,$A^{m},\cdots A^{2m-1}$ 这 $m$ 个要再用前 $m$ 线性表出,由于 $c_i$ 仅有常数个不为 0,可以 $O(m)$ 复杂度把他们写成前 $m$ 个的线性组合。

当 $m>10^4$ 时,我们就有必要写带 NTT 的版本了(NTT 模板可以在我博客中找到)有需求的时候再写吧。

$\sum_{i \equiv r \mod m} \binom{n}{i} \mod M$

数据范围:$1 \leq n \leq 10^{18},\; 2 \leq m \leq 2000,\; 0 \leq r < m,\; 10^8 < M < 10^9$

记 $w$ 满足 $w^m = 1,w^n \neq 1, 0 < n < m$,则答案为 $\frac{1}{m}\sum_{i=0}^{m-1}F(w^i)$,其中 $F(x) = x^{m-r}(1+x)^n \mod x^{m}-1$

假设 $F(x) = \sum_{i=0}^{m-1} a_i x^i$,则答案就是 $a_0$,即答案是 $F(0)$

这个跟上面一样本质是一样的,做带模的多项式运算。

参考:zscoder 的博客

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
#include<bits/stdc++.h>
using namespace std;
using LL = long long;

class Poly{
public:
const static int N = 2003;
int d;
LL a[N];
Poly(int _d):d(_d){
memset(a,0,sizeof(a));
}
auto& operator[](int i){
return a[i];
}
Poly operator*(Poly &A){
Poly R(d+A.d);
for(int i=0;i<=d;++i){
for(int j=0;j<=A.d;++j){
R[i+j] += a[i]*A[j];
}
}
return R;
}
};
Poly mulmod(Poly &A,Poly &B,int m,LL mod){
Poly R = A*B;
for(int i=R.d;i>=m;--i){
R[i-m] += R[i];
R[i]=0;
}
R.d = m;
for(int i=0;i<=m;++i){
R[i]%=mod;
}
return R;
}
template<typename T,typename U>
T powmod(T x,U n,int m,LL mod){
T r(0);
r[0]=1;
while(n){
if(n&1) r=mulmod(r,x,m,mod);
n>>=1; x=mulmod(x,x,m,mod);
}
return r;
}

int main(){
//freopen("in","r",stdin);
std::ios::sync_with_stdio(false);std::cin.tie(nullptr);
LL n,r,mod;
int m;
while(cin>>n>>m>>r>>mod){
Poly A(m-1);
A[0]=A[1]=1;
Poly R = powmod(A,n,m,mod);
cout<<(R[r%m]*m%mod)<<endl;
}
return 0;
}

组合数倒数交错和 $\sum_{i=0}^n \frac{(-1)^i}{{n \choose i}}$

考虑积分 $F(n,m) = \int _0^1 x^m (1-x)^{n-m} dx$,显然 $F(n,m)=F(n,n-m),F(n,0) = \frac{1}{n+1}$
且分部积分即可知道 $F(n,m) = \frac{m}{n-m+1} F(n,m-1) = \frac{1}{(n+1){n \choose m}}$,再等比数列求和就有

To be continued