1: /* taylor_sin.c */
 2: /* taylor 展開 */
 3: /* sin t */
 4: /* 1, 3, 5, ... k_max 次までの展開の誤差を比較 */
 5: /* t のデータ点は 0 - 2 pi を sample 等分 */
 6: #include <stdio.h>
 7: #include <math.h> /* sin 関数を使うのに必要 */
 8: 
 9: /* プロトタイプ宣言 taylor_sin */
10: double taylor_sin( int, double );
11: 
12: 
13: int main( void ){
14: 
15:   /* k_max 次までの taylor 展開を計算 */
16:   int k_max; k_max = 15;
17: 
18:   /* t のデータ点数は sample + 1 */
19:   int sample; sample = 8;
20:   
21:   /* taylor 展開のデータを保存する配列 */
22:   /* 0 行目は時間データ */
23:   /* i 行目は i 次近似 */
24:   /* i + 1 行目は i 次近似の誤差 */
25:   int n; n = k_max + 2;
26:   int m; m = sample + 1;
27:   double data[ n ][ m ];
28: 
29:   /* sin の値を保存する配列 */
30:   double data_sin[ m ];
31:   
32:   /* カウンタ用変数の宣言 */
33:   int i, j;
34: 
35:   /* 円周率 */
36:   double pi; pi = 3.141592;
37:   
38:   /* data[ 0 ][ * ] は t の値 */
39:   /* td は t の刻み幅 */
40:   /* sample は int なので double へキャスト */
41:   double td; td = ( ( 2.0 * pi ) / (double) sample );
42:   for ( j = 0; j < m; j = j + 1 ) {
43: 
44:     data[ 0 ][ j ] = td * j;
45: 
46:   }
47: 
48:   /* sin の値を計算 */
49:   for ( j = 0; j < m; j = j + 1 ) {
50:     
51:     data_sin[ j ] = sin( data[ 0 ][ j ] );
52:     
53:   }
54:   
55:   /* i 次近似 */
56:   /* i は奇数 */
57:   /* data[ i ][ * ] は i 次近似の値 */
58:   /* data[ i + 1 ][ * ] は i 次近似の誤差 */
59:   for ( i = 1; i < k_max + 1; i = i + 2 ) {
60:     for ( j = 0; j < m; j = j + 1 ) {
61:       
62:       data[ i ][ j ] = taylor_sin( i, data[ 0 ][ j ] );
63:       data[ i + 1 ][ j ] = data[ i ][ j ] - data_sin[ j ];
64:     
65:     }
66:   }
67:   
68: 
69:   /* 結果を表示 */
70:   printf ( "       t = ");
71:   for ( j = 0; j < m; j = j + 1 ) {
72: 
73:     printf( "%8.4f", data[ 0 ][ j ] );
74: 
75:   }
76:   printf( "\n\n" );
77:   
78:   for ( i = 1; i < n; i = i + 2 ) {
79: 
80:     printf ( "%2dth sin = ", i );
81:     for ( j = 0; j < m; j = j + 1 ) {
82: 
83:       printf( "%8.4f", data[ i ][ j ] );
84:     
85:     }
86:     printf( "\n" );
87: 
88:     printf( "%2dth err = ", i );
89:     for ( j = 0; j < m; j = j + 1 ) {
90: 
91:       printf( "%8.4f", data[ i + 1 ][ j ] );
92:     
93:     }
94:     printf( "\n\n" );
95: 
96:   }
97:   
98:   return 0;
99: 
100: }
101: 
102: 
103: /* 関数 taylor_sin */
104: /* sin( t ) の i 次 taylor 展開近似を求める */
105: double taylor_sin( int i, double t ){
106: 
107:   /* 出力用変数 */
108:   double z;
109: 
110:   /* i 次近似の第 k 項を保存するための変数 */
111:   double zk;
112: 
113:   /* カウンタ用変数 */
114:   int k;
115:   
116:   /* 1 次近似 */
117:   zk = t;
118:   z = zk;
119: 
120:   /* i >= 3 次近似 */
121:   for ( k = 3; k < i + 1; k = k + 2 ) {
122: 
123:     zk = zk * ( ( ( -1.0 ) / ( ( k - 1 ) * k ) ) * t * t );
124:     z = z + zk;
125:     
126:   }
127:   
128:   return z;
129:   
130: }