Комментарии 21
Не совсем понял, как изложенное в статье объясняет вывод
Я просто высказался, что NumPy — это худший язык для работы с массивами, кроме всех других языков работы с массивами.
Пересаживался с MATLAB на Python и NumPy тогда оказался спасением. На мой взгляд работа с многомерными таблицами - это не та область, где всегда можно обеспечить "очевидность" происходящего.
Когда я читаю какой-то код, чтобы было «очевидно», что он делает.
Вот, кстати, пример кода на MATLAB (кусок численного решения 3D задачи со временем, то есть идет работа с таблицами с 4 измерениями). Тяжело читать? Да, очень. Написать при этом было не очень сложно и сильно пересекается с тем, что было написано на бумаге.
Скрытый текст
fxn = cross(p(:,2:N+1,1:N,:) - p(:,1:N,2:N+1,:), p(:,2:N+1,2:N+1,:) - p(:,1:N,1:N,:),4)/2;
fyn = cross(p(1:N,:,2:N+1,:) - p(2:N+1,:,1:N,:), p(2:N+1,:,2:N+1,:) - p(1:N,:,1:N,:),4)/2;
fzn = cross(p(2:N+1,1:N,:,:) - p(1:N,2:N+1,:,:), p(2:N+1,2:N+1,:,:) - p(1:N,1:N,:,:),4)/2;
% центры граней
fxc = (p(:,1:N,1:N,:) + p(:,1:N,2:N+1,:) + p(:,2:N+1,2:N+1,:) + p(:,2:N+1,1:N,:))/4;
fyc = (p(1:N,:,1:N,:) + p(2:N+1,:,1:N,:) + p(2:N+1,:,2:N+1,:) + p(1:N,:,2:N+1,:))/4;
fzc = (p(1:N,1:N,:,:) + p(1:N,2:N+1,:,:) + p(2:N+1,2:N+1,:,:) + p(2:N+1,1:N,:,:))/4;
% центры ячеек
C = ( p(1:N,1:N,1:N,:) + p(1:N,1:N,2:N+1,:) + p(1:N,2:N+1,1:N,:) + p(1:N,2:N+1,2:N+1,:) + ...
p(2:N+1,1:N,1:N,:) + p(2:N+1,1:N,2:N+1,:) + p(2:N+1,2:N+1,1:N,:) + p(2:N+1,2:N+1,2:N+1,:) )/8;
% объемы ячеек
V = 1/3 * ( - dot(fxn(1:N,:,:,:),fxc(1:N,:,:,:),4) + dot(fxn(2:N+1,:,:,:),fxc(2:N+1,:,:,:),4) ...
- dot(fyn(:,1:N,:,:),fyc(:,1:N,:,:),4) + dot(fyn(:,2:N+1,:,:),fyc(:,2:N+1,:,:),4) ...
- dot(fzn(:,:,1:N,:),fzc(:,:,1:N,:),4) + dot(fzn(:,:,2:N+1,:),fzc(:,:,2:N+1,:),4));
% объемы октаэдров
fxv = 1/3 * ( abs(dot(C(2:N,:,:,:) - p(2:N,1:N,1:N,:),fxn(2:N,:,:,:),4)) + ...
abs(dot(C(1:N-1,:,:,:) - p(2:N,1:N,1:N,:),fxn(2:N,:,:,:),4)) );
fyv = 1/3 * ( abs(dot(C(:,2:N,:,:) - p(1:N,2:N,1:N,:),fyn(:,2:N,:,:),4)) + ...
abs(dot(C(:,1:N-1,:,:) - p(1:N,2:N,1:N,:),fyn(:,2:N,:,:),4)) );
fzv = 1/3 * ( abs(dot(C(:,:,2:N,:) - p(1:N,1:N,2:N,:),fzn(:,:,2:N,:),4)) + ...
abs(dot(C(:,:,1:N-1,:) - p(1:N,1:N,2:N,:),fzn(:,:,2:N,:),4)) );
% коэффициенты alpha
fxa = alpha( x(2:N,1:N,1:N), y(2:N,1:N,1:N), z(2:N,1:N,1:N), ...
x(2:N,1:N,2:N+1), y(2:N,1:N,2:N+1), z(2:N,1:N,2:N+1), ...
x(2:N,2:N+1,2:N+1), y(2:N,2:N+1,2:N+1), z(2:N,2:N+1,2:N+1), ...
x(2:N,2:N+1,1:N), y(2:N,2:N+1,1:N), z(2:N,2:N+1,1:N), ...
C(2:N,:,:,1), C(2:N,:,:,2), C(2:N,:,:,3), ...
C(1:N-1,:,:,1), C(1:N-1,:,:,2), C(1:N-1,:,:,3) );
fya = alpha( x(1:N,2:N,1:N), y(1:N,2:N,1:N), z(1:N,2:N,1:N), ...
x(2:N+1,2:N,1:N), y(2:N+1,2:N,1:N), z(2:N+1,2:N,1:N), ...
x(2:N+1,2:N,2:N+1), y(2:N+1,2:N,2:N+1), z(2:N+1,2:N,2:N+1), ...
x(1:N,2:N,2:N+1), y(1:N,2:N,2:N+1), z(1:N,2:N,2:N+1), ...
C(:,2:N,:,1), C(:,2:N,:,2), C(:,2:N,:,3), ...
C(:,1:N-1,:,1), C(:,1:N-1,:,2), C(:,1:N-1,:,3));
fza = alpha( x(1:N,1:N,2:N), y(1:N,1:N,2:N), z(1:N,1:N,2:N), ...
x(1:N,2:N+1,2:N), y(1:N,2:N+1,2:N), z(1:N,2:N+1,2:N), ...
x(2:N+1,2:N+1,2:N), y(2:N+1,2:N+1,2:N), z(2:N+1,2:N+1,2:N), ...
x(2:N+1,1:N,2:N), y(2:N+1,1:N,2:N), z(2:N+1,1:N,2:N), ...
C(:,:,2:N,1), C(:,:,2:N,2), C(:,:,2:N,3), ...
C(:,:,1:N-1,1), C(:,:,1:N-1,2), C(:,:,1:N-1,3));
можно любить или не любить numpy, но есть ли, чем его заменить? (когда-то давно, для упрощения кода и ускорения вычислений, использовал pdl, но это, простите, perl, и, вероятно не подойдет)
Pandas не лучше?
taichi. Для многомерных массивов ничего лучше не видел. taichi.Vector.field самое гениальное, что можно было придумать, чтобы создавать массивы любых размерностей в каждом элементе которого могут быть другие массивы любых размерностей с полной читаемостью происходящего. Да и cuda из коробки. И киллер фича - конвертация собственных массивов в numpy.
Так и не надо любить numpy.
Надо любить JAX.
Там и нампай, и автодифференцирование и JIT-компиляция хоть под GPU, хоть под TPU, хоть под CPU.
Автоматический батчинг через vmap - просто пишешь код для расчета одного элемента и vmap(foo) возвращает функцию, работающую с целым батчем.
Внезапно оси (axis) превратились в топоры (axes).
Автор, может напишешь свою мини-библиотеку? Возьмёшь ключевые функции и забудешь NumPy как страшный сон. Может ещё зальешь на Git, авось взлетит.
Так делают в других языках, и получается уж как минимум гораздо более внятно. Нативная компиляция в CUDA должна просто происходить отдельным этапом, вот и весь секрет.
Какое Ваше мнение о вот этом? https://cupy.dev/
Понравилась статья.
Пример из области вычматов: с помощью Numpy нельзя запрограммировать численное решение СЛАУ методом Гаусса-Зейделя (метод отчасти учебный).
Всегда удивляло, что на HN (обсуждение там) в таких случаях вытаскивают из шкафа язык APL.
NumPy для массивов как Pandas для табличек в Питоне, все ужасно, но уже поздно что-то делать. Гвидо не обратил внимания на эти пороли боабобов вовремя и вот результат
Я не люблю NumPy