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

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

ただいまの
回答率

88.60%

numpyを使用しないで転置行列の内積を求める方法および行列の内積の内積の求め方の計算方法について

解決済

回答 2

投稿 編集

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

peace0417

score 9

numpyを使用しないで共役勾配法の実装を試みていますが途中の式の計算がうまくできず困っている.

numpyを使用しないで共役勾配法の解を算出するプログラムを作っています.
途中の a = float( np.dot(r0.T,r0) / np.dot(np.dot(p.T, A),p) )という計算をnpを使わないで計算したいのですがうまく変えずに困っております.
自分なりに内包表記やzipを使用して試行錯誤はしたのですが,通りません.
アドバイスいただけたら幸いです.

また,みにくいソースかもしれませんがご容赦ください.

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

ValueError                                Traceback (most recent call last)
<ipython-input-173-78082ced4a70> in <module>
      1 #print(b.shape)
----> 2 ans = cgm(mat,vec,x0)
      3 print(ans)

<ipython-input-172-15436a5d72cf> in cgm(A, b, x_init)
     14 
     15     for i in range(k < k+1):
---> 16         d1 = dot(list(map(list,zip(*r0))),r0)
     17         print(d1)
     18         d2 = dot(dot(list(map(list,zip(*p))),A),p)

<ipython-input-170-9fb70407aafc> in dot(a, b)
      5     for i,v in enumerate(a):
      6         for j,u in enumerate(zip(*b)):
----> 7             tmp[i][j] = sum([x*y for x,y in zip(v,u)])
      8 
      9     return tmp

<ipython-input-170-9fb70407aafc> in <listcomp>(.0)
      5     for i,v in enumerate(a):
      6         for j,u in enumerate(zip(*b)):
----> 7             tmp[i][j] = sum([x*y for x,y in zip(v,u)])
      8 
      9     return tmp

~/.pyenv/versions/3.6.5/lib/python3.6/site-packages/numpy/matrixlib/defmatrix.py in __mul__(self, other)
    218         if isinstance(other, (N.ndarray, list, tuple)) :
    219             # This promotes 1-D vectors to row vectors
--> 220             return N.dot(self, asmatrix(other))
    221         if isscalar(other) or not hasattr(other, '__rmul__') :
    222             return N.dot(self, other)

ValueError: shapes (1,48) and (1,48) not aligned: 48 (dim 1) != 1 (dim 0)

該当のソースコード

import random
from scipy.io import mminfo,mmread


mat = mmread('bcsstm01.mtx').todense()
print(type(mat))
print(mat.shape)

random.seed(1)
vec = [[random.randrange(10) for i in range(48)]]
print(vec)
print(type(vec))

x0 = [[0]*48]
print(type(x0))


def dot(a,b):

    tmp = [[0]*len(b[0]) for i in range(len(a))]

    for i,v in enumerate(a):
        for j,u in enumerate(zip(*b)):
            tmp[i][j] = sum([x*y for x,y in zip(v,u)])

    return tmp


def cgm(A, b, x_init):
    x = x_init

    d = dot(A,x)[0]
    #print(d)
    r0 = [e1 - e2 for e1,e2 in zip(b,d)]
    #r0 = [e1 - e2 for e1,e2 in zip(b,(dot(A,x)[0]))]
    #print(dot(A,x))
    #print(r0)

    p = r0
    k = 0


    for i in range(k < k+1):
        **d1 = dot(list(map(list,zip(*r0))),r0)
        print(d1)
        d2 = dot(dot(list(map(list,zip(*p))),A),p)
        print(d2)
        a = [e3 / e4 for e3,e4 in zip(d1,d2)]**
        #a = dot(r0.T,r0) / dot(dot(p.T, A),p)
        x = x + p*a
        r1 = r0 - dot(A*a, p)
        if np.linalg.norm(r1) < 1.0e-10:
            k = k+1
            return x
        b =  dot(r1.T, r1) / dot(r0.T, r0) 
        p = r1 + b * p
        r0 = r1
    k = k+1
    return x

ans = cgm(mat,vec,x0)
print(ans)

試したこと

a = float( np.dot(r0.T,r0) / np.dot(np.dot(p.T, A),p) )という計算はnumpyを使用しているため下記のようにzipとリスト内包表記を使用したものに変更してみたがエラーは解決できませんでした.
**でくくってあるところの下のコメントa = dot(r0.T,r0) / dot(dot(p.T, A),p)というよな計算を実装したいです.
助言いただけたら幸いです.

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

Python version 3.5.6
jupyter notebook 使用
numpy sympyは使用しないで内積の計算を行いたいためnp.は使用できません.

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

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

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

    クリップを取り消します

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

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

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

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

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

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

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

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

    質問の評価を下げる

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

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

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

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

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

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

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

    詳細な説明はこちら

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

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

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

質問への追記・修正、ベストアンサー選択の依頼

  • meg_

    2019/07/29 20:47

    numpy.dotとarray.Tを使用しないとのことですので、
    ・配列の転置を実行する関数
    ・配列の内積を計算する関数
      と関数を分けて1ステップずつ実装してみてはいかがでしょうか?

    キャンセル

  • peace0417

    2019/07/29 23:19

    質問ありがとうございます.
    仰っている事は理解できているのですが,知識,技術共に乏しく基本的なことしかわからず実装に至らないため,書き方のヒントやどういった関数を使用すれば上手くいくか等のアドバイスをお願いできないでしょうか.

    キャンセル

  • meg_

    2019/07/30 10:45

    内包表記は理解されているということで良いてすか?

    キャンセル

回答 2

checkベストアンサー

0

def tenti(x):
    return [list(j) for j in [i for i in zip(*x)]]



def mydot(x, y):

    y = tenti(y)
    result = []

    for i in x:
        subresult = []
        for j in y:
            subresult.append(sum([k * m for k, m in zip(i, j)]))
        result.append(subresult)
    return result

コードを簡潔にするためエラー処理は一切入れてません。実装する際には適宜入れてください。

a = [[1, 2], [3, 4]]
b = [[5, 6], [7, 8]]

print(mydot(a, b))
>>>[[19, 22], [43, 50]]

a = [[1, 2, 3], [4, 5, 6]]
b = [[1, 4], [2, 5], [3, 6]]

print(mydot(a, b))
>>>[[14, 32], [32, 77]]

投稿

  • 回答の評価を上げる

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

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

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

  • 回答の評価を下げる

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

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

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

  • 2019/08/03 13:25

    返信遅くなり申し訳ございません.
    教えていただいた書き方で無事うまく出力されました.
    ありがとうございました.

    キャンセル

0

二次元のリストを転置する例

x=[[1,2],[3,4]]
[list(j) for j in [i for i in zip(*x)]]

実行結果:[[1, 3], [2, 4]]

※xは要素数が同じ二次元のリストとする

投稿

  • 回答の評価を上げる

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

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

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

  • 回答の評価を下げる

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

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

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

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

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

関連した質問

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