シンプソン 公式。 [C言語] 数値積分(シンプソンの公式)

C++

シンプソン 公式

CONTENTS• スポンサーリンク 数値計算による積分計算の基本~Simpson法~ 積分の計算は高校生までは紙とペンを使用し、様々な積分テクニック(置換積分・部分積分など)を駆使して行っていきます。 大学の初年度に上がると重積分を学習しますが、これも紙とペンを使って問題を解いていき、単位を取得します。 ところが、大学3年生以上になり、専門分野の勉強が進んでいくと、紙とペンでは解けない積分や、紙とペンを解くのが非常に難解だったり、手続きが多くてミスが起こりやすい積分の問題に度々出会うようになります。 そこで、「積分の計算をコンピュータで解く方法」が開発されており、代表例として シンプソン法があります。 今回の資料は、シンプソン法の理論的解説を詳しく行い、 プログラムを解いたり、シンプソン法を使ったプログラムを使用する際に迷いが起こらないようにするための基礎知識を整理するための資料になります。 台形近似より精度の高いシンプソン法 うん! まとめ 今回は、積分をコンピュータで解くのに使用されるシンプソン法について解説しました。 積分をコンピュータで解く初歩的な方法には 台形近似があるのですが、台形近似では非線形性の高い関数の積分を行ったときに精度の甘さが出てきてしまいます。 ここで、 「積分計算を台形の集合体の和ではなく、二次関数の集合体の和に近似しよう。 」というのがシンプソン法の発想です。 ストークも惑わされた 言葉の定義についてですが、積分したい関数を二次関数に近似するだけの公式を単に 「シンプソンの公式」といい、 積分したい関数を二次関数の集合体として取り扱う公式を 「合成シンプソンの公式」と呼びます。 私たちが積分の計算でプログラムを組む時に使われるのは「合成シンプソンの公式」の方になります。 合成シンプソンの公式でググってみるとわかると思いますが、色々な式の形があり、混乱することがあると思います。 奇数だけ足し合わせたり、偶数だけ足し合わせるとか、色々です。 しかし、 基本思想として「2次関数近似の寄せ集め」ということを理解しておけば、自分なりの理解の仕方がマスターでき、それに基づいてプログラムを書けるようになるでしょう。 この資料は私のシンプソン法への自分なりの理解の整理の結果でもあるのです。

次の

C言語のアルゴリズム(シンプソン公式)

シンプソン 公式

