Main Content

A\b のベンチマーク

この例では、クラスターで線形計算の解を求めるベンチマーク方法を示します。A*x = bx を求めるための MATLAB® コードは非常にシンプルです。最も一般的な方法は、行列の左除算 (mldivide またはバックスラッシュ演算子 (\) とも呼ばれる) を使用して x (つまり、x = A\b) を計算するものです。ただし、クラスターでの行列の左除算のパフォーマンスのベンチマークは容易ではありません。

ベンチマークに関する最も重要な課題の 1 つは、システムの全体的なパフォーマンスを 1 つの数値で表そうとする誘惑に負けないようにすることです。性能曲線を確認すると、クラスターにおけるパフォーマンスのボトルネックを特定したり、作成したコードのベンチマークを実行してその結果から有意義な結論を導き出せるようにする方法について理解したりするためにも役立つ可能性があります。

関連する例:

この例のコードは以下の関数に含まれています。

function results = paralleldemo_backslash_bench(memoryPerWorker)

クラスターに適した行列サイズを選択することが非常に重要です。そのための方法として、各ワーカーで使用可能なシステム メモリの量を GB 単位で引数としてこの例の関数に指定します。ここでの既定値は極めて控えめなものであるため、使用するシステムに応じた値を指定してください。

if nargin == 0
    memoryPerWorker = 8.00; % In GB
%    warning('pctexample:backslashbench:BackslashBenchUsingDefaultMemory', ...
%            ['Amount of system memory available to each worker is ', ...
%             'not specified.  Using the conservative default value ', ...
%             'of %.2f gigabytes per worker.'], memoryPerWorker);
end

オーバーヘッドの回避

線形計算の求解の性能に関する正確な測定値を取得するため、オーバーヘッドの要因として考えられるすべてのものを排除する必要があります。そのための作業の一環として、現在の並列プールを取得し、デッドロック検出機能を一時的に無効にします。

p = gcp;
if isempty(p)
    error('pctexample:backslashbench:poolClosed', ...
        ['This example requires a parallel pool. ' ...
         'Manually start a pool using the parpool command or set ' ...
         'your parallel preferences to automatically start a pool.']);
end
poolSize = p.NumWorkers;
pctRunOnAll 'mpiSettings(''DeadlockDetection'', ''off'');'
Starting parallel pool (parpool) using the 'bigMJS' profile ... connected to 12 workers.

ベンチマーク関数

ここでは、行列の左除算 (\) のベンチマークを実行します。spmd ブロックに入るときのコスト、行列作成時間、またはその他のパラメーターのベンチマークは実行しません。そのため、データ生成を線形計算の求解と切り離してから、求解のみにかかる時間を測定します。入力データは 2 次元ブロックサイクリック対話型分散オブジェクトを使用して生成します。これが線形計算の求解のための最も効果的な分散スキームであるためです。このベンチマークでは、すべてのワーカーが線形方程式 A*x = b の求解を完了するまでにかかる時間を測定します。ここでも、考えられるオーバーヘッドの要因が排除されるようにします。

function [A, b] = getData(n)
    fprintf('Creating a matrix of size %d-by-%d.\n', n, n);
    spmd
        % Use the codistributor that usually gives the best performance
        % for solving linear systems.
        codistr = codistributor2dbc(codistributor2dbc.defaultLabGrid, ...
                                    codistributor2dbc.defaultBlockSize, ...
                                    'col');
        A = codistributed.rand(n, n, codistr);
        b = codistributed.rand(n, 1, codistr);
    end
end

function time = timeSolve(A, b)
    spmd
        tic;
        x = A\b; %#ok<NASGU> We don't need the value of x.
        time = gop(@max, toc); % Time for all to complete.
    end
    time = time{1};
end

問題の規模の選択

他の多数の並列アルゴリズムと同様に、並列での線形計算の求解のパフォーマンスは行列のサイズによって大きく異なります。したがって、計算に関する "事前の" 予想は以下のようになります。

  • 小さな行列の場合は、やや非効率的

  • 大きな行列の場合は極めて効率的

  • 行列がシステム メモリのサイズに対して大きすぎるため、オペレーティング システムでディスクへのメモリのスワッピングが開始される場合は非効率的

