2023 年度「計算数学a」 2023-04-21

§3.1 `double` 型

(3.1.1) C 言語には小数点以下のある数を扱える型がいくつかある。この授業では double 型のみを紹介する。

(3.1.2)

#include <stdio.h>

int main() {
  double x;     /* double 型変数 x の宣言 */

  scanf("%lf", &x);
  printf("%f\n", x);
}

double 型の入出力の際には上のように書く。 %lf（パーセント、エル、エフ）と %f との違いに注意せよ。出力は小数点以下 6 桁となる。 20 桁ぐらいの大きな数や、 0.00000001 のような小さな数を入力してみよ。

  printf("%15.8f\n", x);

とすると「全体で少なくとも 15 桁の幅、小数点以下 8 桁」で表示される。

  printf("%15f\n", x);

なら「全体で少なくとも 15 桁の幅」。

  printf("%.8f\n", x);

なら「小数点以下 8 桁」。

(3.1.3) double 型には、十進で言えば 1.234567×10⁸⁹ のような形で値が格納されているので、桁数が大きい整数も代入できる。しかし、近似値になってしまうこともある。

(3.1.4) 自然対数の底 e の値を、e = 1 + 1/1! + 1/2! + 1/3! + ... として計算する。（変数 factorial には n の階乗がはいる。）

#include <stdio.h>

int main() {
  int n;
  double factorial, sum;

  factorial = 1;
  sum = 1;
  for (n = 1; n <= 25; n++) {
    factorial = factorial * n;
    sum = sum + 1/factorial;
    printf("%2d %.40f\n", n, sum);
  }
  printf("   2.718281828459045235360287471352... ← 数表\による値\n");
}

注：「数表\による値」と、表 の後ろに \ がはいっている理由は、 (1.11.12) を見よ。

§3.2 数学関数

(3.2.1) 妙な名前だと思うかもしれないが、 C コンパイラに備え付けの関数で数学に関連するものを「数学関数」と呼ぶ。数学関数の一部をあげる。 K&R2 では §B4 を参照。数学関数を使うには #include <math.h> が必要である。つまり、プログラムの冒頭は次のようになる。

#include <stdio.h>
#include <math.h>

