表5の出力に用いたCプログラム(文字コードはUTF-8)
/* The rho-algorithm */
#include <stdio.h>
#include <math.h>
#define NMAX 10
#define PI 3.141592653589793238462
int main(void){
int n,k,kend;
double s[NMAX+1],tv;
double rho[NMAX+1][NMAX+1]; /* rho[n][k] = \rho_{k}^{(n)} */
tv = PI * PI/6.0;
s[0] = 0.0;
for(n = 1;n <= NMAX;n++){
s[n] = s[n-1] + 1.0 / ((double) n * n);
}
for(n = 1;n <= NMAX;n++){
rho[n][0] = s[n];
if (n > 1) {
rho[n-1][1] = 1.0/(s[n] - s[n-1]);
}
for(k = 2;k < n;k++){
rho[n-k][k] = rho[n-k+1][k-2] + ((double) k)
/ (rho[n-k+1][k-1] - rho[n-k][k-1]);
}
}
for(n = 1;n <= NMAX;n++){
printf("%2d %17.12f",n,rho[n][0]);
if (n <= 2){
printf("\n");
continue;
}
kend = (n-1)/2;
for(k = 1;k <= kend;k++){
printf("%20.15f",rho[n-2*k][2*k]);
}
printf("\n");
}
for(n = 1;n <= NMAX;n++){
printf("%2d %10.2e",n,rho[n][0]-tv);
if (n <= 2){
printf("\n");
continue;
}
kend = (n-1)/2;
for(k = 1;k <= kend;k++){
// printf("%10.2e",rho[n-2*k][2*k]-tv);
printf("%7.2f",-log10(fabs((rho[n-2*k][2*k]-tv)/tv)));
}
printf("\n");
}
return 0;
}
お気づきのことがありましたら、下記のメールアドレスにご一報ください。
メールアドレス 長田直樹 <osada "at" cis.twcu.ac.jp> (スパム防止のため@を"at"で置き換えています。)