このページの内容は最新ではありません。最新版の英語を参照するには、ここをクリックします。

parfeval を使用したパラメーター スイープ中のプロット

この例では、並列パラメーター スイープを parfeval により実行し、その結果を計算中に DataQueue オブジェクトによって戻す方法を示します。

parfeval は MATLAB をブロックしないため、計算の実行中に作業を続行できます。

この例では、ローレンツ常微分方程式系のパラメーター $\sigma$ および $\rho$ に対してパラメーター スイープを実行し、この系のカオス的性質を説明します。

$$\begin{array}{l}\frac{\mathrm{d}}{\mathrm{d}\mathit{t}}\mathit{x}=\sigma \left(\mathit{y}-\mathit{z}\right)\\ \frac{\mathrm{d}}{\mathrm{d}\mathit{t}}\mathit{y}=\mathit{x}\left(\rho -\mathit{z}\right)-\mathit{y}\\ \frac{\mathrm{d}}{\mathrm{d}\mathit{t}}\mathit{z}=\mathit{xy}-\beta \mathit{x}\end{array}$$

パラメーター グリッドの作成

パラメーター スイープで調べるパラメーターの範囲を定義します。

gridSize = 40;
sigma = linspace(5, 45, gridSize);
rho = linspace(50, 100, gridSize);
beta = 8/3;

関数 meshgrid を使用して、パラメーターの 2 次元グリッドを作成します。

[rho,sigma] = meshgrid(rho,sigma);

figure オブジェクトを作成し、'Visible' を true に設定すると、このオブジェクトがライブ スクリプトの外側の新しいウィンドウで開きます。パラメーター スイープの結果を可視化するには、表面プロットを作成します。表面の Z 要素を NaN で初期化すると、空のプロットが作成されることに注意してください。

figure('Visible',true);
surface = surf(rho,sigma,NaN(size(sigma)));
xlabel('\rho','Interpreter','Tex')
ylabel('\sigma','Interpreter','Tex')

並列環境の設定

関数 parpool を使用して並列ワーカーのプールを作成します。

parpool;

Starting parallel pool (parpool) using the 'Processes' profile ...
Connected to the parallel pool (number of workers: 6).

ワーカーからデータを送信するには、DataQueue オブジェクトを作成します。関数 afterEach を使用して、ワーカーがデータを送信するたびに表面プロットを更新する関数を設定します。関数 updatePlot は、この例の最後で定義するサポート関数です。

Q = parallel.pool.DataQueue;
afterEach(Q,@(data) updatePlot(surface,data));

並列パラメーター スイープの実行

パラメーターを定義した後、並列パラメーター スイープを実行できます。

作業負荷を分散すると、parfeval の効率性が向上します。作業負荷を分散するには、調べるパラメーターをグループ化して分割します。この例では、コロン演算子 (:) を使用して、サイズを step として均等に分割します。この結果得られる配列 partitions には、分割の境界が含まれます。最後の分割の終点を追加しなければならないことに注意してください。

step = 100;
partitions = [1:step:numel(sigma), numel(sigma)+1]

partitions = 1×17

           1         101         201         301         401         501         601         701         801         901        1001        1101        1201        1301        1401        1501        1601

最良のパフォーマンスを得るには、次のように分割するようにします。

• 分割のスケジューリングのオーバーヘッドよりも計算時間が長くなる程度に大きい

• すべてのワーカーをビジー状態に維持するために十分な分割数が存在する程度に小さい

並列ワーカーでの関数の実行を表し、その結果を保持するには、future オブジェクトを使用します。

f(1:numel(partitions)-1) = parallel.FevalFuture;

関数 parfeval を使用して、計算を並列ワーカーにオフロードします。parameterSweep はこのスクリプトの最後で定義されている補助関数で、調べるパラメーターの分割についてのローレンツ系を解きます。これには出力引数が 1 つあるため、parfeval では出力の数として 1 を指定しなければなりません。

for ii = 1:numel(partitions)-1
    f(ii) = parfeval(@parameterSweep,1,partitions(ii),partitions(ii+1),sigma,rho,beta,Q);
end

parfeval は MATLAB をブロックしないため、計算の実行中に作業を続行できます。ワーカーは並列で計算を実行し、中間結果を使用できるようになったら DataQueue によって送信します。

parfeval が完了するまで MATLAB をブロックする場合は、future オブジェクトに対して関数 wait を使用します。この後のコードが parfeval の完了に依存する場合、関数 wait を使用すると便利です。

wait(f);

parfeval が計算を完了すると、wait は終了し、さらにコードを実行できるようになります。たとえば、結果として生成される表面の輪郭をプロットします。関数 fetchOutputs を使用して、future オブジェクトに保存された結果を取得します。

results = reshape(fetchOutputs(f),gridSize,[]);
contourf(rho,sigma,results)
xlabel('\rho','Interpreter','Tex')
ylabel('\sigma','Interpreter','Tex')

パラメーター スイープで計算リソースが多く必要で、クラスターにアクセス可能な場合は、parfeval の計算をスケール アップできます。詳細については、デスクトップからクラスターへのスケール アップを参照してください。

補助関数の定義

調べるパラメーターの分割のローレンツ系を解く補助関数を定義します。DataQueue オブジェクトで関数 send を使用して、中間結果を MATLAB クライアントに送信します。

function results = parameterSweep(first,last,sigma,rho,beta,Q)
    results = zeros(last-first,1);
    for ii = first:last-1
        lorenzSystem = @(t,a) [sigma(ii)*(a(2) - a(1)); a(1)*(rho(ii) - a(3)) - a(2); a(1)*a(2) - beta*a(3)];
        [t,a] = ode45(lorenzSystem,[0 100],[1 1 1]);
        result = a(end,3);
        send(Q,[ii,result]);
        results(ii-first+1) = result;
    end
end

新しいデータを受け取った時点で表面プロットを更新する別の補助関数を定義します。

function updatePlot(surface,data)
    surface.ZData(data(1)) = data(2);
    drawnow('limitrate');
end

参考

parpool | parallel.pool.DataQueue | afterEach | parfeval

関連するトピック

• parfor を使用したパラメーター スイープ中のプロット
• デスクトップからクラスターへのスケール アップ