「移動する」球内の複数の配列間でデータを検索するコードのランタイムを高速化する方法

WnGatRC456

時間ステップのサブセットにわたってCFDデータ(スカラーN x M x P配列の形式です。NはY、Mからx、およびPからzに対応します)を平均化しようとしています。以下の希望する平均化プロセスの説明を簡略化しようとしました。

  1. 各タイムステップでグリッドを指定された角度で回転させます(これは、フローが各タイムステップで回転して形状/サイズを変更するコヒーレント構造を持っているためです。それらを重ね合わせて、次のような構造の時間平均形式を見つけたいためです。時間の経過に伴う形状/サイズの変化を説明します
  2. 元の回転していないグリッドを中心とする球を描画する
  3. 球内にあるすべての回転グリッドからグリッドポイントを特定する
  4. 回転した各グリッドのグリッドポイントのインデックスを特定します
  5. インデックスを使用して、球内の回転したグリッドポイントでスカラーデータを見つけます
  6. 球内の値の平均を取る
  7. その新しい平均値を回転していないグリッド上の場所に配置します

自分がやりたいことを正しく実行しているように見えるコードがありますが、計算を完了するのに時間がかかりすぎます。実行速度を上げたいので、必要に応じてコードを変更できます。以下は、小さいバージョンのデータで機能する私のコードのバージョンです。

x = -5:5;       % x - position data
y = -2:.5:5;    % y - position data
z = -5:5;       % z - position data
% my grid is much bigger actually

[X,Y,Z] = meshgrid(x,y,z);    % mesh for plotting data

dX = diff(x)';   dX(end+1) = dX(end);    % x grid intervals
dY = diff(y)';   dY(end+1) = dY(end);    % y grid intervals
dZ = diff(z)';   dZ(end+1) = dZ(end);    % z grid intervals

TestPoints = combvec(x,y,z)';    % you need the Matlab Neural Network Toolbox to run this
dXYZ = combvec(dX',dY',dZ')';

% TestPoints is the unrotated grid

M = length(x);   % size of grid x - direction
N = length(y);   % size of grid y - direction
P = length(z);   % size of grid z - direction

D = randi([-10,10],N,M,P,3);    % placeholder for data for 3 time steps (I have more than 3, this is a subset)
D2{3,M*N*P} = [];
PosAll{3,2} = [];

[xSph,ySph,zSph] = sphere(50);

c = 0.01;   % 1 cm
nu = 8e-6;  % 8 cSt
s = 3*c;    % span for Aspect Ratio 3
r_g = s/sqrt(3);
U_g = 110*nu/c; % velocity for Reynolds number 110
Omega = U_g/r_g;    % angular velocity
T = (2*pi)/Omega;   % period
dt = 40*T/1920;     % time interval
DeltaRotAngle = ((2*pi)/T)*dt;  % angle interval

timesteps = 121:123;   % time steps 121, 122, and 123
for ti=timesteps
    tj = find(ti==timesteps);
    Theta = ti*DeltaRotAngle;
    Rotate = [cos(Theta),0,sin(Theta);...
        0,1,0;...
        -sin(Theta),0,cos(Theta)];
    PosAll{tj,1} = (Rotate*TestPoints')';
end

for i=1:M*N*P
    aa = TestPoints(i,1);
    bb = TestPoints(i,2);
    cc = TestPoints(i,3);
    rs = 0.8*sqrt(dXYZ(i,1)^2 + dXYZ(i,2)^2 + dXYZ(i,3)^2);

    handles.H = figure;
    hs = surf(xSph*rs+aa,ySph*rs+bb,zSph*rs+cc);
    [Fs,Vs,~] = surf2patch(hs,'triangle');
    close(handles.H)

    for ti=timesteps
        tj = find(timesteps==ti);
        f = inpolyhedron(Fs,Vs,PosAll{tj,1},'FlipNormals',false);

        TestPointsR_ti = PosAll{tj,1};
        PointsInSphere = TestPointsR_ti(f,:);

        p1 = [aa,bb,cc];
        p2 = [PointsInSphere(:,1),...
            PointsInSphere(:,2),...
            PointsInSphere(:,3)];

        w = 1./sqrt(sum(...
            (p2-repmat(p1,size(PointsInSphere,1),1))...
            .^2,2));

        D_ti = reshape(D(:,:,:,tj),M*N*P,1);
        D2{tj,i} = [D_ti(f),w];
    end
end

D3{M*N*P,1} = [];
for i=1:M*N*P
    D3{i} = vertcat(D2{:,i});
end

D4 = zeros(M*N*P,1);
for i=1:M*N*P
    D4(i) = sum(D3{i}(:,1).*D3{i}(:,2))/...
                sum(D3{i}(:,2));
end
D_ta = reshape(D4,N,M,P);

各インデックスが、回転していないグリッド内の特定の位置にあるすべてのタイムステップをカバーするすべてのポイントの加重平均であるN x M xP配列を取得することを期待しています。ご覧のとおり、これはまさに私が得たものです。ただし、主な問題は、「実際の」データのより大きなセットを使用する場合に、そうするのにかかる時間の長さです。上記のコードの実行には数分しかかかりませんが、M = 120、N = 24、およびP = 120で、タイムステップ数が24の場合、これにはさらに時間がかかる可能性があります。私の見積もりによると、計算全体が完了するまでに約25日以上かかります。

ダーキー

気にしないでください、私は数学であなたを助けることができます。ここでやろうとしているのは、球の中にあるものを見つけることです。あなたは明確に定義された球を持っているので、これは物事を簡単にします。中心点からすべての点の距離を見つけるだけです。多面体をプロットしたり使用したりする必要はありません。球の中心点で点を変更し、これらの点の距離を計算して、球の半径と比較する行66に注意してください。

% x = -5:2:5;       % x - position data
x = linspace(-5,5,120);
% y = -2:5;    % y - position data
y = linspace(-2,5,24);
% z = -5:2:5;       % z - position data
z = linspace(-5,5,120);
% my grid is much bigger actually

[X,Y,Z] = meshgrid(x,y,z);    % mesh for plotting data

dX = diff(x)';   dX(end+1) = dX(end);    % x grid intervals
dY = diff(y)';   dY(end+1) = dY(end);    % y grid intervals
dZ = diff(z)';   dZ(end+1) = dZ(end);    % z grid intervals

TestPoints = combvec(x,y,z)';    % you need the Matlab Neural Network Toolbox to run this
dXYZ = combvec(dX',dY',dZ')';

% TestPoints is the unrotated grid

M = length(x);   % size of grid x - direction
N = length(y);   % size of grid y - direction
P = length(z);   % size of grid z - direction

D = randi([-10,10],N,M,P,3);    % placeholder for data for 3 time steps (I have more than 3, this is a subset)
D2{3,M*N*P} = [];
PosAll{3,2} = [];

[xSph,ySph,zSph] = sphere(50);

c = 0.01;   % 1 cm
nu = 8e-6;  % 8 cSt
s = 3*c;    % span for Aspect Ratio 3
r_g = s/sqrt(3);
U_g = 110*nu/c; % velocity for Reynolds number 110
Omega = U_g/r_g;    % angular velocity
T = (2*pi)/Omega;   % period
dt = 40*T/1920;     % time interval
DeltaRotAngle = ((2*pi)/T)*dt;  % angle interval

timesteps = 121:123;   % time steps 121, 122, and 123
for ti=timesteps
    tj = find(ti==timesteps);
    Theta = ti*DeltaRotAngle;
    Rotate = [cos(Theta),0,sin(Theta);...
        0,1,0;...
        -sin(Theta),0,cos(Theta)];
    PosAll{tj,1} = (Rotate*TestPoints')';
end
tic
for i=1:M*N*P

    aa = TestPoints(i,1);
    bb = TestPoints(i,2);
    cc = TestPoints(i,3);
    rs = 0.8*sqrt(dXYZ(i,1)^2 + dXYZ(i,2)^2 + dXYZ(i,3)^2);

%     handles.H = figure;
%     hs = surf(xSph*rs+aa,ySph*rs+bb,zSph*rs+cc);
%     [Fs,Vs,~] = surf2patch(hs,'triangle');
%     close(handles.H)

    for ti=timesteps
        tj = find(timesteps==ti);

%         f = inpolyhedron(Fs,Vs,PosAll{tj,1},'FlipNormals',false);
        f = sqrt(sum((PosAll{tj,1}-[aa,bb,cc]).^2,2))<rs;
        TestPointsR_ti = PosAll{tj,1};
        PointsInSphere = TestPointsR_ti(f,:);

        p1 = [aa,bb,cc];
        p2 = [PointsInSphere(:,1),...
            PointsInSphere(:,2),...
            PointsInSphere(:,3)];

        w = 1./sqrt(sum(...
            (p2-repmat(p1,size(PointsInSphere,1),1))...
            .^2,2));

        D_ti = reshape(D(:,:,:,tj),M*N*P,1);
        D2{tj,i} = [D_ti(f),w];
    end
    if ~mod(i,10)
        toc
    end
end

D3{M*N*P,1} = [];
for i=1:M*N*P
    D3{i} = vertcat(D2{:,i});
end

D4 = zeros(M*N*P,1);
for i=1:M*N*P
    D4(i) = sum(D3{i}(:,1).*D3{i}(:,2))/...
                sum(D3{i}(:,2));
end
D_ta = reshape(D4,N,M,P);

実行時間に関しては、私のコンピューターでは、古いコードの実行に57時間かかります。新しいコードには2時間かかります。この時点での主な計算は距離なので、もっと良くなるとは思えません。

この記事はインターネットから収集されたものであり、転載の際にはソースを示してください。

侵害の場合は、連絡してくださいdebugcn@gmail.com

編集
0

コメントを追加

0

関連記事

分類Dev

コードを高速化する方法-データフレームの検索には数時間かかります

分類Dev

パンダのデータフレーム検索を高速化する方法を探しています

分類Dev

iOS / MacOSのメタルコードを高速化する方法

分類Dev

postgresqlのタイムスタンプ列にインデックスを付けることで検索を高速化しますか?

分類Dev

1つの検索バーでデータベース内の複数の列を検索する方法

分類Dev

numpy配列/データフレームの反復プロセスを高速化する方法

分類Dev

PHPで検索するために配列で複数の検索パラメーターを渡す方法

分類Dev

検索に近いジオコーダーを高速化するためのインデックスを追加

分類Dev

Pythonでランダム配列の生成を高速化する方法

分類Dev

複数のリストのインデックスに対するIF / ELIF検索を高速化しますか?

分類Dev

複数のAngularコンポーネント間で静的データを検索する方法は?

分類Dev

Pythonで複数の内積を高速化する方法

分類Dev

ループ内でこの MATLAB コードを高速化する方法は?

分類Dev

コードファイターのこのコードを高速化するjavascriptfirstDuplicate()関数

分類Dev

複数のPythonパンダデータフレーム全体でレコードを検索する

分類Dev

Juliaで複数のブロードキャストを高速化する方法

分類Dev

MySQLデータベースからの複数の結合を高速化するためのヒント

分類Dev

Google Apps Scriptのランタイムを高速化する方法は?

分類Dev

行列計算ランタイムのsympy行列を高速化する方法

分類Dev

1つの文字列内の2つのパターン間で複数の結果を検索する

分類Dev

typeofをループの外に移動することでコードを高速化できますか?

分類Dev

クラスタインデックスの再作成を高速化する方法

分類Dev

Androidアプリ構築の高速化-複数のCPUコアを使用するJavaコンパイラー

分類Dev

NonUniqueフィールドでの検索を高速化するUniquevs NonUnique Clustered Index

分類Dev

複数の IF ステートメントを使用してコードを高速化する方法

分類Dev

データフレーム全体で複数の文字列を検索する

分類Dev

現在のrコードは、worldlcimから気候データを抽出するのに長い時間がかかります。これを高速化する方法

分類Dev

別のデータフレームを参照するパンダのローリングを高速化する

分類Dev

コードを高速化するためのオブジェクトのインスタンス化

Related 関連記事

  1. 1

    コードを高速化する方法-データフレームの検索には数時間かかります

  2. 2

    パンダのデータフレーム検索を高速化する方法を探しています

  3. 3

    iOS / MacOSのメタルコードを高速化する方法

  4. 4

    postgresqlのタイムスタンプ列にインデックスを付けることで検索を高速化しますか?

  5. 5

    1つの検索バーでデータベース内の複数の列を検索する方法

  6. 6

    numpy配列/データフレームの反復プロセスを高速化する方法

  7. 7

    PHPで検索するために配列で複数の検索パラメーターを渡す方法

  8. 8

    検索に近いジオコーダーを高速化するためのインデックスを追加

  9. 9

    Pythonでランダム配列の生成を高速化する方法

  10. 10

    複数のリストのインデックスに対するIF / ELIF検索を高速化しますか?

  11. 11

    複数のAngularコンポーネント間で静的データを検索する方法は?

  12. 12

    Pythonで複数の内積を高速化する方法

  13. 13

    ループ内でこの MATLAB コードを高速化する方法は?

  14. 14

    コードファイターのこのコードを高速化するjavascriptfirstDuplicate()関数

  15. 15

    複数のPythonパンダデータフレーム全体でレコードを検索する

  16. 16

    Juliaで複数のブロードキャストを高速化する方法

  17. 17

    MySQLデータベースからの複数の結合を高速化するためのヒント

  18. 18

    Google Apps Scriptのランタイムを高速化する方法は?

  19. 19

    行列計算ランタイムのsympy行列を高速化する方法

  20. 20

    1つの文字列内の2つのパターン間で複数の結果を検索する

  21. 21

    typeofをループの外に移動することでコードを高速化できますか?

  22. 22

    クラスタインデックスの再作成を高速化する方法

  23. 23

    Androidアプリ構築の高速化-複数のCPUコアを使用するJavaコンパイラー

  24. 24

    NonUniqueフィールドでの検索を高速化するUniquevs NonUnique Clustered Index

  25. 25

    複数の IF ステートメントを使用してコードを高速化する方法

  26. 26

    データフレーム全体で複数の文字列を検索する

  27. 27

    現在のrコードは、worldlcimから気候データを抽出するのに長い時間がかかります。これを高速化する方法

  28. 28

    別のデータフレームを参照するパンダのローリングを高速化する

  29. 29

    コードを高速化するためのオブジェクトのインスタンス化

ホットタグ

アーカイブ