シンプソン則について

質問

大学の課題で　∫0～5　e^(-3t) * cos^2 *(2t) をシンプソン則を用いて求めるプログラムを作成せよ。積分の刻みは0.01毎　　という課題がでました。
私は授業で作ったプログラムを元に、以下のプログラムを作ったのですが、答えが大きすぎて（5より大きいと間違いと言っていました）受け取ってもらえません。どこがいけないのでしょうか？
あと、ｎには　どのような数字をいれればよいのでしょうか？ｎは使わないのでしょうか。
お願いします。

#include "stdio.h"
#include "math.h"

void main()
{
int n;
double i, a, b, h, S, S1, S2, S3;
double f (double);
a = 0.;
b = 5.;

printf ("N=");
scanf ("%d", &n);
h = (b - a) / (double)n;
S1 = 0.;
S2 = 0.;
S3 = f (a) + f (b); /* 両端の値*/
/* 奇数番目の値*/
for (i = 0.01; i < 5; i = i + 0.02)
{
S1 = S1 + 4. * (f (a + (double)i * h));
}
/* 偶数番目の値*/
for (i = 0.0; i < 5; i = i + 0.02)
{
S2 = S2 + 2. * (f (a + (double)i * h));
}
S = (h / 3.) * (S1 + S2 + S3);
printf ("%d %lf\n", n, S);
}
double f (double t)
{
double y;

y = (exp(-3. * t) / 2.) + ((exp(-3.* t) * cos(4. * t)) / 2.);
return y;
}

回答No.10

aquarius_hiro

A No.8 です。

アルゴリズムはそれで合っていると思いますよ。
（細かいところは見ていないですが…。）

ただ i はもともと double で宣言しているので、(double)i としなくて良いと思います。そのまま a + i で良いです。

回答No.9

asuncion

こうだと思います。

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

double f(double x)
{
　　return exp(-3 * x) * cos(2 * x) * cos(2 * x);
}

int main(void)
{
　　double a, b, h, s;
　　int n, i;

　　a = 0, b = 5;
　　h = 0.01;
　　n = (b - a) / h + 0.5;　　　　// +0.5 は丸め誤差への対応
　　for (s = 0, i = 1; i < n; i++)
　　　　s += ((i % 2 == 0) ? 2 : 4) * f(a + i * h);
　　s += f(a) + f(b);
　　s *= h / 3;
　　printf("面積=%f\n", s);
　　return 0;
}

（注）インデントのため、全角空白を使っています。

回答No.8

aquarius_hiro

だんだんまどろっこしくなってきたので…
もう言ってしまいますね。。(^^;)

i は 0.02 ずつ増えています。

a + (double)i*h は、例えば偶数のほうでは、

a, a + 0.02*h, a + 0.04*h, .....

と増えていきますね。

しかし、ほんとうは、a, a + 0.02, a + 0.04, ... でよいのですから、*h が余分ですよね。

a + i だけで良いです。

奇数のほうでも同様です。

もし a + (double)i * h の形にしたければ、

h = 0.01 刻みで増えていくようにするには、i は 2 刻みで増やさないといけませんね。（i=1,3,5,... と i=0,2,4,... ）

お礼

なかなか理解出来ずにすみません。

このプログラムでよいのでしょうか？
本当にありがとうございました！
#include "stdio.h"
#include "math.h"

void main()
{
int n;
double i, a, b, h, S, S1, S2, S3;
double f (double);
a = 0.;
b = 5.;

printf ("N=");
scanf ("%d", &n);
h = (b - a) / (double)n;
S1 = 0.;
S2 = 0.;
S3 = f (a) + f (b); /* 両端の値*/
/* 奇数番目の値*/
for (i = 0.01; i < 5; i = i + 0.02)
{
S1 = S1 + 4. * (f (a + (double)i));
}
/* 偶数番目の値*/
for (i = 0.0; i < 5; i = i + 0.02)
{
S2 = S2 + 2. * (f (a + (double)i));
}
S = (h / 3.) * (S1 + S2 + S3);
printf ("%d %lf\n", n, S);
}
double f (double t)
{
double y;

y = (exp(-3. * t) / 2.) + ((exp(-3.* t) * cos(4. * t)) / 2.);
return y;
}

