import std.stdio;
import std.math;
void main(string[] args)
{
writefln("%d", 3 + 5);
writefln("%d", 3 - 5);
writefln("%d", 3 * 5);
writefln("%d", pow(3, 5));
writefln("%d", 5 / 3);
writefln("%f", 5.0 / 3);
writefln("%f", 5 / 3.0);
writefln("%d", 5 % 3);
writef("%d\n", 5 * 3);
writefln("%3d", 3 * 5);
writefln("%23.20f", 5 / 3.0);
}
import std.stdio;
void main(string[] args)
{
int i = 3 * 5;
writef("3 * 5 = %d\n", i);
}
import std.stdio;
import std.math;
void main(string[] args)
{
const double a = 0;
const double b = 1;
// 蜿ｰ蠖｢蜑・〒遨榊・
int n = 2;
for (int j = 1; j <= 10; j++)
{
double h = (b - a) / n;
double s = 0;
double x = a;
for (int i = 1; i <= n - 1; i++)
{
x += h;
s += f(x);
}
s = h * ((f(a) + f(b)) / 2 + s);
n *= 2;
// 邨先棡繧・ﾏ 縺ｨ豈碑ｼ・
writefln("%2d : %13.10f, %13.10f", j, s, s - PI);
}
}
double f(double x)
{
return 4 / (1 + x * x);
}
import std.stdio;
import std.math;
void main(string[] args)
{
const double a = 0;
const double b = 1;
// 荳ｭ轤ｹ蜑・〒遨榊・
int n = 2;
for (int j = 1; j <= 10; j++)
{
double h = (b - a) / n;
double s = 0;
double x = a + (h / 2);
for (int i = 1; i <= n; i++)
{
s += f(x);
x += h;
}
s *= h;
n *= 2;
// 邨先棡繧・ﾏ 縺ｨ豈碑ｼ・
writefln("%2d : %13.10f, %13.10f", j, s, s - PI);
}
}
double f(double x)
{
return 4 / (1 + x * x);
}
import std.stdio;
import std.math;
void main(string[] args)
{
const double a = 0;
const double b = 1;
// Simpson蜑・〒遨榊・
int n = 2;
for (int j = 1; j <= 5; j++)
{
double h = (b - a) / n;
double s2 = 0;
double s4 = 0;
double x = a + h;
for (int i = 1; i <= n / 2; i++)
{
s4 += f(x);
x += h;
s2 += f(x);
x += h;
}
s2 = (s2 - f(b)) * 2 + f(a) + f(b);
s4 *= 4;
double s = (s2 + s4) * h / 3;
n *= 2;
// 邨先棡繧・ﾏ 縺ｨ豈碑ｼ・
writefln("%2d : %13.10f, %13.10f", j, s, s - PI);
}
}
double f(double x)
{
return 4 / (1 + x * x);
}
import std.stdio;
import std.math;
void main(string[] args)
{
const double a = 0;
const double b = 1;
double t[7][7];
// 蜿ｰ蠖｢蜑・〒遨榊・
int n = 2;
for (int i = 1; i <= 6; i++)
{
double h = (b - a) / n;
double s = 0;
double x = a;
for (int j = 1; j <= n - 1; j++)
{
x += h;
s += f(x);
}
// 邨先棡繧剃ｿ晏ｭ・
t[i][1] = h * ((f(a) + f(b)) / 2 + s);
n *= 2;
}
// Richardson縺ｮ陬懷､匁ｳ・
n = 4;
for (int j = 2; j <= 6; j++)
{
for (int i = j; i <= 6; i++)
{
t[i][j] = t[i][j - 1] + (t[i][j - 1] - t[i - 1][j - 1]) / (n - 1);
if (i == j)
{
// 邨先棡繧・ﾏ 縺ｨ豈碑ｼ・
writefln("%2d : %13.10f, %13.10f", j, t[i][j], t[i][j] - PI);
}
}
n *= 4;
}
}
double f(double x)
{
return 4 / (1 + x * x);
}
import std.stdio;
import std.math;
// 繝・・繧ｿ轤ｹ縺ｮ謨ｰ
const int N = 7;
void main(string[] args)
{
double x[N];
double y[N];
// 1.5蛻ｻ縺ｿ縺ｧ -4.5・・.5 縺ｾ縺ｧ, ・礼せ縺縺大､繧偵そ繝・ヨ
for (int i = 0; i < N; i++)
{
double d = i * 1.5 - 4.5;
x[i] = d;
y[i] = f(d);
}
// 0.5蛻ｻ縺ｿ縺ｧ 荳弱∴繧峨ｌ縺ｦ縺・↑縺・､繧定｣憺俣
for (int i = 0; i <= 18; i++)
{
double d = i * 0.5 - 4.5;
double d1 = f(d);
double d2 = lagrange(d, x, y);
// 蜈・・髢｢謨ｰ縺ｨ豈碑ｼ・
writefln("%5.2f\t%8.5f\t%8.5f\t%8.5f", d, d1, d2, d1 - d2);
}
}
// 蜈・・髢｢謨ｰ
double f(double x)
{
return x - pow(x,3) / (3 * 2) + pow(x,5) / (5 * 4 * 3 * 2);
}
// Lagrange (繝ｩ繧ｰ繝ｩ繝ｳ繧ｸ繝･) 陬憺俣
double lagrange(double d, double x[], double y[])
{
double sum = 0;
for (int i = 0; i < N; i++)
{
double prod = y[i];
for (int j = 0; j < N; j++)
{
if (j != i)
prod *= (d - x[j]) / (x[i] - x[j]);
}
sum += prod;
}
return sum;
}
import std.stdio;
import std.math;
// 繝・・繧ｿ轤ｹ縺ｮ謨ｰ
const int N = 7;
void main(string[] args)
{
double x[N];
double y[N];
// 1.5蛻ｻ縺ｿ縺ｧ -4.5・・.5 縺ｾ縺ｧ, ・礼せ縺縺大､繧偵そ繝・ヨ
for (int i = 0; i < N; i++)
{
double d = i * 1.5 - 4.5;
x[i] = d;
y[i] = f(d);
}
// 0.5蛻ｻ縺ｿ縺ｧ 荳弱∴繧峨ｌ縺ｦ縺・↑縺・､繧定｣憺俣
for (int i = 0; i <= 18; i++)
{
double d = i * 0.5 - 4.5;
double d1 = f(d);
double d2 = neville(d, x, y);
// 蜈・・髢｢謨ｰ縺ｨ豈碑ｼ・
writefln("%5.2f\t%8.5f\t%8.5f\t%8.5f", d, d1, d2, d1 - d2);
}
}
// 蜈・・髢｢謨ｰ
double f(double x)
{
return x - pow(x,3) / (3 * 2) + pow(x,5) / (5 * 4 * 3 * 2);
}
// Neville (繝阪Χ繧｣繝ｫ) 陬憺俣
double neville(double d, double x[], double y[])
{
double w[N][N];
for (int i = 0; i < N; i++)
w[0][i] = y[i];
for (int j = 1; j < N; j++)
{
for (int i = 0; i < N - j; i++)
w[j][i] = w[j-1][i+1] + (w[j-1][i+1] - w[j-1][i]) * (d - x[i+j]) / (x[i+j] - x[i]);
}
return w[N-1][0];
}
import std.stdio;
import std.math;
// 繝・・繧ｿ轤ｹ縺ｮ謨ｰ
const int N = 7;
void main(string[] args)
{
double x[N];
double y[N];
// 1.5蛻ｻ縺ｿ縺ｧ -4.5・・.5 縺ｾ縺ｧ, ・礼せ縺縺大､繧偵そ繝・ヨ
for (int i = 0; i < N; i++)
{
double d = i * 1.5 - 4.5;
x[i] = d;
y[i] = f(d);
}
// 蟾ｮ蛻・膚縺ｮ陦ｨ繧剃ｽ懊ｋ
double d[N][N];
for (int j = 0; j < N; j++)
d[0][j] = y[j];
for (int i = 1; i < N; i++)
{
for (int j = 0; j < N - i; j++)
d[i][j] = (d[i-1][j+1] - d[i-1][j]) / (x[j+i] - x[j]);
}
// ・朱嚴蟾ｮ蛻・膚
double a[N];
for (int j = 0; j < N; j++)
a[j] = d[j][0];
// 0.5蛻ｻ縺ｿ縺ｧ 荳弱∴繧峨ｌ縺ｦ縺・↑縺・､繧定｣憺俣
for (int i = 0; i <= 18; i++)
{
double d1 = i * 0.5 - 4.5;
double d2 = f(d1);
double d3 = newton(d1, x, a);
// 蜈・・髢｢謨ｰ縺ｨ豈碑ｼ・
writefln("%5.2f\t%8.5f\t%8.5f\t%8.5f", d1, d2, d3, d2 - d3);
}
}
// 蜈・・髢｢謨ｰ
double f(double x)
{
return x - pow(x,3) / (3 * 2) + pow(x,5) / (5 * 4 * 3 * 2);
}
// Newton (繝九Η繝ｼ繝医Φ) 陬憺俣
double newton(double d, double x[], double a[])
{
double sum = a[0];
for (int i = 1; i < N; i++)
{
double prod = a[i];
for (int j = 0; j < i; j++)
prod *= (d - x[j]);
sum += prod;
}
return sum;
}
import std.stdio;
import std.math;
// 繝・・繧ｿ轤ｹ縺ｮ謨ｰ
const int N = 7;
const int Nx2 = 14;
void main(string[] args)
{
double x[N];
double y[N];
double yd[N];
// 1.5蛻ｻ縺ｿ縺ｧ -4.5・・.5 縺ｾ縺ｧ, ・礼せ縺縺大､繧偵そ繝・ヨ
for (int i = 0; i < N; i++)
{
double d = i * 1.5 - 4.5;
x[i] = d;
y[i] = f(d);
yd[i] = fd(d);
}
// 蟾ｮ蛻・膚縺ｮ陦ｨ繧剃ｽ懊ｋ
double z[Nx2];
double d[Nx2][Nx2];
for (int i = 0; i < Nx2; i++)
{
int j = i / 2;
z[i] = x[j];
d[0][i] = y[j];
}
for (int i = 1; i < Nx2; i++)
{
for (int j = 0; j < Nx2 - i; j++)
{
if (i == 1 && j % 2 == 0)
d[i][j] = yd[j / 2];
else
d[i][j] = (d[i-1][j+1] - d[i-1][j]) / (z[j+i] - z[j]);
}
}
// ・朱嚴蟾ｮ蛻・膚
double a[Nx2];
for (int j = 0; j < Nx2; j++)
a[j] = d[j][0];
// 0.5蛻ｻ縺ｿ縺ｧ 荳弱∴繧峨ｌ縺ｦ縺・↑縺・､繧定｣憺俣
for (int i = 0; i <= 18; i++)
{
double d1 = i * 0.5 - 4.5;
double d2 = f(d1);
double d3 = hermite(d1, z, a);
// 蜈・・髢｢謨ｰ縺ｨ豈碑ｼ・
writefln("%5.2f\t%8.5f\t%8.5f\t%8.5f", d1, d2, d3, d2 - d3);
}
}
// 蜈・・髢｢謨ｰ
double f(double x)
{
return x - pow(x,3) / (3 * 2) + pow(x,5) / (5 * 4 * 3 * 2);
}
// 蟆朱未謨ｰ
double fd(double x)
{
return 1 - pow(x,2) / 2 + pow(x,4) / (4 * 3 * 2);
}
// Hermite (繧ｨ繝ｫ繝溘・繝・ 陬憺俣
double hermite(double d, double z[], double a[])
{
double sum = a[0];
for (int i = 1; i < Nx2; i++)
{
double prod = a[i];
for (int j = 0; j < i; j++)
prod *= (d - z[j]);
sum += prod;
}
return sum;
}
import std.stdio;
import std.math;
// 繝・・繧ｿ轤ｹ縺ｮ謨ｰ
const int N = 7;
void main(string[] args)
{
double x[N];
double y[N];
// 1.5蛻ｻ縺ｿ縺ｧ -4.5・・.5 縺ｾ縺ｧ, ・礼せ縺縺大､繧偵そ繝・ヨ
for (int i = 0; i < N; i++)
{
double d = i * 1.5 - 4.5;
x[i] = d;
y[i] = f(d);
}
// ・馴・婿遞句ｼ上・菫よ焚縺ｮ陦ｨ繧剃ｽ懊ｋ
double a[N];
double b[N];
double c[N];
double d[N];
for (int i = 1; i < N - 1; i++)
{
a[i] = x[i] - x[i-1];
b[i] = 2.0 * (x[i+1] - x[i-1]);
c[i] = x[i+1] - x[i];
d[i] = 6.0 * ((y[i+1] - y[i]) / (x[i+1] - x[i]) - (y[i] - y[i-1]) / (x[i] - x[i-1]));
}
// ・馴・婿遞句ｼ上ｒ隗｣縺・(繝茨ｼ阪・繧ｹ豕・
double g[N];
double s[N];
g[1] = b[1];
s[1] = d[1];
for (int i = 2; i < N - 1; i++)
{
g[i] = b[i] - a[i] * c[i-1] / g[i-1];
s[i] = d[i] - a[i] * s[i-1] / g[i-1];
}
double z[N];
z[0] = 0;
z[N-1] = 0;
z[N-2] = s[N-2] / g[N-2];
for (int i = N - 3; i >= 1; i--)
{
z[i] = (s[i] - c[i] * z[i+1]) / g[i];
}
// 0.5蛻ｻ縺ｿ縺ｧ 荳弱∴繧峨ｌ縺ｦ縺・↑縺・､繧定｣憺俣
for (int i = 0; i <= 18; i++)
{
double d1 = i * 0.5 - 4.5;
double d2 = f(d1);
double d3 = spline(d1, x, y, z);
// 蜈・・髢｢謨ｰ縺ｨ豈碑ｼ・
writefln("%5.2f\t%8.5f\t%8.5f\t%8.5f", d1, d2, d3, d2 - d3);
}
}
// 蜈・・髢｢謨ｰ
double f(double x)
{
return x - pow(x,3) / (3 * 2) + pow(x,5) / (5 * 4 * 3 * 2);
}
// Spline (繧ｹ繝励Λ繧､繝ｳ) 陬憺俣
double spline(double d, double x[], double y[], double z[])
{
// 陬憺俣髢｢謨ｰ蛟､縺後←縺ｮ蛹ｺ髢薙↓縺ゅｋ縺・
int k = -1;
for (int i = 1; i < N; i++)
{
if (d <= x[i])
{
k = i - 1;
break;
}
}
if (k < 0) k = N - 1;
double d1 = x[k+1] - d;
double d2 = d - x[k];
double d3 = x[k+1] - x[k];
return (z[k] * pow(d1,3) + z[k+1] * pow(d2,3)) / (6.0 * d3)
+ (y[k] / d3 - z[k] * d3 / 6.0) * d1
+ (y[k+1] / d3 - z[k+1] * d3 / 6.0) * d2;
}
import std.stdio;
import std.math;
// 驥榊鴨蜉騾溷ｺｦ
const double g = -9.8;
// 遨ｺ豌玲慣謚嶺ｿよ焚
const double k = -0.01;
// 譎る俣髢馴囈(遘・
const double h = 0.01;
void main(string[] args)
{
// 隗貞ｺｦ
double degree = 45;
double radian = degree * PI / 180.0;
// 蛻晞・250 km/h -> 遘帝溘↓螟画鋤
double v = 250 * 1000 / 3600;
// 豌ｴ蟷ｳ譁ｹ蜷代・騾溷ｺｦ
double vx[2];
vx[0] = v * cos(radian);
// 驩帷峩譁ｹ蜷代・騾溷ｺｦ
double vy[2];
vy[0] = v * sin(radian);
// 邨碁℃遘呈焚
double t = 0.0;
// 菴咲ｽｮ
double x = 0.0;
double y = 0.0;
// Euler豕・
for (int i = 1; y >= 0.0; i++)
{
// 邨碁℃遘呈焚
t = i * h;
// 菴咲ｽｮ
x += h * vx[0];
y += h * vy[0];
writefln("%4.2f\t%8.5f\t%9.5f\t%9.5f\t%8.5f", t, vx[0], vy[0], x, y);
// 騾溷ｺｦ
vx[1] = vx[0] + h * fx(vx[0], vy[0]);
vy[1] = vy[0] + h * fy(vx[0], vy[0]);
vx[0] = vx[1];
vy[0] = vy[1];
}
}
// 遨ｺ豌玲慣謚励↓繧医ｋ豌ｴ蟷ｳ譁ｹ蜷代・貂幃溷・
double fx(double vx, double vy)
{
return k * sqrt(vx * vx + vy * vy) * vx;
}
// 驥榊鴨縺ｨ遨ｺ豌玲慣謚励↓繧医ｋ驩帷峩譁ｹ蜷代・貂幃溷・
double fy(double vx, double vy)
{
return g + (k * sqrt(vx * vx + vy * vy) * vy);
}
import std.stdio;
import std.math;
// 驥榊鴨蜉騾溷ｺｦ
const double g = -9.8;
// 遨ｺ豌玲慣謚嶺ｿよ焚
const double k = -0.01;
// 譎る俣髢馴囈(遘・
const double h = 0.01;
void main(string[] args)
{
// 隗貞ｺｦ
double degree = 45;
double radian = degree * PI / 180.0;
// 蛻晞・250 km/h -> 遘帝溘↓螟画鋤
double v = 250 * 1000 / 3600;
// 豌ｴ蟷ｳ譁ｹ蜷代・騾溷ｺｦ
double vx[3];
vx[0] = v * cos(radian);
// 驩帷峩譁ｹ蜷代・騾溷ｺｦ
double vy[3];
vy[0] = v * sin(radian);
// 邨碁℃遘呈焚
double t = 0.0;
// 菴咲ｽｮ
double x[3];
x[0] = 0.0;
double y[3];
y[0] = 0.0;
// Heun豕・
for (int i = 1; y[0] >= 0.0; i++)
{
// 邨碁℃遘呈焚
t = i * h;
// 菴咲ｽｮ繝ｻ騾溷ｺｦ
x[1] = x[0] + h * vx[0];
y[1] = y[0] + h * vy[0];
vx[1] = vx[0] + h * fx(vx[0], vy[0]);
vy[1] = vy[0] + h * fy(vx[0], vy[0]);
x[2] = x[0] + h * ( vx[0] + vx[1] ) / 2;
y[2] = y[0] + h * ( vy[0] + vy[1] ) / 2;
vx[2] = vx[0] + h * (fx(vx[0], vy[0]) + fx(vx[1], vy[1])) / 2;
vy[2] = vy[0] + h * (fy(vx[0], vy[0]) + fy(vx[1], vy[1])) / 2;
x[0] = x[2];
y[0] = y[2];
vx[0] = vx[2];
vy[0] = vy[2];
writefln("%4.2f\t%8.5f\t%9.5f\t%9.5f\t%9.5f", t, vx[0], vy[0], x[0], y[0]);
}
}
// 遨ｺ豌玲慣謚励↓繧医ｋ豌ｴ蟷ｳ譁ｹ蜷代・貂幃溷・
double fx(double vx, double vy)
{
return k * sqrt(vx * vx + vy * vy) * vx;
}
// 驥榊鴨縺ｨ遨ｺ豌玲慣謚励↓繧医ｋ驩帷峩譁ｹ蜷代・貂幃溷・
double fy(double vx, double vy)
{
return g + (k * sqrt(vx * vx + vy * vy) * vy);
}
import std.stdio;
import std.math;
// 驥榊鴨蜉騾溷ｺｦ
const double g = -9.8;
// 遨ｺ豌玲慣謚嶺ｿよ焚
const double k = -0.01;
// 譎る俣髢馴囈(遘・
const double h = 0.01;
void main(string[] args)
{
// 隗貞ｺｦ
double degree = 45;
double radian = degree * PI / 180.0;
// 蛻晞・250 km/h -> 遘帝溘↓螟画鋤
double v = 250 * 1000 / 3600;
// 豌ｴ蟷ｳ譁ｹ蜷代・騾溷ｺｦ
double vx[2];
vx[0] = v * cos(radian);
// 驩帷峩譁ｹ蜷代・騾溷ｺｦ
double vy[2];
vy[0] = v * sin(radian);
// 邨碁℃遘呈焚
double t = 0.0;
// 菴咲ｽｮ
double x[2];
x[0] = 0.0;
double y[2];
y[0] = 0.0;
// 荳ｭ轤ｹ豕・
for (int i = 1; y[0] >= 0.0; i++)
{
// 邨碁℃遘呈焚
t = i * h;
// 菴咲ｽｮ繝ｻ騾溷ｺｦ
vx[1] = h * fx(vx[0], vy[0]);
vy[1] = h * fy(vx[0], vy[0]);
double wx = vx[0] + vx[1] / 2;
double wy = vy[0] + vy[1] / 2;
vx[0] = vx[0] + h * fx(wx, wy);
vy[0] = vy[0] + h * fy(wx, wy);
x[0] = x[0] + h * wx;
y[0] = y[0] + h * wy;
writefln("%4.2f\t%8.5f\t%9.5f\t%9.5f\t%9.5f", t, vx[0], vy[0], x[0], y[0]);
}
}
// 遨ｺ豌玲慣謚励↓繧医ｋ豌ｴ蟷ｳ譁ｹ蜷代・貂幃溷・
double fx(double vx, double vy)
{
return k * sqrt(vx * vx + vy * vy) * vx;
}
// 驥榊鴨縺ｨ遨ｺ豌玲慣謚励↓繧医ｋ驩帷峩譁ｹ蜷代・貂幃溷・
double fy(double vx, double vy)
{
return g + (k * sqrt(vx * vx + vy * vy) * vy);
}
import std.stdio;
import std.math;
// 驥榊鴨蜉騾溷ｺｦ
const double g = -9.8;
// 遨ｺ豌玲慣謚嶺ｿよ焚
const double k = -0.01;
// 譎る俣髢馴囈(遘・
const double h = 0.01;
void main(string[] args)
{
// 隗貞ｺｦ
double degree = 45;
double radian = degree * PI / 180.0;
// 蛻晞・250 km/h -> 遘帝溘↓螟画鋤
double v = 250 * 1000 / 3600;
// 豌ｴ蟷ｳ譁ｹ蜷代・騾溷ｺｦ
double vx[5];
vx[0] = v * cos(radian);
// 驩帷峩譁ｹ蜷代・騾溷ｺｦ
double vy[5];
vy[0] = v * sin(radian);
// 邨碁℃遘呈焚
double t = 0.0;
// 菴咲ｽｮ
double x[5];
x[0] = 0.0;
double y[5];
y[0] = 0.0;
// Runge-Kutta豕・
for (int i = 1; y[0] >= 0.0; i++)
{
// 邨碁℃遘呈焚
t = i * h;
// 菴咲ｽｮ繝ｻ騾溷ｺｦ
x[1] = h * vx[0];
y[1] = h * vy[0];
vx[1] = h * fx(vx[0], vy[0]);
vy[1] = h * fy(vx[0], vy[0]);
double wx = vx[0] + vx[1] / 2;
double wy = vy[0] + vy[1] / 2;
x[2] = h * wx;
y[2] = h * wy;
vx[2] = h * fx(wx, wy);
vy[2] = h * fy(wx, wy);
wx = vx[0] + vx[2] / 2;
wy = vy[0] + vy[2] / 2;
x[3] = h * wx;
y[3] = h * wy;
vx[3] = h * fx(wx, wy);
vy[3] = h * fy(wx, wy);
wx = vx[0] + vx[3];
wy = vy[0] + vy[3];
x[4] = h * wx;
y[4] = h * wy;
vx[4] = h * fx(wx, wy);
vy[4] = h * fy(wx, wy);
x[0] += ( x[1] + x[2] * 2 + x[3] * 2 + x[4]) / 6;
y[0] += ( y[1] + y[2] * 2 + y[3] * 2 + y[4]) / 6;
vx[0] += (vx[1] + vx[2] * 2 + vx[3] * 2 + vx[4]) / 6;
vy[0] += (vy[1] + vy[2] * 2 + vy[3] * 2 + vy[4]) / 6;
writefln("%4.2f\t%8.5f\t%9.5f\t%9.5f\t%9.5f", t, vx[0], vy[0], x[0], y[0]);
}
}
// 遨ｺ豌玲慣謚励↓繧医ｋ豌ｴ蟷ｳ譁ｹ蜷代・貂幃溷・
double fx(double vx, double vy)
{
return k * sqrt(vx * vx + vy * vy) * vx;
}
// 驥榊鴨縺ｨ遨ｺ豌玲慣謚励↓繧医ｋ驩帷峩譁ｹ蜷代・貂幃溷・
double fy(double vx, double vy)
{
return g + (k * sqrt(vx * vx + vy * vy) * vy);
}
import std.stdio;
import std.math;
// 驥榊鴨蜉騾溷ｺｦ
const double g = -9.8;
// 遨ｺ豌玲慣謚嶺ｿよ焚
const double k = -0.01;
// 譎る俣髢馴囈(遘・
const double h = 0.01;
void main(string[] args)
{
// 隗貞ｺｦ
double degree = 45;
double radian = degree * PI / 180.0;
// 蛻晞・250 km/h -> 遘帝溘↓螟画鋤
double v = 250 * 1000 / 3600;
// 豌ｴ蟷ｳ譁ｹ蜷代・騾溷ｺｦ
double vx[5];
vx[0] = v * cos(radian);
// 驩帷峩譁ｹ蜷代・騾溷ｺｦ
double vy[5];
vy[0] = v * sin(radian);
// 邨碁℃遘呈焚
double t = 0.0;
// 菴咲ｽｮ
double x[5];
x[0] = 0.0;
double y[5];
y[0] = 0.0;
// Runge-Kutta-Gill豕・
for (int i = 1; y[0] >= 0.0; i++)
{
// 邨碁℃遘呈焚
t = i * h;
// 騾溷ｺｦ
x[1] = h * vx[0];
y[1] = h * vy[0];
vx[1] = h * fx(vx[0], vy[0]);
vy[1] = h * fy(vx[0], vy[0]);
double wx = vx[0] + vx[1] / 2;
double wy = vy[0] + vy[1] / 2;
x[2] = h * wx;
y[2] = h * wy;
vx[2] = h * fx(wx, wy);
vy[2] = h * fy(wx, wy);
wx = vx[0] + vx[1] * ((sqrt(2.0) - 1) / 2) + vx[2] * (1 - 1 / sqrt(2.0));
wy = vy[0] + vy[1] * ((sqrt(2.0) - 1) / 2) + vy[2] * (1 - 1 / sqrt(2.0));
x[3] = h * wx;
y[3] = h * wy;
vx[3] = h * fx(wx, wy);
vy[3] = h * fy(wx, wy);
wx = vx[0] - vx[2] / sqrt(2.0) + vx[3] * (1 + 1 / sqrt(2.0));
wy = vy[0] - vy[2] / sqrt(2.0) + vy[3] * (1 + 1 / sqrt(2.0));
x[4] = h * wx;
y[4] = h * wy;
vx[4] = h * fx(wx, wy);
vy[4] = h * fy(wx, wy);
x[0] += ( x[1] + x[2] * (2 - sqrt(2.0)) + x[3] * (2 + sqrt(2.0)) + x[4]) / 6;
y[0] += ( y[1] + y[2] * (2 - sqrt(2.0)) + y[3] * (2 + sqrt(2.0)) + y[4]) / 6;
vx[0] += (vx[1] + vx[2] * (2 - sqrt(2.0)) + vx[3] * (2 + sqrt(2.0)) + vx[4]) / 6;
vy[0] += (vy[1] + vy[2] * (2 - sqrt(2.0)) + vy[3] * (2 + sqrt(2.0)) + vy[4]) / 6;
writefln("%4.2f\t%8.5f\t%9.5f\t%9.5f\t%9.5f", t, vx[0], vy[0], x[0], y[0]);
}
}
// 遨ｺ豌玲慣謚励↓繧医ｋ豌ｴ蟷ｳ譁ｹ蜷代・貂幃溷・
double fx(double vx, double vy)
{
return k * sqrt(vx * vx + vy * vy) * vx;
}
// 驥榊鴨縺ｨ遨ｺ豌玲慣謚励↓繧医ｋ驩帷峩譁ｹ蜷代・貂幃溷・
double fy(double vx, double vy)
{
return g + (k * sqrt(vx * vx + vy * vy) * vy);
}
import std.stdio;
import std.math;
void main(string[] args)
{
double a = 1;
double b = 2;
writefln("%12.10f", bisection(a, b));
}
double bisection(double a, double b)
{
double c;
while (true)
{
// 蛹ｺ髢・(a, b) 縺ｮ荳ｭ轤ｹ c = (a + b) / 2
c = (a + b) / 2;
writefln("%12.10f\t%13.10f", c, c - sqrt(2.0));
double fc = f(c);
if (fabs(fc) < 0.0000000001) break;
if (fc < 0)
{
// f(c) < 0 縺ｧ縺ゅｌ縺ｰ, 隗｣縺ｯ蛹ｺ髢・(c, b) 縺ｮ荳ｭ縺ｫ蟄伜惠
a = c;
}
else
{
// f(c) > 0 縺ｧ縺ゅｌ縺ｰ, 隗｣縺ｯ蛹ｺ髢・(a, c) 縺ｮ荳ｭ縺ｫ蟄伜惠
b = c;
}
}
return c;
}
double f(double x)
{
return x * x - 2;
}
import std.stdio;
import std.math;
void main(string[] args)
{
double a = 1;
double b = 2;
writefln("%12.10f", falseposition(a, b));
}
double falseposition(double a, double b)
{
double c;
while (true)
{
// 轤ｹ (a,f(a)) 縺ｨ 轤ｹ (b,f(b)) 繧堤ｵ舌・逶ｴ邱壹→ x霆ｸ縺ｮ莠､轤ｹ
c = (a * f(b) - b * f(a)) / (f(b) - f(a));
writefln("%12.10f\t%13.10f", c, c - sqrt(2.0));
double fc = f(c);
if (fabs(fc) < 0.0000000001) break;
if (fc < 0)
{
// f(c) < 0 縺ｧ縺ゅｌ縺ｰ, 隗｣縺ｯ蛹ｺ髢・(c, b) 縺ｮ荳ｭ縺ｫ蟄伜惠
a = c;
}
else
{
// f(c) > 0 縺ｧ縺ゅｌ縺ｰ, 隗｣縺ｯ蛹ｺ髢・(a, c) 縺ｮ荳ｭ縺ｫ蟄伜惠
b = c;
}
}
return c;
}
double f(double x)
{
return x * x - 2;
}
import std.stdio;
import std.math;
void main(string[] args)
{
double x = 1;
writefln("%12.10f", iterative(x));
}
double iterative(double x0)
{
double x1;
while (true)
{
x1 = g(x0);
writefln("%12.10f\t%13.10f", x1, x1 - sqrt(2.0));
if (fabs(x1 - x0) < 0.0000000001) break;
x0 = x1;
}
return x1;
}
double g(double x)
{
return (x / 2) + (1 / x);
}
import std.stdio;
import std.math;
void main(string[] args)
{
double x = 2;
writefln("%12.10f", newton(x));
}
double newton(double x0)
{
double x1;
while (true)
{
x1 = x0 - (f0(x0) / f1(x0));
writefln("%12.10f\t%13.10f", x1, x1 - sqrt(2.0));
if (fabs(x1 - x0) < 0.0000000001) break;
x0 = x1;
}
return x1;
}
double f0(double x)
{
return x * x - 2;
}
double f1(double x)
{
return 2 * x;
}
import std.stdio;
import std.math;
void main(string[] args)
{
double x = 2;
writefln("%12.10f", bailey(x));
}
double bailey(double x0)
{
double x1;
while (true)
{
x1 = x0 - (f0(x0) / (f1(x0) - (f0(x0) * f2(x0) / (2 * f1(x0)))));
writefln("%12.10f\t%13.10f", x1, x1 - sqrt(2.0));
if (fabs(x1 - x0) < 0.0000000001) break;
x0 = x1;
}
return x1;
}
double f0(double x)
{
return x * x - 2;
}
double f1(double x)
{
return 2 * x;
}
double f2(double x)
{
return 2;
}
import std.stdio;
import std.math;
void main(string[] args) {
double x0 = 1;
double x1 = 2;
writefln("%12.10f", secant(x0, x1));
}
double secant(double x0, double x1)
{
double x2;
while (true)
{
x2 = x1 - f(x1) * (x1 - x0) / (f(x1) - f(x0));
writefln("%12.10f\t%13.10f", x1, x1 - sqrt(2.0));
if (fabs(x2 - x1) < 0.0000000001) break;
x0 = x1;
x1 = x2;
}
return x2;
}
double f(double x)
{
return x * x - 2;
}
import std.stdio;
import std.math;
const int N = 4;
void main(string[] args)
{
double[N][N] a = [[9,2,1,1],[2,8,-2,1],[-1,-2,7,-2],[1,-1,-2,6]];
double[N] b = [20,16,8,17];
double[N] c = [0,0,0,0];
jacobi(a,b,c);
writefln("X");
for (int i = 0; i < N; i++)
writef("%14.10f\t", c[i]);
writefln("");
}
void jacobi(double[N][N] a, double[N] b, ref double[N] x0)
{
while (true)
{
double[N] x1;
bool finish = true;
for (int i = 0; i < N; i++)
{
x1[i] = b[i];
for (int j = 0; j < N; j++)
if (j != i)
x1[i] -= a[i][j] * x0[j];
x1[i] /= a[i][i];
if (fabs(x1[i] - x0[i]) > 0.0000000001) finish = false;
}
for (int i = 0; i < N; i++)
{
x0[i] = x1[i];
writef("%14.10f\t", x0[i]);
}
writefln("");
if (finish) return;
}
}
import std.stdio;
import std.math;
const int N = 4;
void main(string[] args)
{
double[N][N] a = [[9,2,1,1],[2,8,-2,1],[-1,-2,7,-2],[1,-1,-2,6]];
double[N] b = [20,16,8,17];
double[N] c = [0,0,0,0];
gauss(a,b,c);
writefln("X");
for (int i = 0; i < N; i++)
writef("%14.10f\t", c[i]);
writefln("");
}
void gauss(double[N][N] a, double[N] b, ref double[N] x0)
{
while (true) {
double x1;
bool finish = true;
for (int i = 0; i < N; i++)
{
x1 = 0;
for (int j = 0; j < N; j++)
if (j != i)
x1 += a[i][j] * x0[j];
x1 = (b[i] - x1) / a[i][i];
if (fabs(x1 - x0[i]) > 0.0000000001) finish = false;
x0[i] = x1;
}
for (int i = 0; i < N; i++)
writef("%14.10f\t", x0[i]);
writefln("");
if (finish) return;
}
}
import std.stdio;
import std.math;
const int N = 4;
void main(string[] args)
{
double[N][N] a = [[-1,-2,7,-2],[1,-1,-2,6],[9,2,1,1],[2,8,-2,1]];
double[N] b = [8,17,20,16];
pivoting(a,b);
writefln("pivoting");
disp_progress(a,b);
forward_elimination(a,b);
writefln("forward_elimination");
disp_progress(a,b);
backward_substitution(a,b);
writefln("X");
for (int i = 0; i < N; i++)
writef("%14.10f\t", b[i]);
writefln("");
}
void forward_elimination(ref double[N][N] a, ref double[N] b)
{
for (int pivot = 0; pivot < N - 1; pivot++)
{
for (int row = pivot + 1; row < N; row++)
{
double s = a[row][pivot] / a[pivot][pivot];
for (int col = pivot; col < N; col++)
a[row][col] -= a[pivot][col] * s;
b[row] -= b[pivot] * s;
}
}
}
void backward_substitution(ref double[N][N] a, ref double[N] b)
{
for (int row = N - 1; row >= 0; row--)
{
for (int col = N - 1; col > row; col--)
b[row] -= a[row][col] * b[col];
b[row] /= a[row][row];
}
}
void pivoting(ref double[N][N] a, ref double[N] b)
{
for(int pivot = 0; pivot < N; pivot++)
{
int max_row = pivot;
double max_val = 0;
for (int row = pivot; row < N; row++)
{
if (fabs(a[row][pivot]) > max_val)
{
max_val = fabs(a[row][pivot]);
max_row = row;
}
}
if (max_row != pivot)
{
double tmp;
for (int col = 0; col < N; col++) {
tmp = a[max_row][col];
a[max_row][col] = a[pivot][col];
a[pivot][col] = tmp;
}
tmp = b[max_row];
b[max_row] = b[pivot];
b[pivot] = tmp;
}
}
}
void disp_progress(double[N][N] a, double[N] b)
{
for (int i = 0; i < N; i++)
{
for (int j = 0; j < N; j++)
writef("%14.10f\t", a[i][j]);
writefln("%14.10f\t", b[i]);
}
writeln();
}
import std.stdio;
import std.math;
const int N = 4;
void main(string[] args)
{
double[N][N] a = [[-1,-2,7,-2],[1,-1,-2,6],[9,2,1,1],[2,8,-2,1]];
double[N] b = [8,17,20,16];
pivoting(a,b);
writefln("pivoting");
disp_progress(a,b);
forward_elimination(a,b);
writefln("forward_elimination");
disp_progress(a,b);
backward_substitution(a,b);
writefln("X");
for (int i = 0; i < N; i++)
writef("%14.10f\t", b[i]);
writeln();
}
void forward_elimination(ref double[N][N] a, ref double[N] b)
{
for (int pivot = 0; pivot < N; pivot++)
{
double s = a[pivot][pivot];
for (int col = 0; col < N; col++)
a[pivot][col] /= s;
b[pivot] /= s;
}
writeln("forward_elimination");
disp_progress(a,b);
for (int pivot = 0; pivot < N; pivot++)
{
for (int row = 0; row < N; row++)
{
if (row == pivot) continue;
double s = a[row][pivot] / a[pivot][pivot];
for (int col = pivot; col < N; col++)
a[row][col] -= a[pivot][col] * s;
b[row] -= b[pivot] * s;
}
}
}
void backward_substitution(double[N][N] a, ref double[N] b)
{
for (int pivot = 0; pivot < N; pivot++)
b[pivot] /= a[pivot][pivot];
}
void pivoting(ref double[N][N] a, ref double[N] b)
{
for(int pivot = 0; pivot < N; pivot++)
{
int max_row = pivot;
double max_val = 0;
for (int row = pivot; row < N; row++)
{
if (fabs(a[row][pivot]) > max_val)
{
max_val = fabs(a[row][pivot]);
max_row = row;
}
}
if (max_row != pivot)
{
double tmp;
for (int col = 0; col < N; col++) {
tmp = a[max_row][col];
a[max_row][col] = a[pivot][col];
a[pivot][col] = tmp;
}
tmp = b[max_row];
b[max_row] = b[pivot];
b[pivot] = tmp;
}
}
}
void disp_progress(double[N][N] a, double[N] b)
{
for (int i = 0; i < N; i++)
{
for (int j = 0; j < N; j++)
writef("%14.10f\t", a[i][j]);
writefln("%14.10f\t", b[i]);
}
writeln();
}
import std.stdio;
import std.math;
const int N = 4;
void main(string[] args)
{
double[N][N] a = [[-1,-2,7,-2],[1,-1,-2,6],[9,2,1,1],[2,8,-2,1]];
double[N] b = [8,17,20,16];
pivoting(a,b);
writefln("pivoting");
disp_progress(a);
// LU
double[N] x = [0,0,0,0];
decomp(a,b,x);
writeln("X");
for (int i = 0; i < N; i++)
writef("%14.10f\t", x[i]);
writeln();
}
// LU
void decomp(ref double[N][N] a, ref double[N] b, ref double[N] x)
{
//
forward_elimination(a,b);
//
writeln("L");
for (int i = 0; i < N; i++)
{
for (int j = 0; j < N; j++)
if (i > j)
writef("%14.10f\t", a[i][j]);
else if (i == j)
writef("%14.10f\t", 1.0);
else
writef("%14.10f\t", 0.0);
writeln();
}
writeln();
//
writeln("U");
for (int i = 0; i < N; i++)
{
for (int j = 0; j < N; j++)
if (i <= j)
writef("%14.10f\t", a[i][j]);
else
writef("%14.10f\t", 0.0);
writeln();
}
writeln();
// Ly=b
double y[N] = [0,0,0,0];
for (int row = 0; row < N; row++)
{
for (int col = 0; col < row; col++)
b[row] -= a[row][col] * y[col];
y[row] = b[row];
}
// y
writeln("Y");
for (int i = 0; i < N; i++)
writef("%14.10f\t", y[i]);
writeln("\n");
// Ux=y
for (int row = N - 1; row >= 0; row--)
{
for (int col = N - 1; col > row; col--)
y[row] -= a[row][col] * x[col];
x[row] = y[row] / a[row][row];
}
}
//
void forward_elimination(ref double[N][N] a, double[N] b)
{
for (int pivot = 0; pivot < N - 1; pivot++)
{
for (int row = pivot + 1; row < N; row++)
{
double s = a[row][pivot] / a[pivot][pivot];
for (int col = pivot; col < N; col++)
a[row][col] -= a[pivot][col] * s; //
a[row][pivot] = s; //
// b[row] -= b[pivot] * s; //
}
}
}
void pivoting(ref double[N][N] a, ref double[N] b)
{
for(int pivot = 0; pivot < N; pivot++)
{
int max_row = pivot;
double max_val = 0;
for (int row = pivot; row < N; row++)
{
if (fabs(a[row][pivot]) > max_val)
{
max_val = fabs(a[row][pivot]);
max_row = row;
}
}
if (max_row != pivot)
{
double tmp;
for (int col = 0; col < N; col++) {
tmp = a[max_row][col];
a[max_row][col] = a[pivot][col];
a[pivot][col] = tmp;
}
tmp = b[max_row];
b[max_row] = b[pivot];
b[pivot] = tmp;
}
}
}
void disp_progress(double[N][N] a)
{
for (int i = 0; i < N; i++)
{
for (int j = 0; j < N; j++)
writef("%14.10f\t", a[i][j]);
writeln();
}
}
import std.stdio;
import std.math;
const int N = 4;
void main(string[] args)
{
double[N][N] a = [[5,2,3,4],[2,10,6,7],[3,6,15,9],[4,7,9,20]];
double[N] b = [34,68,96,125];
writefln("A");
disp_progress(a);
// LU
double[N] x = [0,0,0,0];
decomp(a,b,x);
writeln("X");
for (int i = 0; i < N; i++)
writef("%14.10f\t", x[i]);
writeln();
}
// LU
void decomp(ref double[N][N] a, ref double[N] b, ref double[N] x)
{
//
forward_elimination(a,b);
//
writeln("L");
for (int i = 0; i < N; i++)
{
for (int j = 0; j < N; j++)
if (j <= i)
writef("%14.10f\t", a[i][j]);
else
writef("%14.10f\t", a[j][i]);
writeln();
}
writeln();
//
// Ly=b
double y[N] = [0,0,0,0];
for (int row = 0; row < N; row++)
{
for (int col = 0; col < row; col++)
b[row] -= a[row][col] * y[col];
y[row] = b[row] / a[row][row];
}
// y
writeln("Y");
for (int i = 0; i < N; i++)
writef("%14.10f\t", y[i]);
writeln("\n");
// Ux=y
for (int row = N - 1; row >= 0; row--)
{
for (int col = N - 1; col > row; col--)
y[row] -= a[col][row] * x[col];
x[row] = y[row] / a[row][row];
}
}
//
void forward_elimination(ref double[N][N] a, double[N] b)
{
for (int pivot = 0; pivot < N; pivot++)
{
double s = 0;
for (int col = 0; col < pivot; col++)
s += a[pivot][col] * a[pivot][col];
a[pivot][pivot] = sqrt(a[pivot][pivot] - s);
for (int row = pivot + 1; row < N; row++)
{
s = a[row][pivot] / a[pivot][pivot];
s = 0;
for (int col = 0; col < pivot; col++)
s += a[row][col] * a[pivot][col]; //
a[row][pivot] = (a[row][pivot] - s) / a[pivot][pivot];
}
}
}
void disp_progress(double[N][N] a)
{
for (int i = 0; i < N; i++)
{
for (int j = 0; j < N; j++)
writef("%14.10f\t", a[i][j]);
writeln();
}
}
import std.stdio;
import std.math;
const int N = 4;
void main(string[] args)
{
double[N][N] a = [[5,2,3,4],[2,10,6,7],[3,6,15,9],[4,7,9,20]];
double[N] b = [34,68,96,125];
writefln("A");
disp_progress(a);
// LU
double[N] x = [0,0,0,0];
decomp(a,b,x);
writeln("X");
for (int i = 0; i < N; i++)
writef("%14.10f\t", x[i]);
writeln();
}
// LU
void decomp(ref double[N][N] a, ref double[N] b, ref double[N] x)
{
//
forward_elimination(a,b);
//
writeln("L");
for (int i = 0; i < N; i++)
{
for (int j = 0; j < N; j++)
if (j <= i)
writef("%14.10f\t", a[i][j]);
else
writef("%14.10f\t", a[j][i]);
writeln();
}
writeln();
//
// Ly=b
double y[N] = [0,0,0,0];
for (int row = 0; row < N; row++)
{
for (int col = 0; col < row; col++)
b[row] -= a[row][col] * y[col];
y[row] = b[row];
}
// y
writeln("Y");
for (int i = 0; i < N; i++)
writef("%14.10f\t", y[i]);
writeln("\n");
// Ux=y
for (int row = N - 1; row >= 0; row--)
{
for (int col = N - 1; col > row; col--)
y[row] -= a[col][row] * a[row][row] * x[col];
x[row] = y[row] / a[row][row];
}
}
//
void forward_elimination(ref double[N][N] a, double[N] b)
{
for (int pivot = 0; pivot < N; pivot++)
{
double s;
// pivot < k
for (int col = 0; col < pivot; col++)
{
s = a[pivot][col];
for (int k = 0; k < col; k++)
s -= a[pivot][k] * a[col][k] * a[k][k];
a[pivot][col] = s / a[col][col];
}
// pivot == k
s = a[pivot][pivot];
for (int k = 0; k < pivot; k++)
s -= a[pivot][k] * a[pivot][k] * a[k][k];
a[pivot][pivot] = s;
}
}
void disp_progress(double[N][N] a)
{
for (int i = 0; i < N; i++)
{
for (int j = 0; j < N; j++)
writef("%14.10f\t", a[i][j]);
writeln();
}
}
|