Частные производные первого и второго порядка функции с двумя переменными в MATLAB - PullRequest
Новая
       61

Частные производные первого и второго порядка функции с двумя переменными в MATLAB

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

Я пытаюсь реализовать различные численные методы в MATLAB без использования встроенной функции, такой как градиент или del2.Пока это мой код: 

clear all
close all
x = [-1:0.1:1];
y = [-2:0.1:2];
vel = @(x,y) x+exp(-((x-x(1)).^2+(y-y(1)).^2));
nx = length(x);
ny = length(y);
derivx = zeros(nx-1,ny-1)
% The partial derivative with respect to x
for ii = 1:nx-1
    for jj = 1:ny-1
        derivx(ii,jj) = (vel(ii+1,jj) - vel(ii,jj))./(x(jj+1,ii)-x(jj,ii));
    end
end
% The partial with respect to y
derivy = zeros(ny-1,nx-1)
for ii = 1:ny-1
    for jj = 1:nx-1
    derivy(ii,jj) = (vel(ii+1,jj) - vel(ii,jj))./(y(jj+1,ii)-y(jj,ii));
    end
end

Этот код не работает с сообщением об ошибке, в котором говорится, что индексы матрицы превышены.

Index in position 1 exceeds array bounds (must not exceed 1).
Error in untitled6 (line 13)
    derivx(ii,jj) = (vel(ii+1,jj) - vel(ii,jj))./(x(jj+1,ii)-x(jj,ii));

И как мне рассчитать партиалы второго порядка с повторением по x и y (не смешанным)?

Заранее спасибо за помощь!

1 Ответ

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

Чтобы ответить на вопрос, который вы задали, проблема заключается в следующем: (x(jj+1,ii)-x(jj,ii).x - это вектор, но вы рассматриваете его как матрицу.Тем не менее, я думаю, что в вашем коде есть более глубокие проблемы.Например, то, как вы относитесь к vel, довольно необычно.Вы написали как функцию x и y, где x и y - предположительно векторы (или матрицы), но вы когда-либо вызывали его только со скалярами.Если бы мне нужно было сделать предположение, я бы предположил, что вы хотите написать vel как: 

x = [-1:0.1:1];
y = [-2:0.1:2]'; % Note the transpose here
vel = x+exp(-((x-x(1)).^2+(y-y(1)).^2));

. Это создаст vel как 2D матрицу, где каждый элемент (a, b) являетсязначение vel оценивается при x=a и y=b.После того, как вы это сделали, вы можете покончить с дважды вложенными циклами for (которые почти никогда не являются хорошей идеей в MATLAB): 

derivx = (vel(2:end,1:end-1) - vel(1:end-1,1:end-1)./(x(2:end)-x(1:end-1));
...