2018年11月26日

【Javaで数値計算入門】その2:ラグランジェ補間

(n+1)個の既知の点 [x,f(x)] を n次多項式で近似し
求めたい点x=pのf(p)を推定する方法。
public class LagrangeInterpolation {
public static void main(String... args) {
int n_points = 4;

double p = 20.; // x-value for f(p)
double ans=0.;

double L_bunshi1;
double L_bunbo1;
double L_bunshi2;
double L_bunbo2;

double lj;

double[] x = new double[n_points];
double[] f = new double[n_points];
//double[] l = new double[n_points];

x[0]=19.9466;
f[0]=2.712;
x[1]=19.9907;
f[1]=2.714;
x[2]=20.0349;
f[2]=2.716;
x[3]=20.0793;
f[3]=2.718;

for (int j = 0; j < n_points; j++) {
L_bunshi1=1.;
L_bunbo1=1.;
L_bunshi2=1.;
L_bunbo2=1.;

for (int jj = 0; jj <= j - 1; jj++) {
L_bunshi1=L_bunshi1*(p-x[jj]);
L_bunbo1=L_bunbo1*(x[j]-x[jj]);
}
for (int jj = j+1; jj < n_points; jj++) {
L_bunshi2=L_bunshi2*(p-x[jj]);
L_bunbo2=L_bunbo2*(x[j]-x[jj]);
}
lj=(L_bunshi1*L_bunshi2)/(L_bunbo1*L_bunbo2);
ans=ans+f[j]*lj;
}
System.out.printf("x = %8.6f, f(x) = %11.9f" , p,ans); //True f(20) = 2.7144176
}
}

以上を実行してみると
x = 20.000000, f(x) = 2.714421340
が得られた(環境により変化する)

※参考:「理工学のための数値計算法」p.24
    「数値計算入門 C言語版」p.97-98
  


Posted by 和歌山サイエンスカフェインフィニティ at 22:31Comments(0)javaで数値計算入門

2018年11月25日

【Javaで数値計算入門】その1:エイトケン加速

級数の和の加速法の一つエイトケン加速

円周率を求める級数の一つを加速させて計算する。
π/4なので、4倍すると円周率になる。

public class AitkenProcess {

public static void main(String... args) {
double[][] s = new double[10][10];

for (int n = 0; n < 10; n++) {
for (int k = 0; k < n + 1; k++) {
s[n][0] += Math.pow(-1., k) / (2 * k + 1);
}
System.out.printf("n = %3d, k= 0, s = %10.8f \n", n, s[n][0]);
}

for (int n = 1; n < 9; n++) {
for (int k = 0; k < n; k++) {
s[n][k + 1] = (s[n + 1][k] * s[n - 1][k] - (s[n][k] * s[n][k]))
/ (s[n + 1][k] + s[n - 1][k] - 2 * s[n][k]);

System.out.printf("n = %3d, k= %3d, s = %10.8f \n", n, k + 1, s[n][k + 1]);
}
}
}
}


※参考:「理工学のための数値計算法」p.9  


Posted by 和歌山サイエンスカフェインフィニティ at 21:15Comments(0)javaで数値計算入門