Project_Euler-12 题解

Project_Euler-12 题解

标题

题目

1
2

思路

我们可以从小到大枚举每一个三角形数,然后计算他们的约数个数,从而得到结果。

代码

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <string.h>
#include <time.h>
//#include "./DBG_text.h"

#define MAX_N 100

typedef long long ll;

ll get_factory(ll n) {
    ll cnt = 0;
    for (ll i = 1, I = sqrt(n); i <= I; i++) {
        if (n % i) continue;
        cnt += 1;
        cnt += (i != n / i);
    }
    return cnt;
}

ll triangle(ll n) {
    return n * (n + 1) >> 1;
}

ll main() {
    ll n = 0;
    while (1) {
        n++;
        // 求第n个三角形数
        ll val = triangle(n);
        // 计算因子个数
        ll cnt = get_factory(val);
        //DBG("cnt = %lld\n", cnt);
        if (cnt <= 500) continue;
        printf("val = %lld, cnt = %lld\n", val, cnt);
        break;
    }
    return 0;
}

优化思路

这里不得不提到 一个定理了:


约数个数定理

约数个数定理:如果一个大于1的正整数 n 可以分解为质因数的形式:n = p1^a1 * p2^a2 * ... * pk^ak,其中 p1, p2, ..., pkn 的不同质因数,a1, a2, ..., ak 是对应的幂次,那么 n 的约数个数就是 (a1+1) * (a2+1) * ... * (ak+1)

举个例子,假设 n = 2^3 * 3^2 * 5^1,那么 n 的约数个数就是 (3+1) * (2+1) * (1+1) = 24


有了这个定理可以得到什么呢?当一个数字n被确定时,我们通过素数快速计算出它的约数个数。

但是前提,我们需要先知道它的素因子是什么,这就需要我们使用线性筛来求解。

不管先不着急,因为我们要求的是三角形数字的因数个数,而不是n的,因此我们需要进一步处理:

具体是这样的做的:

我们设 F ( n ) F(n) F(n) 表示数字 n n n 的因数个数, f ( n ) f(n) f(n) 表示第 n n n 个三角形数字,那么我们可以得到下面的关系:

F ( f ( n ) ) = F ( n ∗ ( n + 1 ) 2 ) F(f(n)) = F(\frac {n * (n + 1)} {2}) F(f(n))=F(2n(n+1))

此时我们发现 n n n n + 1 n + 1 n+1 一定是一奇一偶的关系,并且他们一定互素。


两个互素的数相乘得到的数的因数 ( F ( n ∗ ( n + 1 ) ) F(n * (n + 1)) F(n(n+1))) 的个数等于什么?

F ( n ∗ ( n + 1 ) ) = F ( n ) ∗ F ( n + 1 ) F(n*(n + 1)) = F(n) * F(n + 1) F(n(n+1))=F(n)F(n+1)

为什么呢?因为两个数的最大公因数为1,所以他们没有相同的因数,因此他们两个相乘得到的因数的个数就是他们自己拆解后的因数不断排列组合得到的因数的个数。

举个例子: 5 5 5 6 6 6 5 5 5 的因数个数为 2 2 2 6 6 6 的因数个数为 4 4 4 ,那么他们的乘积就是 30 30 30 30 30 30 的因数个数就是 2 ∗ 4 = 8 2 * 4 = 8 24=8 个。


那么回到上面的式子:

F ( f ( n ) ) = F ( n ∗ ( n + 1 ) 2 ) F(f(n)) = F(\frac {n * (n + 1)} {2}) F(f(n))=F(2n(n+1))

就可以化简为:

F ( f ( n ) ) = F ( n ) ∗ F ( n + 1 ) 2 F(f(n)) = \frac{F(n) * F(n + 1)}{2} F(f(n))=2F(n)F(n+1)

对n分奇偶情况讨论,可以将2化简掉:

