Aitken's delta-squared process

久しぶりの更新です。

参照
http://en.wikipedia.org/wiki/Aitken%27s_delta-squared_process

つまり少ない計算量で数列の収束を速めるためのアルゴリズムの一つである。(wiki参照)

#include<cstdio>

#define N 50

void delta(double array[]){
   double p, q;

   for(int i = 0;i < N;i++){
      for (int j = i; j >= 2; j -= 2) {
         q = array[j] - 2 * array[j - 1] + array[j - 2];
         if (q == 0) {
            printf("convergence\n");
            return;
         }
         p = array[j] - array[j - 1];
         array[j - 2] = array[j] - p * p / q;
      }
      printf("%d % .20e\n", i, array[i]);
   }
   
}

int main(void){
   static double array[N];
   double t = 1;

   array[0] = 1;

   printf("/////////////////////////////////////////////////\n");
   
   for(int i = 1;i < N;i++){
      t *= (-0.1) * i;
      //std :: cout << t << std :: endl;
      array[i] = array[i - 1] + t;
      printf("%d % .20e\n", i, array[i]);
   }
   
   printf("/////////////////////////////////////////////////\n\n");
   
   delta(array);
   
   return 0;
}

result

/////////////////////////////////////////////////
1 9.00000000000000022204e-01
2 9.20000000000000039968e-01
3 9.14000000000000034639e-01
4 9.16399999999999992362e-01
5 9.15200000000000013500e-01
6 9.15920000000000067431e-01
7 9.15416000000000118497e-01
8 9.15819200000000166462e-01
9 9.15456320000000212112e-01
10 9.15819200000000166462e-01
11 9.15420032000000216676e-01
12 9.15899033600000245237e-01
13 9.15276331520000208108e-01
14 9.16148114432000171270e-01
15 9.14840440064000115505e-01
16 9.16932719052800160320e-01
17 9.13375844771840195158e-01
18 9.19778218477568243472e-01
19 9.07613708436684984981e-01
20 9.31942728518451390940e-01
21 8.80851786346741882916e-01
22 9.93251859124502867182e-01
23 7.34731691735652558961e-01
24 1.35518009346889334310e+00
25 -1.95940910864208284181e-01
26 3.83697370040185603557e+00
27 -7.05189575001651824948e+00
28 2.34369387111549301039e+01
29 -6.49806812262422965887e+01
30 2.00272178585949376384e+02
31 -6.22011686831844826884e+02
32 2.00929668250509689642e+03
33 -6.67402093630681156355e+03
34 2.28492589676536808838e+04
35 -8.04822206962080381345e+04
36 2.91511106093694223091e+05
37 -1.08486420302894408815e+06
38 4.14536197163708135486e+06
39 -1.62525201095604188740e+07
40 6.53390082152295857668e+07
41 -2.69186257916409492493e+08
42 1.13581985983647465706e+09
43 -4.90570644650092697144e+09
44 2.16770093013836402893e+10
45 -9.79452115640969238281e+10
46 4.52317004417113769531e+11
47 -2.13391541069457666016e+12
48 1.02800001818415390625e+13
49 -5.05481862215854296875e+13
/////////////////////////////////////////////////

0 1.00000000000000000000e+00
1 9.00000000000000022204e-01
2 9.20000000000000039968e-01
3 9.14000000000000034639e-01
4 9.16399999999999992362e-01
5 9.15200000000000013500e-01
6 9.15920000000000067431e-01
7 9.15416000000000118497e-01
8 9.15819200000000166462e-01
9 9.15456320000000212112e-01
10 9.15819200000000166462e-01
11 9.15420032000000216676e-01
12 9.15899033600000245237e-01
13 9.15276331520000208108e-01
14 9.16148114432000171270e-01
15 9.14840440064000115505e-01
16 9.16932719052800160320e-01
17 9.13375844771840195158e-01
18 9.19778218477568243472e-01
19 9.07613708436684984981e-01
convergence

多分合ってると思う。