数列 a(1)=1, a(n+1)=a(n)[1-a(n)/3], 求证 5/2<100a(100)<3

TSC999

问题 ①： 2022 年浙江高考数学卷第 10 题：



问题 ②： 如何用软件 mathematica 在短时间内求出该数列第 100 项的近似数值？ 有一个指令可以使用，但是时间太长，例如计算第 25 项，需要用 18 秒。
想算到第 100 项可能需要几十个小时甚至几个月？

1. Timing[
   
2. N@RecurrenceTable[{a[n + 1] == a[n] (1 - a[n]/3), a[1] == 1},
   
3.    a, {n, 25, 25}]
   
4. ]

复制代码

计算结果是 a(25)≈0.101709。

还有一个程序计算稍快一点，计算第 25 项用时 12 秒，虽有改进，也不咋的：

1. Timing[a[1] = 1; a[2] = 2/3;
   
2. a[n_] := a[n] = a[n - 1] - a[n - 1] a[n - 1]/3;
   
3. N@a[25]
   
4. ]

复制代码

上面这两个程序能不能改进到用时大量缩短，达到计算第 100 项不超过 10 分钟？

wayne

数学归纳法证明$a_n<=\frac{3}{n+2}$.对于 $n=1$成立。 对于一般情况，因为$f(x) =x -1/3x^2$在$(0,3/2)$内单调递增。 所以如果 $a_n< \frac{3}{n+2}$,那么，$a_{n+1}=f(a_n)<f(\frac{3}{n+2}) <3/{n+3}$ 确实恒成立,所以，对于任意的$n>1$,$a_n<\frac{3}{n+2}$

同样，可以证明$na_n$单调递增。因为$(n+1)a_{n+1}-na_n =a_n(1-\frac{n+1}{3}a_n)>=\frac{1}{n+2}a_n >0$

