Векторизация Matlab с циклами if и for, использующими ODE45 для интеграции - PullRequest
Новая
       47

Векторизация Matlab с циклами if и for, использующими ODE45 для интеграции

0 голосов
/ 06 января 2020

Я заинтересован в оптимизации скорости моего кода, и при использовании «Run and Time» - это функция в моем коде, которая сильно влияет на скорость, но мне трудно понять, как правильно векторизовать эту функцию, как я обычно просто сделать зацикливание, в попытке я наткнулся на ошибку, так как она также используется в интеграции, моя первоначальная функция выглядит следующим образом и не приводит к ошибке 

function [dotStates] = ODEFunc(t,states,params)
%ODE function
% Loading in and assigning the variables from parameters
K = params(1);                                                             
N = params(2);
nn = params(3);
% magnitude of the coupling based on the number of neighbours
kn = K/nn;
w = params(4:end);
dotStates=states;
    % For each oscillator
    for i=1:N
        % Use the oscillators natural frequency
        dotStates(i) = w(i);
        % For j number of neighbours
        for j=(i-nn):(i+nn)
            % neighbour number is positive and shorter than # of oscilators
            if (j > 0) && (j < length(dotStates))
                dotStates(i) =  dotStates(i) + (kn * sin( states(j)-states(i) ));
            end
        end
    end
end

Я пытался следуя руководству по векторизации mathworks: https://se.mathworks.com/help/matlab/matlab_prog/vectorization.html

До сих пор я пытался использовать некоторые из входных данных, которые они используют, такие как использование маски, и сгенерировал следующий код 

function [dotStates] = ODEFunc(t,states,params)
%ODE function
% Loading in and assigning the variables from parameters
K = params(1);                                                             
N = params(2);
nn = params(3);
% magnitude of the coupling based on the number of neighbours
kn = K/nn;
w = params(4:end);
dotStates=states;
% Use the oscillators natural frequency
dotStates = w';
% Mask of j states
j = (i-nn):(i+nn);
% neighbours cannot exceed boundaries
j = j(j>0 & j <=length(dotStates));
jstate = states(j);
jstate(numel(states)) = 0;
dotStates =  dotStates + (kn * sin( jstate'-states ));
end

Я получил вектор, который короче того, что записывается, и я решил просто добавить несколько нулей в переменную "jstate", чтобы компенсировать разницу, но это не так чувствую себя как правильная векторизация и когда я запускаю код, я получаю следующую ошибку, которая связана с шагом интеграции

Wa rning: операнды двоеточия должны быть настоящими скалярами

В RK_ODE_2411> ODEFun c (строка 99) В RK_ODE_2411> @ (t, состояния) ODEFun c (t, состояния, параметры) В ode45 (строка 324) В RK_ODE_2411 (строка 58)

функция в свою очередь используется в следующем сегменте для интеграции с использованием ODE45 

%% Integration via ODE45
for K = 0:.1:Klen
    params(1) = K;
    K_count = K_count+1;
    nn_count = 0;
    for nn = nnlen:nnlen
        params(3) = nn;
        % index counter
        nn_count = nn_count+1;
        % 6th order runge kutta

        sol(K_count,nn_count) = ode45(@(t,states) ODEFunc(t,states,params),tSpan,init,options);
    end
end

, где строка 58 равна 

sol(K_count,nn_count) = ode45(@(t,states) ODEFunc(t,states,params),tSpan,init,options);

EDIT: строка 99 in ODEFun c is 

j = (i-nn):(i+nn);

1 Ответ

1 голос
/ 07 января 2020

Попробуйте этот фрагмент 

% For each oscillator
for i=1:N
    % For j number of neighbours
    j=(i-nn):(i+nn);

    % neighbour number is positive and shorter than # of oscilators
    lg = (j > 0) & (j < length(dotStates));
    dotStates(i) =  w(i) + sum(kn * sin( states(lg)-states(i) ));
end

, но самое важное, что dotStates не будет больше, чем stats, так как это заставит matlab переставить свою память, что значительно замедлит код .

...