5 数値積分(定積分) ここでは,1変数の連続関数の定積分を数値的に求めることを扱う。 一般に微分方程式の解を数値的に求めることも数値積分ということになるが,これは後で扱う。 不定積分が代数関数や超越関数(三角関数,対数関数,指数関数など)で与えられる場合は,定積分を求めるのにわざわざ数値積分する必要はない。 そうでない場合でも,よく出てくるものについては定積分が特殊関数として定義され(誤差関数,楕円関数など),市販の数学ライブラリに納められていたり,数値計算のデータブックで優秀な近似式が見つかることも多い。 しかしながら,比較的性質の良い関数なら,変形してライブラリを探したりするよりは,自分でプログラムを組んでしまう方が手っ取り早く,計算時間をそれほど気にしなくても,実用的には十分な精度で数値解が得られるものである。 しかしながら階段ではなく,各区間で関数を右図のように直線で近似して台形の面積の和を求めることにすれば,プログラム上は以下のようにほんのわずかな修正で精度は格段に上がる。 右端(または左端)は加えない。 理由を考えてみよう。 階段で近似する場合の誤差は,図のような単調関数の場合,明らかに両端の柱部分である。 このことから一般に N 等分した場合の近似値 s N と真の値 s との差は s N - s < f max - f min h f max,f min は, a,b での最大値と最小値 であると考えてよい。 a である。 である。 の無い形で表すことが可能である。 問題は 「ノートNo. 2」 で述べた「丸め誤差の蓄積」である。 7, F12. 14159265 END DO END! 1389885 0. 999171 16 3. 1409416 0. 999793 32 3. 1414294 0. 999948 64 3. 1415517 0. 999987 128 3. 1415832 0. 999997 256 3. 1415906 0. 999999 512 3. 1415911 0. 999999 1024 3. 1415923 1. 000000 2048 3. 1415927 1. 000000 4096 3. 1415932 1. 000000 8192 3. 1415925 1. 000000 16384 3. 1415939 1. 000000 32768 3. 1415918 1. 000000 65536 3. 1415901 0. 999999 131072 3. 1415541 0. 999988 262144 3. 1415567 0. 999989 524288 3. 1412396 0. 999888 1048576 3. 1404257 0. 999628 区分数をむやみに増やしても,何の役にもたたず,かえって誤差が増えていることがわかるだろう。 この区分数を2のn乗の形にとることは,次のシンプソン公式との関係で意味がある。 (2)シンプソン公式 連続する3点を放物線で近似することにより,さらに精度はよくなる。 ] となる。 プログラムはやや複雑になる。 (関数の定義部分は前出のプログラムと同じだから省略した。 ) INTEGER :: j, k, N 1 4 1 1 4 1 1 4 1... REAL :: h, s, x, f 1 4 1 1 4 1..... PRINT " I10, F12. 7, F12. 1415927 1. 000000 16 3. 1415927 1. 000000 32 3. 1415927 1. 000000 64 3. 1415932 1. 000000 128 3. 1415930 1. 000000 256 3. 1415913 1. 000000 512 3. 1415930 1. 000000 1024 3. 1415930 1. 000000 2048 3. 1415913 1. 000000 4096 3. 1415951 1. 000001 8192 3. 1415949 1. 000001 16384 3. 1415930 1. 000000 32768 3. 1415961 1. 000001 65536 3. 1416073 1. 000005 131072 3. 1416357 1. 000014 262144 3. 1414626 0. 999959 524288 3. 1418412 1. 000079 1048576 3. 1398065 0. 999431 結果を見てわかるように,このような滑らかな関数であれば,有効数字6桁の範囲でなら,実に8等分でも十分なのである。 単一の定積分を求めるだけなら,台形公式でも十分役立つ。 アルゴリズムを暗記できるくらいにプログラムが簡単だからである。 しかしながら大きなプログラムの中で定積分を繰り返し何度も必要とする場合にはシンプソン公式が有効である。 ただし上の例ではプログラムを分かりやすくするため,これを使ってはいない。 こうして計算して途中の結果も出力すれば,結果的には同じ計算回数でもって近似が上がっていく(あるいは次数が高すぎると悪くなっていく)経過を見ることができる。 関数は,ここで用いたものでもよいし,そのほか予め定積分の分かっているものがよい。 [Visual Fortran を使用する場合は,コンパイルの際に実数処理最適化を無効にするオプション「 -Op 」 をつけてコンパイルすること。 上の方法で中間点の値の和を求める時にも利用できる。 これはあくまでも計算機の誤差の特徴や実数計算の最適化を理解してもらうためであって,実際にはこんな計算を行うことは愚かであり, 結果としてほしい精度を得るのに必要最小限の区分数を用いればよい。 (4)ロンバーグ公式 倍精度(有効数字15桁)で結果を得たいときには,原理的にはシンプソン公式でも1000〜5000等分くらいの区分数が必要である。 そこで,さらに高次の多項式で近似することが考えられよう。 台形公式とシンプソン公式の関係は次のように書くこともできる。 と展開される( オイラー-マクローリンの公式)。 ここでは右辺の展開係数の具体的な形を気にすることはない。 となり近似の次数が向上する。 区分幅が2倍のとき,先頭の項は今度は 2 4 倍になる。 高次の多項式を用いて補間近似を行う公式の最も分かりやすい例であろう。 , 2 p 区分の台形公式による値, S 0, 0 , S 0, 1 ,... , S 0, p-1 , S 0, p を求めておき,せっかくそこまで計算したのであれば,次に上の漸化式を用いて S p, p まで求めれば,計算操作をさほど増やすことなく与えられた区分数において原理的には最良の近似を得ることができる。 (5)無限区間の積分 積分区間が無限になる場合,あるいは区間は有限でも端(または内部のどこか)で関数が発散している場合や,関数そのものは発散しなくても微分が発散しているような場合,以上のようにはいかない。 思うようにはいかないであろう。 ) このような場合には,一般には変数変換によって特異性を除くことができれば,ここまで述べてきた方法を適用することができる。 あるいは,予め関数の漸近形が評価できる場合には,発散点を挟む狭い領域を除くとか積分区間を有限区間で近似した場合の誤差を評価することができるから,素朴にはこのようにして積分値を求めることも考えられよう。 最後におそらく意外な事実をあげておこう。 これは(4)で述べたを見れば予想できよう。 この場合の誤差は,区分幅ではなく打ち切りから来るものの方が本質的である。 関数の収束が十分に速い場合には,信じられないほどの粗さで驚くべき精度の結果が得られる。 0625, 0. 125, 0. 25, 0. 5 でである。 この事実を利用すれば,有限区間の積分であっても,変数変換により収束の速い関数の無限区間の積分に変換することによって,計算効率を上げることが考えられる(二重指数関数法)が,ここでは割愛する。 単発の積分なら,ここで述べた拡張定積分も含めてやってのけてくれる。 「プログラム」「Maple V」を立ち上げてから以下を実行してみよ。 (>はプロンプトで,行入力は ;(セミコロン)で閉じさせること。 b と指定する。