回答No.7

aquarius_hiro

そうですね、i は 0.02 ずつ増えるわけですね。

そのとき、

f() の引数

a + (double)i * h

は、a = 0 （の近く）から、どんなふうに増えていくと思いますか？

具体的に数字を順に書いてみたらどうでしょうね？

b = 5 （の近く）まで変化していかないといけないわけですが…。

お礼

a + (double)i * h の aも0.02ずつ値が大きくなるということでしょうか？つまり、　i + (double)i * h でしょうか？

回答No.6

asuncion

> 半角の公式や余弦定理などを使い　分解してみました！

シンプソンの公式にはそういう変形が必要だとは書いていないです。

回答No.5

asuncion

> ∫0～5　e^(-3t) * cos^2 *(2t)

積分したい関数と

> y = (exp(-3. * t) / 2.) + ((exp(-3.* t) * cos(4. * t)) / 2.);

実際のコードに書かれたこの式との関係がわかりません。

お礼

半角の公式や余弦定理などを使い　分解してみました！

回答No.4

aquarius_hiro

できるだけ自分で考えてほしいので、まどろっこしい書き方になってしまいますが…。

> double fの中での　t にうまく代入されていないということですか？

そういうことではなく…

for (i = 0.01; i < 5; i = i + 0.02)

では、i はどのように増えていっていますか？

そして、

f(a + (double)i * h) と書かれていますが、

この変数は意図したように値が増えていっていますか？

よ～く、考えてみてください。

h = ?

もしこの式なら、i はどう増えていくべきですか？

修正の仕方は二通りあります。

お礼

for (i = 0.01; i < 5; i = i + 0.02)は奇数だけを考えるようにしたので0.02ずつ増えていくようにしました。その代わり偶数を考えるfor文も作り、こぼれがないようにしたつもりです。

f(a + (double)i * h) の式の意味が分からなくなってしまいました。今回の問題には使えないような気がしてきました・・・。

回答No.3

Tacosan

あ, h はその意味なんですね. それなら 3 ですね. 了解です.
あと, ちょっとそのプログラムを動かしてみたんだけど危険なところがいくつかありますね. 簡単に列挙すると
・stdio.h や math.h は #include <～> でインクルードすべし
・f の宣言は main の中でしない方が安全
・やたらと無駄なキャストをし過ぎ
・for で double の変数を動かすのは危険
くらい?
ついでにいうと, #2 で挙げられたコメントとは別件で, 多分もう 1つ間違ってますね. 「偶数番目の値」のコメントのあとの for 文があやしい感じ.

お礼

y = (exp(-3. * t) / 2.) + ((exp(-3.* t) * cos(4. * t)) / 2.)　のtに値を入れていくのはdouble fのなかで　forで回せばいいのですか？

for (i = 0.02; i < 5; i = i + 0.02)にすればよいでしょうか？
質問ばかりですみません！

回答No.2

aquarius_hiro

どこが間違っているのか指摘したら、getherさんの勉強にならないと思うので、ヒントをお教えします。

きざみ幅が0.01で、区間の幅が5なら、きざみの数 n はどうなるか考えてみてください。また、h = (b-a)/n と書いた h は何を表しているのかも考えてみてください。0.01とhの両方がプログラムの中に書かれていますが…。

そして、f(x) の x は loop の中で思ったように増えていっているか、よく考えてみてください。

お礼

ｎ＝500ということでしょうか？
ｈは0.01刻みごとの幅だと考えています。

double fの中での　t にうまく代入されていないということですか？
回答ありがとうございました！！

回答No.1

Tacosan

あれ? シンプソン則って, 3 で割るんだっけ? なんか, 6 だったような気がするんだけど....

お礼

ネットで調べたところ　3で正しいようです。
迅速な回答ありがとうございます。