1. With[{m = 5000}, RecurrenceTable[{a[n + 1] == a[n] (1 - a[n]/3), a[1] == 1`100}, a, {n, m - 10, m}]] // Reverse

复制代码

更进一步的，得到$a_n <=\frac{3}{n+2}-\frac{275 (n-1)}{54 (n+2)^3}$ 等号在$n=1,3$的时候取得。$n=12$取得最大残差。

。

lsr314

需要证明$na_n$是递增的

mathe

查看:
https://bbs.emath.ac.cn/forum.ph ... =8790&pid=61918
先定义$b_n=\frac1{a_n}$,
记\(h(x)=\frac{x^2}{x-\frac13}\)得到\(b_{n+1}=h(b_n)\)

wayne

受楼上启发，设$c_n=b_n-1/3 = 1/a_n-1/3$，那么，$c_1=2/3, c_n =\frac{n+1}{3} +1/9\sum_{k=1}^{n-1}1/c_k$

TSC999

mathe

https://bbs.emath.ac.cn/forum.ph ... 2433&fromuid=20
$v(x)=x^2/(x-1/3)=x+1/3+1/{9x}+1/{27x^2}+1/{81x^3}+\cdots$
于是我们通过$a_1=1, a_{n+1}=v(a_n)$递推式得出来的数列和本题中数列互为倒数。
对比链接中记号有$b_0=1/3,b_1=1/9$

于是$1/3 \ln(m)-1/3 m = 1/3 \ln(m')-1/3 (m'-1)$
即$\ln(m')-m'+1=\ln(m)-m$, 即$ m' =m +1 +\ln({m'}/m)$
于是$m'=J(m)=m+1 + \ln(\frac{J(m)}m)$
使用迭代计算

1. ? x+O(x^2)
   
2. %8 = x + O(x^2)
   
3. ? log(1+x+x*%)
   
4. %9 = x + 1/2*x^2 + O(x^3)
   
5. ? log(1+x+x*%)
   
6. %10 = x + 1/2*x^2 - 1/6*x^3 + O(x^4)
   
7. ? log(1+x+x*%)
   
8. %11 = x + 1/2*x^2 - 1/6*x^3 - 5/12*x^4 + O(x^5)
   
9. ? log(1+x+x*%)
   
10. %12 = x + 1/2*x^2 - 1/6*x^3 - 5/12*x^4 - 1/20*x^5 + O(x^6)
    
11. ? log(1+x+x*%)
    
12. %13 = x + 1/2*x^2 - 1/6*x^3 - 5/12*x^4 - 1/20*x^5 + 49/120*x^6 + O(x^7)
    
13. ? log(1+x+x*%)
    
14. %14 = x + 1/2*x^2 - 1/6*x^3 - 5/12*x^4 - 1/20*x^5 + 49/120*x^6 + 15/56*x^7 + O(x^8)
    
15. ? log(1+x+x*%)
    
16. %15 = x + 1/2*x^2 - 1/6*x^3 - 5/12*x^4 - 1/20*x^5 + 49/120*x^6 + 15/56*x^7 - 457/1260*x^8 + O(x^9)
    
17. ? log(1+x+x*%)
    
18. %16 = x + 1/2*x^2 - 1/6*x^3 - 5/12*x^4 - 1/20*x^5 + 49/120*x^6 + 15/56*x^7 - 457/1260*x^8 - 49/90*x^9 + O(x^10)
    
19. ? log(1+x+x*%)
    
20. %17 = x + 1/2*x^2 - 1/6*x^3 - 5/12*x^4 - 1/20*x^5 + 49/120*x^6 + 15/56*x^7 - 457/1260*x^8 - 49/90*x^9 + 391/2016*x^10 + O(x^11)

复制代码

得到$J(1/x)=1/x+1+x + 1/2*x^2 - 1/6*x^3 - 5/12*x^4 - 1/20*x^5 + 49/120*x^6 + 15/56*x^7 - 457/1260*x^8 - 49/90*x^9 + 391/2016*x^10 + O(x^11)$
然后得到渐进式$s(x)\approx 1/{3x}-(x/6+{5x^2}/18+{13x^3}/27+{2009x^4}/2160+{6973x^5}/3600+{1923143x^6}/453600+{30503779x^7}/3175200+{850720483x^8}/38102400+{1209787591x^9}/22861440+{3204878097563x^10}/25147584000+{43027959566971x^11}/138311712000+\cdots)$

\(s(\frac1{J(m)}) = v(s(\frac1m))\)， 而且\(s(x)=1/{3x}+P_1 x+P_2x^2+...\)
根据现有展开结果，我们可以先计算$a_{1000}=0.0029760594191188835248176874067736714838$, 计算\(\frac1m=s^{-1}(1/{a_{1000}})=0.00099201931744002304343664937365698435425\)。
所以对应\(m=1008.0448862433159869955\), 那么对应极限应该约等于$c_{a_0} =0.37637275522427338461471882180513818568$
但是代回n=2000情况，需要计算对应的m为$1/3 \ln(m)+c_{a_0}-1/3(m-n)=0$，对应m=2008.7343784059206482994905154424269013
然后计算对应数列倒数值为\(s(1/m)=669.57804309542223061154167143880331827\),得出原始数列第2000项为\(1/{s(1/m)}=0.0014934778855307969063639720682167818331\).
直接迭代计算的结果为0.0014934778855307969063639720682167818339，差别应该是舍入误差引起的。

于是对于本题中数列中较大的n,可以通过如下方法计算第n项:
i) 解关于m的方程 $\ln(m)-(m-n)+3*0.37637275522427338461471882180513818568=0$,得到实数m (m略大于n,可以在(n,2n)中使用二分法找m)
ii) 计算$1/{s(1/m)}$, 即为原数列第n项。

mathe

如果我们尝试把上面的s(x)的渐进展开式试着转化为连分式形式，会得到
$s(x)=1/{3x}+[-x/6, -{5x}/3, -x/15, -{619x}/210, {419x}/120,  {577191x}/9077635, - {4163475702921x}/1746509741095,
  {39018927590570016701x}/84109224655563621885, {8320244433144139792056998559830x}/8142103810321752087542998506339,
- {146451852647634675917906267906355624753965x}/57252060088539150878523070380820686585512]$
比较有意思的是不像前面渐进展开式，系数绝对值越来越大，连分式展开式中每项的值都不是特别大。

mathe

47#另外一种推导方法有\(C(x)\approx 3x-\ln(x)+\alpha_0 + \sum_{h=1}^{\infty} \frac{\alpha_h}{x^h}\), 而且$C(v(x))=C(x)+1$
设$x=1/t, C(x)=3x-\ln(x)+\alpha_0 +h(\frac1x)$,得到\(h(t)-h(t-\frac{t^2}3)=\frac t{3-t}-\ln(\frac3{3-t})\)
由此得到\(h(t)=\frac t6+\frac{t^2}{27}+\frac{13t^3}{972}+\frac{113t^4}{19440}+\frac{1187t^5}{437400}+\frac{877t^6}{688905}+\frac{14569t^7}{25719120}+\frac{176017t^8}{793618560}+\frac{1745717t^9}{26784626400}+\frac{88217t^{10}}{15345358875}-\frac{147635381t^{11}}{19445638766400}-\frac{3238110769t^{12}}{827323540243200}+\frac{63045343657t^{13}}{37643221081065600}+\frac{24855467017t^{14}}{7125970336860375}+O(t^{15})\)，看上去这是一个收敛的级数。
而我们取倒数以后数列第100项为：35.261925907006479247389066690274451303，
计算$C(35.261925907006479247389066690274451303)=102.22773055434092984523940170232737503+\alpha_0$，计算中迭代式最后一项值为8E-28,说明至少精确到小数点后27位。所以我们得到$C(a_{1000})=900+102.22773055434092984523940170232737503+\alpha_0$,二分法（牛顿法也应该可行）求解得到$a_{1000}=336.01479647071971871965605047571202386$,
如果回到原始数列还需要求一次倒数，就是0.0029760594191188835248176874067737027128精度超过了27位。
可以看出上面计算过程中$\alpha_0$的真正取值我们不需要关心。

现在我们再查看\(h(t)-h(t-\frac{t^2}3)=\frac t{3-t}-\ln(\frac3{3-t})\)，设\(h(t)=\sum_{k=1}^{\infty} h_k t^k\)
于是我们有
\(h_k - \sum_{s=0}^{\lfloor \frac k2\rfloor} {k-s\choose s}(-\frac13)^s h_{k-s}=\frac1{3^k}-\frac 1{k 3^k}\)
即
\(\sum_{s=1}^{\lfloor \frac k2\rfloor} {k-s\choose s}(-1)^{s-1}3^{k-s} h_{k-s}=\frac{k-1}{k}\)
由此可以计算出$h_k$,比较有意思，在k不超过48时，$h_k$的绝对值都不大，但是超过48以后，会快速增加，比如50以内的近似值如下:
[0.16666666666666666666666666666666666667, 0.037037037037037037037037037037037037037, 0.013374485596707818930041152263374485597, 0.0058127572016460905349794238683127572016, 0.0027137631458619112940100594421582075903, 0.0012730347435422881239067795995093663132, 0.00056646572666560908771373204059858968736, 0.00022179042788515429881075361947180267558, 6.5176081754121461257342756888332032139 E-5, 5.7487739920973337288600883242621460034 E-6, -7.5922104063301982988954141657144179788 E-6, -3.9139594263788566532322554993828840680 E-6, 1.6748126713500501434084340204178928013 E-6, 3.4880115748490539382016935213242119700 E-6, 1.7363494737148565213694432284329098707 E-6, -7.9799420354630179737419381898282137015 E-7, -1.6663666605547081554496798272929129782 E-6, -4.8980993654953182067160837188121197326 E-7, 1.1105335135062783470551270653313374172 E-6, 1.2155150887944245256861254574297638516 E-6, -3.7262538223826892478950604364811930886 E-7, -1.6859236599273004152507711520678032307 E-6, -6.4634511231469059808323594800515740982 E-7, 1.9668726787931563176411535184133687908 E-6, 2.3249532336163214280369647854475855145 E-6, -1.7590830472271215203914538899448131193 E-6, -5.4348511865852951752738608483990514820 E-6, -1.7776362697917386439362668349121145023 E-7, 1.1443113858678654387830285865536416390 E-5, 8.0991151110068231025564920218564512583 E-6, -2.2571464149286119366016217567932791543 E-5, -3.6372648930010308475889577225317185439 E-5, 3.8558160191635557768001768419753515167 E-5, 0.00013401763051510425204805030575272547318, -3.1649076903916762404526823514940657499 E-5, -0.00047010120902109810550717849867535525221, -0.00020380104534926813896542143755473092232, 0.0016252503737284226898849132052924148248, 0.0018988820328237181903767321288950612215, -0.0055023715562717999535623302581466140326, -0.012034843927255075439229442533738015443, 0.017308691299680244760881895244865116243, 0.069955901283703582133372768457402055335, -0.040971229207006948112488694223475275324, -0.40002895638879882523958547281664757706, -0.040430165593565901058044155533920062087, 2.3033291921212683769970128889917119916, 1.7931348656167185021625926103215511738, -13.420200122406705514524236511193599842, -20.882362404153083509726322768841702422]
但是51项到60项如下：
[78.513621003297250420074188707102682220, 202.58713829274578879553773336619242812, -449.99520262010819218204869756666128052, -1873.4863539770463771221218334858172231, 2366.9955239364074357260471634171465472, 17249.065560767631956266214434624779288, -8940.7592181048466843451670448875449069, -160960.78661872319570785726231645359908, -24829.450164495412949086768442609647625, 1533556.2377505922967945825482291368540]

ShuXueZhenMiHu

IMG_20220624_0001_compressed.pdf (110.59 KB, 下载次数: 1)

2022-6-25 02:42 上传

点击文件名下载附件