whileの条件文が真になってもループを抜けない

Question

前提・実現したいこと

google cloud shell ,cを用いて、0.5 + 25x - 200xx + 675xxx - 900xxxx + 400xxxx*x　という関数を0～0．8までシンプソン1／3則を用いて、厳密な積分値1．880533になるまで分割数ｎを変化させるというプログラムを作っています。

発生している問題・エラーメッセージ

whileを使っているのですが、厳密な積分値になってもループが止まりません。

該当のソースコード

ｃ

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

int main(){
double a, b, h, x, f, f0, fl;
double A, B, I;
int i, n;

a = 0, b = 0.8, n = 0;

while(I != 1.880533){
    h = (b - a) / (double)n;

    A = 0, B = 0;

    for(i = 2; i <= n - 2; i += 2){/*偶数*/
        x = a + h * (double)i;
        f = 0.5 + 25 * x - 200 * pow(x,2) + 675 * pow(x,3) - 900 * pow(x,4) + 400 * pow(x,5);

        A += f;    
    }

    for(i = 1; i <= n - 1; i += 2){/*奇数*/
        x = a + h * (double)i;
        f = 0.5 + 25 * x - 200 * pow(x,2) + 675 * pow(x,3) - 900 * pow(x,4) + 400 * pow(x,5);

        B += f;    
    }

    x = 0;
    f0 = 0.5 + 25 * x - 200 * pow(x,2) + 675 * pow(x,3) - 900 * pow(x,4) + 400 * pow(x,5);

    x = a + h * (double)n;
    fl = 0.5 + 25 * x - 200 * pow(x,2) + 675 * pow(x,3) - 900 * pow(x,4) + 400 * pow(x,5);

    I = h * (f0 + 4 * B + 2 * A + fl) / 3;

    printf("n=%d I=%lf\n", n, I);    

    n += 2;

}

return 0;

}

Comments

[jimbe]

> whileの条件文が真になってもループを抜けない
while は「条件文が真の間ループ」しますので、タイトルは逆では無いでしょうか。

Answer 1

仕組み上、小数を厳密に比較することは難しいです。
目標値との差分を取って判定するのが一般的かと思います。

while(fabs(I - 1.880533) <= 0.0001) {

0.0001は適当に決めた精度です。

参考： 浮動小数点数の誤差を考慮した比較【double/float型の正しい比較方法】 | MaryCore

Answer 2

while(I != 1.880533)ですが、浮動小数点数は厳密な比較には適しません（ごく近いけど別な値に収束していることが考えられます）。

Answer 3

コードの投稿の仕方が間違っています。teratail のヘルプを見て質問を書き直してください。

変数 I が初期化されないまま、1.88053 と比較されています。
さらに　n の初期値が 0 なので (b - a) / (double)n がとんでもないことになります。
n は 2 から始めないといけないでしょう。

I が 1.880533 になってほしいということは、
I が 1.8805325以上 1.8805335未満になればいいということですから、
while(I != 1.880533){ を
while (fabs(I - 1.880533) > 0.0000005) { に変えてみてください。

同じ計算が何度も出てくるのは見にくいものです。
関数を使って見やすくなるようにしましょう。

#include <stdio.h>  // printf
#include <math.h>   // fabs

double f(double x)
{
	return ((((400*x - 900)*x + 675)*x - 200)*x + 25)*x + 0.5;
}

double simpson(double a, double b, int n)
{
	double h = (b - a) / n;
	double o = 0, e = 0;
	for (int i = 1; i < n; i += 2) o += f(a + h*i);
	for (int i = 2; i < n; i += 2) e += f(a + h*i);
	return (f(a) + 4*o + 2*e + f(b)) * h / 3;
}

int main(void)
{
	int n = 2;
	double v = 0;
	while (fabs(v - 1.880533) >= 0.0000005) {
		v = simpson(0, 0.8, n);
		printf("%4d %f\n", n, v);
		n += 2;
	}
}