したがって、サイズが異なる多数の行列の計算時間を測定し、この場合において行列が "小さい"、"大きい"、"大きすぎる" ということがそれぞれどのようなことを意味するか理解することが重要です。これまでの検証に基づき、以下のように想定されます。

  • "小さすぎる" 行列のサイズは 1,000 行 1,000 列

  • "大きな" 行列は各ワーカーで使用可能なメモリの 45% 弱を占有

  • "大きすぎる" 行列は各ワーカーで使用可能なシステム メモリの 50% 以上を占有

これらは、経験則に基づく値であるため、正確な値はリリースによって異なる場合があります。したがって、この範囲全体にわたる行列サイズを使用して、予測されるパフォーマンスを確認することが重要です。

ここでは、ワーカーの数に応じて問題の規模を変更するので、ウィーク スケーリングを使用します。ブラックジャックを使用した parfor のシンプルなベンチマーククラスターにおける独立ジョブのベンチマークなど、他のベンチマークの例でもウィーク スケーリングを使用しています。それらの例ではタスク並列計算のベンチマークを実行するので、反復の数がワーカー数に比例するようなウィーク スケーリングとなっています。一方、この例では、データ並列計算のベンチマークを実行するので、行列の上限サイズとワーカー数の関連性を確認します。

% Declare the matrix sizes ranging from 1000-by-1000 up to 45% of system
% memory available to each worker.
maxMemUsagePerWorker = 0.45*memoryPerWorker*1024^3; % In bytes.
maxMatSize = round(sqrt(maxMemUsagePerWorker*poolSize/8));
matSize = round(linspace(1000, maxMatSize, 5));

パフォーマンスの比較: ギガフロップス

パフォーマンスを測定するための値として 1 秒当たりの浮動小数点演算の回数 (フロップス) を使用します。そのようにすると、さまざまな行列サイズとさまざまなワーカー数についてアルゴリズムのパフォーマンスを比較できるようになります。十分に広い範囲の行列サイズについて行列の左除算のパフォーマンス テストが正常に完了した場合、パフォーマンスのグラフは以下のようになります。

このようなグラフを生成することによって、以下の疑問を解消できます。

  • 最も小さい行列はパフォーマンスが低下するほど小さいものなのか

  • 行列が大きすぎて合計システム メモリの 45% を占有する場合にパフォーマンスの低下が見られるか

  • ある数のワーカーについて実現できる最大パフォーマンスはどの程度か

  • ワーカー数を 8 ではなく 16 にした場合にはどの行列サイズでパフォーマンスが向上するか

  • システム メモリによってピーク パフォーマンスが制限されるか

行列サイズを指定すると、ベンチマーク関数によって行列 A と右辺 b が 1 回作成されてから、A\b の求解が複数回実行されて、所要時間の正確な測定値が取得されます。ここでは、HPC チャレンジの浮動小数点演算の回数を使用します。したがって、n 行 n 列の行列の場合、浮動小数点演算の回数は 2/3*n^3 + 3/2*n^2 になります。

function gflops = benchFcn(n)
    numReps = 3;
    [A, b] = getData(n);
    time = inf;
    % We solve the linear system a few times and calculate the Gigaflops
    % based on the best time.
    for itr = 1:numReps
        tcurr = timeSolve(A, b);
        if itr == 1
            fprintf('Execution times: %f', tcurr);
        else
            fprintf(', %f', tcurr);
        end
        time = min(tcurr, time);
    end
    fprintf('\n');
    flop = 2/3*n^3 + 3/2*n^2;
    gflops = flop/time/1e9;
end

ベンチマークの実行

これでセットアップがすべて完了したので、ベンチマークを容易に実行できるようになりました。ただし、計算が完了するまでに時間がかかることがあるため、それぞれの行列サイズのベンチマークが完了するたびに中間ステータス情報を出力するものとします。

fprintf(['Starting benchmarks with %d different matrix sizes ranging\n' ...
         'from %d-by-%d to %d-by-%d.\n'], ...
        length(matSize), matSize(1), matSize(1), matSize(end), ...
        matSize(end));
gflops = zeros(size(matSize));
for i = 1:length(matSize)
    gflops(i) = benchFcn(matSize(i));
    fprintf('Gigaflops: %f\n\n', gflops(i));
end
results.matSize = matSize;
results.gflops = gflops;
Starting benchmarks with 5 different matrix sizes ranging
from 1000-by-1000 to 76146-by-76146.
Creating a matrix of size 1000-by-1000.
Analyzing and transferring files to the workers ...done.
Execution times: 1.038931, 0.592114, 0.575135
Gigaflops: 1.161756

