[数学] 沉淀三角中求某奇数的倍数的个数

有一个正整数M（假设取10），构造如下的倒三角数列（因下一行每一个数字都是上一行两个数据之和，因而我给它简单命名为沉淀三角，希望有人知道学名告诉我）；

第一行：1 2 3 4 5 ...... M-2 M-1 M

第二行： 3 5 7 9 ...... 2M-3 2M-1

第三行：  8 12 16 ...... 4M-4

......

第M行：            ...

在这个数列中求可整除奇数X（1除外）的数字的个数；

 

好了，来解这道数学题：

1. 首先，很容易看出，每一行是一个等差数列，差为2^n-1；

2. 然后，求解每行的首数，计算方法为：First_n = ∑Cⁱ_n-1*(i-1) = C⁰_n-1 *1 + C¹_n-1*2 + ... + Cⁱ_n-1*(i-1) + ... + C^n-1_n-1*n

    则，2 * First_n = ∑Cⁱ_n-1*(n+1) = 2^n-1*(n+1)；

    则，First_n = 2^n-2*(n+1), (2<=n<=M)

3. 这时可知每一行的数据为：2^n-2*(n+1+2*y), (0<=y<=M-n), (2<=n<=M)，在这个范围内，求每一个数据元素，检查是否满足

    2^n-2*(n+1+2*y) ≡ 0 (mod X), (0<=y<=M-n), (2<=n<=M)

4. 因为数X（1以外的奇数）肯定与2^n-2互质，则以上同余式化为

    (n+1+2*y) ≡ 0 (mod X), (0<=y<=M-n), (2<=n<=M) （解读：从n+1到2M+1-n，以2为公差的等差数列对X求模）

5. 同理，因为X与2互质，所以，从第一个出现可整除X的数开始，每增加2X出现一次；

6. 第一个出现的可整除X的数为

    first_match = ceil( (n+1) / X ) * X（如果注意，这样求出的理论第一个数并不一定在这一行的数列中，还需要保证与n+1同奇偶；当first_match > 2M+1-n，这个数不存在）；

    之后出现的可整除X的数的个数为 floor( (2M+1-n - first_match) / 2X )；

7. 以上公式对n的范围要求是2~M，所以第一行要单独计算一下；

 

先简单给个构造数列的代码，把数列输出到屏幕，然后使用上面求解的公式算一下结果：

# -*- coding: utf-8 -*-  

import math

# 构造triangle数列
def construct(triangle, M):
    last = []
    for m in range(M, 0, -1):
        if len(last) == 0: # no last exist
            for i in range(1, m+1):
                last.append(i)
        else: # last exist
            tmp = []
            for i in range(m):
                tmp.append(last[i] + last[i+1])
            last = tmp

        triangle.append(tuple(last))

# 将triangle输出到屏幕，M很大的时候不要用
def show_triangle(triangle):
    print triangle
    M = len(triangle)
    pre = ""
    sep = "       "
    for line in triangle:
        print pre,
        for e in line:
            print "%4d" % e, sep,
        print
        pre += sep

# 计算X倍数的个数
def num_multiple(triangle, X):
    result = 0
    M = len(triangle)
    # 第一行不符合计算公式，单独处理
    if X < M:
        result += 1
    # 2~M行按公式计算
    for n in range(2, M+1):
        # 划定范围
        low = n+1
        high = 2*M+1-n
        # 找到第一个X的倍数
        first_match = math.ceil(float(low)/X) * X
        while (first_match % 2) != (low % 2):
            first_match += X 
        if first_match > high:
            continue
        result += int(1 + math.floor((high-first_match) / (2*X)))
    return result

if __name__ == "__main__":
    triangle = []
    construct(triangle, 10)
    show_triangle(triangle)
    ret = num_multiple(tuple(triangle), 3)
    print ret