int main() {

(3.2.2) x, y は double 型、答えもすべて double 型である。角度の単位はラジアンである。

sin(x), cos(x), tan(x) ... 三角関数
asin(x), acos(x), atan(x) ... 逆三角関数
atan2(x, y) ... tan^-1(y/x)
sinh(x), cosh(x), tanh(x) ... 双曲線関数
exp(x) ... exp(x)（e^x のこと）
log(x) ... log_e(x)
log10(x) ... log₁₀(x)
sqrt(x) ... ルート x

(3.2.3) 三角関数表の一部を出力するプログラム。 PI は円周率の近似値である。 πの呼び名にちなむ。なお、これを「パイ」と読むのは英語式である。

#include <stdio.h>
#include <math.h>

#define PI 3.141592653589793238462643383279

int main() {
  double x;

  for (x = 0; x <= 90; x = x + 5) {
    printf("%9f %6f\n", x, sin(x / 180 * PI));
  }
}

§3.3 コマンドライン引数の使い方（飛ばしてもよい）

(3.3.1) 「二次方程式を解くプログラムを書け」と言われたら、次の二つのタイプが考えられる。

(3.3.2) 第一は、プログラムの中に数値を書き込んでしまう方法。第二は、プログラムが動いてから、ユーザに数値を入力させる方法。

(3.3.3) 第三の方式は、プロンプトの出ているところに「a 3 4.56 7.8」または「a.exe 3 4.56 7.8」のように打ち込んで動かし、プログラムの側でそれらの値を利用する方法である。この「3 4.56 7.8」をコマンドライン引数という。（初回のプリントの「付」でも述べた方法である。）

#include <stdio.h>
#include <stdlib.h> /* atoi(), atof() を使うのに必要 */

int main(int argc, char *argv[ ]) {     /* 小かっこの中は、常にこう書く。きょうは説明しない */
  int x;
  double y, z;

  if (argc != 4) {                      /* argc は、コマンドラインに打ち込んだものの個数 */
    printf("使い方が違います.\n");
  } else {
    x = atoi(argv[1]);                  /* int 型に収めるときは atoi() */
    y = atof(argv[2]);                  /* double 型に収めるときは atof() */
    z = atof(argv[3]);
    printf("x = %d, y = %f, z = %f\n", x, y, z);
                                        /* きちんと収まったか、確かめる（だけ）*/
  }
}

§3.4 練習問題（必須）

(3.4.1) 実数を係数とする二次方程式を解くプログラムを書け。解が虚数になる場合は、「実数解を持ちません」のように出力するだけでもよい。上で述べたうちのどのスタイルで書いてもよい。よくあるミスを見抜くために、2x² - 3x + 1 = 0 も解かせてみること。変数はすべて double 型にするのがよかろう。

§3.5 練習問題（飛ばしてもよい）

(3.5.1) 未知数が二つの、連立一次方程式を解くプログラムを書け。

(3.5.2) フィボナッチ数列を、double 型を用いて計算し、隣り合った項の間の比を表示させてみよ。

§3.6 `double` 型についての補足

(3.6.1) 階乗を計算するプログラム。あるところから先は“無限大”になる。（コマンドプロンプトを最大化して試すとよい。）

#include <stdio.h>

int main() {
  int n;
  double factorial;

  factorial = 1;
  for (n = 1; n <= 200; n++) {
    factorial = factorial * n;
    printf("%2d %f\n", n, factorial);
  }
}

(3.6.2) int 型と double 型を混ぜて使う際の注意。 int 型を double 型に代入すると、（たぶん）そのまま格納される。 double 型を int 型に代入すると、（確か）小数点以下が切り捨てられる。

やや間違えやすいケースをプログラムで示す。

一番めの x = 1/3; では、 int 型の 1 を int 型の 3 で割ってから double 型変数の x に代入するので、 int 型の割り算の規則により、切り捨てた値 0 が x に代入される。
二番めでは、1.0/3 と書くことで、double 型の 1 を int 型の 3 で割ることになり、 x には 0.333... が代入される。
三番めでは、double 型に代入された 1 を int 型の 3 で割るので、x は 0.333... となる。

#include <stdio.h>

int main() {
  double x;

  x = 1/3;
  printf("%f\n", x);

  x = 1.0/3;
  printf("%f\n", x);

  x = 1;
  x = x/3;
  printf("%f\n", x);
}

これらは、C 言語のやや使いにくい点である。（ほかにもまだまだいろいろな型があるので、複雑な規則がある。）

(3.6.3) double 型は二進で数値を保存しており、十進に変換して出力される。十進の 0.1 は二進では 1/1010 なので循環小数となり、十進の 1.333 のような近似値として扱われる。（以下のプログラムでそのことがわかる。）

(3.6.4) （注意）二進の有限小数は十進でも有限小数である。

(3.6.5) 次のプログラムは、0.9 までが出力されるように思えるが、 gcc ではそうはならない。

#include <stdio.h>

int main() {
  double x;

  for (x = 0; x < 1; x = x + 0.1) {
    printf("%f\n", x);
  }
}

1 も出力される。それは、0.1 の十倍が、誤差のため、1 より少し小さい値になるからである。より多くの桁数を出力させ、確かめてみよ。

(3.6.6) 以下のプログラムを動かしてみよ。 double 型の誤差がどのように出るか、わかるであろう。（100 とか 150 といった定数に特別の意味はない。適宜、変えてみてもよい。）

#include <stdio.h>

int main() {
  int n;
  double x;

  x = 1;
  for (n = 0; n < 100; n++) {
    printf("%f\n", x);
    x = x * 10;                         /* 次々と 10 倍する */
  }
}

(3.6.7)

#include <stdio.h>

int main() {
  int n;
  double x;

  x = 1;
  for (n = 0; n < 150; n++) {
    printf("%.150f\n", x);
    x = x / 10;                         /* 次々と 10 で割る */
  }
}

(3.6.8)

#include <stdio.h>

int main() {
  int n;
  double x;

  x = 1;
  for (n = 0; n < 450; n++) {
    printf("%150f\n", x);
    x = x * 2;                          /* 次々と 2 倍する */
  }
}

(3.6.9)

#include <stdio.h>

int main() {
  int n;
  double x;

  x = 1;
  for (n = 0; n < 450; n++) {
    printf("%.150f\n", x);
    x = x / 2;                          /* 次々と 2 で割る */
  }
}

(3.6.10)

  printf("%g\n", x);

とすると、1.2345×10⁶⁷ のように表示される。実際の出力は 1.2345e+067 となる。

#include <stdio.h>

int main() {
  int n;
  double x;

  x = 1;
  for (n = 0; n < 400; n++) {
    printf("%.50g\n", x);
    x = x / 10;                         /* 次々と 10 で割る */
  }
}

どこかから先は 0 になる。

付：無謀なプログラムの例

1/1 + 1/2 + 1/3 + 1/4 + 1/5 + ... が ∞ に発散することを見ようとしている。

#include <stdio.h>

int main() {
  int n;
  double sum;

  sum = 0;
  for (n = 1;     ; n++) {      /* 永遠に繰り返せ、の意味になる */
    printf("%d %.25lf\n", n, sum);
    sum = sum + 1.0 / n;
  }
}

次は、上のプログラムに細工をして、 n が int 型で表せる最大数の付近のところだけを出力させたもの。

2147483638 22.0647782573692480000000000
2147483639 22.0647782578349090000000000
2147483640 22.0647782583005710000000000
2147483641 22.0647782587662320000000000
2147483642 22.0647782592318930000000000
2147483643 22.0647782596975550000000000
2147483644 22.0647782601632160000000000
2147483645 22.0647782606288770000000000
2147483646 22.0647782610945380000000000
2147483647 22.0647782615602000000000000

岩瀬順一