当 n 为偶数时: F ( f ( n ) ) = F ( n 2 ) ∗ F ( n + 1 ) 当 n 为奇数时: F ( f ( n ) ) = F ( n ) ∗ F ( n + 1 2 ) 当n为偶数时: F(f(n)) = F(\frac {n}{2}) * F(n + 1) \\ 当n为奇数时: F(f(n)) = F(n) * F(\frac {n + 1} {2}) n为偶数时:F(f(n))=F(2n)F(n+1)n为奇数时:F(f(n))=F(n)F(2n+1)

得到这个公式,我们就可以通过 n n n 来表达第 n n n 个三角形数的因数个数了。

程序的主框架就完成了:

int main() {
    // ... init();
    ll n = 0; 
    while (1) {
        n++;
        ll val;
        val = (n & 1) ? f[n] * f[(n + 1) >> 1] : f[n >> 1] * f[n + 1];
        if (val <= 500) continue;
        printf("n = %lld\n", n * (n + 1) >> 1);
        break;
    }
    return 0;
}

接下来就是怎么确定 F ( n ) F(n) F(n) ,我们必须保证,当我们拿出一个 n n n 的时候,要得到 n n n 对应的因数的个数。

我们可以将 n n n 对应的因数个数存放在数组里。

我们需要找到 n n n 的最小素因子(prime[j]),这个可以通过遍历素数表来实现,但是对于 n n n 本身,我们需要在线性筛中实现( n = i ∗ p r i m e [ j ] n = i * prime[j] n=iprime[j]), i i i 是线性筛的遍历值,在筛中有以下几种情况讨论:

  • n n n为素数时:其因数个数只有1和它本身,因此就是2
  • i i i与某个素数 j j j互素的时候,f[prime[j] * i] = f[prime[j]] * f[i];
  • i i i与某个素数 j j j非互素的时候 f[prime[j] * i] = (f[i] / (cnt[i] + 1)) * (cnt[i] + 2);

优化代码

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <string.h>
#include <time.h>
#include "./DBG_text.h"
#define MAX_N 100000

typedef long long ll;

int prime[MAX_N + 5];
int f[MAX_N + 5];
int cnt[MAX_N + 5]; // 表示的是数字i中最小素因子的幂次
void init() {
    for (int i = 2; i <= MAX_N; i++) {
        if (!prime[i]) {
            prime[++prime[0]] = i;
            f[i] = 2;
            cnt[i] = 1;
        }
        for (int j = 1; j <= prime[0]; j++) {
            if (i * prime[j] > MAX_N) break;
            prime[i * prime[j]] = 1;
            // == 0 说明不互素并且prime[j]是i最小的素因数,否则说明互素
            if (i % prime[j] == 0) {
            	// 这里的prime[j]就是n的最小素因数,cnt[i]中记录的就是它的幂数
            	// prime[j] * i说明n这个数相对于i这个数的最小素因数增加了1幂,此时需要先让f[i]先除掉原来的最小幂
            	// 然后再乘以新的幂数也就是+2
            	// 这个2中有一个1表示prime[j] * i之后,增加的1幂
            	// 另一个1表示根据约数个数定理而必须+1
                f[prime[j] * i] = (f[i] / (cnt[i] + 1)) * (cnt[i] + 2); // 这句话较难理解
                cnt[prime[j] * i] = cnt[i] + 1;
                break;
            } else {
                f[prime[j] * i] = f[prime[j]] * f[i];
                cnt[prime[j] * i] = 1;
            }
        }
    }
    return ;
}

int main() {
    init();
    ll n = 0; 
    while (1) {
        n++;
        ll val;
        val = (n & 1) ? f[n] * f[(n + 1) >> 1] : f[n >> 1] * f[n + 1];
        if (val <= 500) continue;
        printf("n = %lld\n", n * (n + 1) >> 1);
        break;
    }
    return 0;
}

补充

对于这段代码:
f[prime[j] * i] = (f[i] / (cnt[i] + 1)) * (cnt[i] + 2);

首先要清楚, p r i m e [ i ] prime[i] prime[i] 一定是 i i i 的最小素因数

之后可以用下面的例子来进行理解:

