今日は原始的なプログラム。
ニュートン法で解を求めてみます。
そのうちこのﾌﾞﾛｸﾞでC言語の書き方をまとめたものなんか
出して行きたいですね。//-------------------------------------------------------------------------C
// 2007/04/20 Akitoshi Takayasu
// ニュートン法の計算 f(x)=x^3-3x^2+x+3
// ** 安定性に欠ける **
//-------------------------------------------------------------------------C
#include
//double x[51];
//double y[51],z[51];
FILE *fp;
double x,y;
int i,N;
void main() {
fp=fopen("Newton02.txt","w");
printf("N (:Number of Iteration) ? \n");
scanf("%d", &N );
//x[0]=1;
//y[0]=1;
//z[0]=1;
x=1;
for(i=1; i うまく載らないな。。なぜ？？
情報数学Ⅰのレポート意味もなくアップします。
もしかしたら困っている人がたどり着くかもだしね。
情報数学Ⅰ前期レポート
1E04J×××× ××××
・回答レベル: 応用
・プログラムリスト、解説:
//-------------------------------------------------------------------C
// 情報数学Ⅰ前期レポート(応用)
// Akitoshi Takayasu 2007/07/09
// File name = IMI.c
//-------------------------------------------------------------------C
// 各種ヘッダファイルを宣言
#include
#include
#include
// グローバル変数宣言
char line[20];
double x, y, pi, eps;
int i, j, count;
//===================================================================C
void main(){
//===================================================================C
// モンテカルロ法による円周率の計算
//-------------------------------------------------------------------C
// pi : πの近似値
// eps : 相対誤差
// count : ランダムに与えられた点のうち扇型に入る点の個数
// i : 反復回数,点の総数
// pi = 4.0*(count/i)
//-------------------------------------------------------------------C
// 初期設定
count = 0;
i = 1;
j = 100;
line[0] = 'y';
// 擬似乱数のシードを指定
srand((unsigned int)time(NULL));
// Main Loop-----------------------------------------------------------C
while((line[0] != 'n')&&(line[0] != 'N')){
// x,yにそれぞれ0以上1未満の実数を代入し、扇型内に含まれるか判定
// そしてπの近似値を計算している
x = (double)rand()/(RAND_MAX + 1.0);
y = (double)rand()/(RAND_MAX + 1.0);
if((x*x) + (y*y) < pi =" 4.0*((double)count/i);" i ="=" eps =" fabs(M_PI" pi="%f" i="=" count =" 0;" i =" 0;" j =" 100;">
・プログラムの実行結果:
(繰り返し回数: 100回) pi=2.800000 (相対誤差:0.108732)
(繰り返し回数: 1000回) pi=3.148000 (相対誤差:0.002040)
(繰り返し回数: 10000回) pi=3.140000 (相対誤差:0.000507)
(繰り返し回数: 100000回) pi=3.144520 (相対誤差:0.000932)
(繰り返し回数: 1000000回) pi=3.144516 (相対誤差:0.000931)
(繰り返し回数:10000000回) pi=3.142208 (相対誤差:0.000196)
もう一度やりますか？(y/n)
y
(繰り返し回数: 100回) pi=2.960000 (相対誤差:0.057803)
(繰り返し回数: 1000回) pi=3.156000 (相対誤差:0.004586)
(繰り返し回数: 10000回) pi=3.117200 (相対誤差:0.007764)
(繰り返し回数: 100000回) pi=3.139640 (相対誤差:0.000622)
(繰り返し回数: 1000000回) pi=3.141672 (相対誤差:0.000025)
(繰り返し回数:10000000回) pi=3.141430 (相対誤差:0.000052)
もう一度やりますか？(y/n)
h
(繰り返し回数: 100回) pi=3.120000 (相対誤差:0.006873)
(繰り返し回数: 1000回) pi=3.220000 (相対誤差:0.024958)
(繰り返し回数: 10000回) pi=3.147200 (相対誤差:0.001785)
(繰り返し回数: 100000回) pi=3.136440 (相対誤差:0.001640)
(繰り返し回数: 1000000回) pi=3.139044 (相対誤差:0.000811)
(繰り返し回数:10000000回) pi=3.141416 (相対誤差:0.000056)
もう一度やりますか？(y/n)
N
久しぶりの投稿になってしまいました。
無事にNAS2007での発表も終わり、
本格的に研究を進めています。
まずはテストプログラムをC言語に翻訳しました。（今日までに終了）
今日はそのプログラムを公開します。
FORTRAN⇔C言語の翻訳は大変でした。
初めての大型プログラムなので、汚いですが
頑張ったので載せます。
//-------------------------------------------------------------------------C
// 2007/06/28 Akitoshi Takayasu
// 係数の算出（mod f(θ)）
// Test program for Function Multiplication and Reduction
// File name = KEISU.c
//-------------------------------------------------------------------------C
#include
#define ND 101
FILE *fp1, *fp2;
int N, i, j, k, NV;
double C[3], A[3], B[3][ND];
double A1, A2, B1, B2;
//=========================================================================C
void KEISU(){
//=========================================================================C
// Function Multiplication and Reduction
// f(X) = A[0]*X*X + A[1]*X + A[2]
// F(X) = (B[1][1]*X + B[2][1])*...*(B[1][N]*X + B[2][N])
// g(X) = C[1]*X + C[2]
// g(X) = F(X) mod f(X)
//-------------------------------------------------------------------------C
/*
// Input
puts("N(=イデアルの個数),A(0)*X**2+A(1)*X+A(2)");
puts("(B(1,1)*X+B(2,1))*...*(B(1,N)*X+B(2,N))を入力してください。");
scanf_s("%d",&N);
for(i=0; i<=2; i++){ scanf_s("%lf",&A[i]); printf("A(%d)=%lf\n",i,A[i]); } for(j=1; j<=N; j++){ scanf_s("%lf",&B[1][j]); scanf_s("%lf",&B[2][j]); printf("B(1)(%d)=%lf\n",j,B[1][j]); printf("B(2)(%d)=%lf\n",j,B[2][j]); } */ // First Multiplication if(N%2 == 1){ C[1]=A[0]*B[1][1]; C[2]=A[0]*B[2][1]; }else{ C[1]=B[1][1]; C[2]=B[2][1]; } //printf("A(1)=%lf,A(2)=%lf\n",A[1],A[2]); // Main Loop for(k=2; k<=N; k++){ A1=C[1]; A2=C[2]; B1=B[1][k]; B2=B[2][k]; C[1]=A[0]*A2*B1 + A1*(A[0]*B2-B1*A[1]); C[2]=A[0]*A2*B2 - A1*B1*A[2]; } //printf("C(1)=%f,C(2)=%f\n",C[1],C[2]); } main(){ // Initial Set for A /*puts("A(0)*X**2+A(1)*X+A(2)の係数を入力"); for(i=0; i<=2; i++){ scanf_s("%lf",&A[i]); printf("A(%d)=%f\n",i,A[i]); }*/ A[0]=2; A[1]=0; A[2]=-27; // File Open and Read fp1 = fopen("INFunc.txt","r"); fp2 = fopen("OutFunc.txt","w"); while(100){ if(fscanf(fp1,"%d",&N)==EOF) {break;} //fprintf(fp2,"%d, \n",N); for(i=1; i<=N; i++){ fscanf(fp1,"%lf",&B[1][i]); fscanf(fp1,"%lf",&B[2][i]); //fprintf(fp2,"%lf, %lf,\n",B[1][i],B[2][i]); } //fprintf(fp2,"\n"); // Computation Function KEISU(); NV = (N-1)/2; // Output fprintf(fp2,"%17.0f%17.0f%5d\n",C[1],C[2],NV); } } これはふるいで得られたデータをもとに modをとり次数を落としていきます。代数平方根の計算のステップ1です。
