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

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

ただいまの
回答率

90.50%

  • Python

    8026questions

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

  • Python 3.x

    6438questions

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

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

解決済

回答 2

投稿

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

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

投稿

  • 回答の評価を上げる

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

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

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

  • 回答の評価を下げる

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

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

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

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

  • ただいまの回答率 90.50%
  • 質問をまとめることで、思考を整理して素早く解決
  • テンプレート機能で、簡単に質問をまとめられる

関連した質問

  • 解決済

    pythonで実行時間を求めたい

    実現したいこと 10種類の疑似乱数をバブルソートで実行。 その実行時間をしりたい。 エラーメッセージ t=timeit.Tiimer("Bublesort(Random

  • 解決済

    ソート結果が正しくない

    前提・実現したいこと csvファイルdata1.csvおよびそれと同じフォーマットのcsvファイルdata2.csv, csvファイルdata3.csvのそれぞれに対して,温度(

  • 解決済

    Python 3.x 辞書のキー値によって変換する場合の高速化

    Pythonにて、辞書(dict({key,value})を使って、list型の全要素をValue値に変換する際の、 高速化が可能かどうかをご教授いただきたいです。 dict1

  • 解決済

    高階関数の構造が理解できません。

    Python初心者です。書籍を購入して勉強していますが、内容を読んでも意味が分からない箇所があったので質問します。 書籍:Pythonプログラミング パーフェクトマスター (秀

  • 解決済

    構造の異なる辞書をソートしたい

    構造の異なる辞書をソートしたい。 辞書dfには df = [{"A":{"A1":{"a1":100,"b1":200,"c1":300}}},{"B":{"B1":{"a2"

  • 解決済

    pythonで_chrの意味がわかりません

    def long_repeat(line): #文字列の中の最大連続数表示の関数 if line == '': return(0) else:

  • 受付中

    Pythonにも『ラムダ式の初期化キャプチャ』のような機能がありますか?

    C++14には所謂『ラムダ式の初期化キャプチャ』という機能があります。 即ち、関数の中の部分パラメーターが環境やコンテクストに依存するものの、値は固定です。 なので、これらの環境依

  • 解決済

    Python 素数の表示

    現在 for文 と if文程度の基礎的な文法のみで 素数を小さいほうから2000番目を表示するプログラムを考えています。 アドバイス お願いいたします。 a = 1000000

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

  • Python

    8026questions

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

  • Python 3.x

    6438questions

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