Creating a matrix of size 19787-by-19787.
Execution times: 119.402579, 118.087116, 119.323904
Gigaflops: 43.741681

Creating a matrix of size 38573-by-38573.
Execution times: 552.256063, 549.088060, 555.753578
Gigaflops: 69.685485

Creating a matrix of size 57360-by-57360.
Execution times: 3580.232186, 3726.588242, 3113.261810
Gigaflops: 40.414533

Creating a matrix of size 76146-by-76146.
Execution times: 9261.720799, 9099.777287, 7968.750495
Gigaflops: 36.937936

パフォーマンスのプロット

次に、結果をプロットし、上記の予測グラフと比較します。

fig = figure;
ax = axes('parent', fig);
plot(ax, matSize/1000, gflops);
lines = ax.Children;
lines.Marker = '+';
ylabel(ax, 'Gigaflops')
xlabel(ax, 'Matrix size in thousands')
titleStr = sprintf(['Solving A\\b for different matrix sizes on ' ...
                    '%d workers'], poolSize);
title(ax, titleStr, 'Interpreter', 'none');

ベンチマークの結果が予測したほど良好なものではない場合は、以下の点を検討してみてください。

  • 基礎となる実装では優れたパフォーマンスで定評がある ScaLAPACK が使用されています。したがって、アルゴリズムまたはライブラリが原因で効率性が低下する可能性は極めて低く、以下の項目で説明するように、どちらかというとその使用方法が原因であると考えられます。

  • 行列がクラスターに対して大きすぎたり、小さすぎたりする場合、最終的なパフォーマンスは低くなります。

  • ネットワーク通信速度が遅い場合はパフォーマンスに著しく影響します。

  • CPU 速度とネットワーク通信速度が両方とも非常に速いにもかかわらずメモリ量が限られている場合は、使用可能な CPU とネットワーク帯域幅を最大限に活用できるほど十分に大きい行列のベンチマークを実行できないことがあります。

  • 最大限のパフォーマンスを実現するには、ネットワークのセットアップに合わせて調整された MPI のバージョンを使用し、共有メモリを介して可能な限り多くの通信が行われるような方法でワーカーが実行されるようにすることが重要です。ただし、これらの問題を特定および解決する方法はこの例の対象外であるため、ここでは説明しません。

さまざまな数のワーカーの比較

ここでは、さまざまな数のワーカーを比較する方法を紹介します。具体的には、ワーカーの数を変更しながらこの例を実行して取得されたデータを確認します。データは上記のものとは異なるクラスターで取得します。

クラスターにおける独立ジョブのベンチマークなどの他の例で、異なるワーカー数について並列アルゴリズムのベンチマークを実行する場合は一般的にウィーク スケーリングが使用されることを説明しました。つまり、ワーカー数の増加に応じて問題の規模を大きくするということです。行列の左除算の場合、除算のパフォーマンスは行列のサイズによって大幅に異なるので、一層注意する必要があります。以下のコードを実行すると、さまざまなワーカー数でテストされたすべての行列サイズについてのパフォーマンス (単位: ギガフロップス) のグラフが作成され、"この特定のクラスター" における行列の左除算のパフォーマンス特性について詳しく確認することができます。

s = load('pctdemo_data_backslash.mat', 'workers4', 'workers8', ...
         'workers16', 'workers32', 'workers64');
fig = figure;
ax = axes('parent', fig);
plot(ax, s.workers4.matSize./1000, s.workers4.gflops, ...
     s.workers8.matSize./1000, s.workers8.gflops, ...
     s.workers16.matSize./1000, s.workers16.gflops, ...
     s.workers32.matSize./1000, s.workers32.gflops, ...
     s.workers64.matSize./1000, s.workers64.gflops);
lines = ax.Children;
set(lines, {'Marker'}, {'+'; 'o'; 'v'; '.'; '*'});
ylabel(ax, 'Gigaflops')
xlabel(ax, 'Matrix size in thousands')
title(ax, ...
      'Comparison data for solving A\\b on different numbers of workers');
legend('4 workers', '8 workers', '16 workers', '32 workers',  ...
       '64 workers', 'location', 'NorthWest');

