ベッセル関数 (J0(x) と J1(x)) の数表 (例えば x は 0~10 の範囲で 0.1 刻み) を出力する C のプログラムを作れ。という課題なのですができません教えて下さい
> itagi_juniorさん
知恵袋過去を探したら
detail.chiebukuro.yahoo.co.jp/qa/question_detail/q1112026261
がありました。
あとGSLと言うライブラリにもベッセル関数含まれてます。
それのソース見れば、答になってるはず。
http://ftp.jaist.ac.jp/pub/GNU/gsl/gsl-1.15.tar.gz
展開してそれらしい名前のファイルが以下
gsl-1.15/specfunc/bessel_J0.c
gsl-1.15/specfunc/bessel_J1.c
答えあわせも使えると思います。
よろしくお願い致します。
****追記1
あと、何重にも同じ投稿されているので、他消していただけないでしょうか。
よろしくお願い致します。
*** 補足分
私も良く見てませんが、テイラー展開(かな?)した項を計算していって、項の値がある一定以下になるまで計算してるってパターンじゃないのでしょうか?
後で過去のプログラムとGS見ときます。何かわかったら追記します。
よろしくお願い致します。
*** 追記2
上のプログラム後半削って、最初の方にあるループの範囲を変更したらJ0だと思う。
数値未確認、GSLと比較したら検算になると思います。
後は以下を頑張ってください。
1)検算
2)検算合ってたら、下のプログラムを改造してJ1を計算する
よろしくお願い致します。
--- em163.c
#include <stdio.h>
#include <math.h>
double Bess0( // J0(x) を計算し、戻り値で返す
double x, //[i] 独立変数
double eps, //[i] 打ち切り誤差
int max, //[i] 級数項数の上限
int* lp) //[o] 計算した項数
{
int k;
double fct, dJ, J0;
k = 0; J0 = 1;
dJ = 1; fct = -pow(x/2,2.);
do {
if ( ++k > max ) break;
dJ *= fct/(k*k);
J0 += dJ;
} while ( fabs(dJ) > eps );
*lp = k;
return J0;
}
int main(){
int loop, maxlp=1000;
double x, J0, eps=1.e-12;
printf("1.0<x<10.0 でのJ0(x)の分布\n");
printf(" x N J0\n");
for(x=1.;x<=10.;x+=0.1) {
J0 = Bess0( x, eps, maxlp, &loop );
printf("%f %3d %9f\n",x,loop,J0);
}
printf("\n");
return 0;
}
---
0 件のコメント:
コメントを投稿