3

私は、3次元空間の規則的な間隔のグリッド位置で、関心のあるさまざまな測定値（温度、密度、圧力など）を与える有限要素プログラムの結果を持っています。

値は各座標に沿って等間隔に配置されますが、この間隔は座標ごとに異なる場合があります。例えば、

x1 = [0, 0.1, 0.2, ..., 1.0]      (a total of NX1 pts) 
x2 = [0, 0.5, 1.0, ..., 20]       (a total of NX2 pts) 
x3 = [0, 0.2, 0.4, ..., 15]       (a total of NX3 pts)

ソフトウェアから出力される結果は、次の形式です。

x1_1, x2_1, x3_1, f_x, g_x, h_x
x1_1, x2_1, x3_2, f_x, g_x, h_x
x1_1, x2_1, x3_3, f_x, g_x, h_x
...
x1_1, x2_2, x3_1, f_x, g_x, h_x
x1_1, x2_2, x3_2, f_x, g_x, h_x
x1_1, x2_2, x3_3, f_x, g_x, h_x
...
x1_2, x2_1, x3_1, f_x, g_x, h_x
x1_2, x2_1, x3_2, f_x, g_x, h_x
x1_2, x2_1, x3_3, f_x, g_x, h_x
...

ここで、f_x、g_x、h_xは、特定のグリッドポイントで対象となるメジャーです。

上記のデータ形式を変換して、f、g、hの（NX1 x NX2 x NX3）numpy配列を取得したいと思います。

一部の結果セットはかなり大きい（80 x 120 x 100）。

この変換を効率的に行うためのヒントはありますか？

4

2 に答える 2

1

data配列全体を形状の配列としてメモリに読み込んだとします(Nx1 * Nx2 * Nx3, 6)。

data = np.loadtxt('data.txt', dtype=float, delimiter=',')

例が示すように、ポイントが辞書式順序で生成される場合は、列を取得して、それらの形状を変更するだけで済みfます。gh

f = data[:, 3].reshape(Nx1, Nx2, Nx3)
g = data[:, 4].reshape(Nx1, Nx2, Nx3)
h = data[:, 5].reshape(Nx1, Nx2, Nx3)

Nx1、Nx2およびが何であるかを理解する必要がある場合Nx3は、次を使用できますnp.unique。

Nx1 = np.unique(data[:, 0]).shape[0]
Nx2 = np.unique(data[:, 1]).shape[0]
Nx3 = np.unique(data[:, 2]).shape[0]

ポイントの順序が保証されていない場合のより堅牢なソリューションはnp.unique、グリッド値へのインデックスを抽出するために使用することです。

Nx1, idx1 = np.unique(data[:, 0], return_inverse=True)
Nx1 = Nx1.shape[0]
Nx2, idx2 = np.unique(data[:, 1], return_inverse=True)
Nx2 = Nx2.shape[1]
Nx3, idx3 = np.unique(data[:, 2], return_inverse=True)
Nx3 = Nx3.shape[0]

f = np.empty((Nx1, Nx2, Nx3))
f[idx1, idx2, idx3] = data[:, 3]
g = np.empty((Nx1, Nx2, Nx3))
g[idx1, idx2, idx3] = data[:, 4]
h = np.empty((Nx1, Nx2, Nx3))
h[idx1, idx2, idx3] = data[:, 5]

これにより、元の配列へのビューではなく、、およびのf新しいg配列が作成されるため、より多くのメモリを使用します。hdata

そしてもちろん、上記の醜いコードですべてを3回繰り返す代わりに、ループまたはリスト内包表記を使用する必要があります。

于 2013-02-27T15:56:27.603 に答える
0

何があってもファイル全体をループする必要があるので、配列を初期化して値を入力してみませんか？

トリッキーな部分は、形状のある配列が必要な(NX1,NX2,NX3)場合でも、x1,x2,x3値がfloatsの場合は、何らかの方法で配列にインデックスを付ける必要があることです。おそらくこれにはデータ構造が存在しますが、次のようなものを使用できます

def xyz_index((x,y,z),(n1,n2,n3)):
    """ return integer indices for x,y,z position
        given a constant step """
    return tuple(map(int,[x/n1,y/n2,z/n3]))

import numpy as np
NX1,NX2,NX3 =  (80, 120, 100)
ns = n1, n2, n3 =   (.1,.5,.2)
x1, x2, x3 = np.arange(0,1+n1,n1), np.arange(0,20+n2,n2), np.arange(0,15+n3,n3),

data = np.empty((NX1,NX2,NX3),dtype=[('f',float),('g',float),('h',float)])
with open(filename,'r') as f:
    for line in f:
        x,y,z,f,g,h = map(float,line.split(', '))
        data[xyz_index((x,y,z),ns)] = (f,g,h)

次に、次のようにデータにアクセスできます。

ポイントのh-valueにはx1,x2,x3 = .2, .5, 0.、

data[xyz_index((.2,.5,0),ns)]['h']

がないと、上記のタプル['h']が返されます。(f,g,h)dtype

(NX1,NX2,NX3)インデックスがないと、すべての値の配列が返されhます。

これを見てみましょn1, n2, n3う。常に同じである場合は、関数内でそれらを定義することをお勧めします。そうすれば、毎回inxyz_indexを渡す必要はありません。ns

data[xyz_index(.2,.5,0)]['h']
于 2013-02-27T15:33:41.290 に答える