上記のグラフを見て最初に気づくのは、ワーカーの数が 64 個の場合はわずか 4 個のワーカーの場合に比べてはるかに大きいサイズの線形方程式の解を求めることが可能であるということです。また、60,000 行 60,000 列の行列を 4 個のワーカーで処理できたとしてもパフォーマンスは 10 ギガフロップス程度にしかならないということもわかります。したがって、4 個のワーカーでそのような大きさの問題の解を求めるのに十分なメモリを使用できる場合でも、やはり、64 個のワーカーを使用した方がはるかに優れたパフォーマンスが実現されることになります。

4 個のワーカーの曲線の変化を見ると、3 番目に大きいサイズから最大サイズまでの行列ではパフォーマンスは大幅には向上しないことがわかります。これを、さまざまなサイズの行列に関する A\b の予測パフォーマンスを示す前掲のグラフと比較すると、4 個のワーカーでは行列サイズが 7,772 行 7,772 列の場合にパフォーマンスがほぼピークに達することがわかります。

ワーカー数が 8 個の場合と 16 個の場合の曲線を見ると、行列サイズが最大のときにパフォーマンスが低下していることがわかります。これは、つまり、使用可能なシステム メモリがすべて消費されそうになっている状態または既に消費された状態であることを示します。一方、2 番目と 3 番目に大きな行列の間ではパフォーマンスがごく緩やかに向上しており、ある程度安定していることがわかります。したがって、8 個または 16 個のワーカーについては、システム メモリを増設して比較的大きな行列サイズについてテストを実行した場合に、ギガフロップスの大幅な上昇は見込めないものと推測されます。

ワーカーが 32 個と 64 個の場合の曲線を見ると、2 番目と 3 番目に大きな行列サイズの間でパフォーマンスが大幅に向上しています。64 個のワーカーでは、2 番目に大きいサイズから最大サイズの行列の間でもパフォーマンスが大幅に向上しています。したがって、ワーカーが 32 個と 64 個の場合にはピーク パフォーマンスに達する前にシステム メモリがすべて消費されるものと推測されます。その推測が正しければ、コンピューターのメモリを増設することによって、より大きな問題が解決でき、また行列のサイズが比較的大きい場合のパフォーマンスを向上させることができます。

高速化

バックスラッシュなどの線形代数アルゴリズムによって実現される高速化を測定するために従来から使用されている方法は、ピーク パフォーマンスの比較です。そこで、それぞれの数のワーカーについてギガフロップスの最大数を求めてみることにします。

peakPerf = [max(s.workers4.gflops), max(s.workers8.gflops), ...
            max(s.workers16.gflops), max(s.workers32.gflops), ...
            max(s.workers64.gflops)];
disp('Peak performance in Gigaflops for 4-64 workers:')
disp(peakPerf)

disp('Speedup when going from 4 workers to 8, 16, 32 and 64 workers:')
disp(peakPerf(2:end)/peakPerf(1))
Peak performance in Gigaflops for 4-64 workers:
   10.9319   23.2508   40.7157   73.5109  147.0693

Speedup when going from 4 workers to 8, 16, 32 and 64 workers:
    2.1269    3.7245    6.7244   13.4532

その結果、ワーカー数を 4 から 64 へと 16 倍に増加させると約 13.5 倍の高速化が実現されることがわかります。既に述べたように、パフォーマンスのグラフから、64 個のワーカーの場合にはクラスター コンピューターのシステム メモリを増設することによってパフォーマンスを向上 (ひいては、さらなる高速化を実現) できるものと推測されます。

使用したクラスター

このデータの生成には、デュアル プロセッサとオクタコアを装備した、64 GB のメモリをもつコンピューター 16 台をギガビット イーサネットで接続して使用しました。4 個のワーカーの場合、すべてのワーカーを 1 台のコンピューターで実行しました。また、8 個のワーカーの場合は 2 台のコンピューター、16 個のワーカーの場合は 4 台のコンピューターというような設定で行いました。

デッドロック検出の再有効化

これでベンチマークが完了したので、現在の並列プールでデッドロック検出を再び有効にしても問題ありません。

pctRunOnAll 'mpiSettings(''DeadlockDetection'', ''on'');'
end
ans = 

  struct with fields:

    matSize: [1000 19787 38573 57360 76146]
     gflops: [1.1618 43.7417 69.6855 40.4145 36.9379]