假设 i = p a i = p^a i=pa,那么 f [ i ] = a + 1 f[i] = a + 1 f[i]=a+1 c n t [ i ] = a cnt[i] = a cnt[i]=a。当 p r i m e [ j ] = p prime[j] = p prime[j]=p 时, p r i m e [ j ] ∗ i = p ( a + 1 ) prime[j] * i = p^{(a+1)} prime[j]i=p(a+1),其约数个数 f [ p r i m e [ j ] ∗ i ] f[prime[j] * i] f[prime[j]i] 应该是 a + 2 a + 2 a+2。这正好等于 f [ i ] / ( c n t [ i ] + 1 ) ∗ ( c n t [ i ] + 2 ) f[i] / (cnt[i] + 1) * (cnt[i] + 2) f[i]/(cnt[i]+1)(cnt[i]+2)

好的,我可以为您讲解一下拉格朗日方程的推导过程,以一维情况为例。 假设一个质点在一维坐标系上运动,其位置用 $q$ 表示,时间用 $t$ 表示。该质点的运动可以用拉格朗日函数 $L(q, \dot{q}, t)$ 描述,其中 $\dot{q}=\frac{dq}{dt}$ 表示质点的速度。 为了推导欧拉-拉格朗日方程,我们首先需要定义一个重要的量,即作用量 $S$。作用量定义为 $$S = \int_{t_1}^{t_2} L(q, \dot{q}, t) dt$$ 其中 $t_1$ 和 $t_2$ 表示质点运动的起点和终点。作用量可以理解为一个积分,它是拉格朗日函数在时间 $t_1$ 到 $t_2$ 内的时间积分。 接下来,我们需要考虑如何对作用量进行变分。变分是一种数学操作,它类似于求导,但是是对函数进行微小的偏移,即对函数进行微小的扰动。对于一个作用量 $S$,其变分可以表示为 $$\delta S = S[q+\delta q] - S[q]$$ 其中 $\delta q$ 表示对 $q$ 进行微小的扰动。 现在我们来考虑如何对作用量进行变分。首先,我们将作用量中的积分拆开,得到 $$\delta S = \int_{t_1}^{t_2} \left( \frac{\partial L}{\partial q} \delta q + \frac{\partial L}{\partial \dot{q}} \delta \dot{q} \right) dt$$ 其中第一个式子中的 $\frac{\partial L}{\partial q}$ 表示对 $L$ 关于 $q$ 的偏导数,第二个式子中的 $\frac{\partial L}{\partial \dot{q}}$ 表示对 $L$ 关于 $\dot{q}$ 的偏导数。 现在我们需要将 $\delta \dot{q}$ 转化为 $\delta q$。由于 $\dot{q}=\frac{dq}{dt}$,我们可以得到 $$\delta \dot{q} = \frac{d}{dt} \delta q$$ 将上式代入到 $\delta S$ 中,得到 $$\delta S = \int_{t_1}^{t_2} \left( \frac{\partial L}{\partial q} \delta q + \frac{\partial L}{\partial \dot{q}} \frac{d}{dt} \delta q \right) dt$$ 接下来,我们需要对第二个式子进行分部积分,得到 $$\delta S = \int_{t_1}^{t_2} \left( \frac{\partial L}{\partial q} - \frac{d}{dt} \frac{\partial L}{\partial \dot{q}} \right) \delta q dt + \left[ \frac{\partial L}{\partial \dot{q}} \delta q \right]_{t_1}^{t_2}$$ 现在我们需要对第二个式子进行讨论。由于质点在运动时,其位置 $q$ 和速度 $\dot{q}$ 在起点和终点上都是确定的,因此 $\delta q$ 在 $t_1$ 和 $t_2$ 处的值都应该为 0。因此,第二个式子等于 0。 最终,我们得到了欧拉-拉格朗日方程: $$\frac{d}{dt} \frac{\partial L}{\partial \dot{q}} - \frac{\partial L}{\partial q} = 0$$ 这个方程描述了质点的运动。如果我们能够求出拉格朗日函数 $L$,那么欧拉-拉格朗日方程就可以帮助我们计算质点的运动。
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

若亦_Royi

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值