質問をすることでしか得られない、回答やアドバイスがある。

15分調べてもわからないことは、質問しよう!

ただいまの
回答率

90.33%

  • Python

    9209questions

    Pythonは、コードの読みやすさが特徴的なプログラミング言語の1つです。 強い型付け、動的型付けに対応しており、後方互換性がないバージョン2系とバージョン3系が使用されています。 商用製品の開発にも無料で使用でき、OSだけでなく仮想環境にも対応。Unicodeによる文字列操作をサポートしているため、日本語処理も標準で可能です。

  • Python 3.x

    7401questions

    Python 3はPythonプログラミング言語の最新バージョンであり、2008年12月3日にリリースされました。

Python ガウス消去法の中に出てくるZeroDivisionErrorが解決できません

解決済

回答 2

投稿

  • 評価
  • クリップ 0
  • VIEW 210

Akiranaba

score 2

 前提・実現したいこと

python3.xで、連立一次方程式をPivot選択付きGauss消去法で解こうとしています。前進消去法の計算過程でエラーが出てしまいます。
Ax=bについては、Aは5次元ヒルベルト行列、x=(1,1,1,1,1)'、bはxの値に伴ってAの列ごとの和を成分に持つ5次元ベクトルです。

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

File "<ipython-input-195-0db97456e826>", line 1, in <module>
    forwardelimination(a,b)

  File "<ipython-input-194-b8398e40a32d>", line 4, in forwardelimination
    aik = a[i][k] / a[k][k]

ZeroDivisionError: float division by zero

 該当のソースコード

N = 5
a = [[0 for i in range(5)] for j in range(5)]
b = [0 for i in range(5)]

#ヒルベルト行列(A)の作成#
def makehilbertmatrix (N, a):
    for i in range(0, N, 1):
        for j in range(0, N, 1):
            a[i][j] = 1 / (i + 1 + j + 1 - 1)

#ベクトルbの作成#
def makeb (a, b):
    for i in range(0, N, 1):
        for j in range (0, N, 1):
            b[i] += a[i][j]

#ピボット選択付き前進消去法#
def forwardelimination (a, b):
    for k in range(0, N, 1):
        for i in range(k, N - 1, 1):
            aik = a[i][k] / a[k][k]
            for j in range(k, N ,1):
                a[i][j] -= a[k][j] * aik
                b[i] -= b[k] * aik
                #ピボット選択#
                for k in range (0, N, 1):
                    pivot_i = k
                    pivot_val = 0.0
                    for i in range(k, N, 1):
                        if (abs(a[i][k]) > pivot_val):
                            pivot_i = i
                            pivot_val = abs(a[i][k])
                            if (pivot_i != k):
                                tmp = 0.0
                                for j in range(0, N, 1):
                                    tmp = a[pivot_i][j]
                                    a[pivot_i][j] = a[k][j]
                                    a[k][j] = tmp
                                    tmp = b[i]
                                    b[i] = b[k]
                                    b[k] = tmp

#関数の実行#
makehilbertmatrix(N, a)
makeb(a, b)
forwardelimination(a, b)

 試したこと

ピボット選択の部分を外して実行しましたが同様のエラーが出ます。
それぞれのa[k][k](エラーが出ているコマンドの割り算の分母の部分)を表示させると

def func (a, b):
    for k in range(0, N, 1):
        for i in range(k, N - 1, 1):
            print(a[k][k])


func(a,b)
0.3333333333333333
0.3333333333333333
0.3333333333333333
0.3333333333333333
0.5
0.5
0.5
0.25
0.25
0.14285714285714285


このようになり、0は入っていません。

 補足情報(FW/ツールのバージョンなど)

