Преобразование явного эйлера в неявный эйлер (с помощью итерации с фиксированной запятой) - PullRequest
0 голосов
/ 08 февраля 2019

Таким образом, у меня есть школьное задание, в котором мне нужно вычислить положение группы автомобилей, следующих друг за другом на дороге (АКА едет по линии, поэтому, если автомобиль 10 [автомобиль, который находится в очереди) тормозит, тоАвтомобиль 9 тормозит, а когда автомобиль 9 тормозит, автомобиль 8 должен тормозить и т. Д.).

Каждый автомобиль следует «правилу трех секунд» (кроме автомобиля, который стоит первым в очереди, он может ехать так же быстро, как и онхочет и все другие машины в очереди регулируют свою скорость соответственно).Скорость каждого автомобиля представлена ​​выражением:

enter image description here

, где «i» - это индекс автомобиля, а «t» - это момент времени (автомобиль с наивысшим индексом - это автомобиль, который является первым в очереди), а функция 'f' представлена ​​этим кодом:

% Input x: Distance to the car in front (of you. NOT the car that is at the
% very front of the line).
% Output: Car velocity
function output = carFunc(x)
    if (x >= 75)
        output = 25;
    else
        output = x/3;
    end
end

Автомобиль, который находится на переднем крае, просто имеет постоянную скорость 'g ':

enter image description here

В любом случае , теперь, когда вы знаете контекст задачи, давайте фактически перейдем к самой задаче.Эта школьная задача состоит из нескольких этапов, и первым шагом является вычисление положения каждого автомобиля с течением времени с использованием прямого / явного Эйлера, что я и сделал.Вот код для этого:

% Euler's Method
% Initial conditions and setup
g = 25; % Velocity of car M (the car that is first in line, ahead of all others) [m/s]
h = 0.01;  % step size
M = 10; % Number of cars
x = 0:h:40;  % the range of x (time)
n = numel(x);  % the number of timesteps
posMat = zeros(M,n); % Matrix that holds positions of each car (car = row, time = column)
for i=1:M
    posMat(i,1) = 10*i;  % the initial positions for the cars
end
for t=1:(n-1)
    % Calculate position of all cars EXCEPT car M (using the first
    % equation)
    for c=1:(M-1)
    f = carFunc(posMat(c+1,t) - posMat(c,t)); % Velocity of car 'c'
    posMat(c,t+1) = posMat(c,t) + h * f; % Explicit Euler
    end
    % Calculate positon of last car M (first car in line) here:
    % x_M^(n+1) = x_M^n + h*g
    posMat(M,t+1) = posMat(M,t) + h*g;
end

Однако проблема (в которой я застрял) - это второй шаг, где мне нужно теперь изменить свой код на первом шаге, чтобы использовать обратный / неявный эйлерчерез итерацию с фиксированной запятой. Я написал функцию, которая выполняет итерацию с фиксированной запятой, но кроме этого я не знаю, что делать.Вот мой код для итерации с фиксированной запятой:

%Fixed-Point Iteration
%Computes approximate solution of g(x)=x
%Input: function handle g, starting guess x0,
% number of iteration steps k
%Output: Approximate solution

function out=fpi(g, x0, k)
    x = zeros(1, k+1);
    x(1)=x0;
        for i=1:k
            x(i+1)=g(x(i));
        end
    out=x(k+1);
end

Любая помощь приветствуется.Извините за длинный текст.Верхняя часть моего поста - это, в основном, «краткое» описание контекста задачи.Это не является абсолютно необходимым (и здесь это не главное), но я все равно добавил его, чтобы вы, ребята, знали, что происходит в моем коде.

Спасибо!

1 Ответ

0 голосов
/ 09 февраля 2019

Проблема в том, что вы написали свою итерацию с фиксированной точкой для скаляров , но у вас есть система дифференциальных уравнений.Если переписать систему в векторном виде, она станет намного понятнее.Я добавил комментарий о том, как будет выглядеть явное и неявное уравнение таким образом, чтобы они действительно имели ту же стандартную форму y(t+1) = y(t) + h * f(y(t),t) (соответственно y(t+1) = y(t) + h * f(y(t+1),t+1)), которую вы можете найти в текстовой книге.Тогда вы можете легко записать обновления:

% Euler's Method
% Initial conditions and setup
g = 25; % Velocity of car M (the car that is first in line, ahead of all others) [m/s]
h = 0.01;  % step size
M = 10; % Number of cars
x = 0:h:40;  % the range of x (time)
n = numel(x);  % the number of timesteps
posMat = zeros(M,n); % Matrix that holds positions of each car (car = row, time = column)
k = 20; % number of fixed point iterations
explicit = true;
for i=1:M
    posMat(i,1) = 10*i;  % the initial positions for the cars
end

for t=1:n-1 
    % Calculate position of all cars EXCEPT car M (using the first
    % equation)
    c=1:M-1;
    if explicit
        f = [carFunc(posMat(c+1,t) - posMat(c,t)); g]; % Velocity of car 'c'
        posMat(:,t+1) = posMat(:,t) + h * f; % Explicit Euler
    else 
        %explicit euler:
        %posMat(:,t+1) = posMat(:,t) + h * [carFunc(posMat(c+1,t) - posMat(c,t)); g];
        %implicit euler: (needs to be solved)
        %posMat(:,t+1) = posMat(:,t) + h * [carFunc(posMat(c+1,t+1) - posMat(c,t+1)); g];

        %fixed point iteration:
        posMat(:,t+1) = posMat(:,t); % initialization
        for m=1:k
            posMat(:,t+1) = posMat(:,t) + h * [carFunc(posMat(c+1,t+1) - posMat(c,t+1)); g]; 
        end
    end
end
plot(posMat');


% Input x: Distance to the car in front (of you. NOT the car that is at the
% very front of the line).
% Output: Car velocity
function output = carFunc(x)
    mask = x >= 75;
    output = zeros(size(x));
    output(mask) = 25;
    output(~mask) = x(~mask)/3;
end
...