7

Pythonで数値コードを試してみる前の簡単な演習として、LDLTアルゴリズムを作成しようとしています。「足を濡らす」ためだけに。

しかし、私はnumpy配列の基本的な理解が不足しているようです。次の例を参照してください。

def ldlt(Matrix):
    import numpy

    (NRow, NCol) = Matrix.shape

    for col in range(NCol):
        Tmp = 1/Matrix[col,col]
        for D in range(col+1, NCol):
            Matrix[col,D] = Matrix[D,col]*Tmp  
            
if __name__ == '__main__':
    import numpy
    A = numpy.array([[2,-1,0],[-1,2,-1],[0,-1,2]])
    ldlt(A)

この例は、私が取り組んでいる完全なコードではありません。ただし、試して実行し、Matrix [col、D]=..にブレークポイントを設定してください。

最初の評価で期待するのは、行0の列1（開始値-1）が= -1 *（1/2）=-0.5に等しく設定されることです。

ただし、コードを実行すると、0に設定されているように見えます。なぜですか？私が本当に理解していない基本的な何かがあるに違いありませんか？

皆さん、私を助けてくれてありがとう。

編集1：

Python Ver .: 3.3 Tmp .: 0.5になります（私のデバッガーによって報告されたように）。

4

2 に答える 2

4

以下は、何が起こっているかを示している可能性があります。

>>> A = np.array([[2,-1,0],[-1,2,-1],[0,-1,2]])
>>> A.dtype
dtype('int32')
>>> A[0, 1]
-1
>>> A[0, 1] * 0.5
-0.5
>>> A[0, 1] *= 0.5
>>> A[0, 1]
0
>>> int(-0.5)
0

配列は32ビット整数しか保持できないため、配列に割り当てようとする浮動小数点値は、int32にキャストされます。つまり切り捨てられます。

同じ価格で、あなたが求めていたものを実行するためのよりnumpythonicな方法があります：forループは、numpyの目的全体を無効にするため、通常は回避する必要があります。

def ldlt_np(arr) :
    rows, cols = arr.shape
    tmp = 1 / np.diag(arr) # this is a float array
    mask = np.tril_indices(cols)
    ret = arr * tmp[:, None] # this will also be a float array
    ret[mask] = arr[mask]

    return ret

>>> A = np.array([[2,-1,0],[-1,2,-1],[0,-1,2]])
>>> ldlt_np(A)
array([[ 2. , -0.5,  0. ],
       [-1. ,  2. , -0.5],
       [ 0. , -1. ,  2. ]])
于 2013-02-22T16:29:34.787 に答える
0

numpy配列の型は固定されています。int後で配列をfloatに変更することはできません。配列をfloatの配列として初期化します。

A = numpy.array([[2, -1, 0], [-1, 2, -1], [0, -1, 2]], numpy.float)
于 2013-02-22T16:22:02.297 に答える