Четырехкратный интеграл с использованием вложенного интеграла 2 в Matlab

Question

Я пытаюсь решить проблему следующего вида:

f=@(x,y,z,w) x.*y.*z.*w;          % A complicated black box function
a=1;b=1;c=1;d=1;                  % Integration limits
I=integral2(@(x,y)integral2(@(z,w)f(x,y,z,w),c,-c,d,-d),a,-a,b,-b);

Используя эту реализацию, я получаю следующую ошибку:

Error using  .* 
Matrix dimensions must agree.

Проблема в том, что x, y, z и w не имеют одинаковый размер. Для оценки первой функции все входные данные имеют одинаковый размер, но затем для оценки второй функции x и y не совпадают с размерами z и w.

Как я могу устранить эту ошибку?

Этот вопрос похож на этот оставшийся без ответа вопрос: ошибка размера входного массива для четырехкратного интегрирования с использованием вложенного интеграла2

================================================== ================================

В ответ на ответ:

I=integral(@(x)integral3(@(y,z,w)f(x,y,z,w),b,-b,c,-c,d,-d),a,-a,'ArrayValued',true);

Это решает проблему, однако, для меня не очевидно, почему это работает. Я действительно видел это решение раньше, но забыл упомянуть его в своем вопросе ( http://www.mathworks.com/matlabcentral/answers/77571-how-to-perform-4d-integral-in-matlab).

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

I=integral(@(x)integral3(@(y,z,w)f(x,y,z,w),b,-b,c,-c,d,-d,'Method','iterated'),a,-a,'ArrayValued',true);

Answer 1

integral2 вызывает его подынтегральное выражение с двумя матричными аргументами одинакового размера. Проблема в том, что вы не можете просто смешать переменные x,y,z,w в вызове функции f(x,y,z,w), как вы делали это в своем вопросе, потому что размерность x и y определяется внешний integral2тогда как размерность z и w определяется внутренним integral2, поэтому не гарантируется, что размеры одинаковы. integralX функции в любом случае не векторизованы, каждый вызов может предоставить только одно выходное значение.

Функция integral предоставляет опцию, которая вызывает его подынтегральное выражение только со скалярными значениями, а со скалярным расширением Matlab это работает вместе с трехмерными массивами того же размера, что и внутренняя функция integral3 обеспечивает.

I=integral(@(x)integral3(@(y,z,w)f(x,y,z,w),b,-b,c,-c,d,-d),a,-a,'ArrayValued',true);

Вы можете добиться того же (вызывая со скалярами), инкапсулируя внутренний вызов integra2 с помощью arrayfun:

I=integral2(@(x,y)arrayfun(@(x,y)integral2(@(z,w)f(x,y,z,w),c1,c2,d1,d2),x,y),a1,a2,b1,b2)

Последний был примерно в шесть раз быстрее в моих экспериментах.

Answer 2

Чтобы добавить ответ Даниэля, я, наконец, сломался и написал представление об обмене файлами в MATLAB, которое автоматически выполняет вложение для вас, и может обрабатывать более утомительные случаи, когда ограничения являются функциями, а не константами. Он называется интегральным N и написан для обработки 4-кратных, 5-кратных и 6-кратных интегралов. Я не думаю, что это отличный способ атаковать многомерную интеграцию, но если он оказывается достаточно быстрым, чтобы получить работу для вас, тогда его можно использовать.