ここにより詳細な情報を記載してください。

  • 気になる質問をクリップする

    クリップした質問は、後からいつでもマイページで確認できます。

    またクリップした質問に回答があった際、通知やメールを受け取ることができます。

    クリップを取り消します

  • 良い質問の評価を上げる

    以下のような質問は評価を上げましょう

    • 質問内容が明確
    • 自分も答えを知りたい
    • 質問者以外のユーザにも役立つ

    評価が高い質問は、TOPページの「注目」タブのフィードに表示されやすくなります。

    質問の評価を上げたことを取り消します

  • 評価を下げられる数の上限に達しました

    評価を下げることができません

    • 1日5回まで評価を下げられます
    • 1日に1ユーザに対して2回まで評価を下げられます

    質問の評価を下げる

    teratailでは下記のような質問を「具体的に困っていることがない質問」、「サイトポリシーに違反する質問」と定義し、推奨していません。

    • プログラミングに関係のない質問
    • やってほしいことだけを記載した丸投げの質問
    • 問題・課題が含まれていない質問
    • 意図的に内容が抹消された質問
    • 広告と受け取られるような投稿

    評価が下がると、TOPページの「アクティブ」「注目」タブのフィードに表示されにくくなります。

    質問の評価を下げたことを取り消します

    この機能は開放されていません

    評価を下げる条件を満たしてません

    評価を下げる理由を選択してください

    詳細な説明はこちら

    上記に当てはまらず、質問内容が明確になっていない質問には「情報の追加・修正依頼」機能からコメントをしてください。

    質問の評価を下げる機能の利用条件

    この機能を利用するためには、以下の事項を行う必要があります。

回答 2

+1

手元で計算過程を出力してみましたが、ガウス消去法で計算している過程で対角成分が 0 になるケースが実際に発生しています。

def forwardelimination (a, b):
    for k in range(0, N, 1):
        for i in range(k, N - 1, 1):
            print(a)
            aik = a[i][k] / a[k][k]

のようにしてみればどんな値になっているのかわかると思います。

ところで、ここまで for や if がネストしていると大変読みにくいので、いくつかの関数に分解するか処理を見直すことをおすすめします。

投稿

  • 回答の評価を上げる

    以下のような回答は評価を上げましょう

    • 正しい回答
    • わかりやすい回答
    • ためになる回答

    評価が高い回答ほどページの上位に表示されます。

  • 回答の評価を下げる

    下記のような回答は推奨されていません。

    • 間違っている回答
    • 質問の回答になっていない投稿
    • スパムや攻撃的な表現を用いた投稿

    評価を下げる際はその理由を明確に伝え、適切な回答に修正してもらいましょう。

  • 2018/06/25 20:10

    ループ2回目でいきなり0ですね。

    キャンセル

check解決した方法

0

前進消去法の関数内の2つ目のfor文の範囲が間違っていました.消去するのはピボット選択された行以降の行でした.下が訂正したものになります

#ピボット選択付き前進消去法#
def forwardelimination (a, b):
    for k in range(0, N, 1):
        for i in range(k + 1, N, 1):
            aik = a[i][k] / a[k][k]
            for j in range(k, N ,1):
                a[i][j] -= a[k][j] * aik
                b[i] -= b[k] * aik
                #ピボット選択#
                for k in range (0, N, 1):
                    pivot_i = k
                    pivot_val = 0.0
                    for i in range(k, N, 1):
                        if (abs(a[i][k]) > pivot_val):
                            pivot_i = i
                            pivot_val = abs(a[i][k])
                            if (pivot_i != k):
                                tmp = 0.0
                                for j in range(0, N, 1):
                                    tmp = a[pivot_i][j]
                                    a[pivot_i][j] = a[k][j]
                                    a[k][j] = tmp
                                    tmp = b[i]
                                    b[i] = b[k]
                                    b[k] = tmp

投稿

  • 回答の評価を上げる

    以下のような回答は評価を上げましょう

    • 正しい回答
    • わかりやすい回答
    • ためになる回答

    評価が高い回答ほどページの上位に表示されます。

  • 回答の評価を下げる

    下記のような回答は推奨されていません。

    • 間違っている回答
    • 質問の回答になっていない投稿
    • スパムや攻撃的な表現を用いた投稿

    評価を下げる際はその理由を明確に伝え、適切な回答に修正してもらいましょう。

同じタグがついた質問を見る

  • Python

    9209questions

    Pythonは、コードの読みやすさが特徴的なプログラミング言語の1つです。 強い型付け、動的型付けに対応しており、後方互換性がないバージョン2系とバージョン3系が使用されています。 商用製品の開発にも無料で使用でき、OSだけでなく仮想環境にも対応。Unicodeによる文字列操作をサポートしているため、日本語処理も標準で可能です。

  • Python 3.x

    7401questions

    Python 3はPythonプログラミング言語の最新バージョンであり、2008年12月3日にリリースされました。