次の

数値計算入門I〜第7回〜

シンプソン 公式

シンプソンの公式マクロの使い方 1)マクロをコピペする Visual Basicのエディタを開き、マクロをコピペして下さい。 マクロの使い方が分からない方は、をご覧ください。 '----------------------------------------------------------------------- ' Summary : シンプソン法により、定積分を計算する ' Comments: 初めは単純に本関数を呼び出す ' : 次からは、関数、下限値、上限値、分割数を変えて呼び出す ' : 他のブックからの呼び出し方法: ' : Application. Run "Excelのファイル名! Simpson" ' : 例:個人用マクロブックにこのマクロを入れた場合は下行になる ' : Application. Run "PERSONAL. XLSB! Add sht. Range "A1". Range "A2". Range "A3". Range "A4". Range "A5". Range "A6". Names. Names. Names. Names. Names. Names. Range "fx". Range "x". Range "a". Range "b". Range "N". Range "D:E". Range "a". Range "b". Range "N". Range "A2". Range "A1" ' 変数xに下限値をセット. Range "x". Value ' xの下限値とその時の関数値をD2とE2にそれぞれ記入. Range "x". Range "x". Range "fx". Range "fx". Value End If ' 変数xと関数f x の値をそれぞれD、E列に出力する. Range "x". Range "fx" Next i ' 変数xに上限値をセット. Range "x". Range "fx". Range "x". Range "fx" ' 積分の結果を表示する. Range "Result". 「シンプソン」というシートができて、そこに勝手に値が書きこまれていきます。 3秒以内で終わると思います。 B2は確率変数xの値で、マクロ終了後には積分範囲の上限値が入っています。 B3は積分範囲の下限a(-3)です。 B4は積分範囲の上限b(3)です。 B5は分割数N(100)です。 B6は積分の結果(厳密には0. 997300192454353)です。 99730020393674」なので、差は「-1. 14824E-08」です。 どちらが正解に近いのか分かりませんが、なかなかの高精度で合っているので、マクロのコードに致命的なミスはなさそうです。 3)所望の関数に書き換える B1~B5を書き換えることで、所望の関数を定積分することができます。 すなわち、 です。 B1セルの関数定義には、xという文字が使えるので、お使い下さい。 B3 a を「0」とします。 積分の結果は「2. 0000000108245」でした。 正解は「2」ですから、近似値としては十分だと思います。 5)グラフを書く D、E列に、変数xとf x の数値表が出来あがりますので、必要に応じてグラフを作成下さい。

次の