Ответ 1
Ниже приведен код, приведенный ниже: http://pastebin.com/4PWWxGhB. Просто скопируйте и вставьте его в блокнот, чтобы проверить его.
Я действительно пытался сделать несколько функциональных способов вычисления матриц, так как я что функциональный способ (который обычно идиоматичен в Mathematica) более эффективен.
В качестве одного примера у меня была эта матрица, состоящая из двух списков:
In: L = 1200;
e = Table[..., {2L}];
f = Table[..., {2L}];
h = Table[0, {2L}, {2L}];
Do[h[[i, i]] = e[[i]], {i, 1, L}];
Do[h[[i, i]] = e[[i-L]], {i, L+1, 2L}];
Do[h[[i, j]] = f[[i]]f[[j-L]], {i, 1, L}, {j, L+1, 2L}];
Do[h[[i, j]] = h[[j, i]], {i, 1, 2 L}, {j, 1, i}];
Моим первым шагом было все время.
In: h = Table[0, {2 L}, {2 L}];
AbsoluteTiming[Do[h[[i, i]] = e[[i]], {i, 1, L}];]
AbsoluteTiming[Do[h[[i, i]] = e[[i - L]], {i, L + 1, 2 L}];]
AbsoluteTiming[
Do[h[[i, j]] = f[[i]] f[[j - L]], {i, 1, L}, {j, L + 1, 2 L}];]
AbsoluteTiming[Do[h[[i, j]] = h[[j, i]], {i, 1, 2 L}, {j, 1, i}];]
Out: {0.0020001, Null}
{0.0030002, Null}
{5.0012861, Null}
{4.0622324, Null}
DiagonalMatrix[...]
был медленнее, чем циклы do, поэтому я решил просто использовать циклы Do
на последнем шаге. Как вы можете видеть, использование Outer[Times, f, f]
в этом случае было намного быстрее.
Затем я написал эквивалент, используя Outer
для блоков в верхнем правом и нижнем левом углу матрицы, и DiagonalMatrix
для диагонали:
AbsoluteTiming[h1 = ArrayPad[Outer[Times, f, f], {{0, L}, {L, 0}}];]
AbsoluteTiming[h1 += Transpose[h1];]
AbsoluteTiming[h1 += DiagonalMatrix[Join[e, e]];]
Out: {0.9960570, Null}
{0.3770216, Null}
{0.0160009, Null}
DiagonalMatrix
был фактически медленнее. Я мог бы заменить его только контурами Do
, но я сохранил его, потому что он выглядел чище.
Текущая цифра составляет 9,06 секунды для наивного цикла Do
и 1,389 секунды для моей следующей версии с использованием Outer
и DiagonalMatrix
. Примерно в 6,5 раз быстрее, не так уж плохо.
Звучит намного быстрее, не так ли? Попробуйте использовать Compile
сейчас.
In: cf = Compile[{{L, _Integer}, {e, _Real, 1}, {f, _Real, 1}},
Module[{h},
h = Table[0.0, {2 L}, {2 L}];
Do[h[[i, i]] = e[[i]], {i, 1, L}];
Do[h[[i, i]] = e[[i - L]], {i, L + 1, 2 L}];
Do[h[[i, j]] = f[[i]] f[[j - L]], {i, 1, L}, {j, L + 1, 2 L}];
Do[h[[i, j]] = h[[j, i]], {i, 1, 2 L}, {j, 1, i}];
h]];
AbsoluteTiming[cf[L, e, f];]
Out: {0.3940225, Null}
Теперь он работает в 3,56 раза быстрее, чем моя последняя версия, и в 23,23 раза быстрее, чем первая. Следующая версия:
In: cf = Compile[{{L, _Integer}, {e, _Real, 1}, {f, _Real, 1}},
Module[{h},
h = Table[0.0, {2 L}, {2 L}];
Do[h[[i, i]] = e[[i]], {i, 1, L}];
Do[h[[i, i]] = e[[i - L]], {i, L + 1, 2 L}];
Do[h[[i, j]] = f[[i]] f[[j - L]], {i, 1, L}, {j, L + 1, 2 L}];
Do[h[[i, j]] = h[[j, i]], {i, 1, 2 L}, {j, 1, i}];
h], CompilationTarget->"C", RuntimeOptions->"Speed"];
AbsoluteTiming[cf[L, e, f];]
Out: {0.1370079, Null}
Большая часть скорости поступала от CompilationTarget->"C"
. Здесь я получил еще 2,84 ускорения по сравнению с самой быстрой версией и в 66,13 раза быстрее, чем первая версия. Но все, что я сделал, это просто скомпилировать его!
Теперь это очень простой пример. Но это реальный код, который я использую для решения проблемы физики конденсированных сред. Поэтому не отвергайте его как возможно "игрушечный пример".
Как насчет другого примера метода, который мы можем использовать? У меня есть относительно простая матрица, которую я должен создать. У меня есть матрица, состоящая из ничего, кроме только от начала до некоторой произвольной точки. Наивный способ может выглядеть примерно так:
In: k = L;
AbsoluteTiming[p = Table[If[i == j && j <= k, 1, 0], {i, 2L}, {j, 2L}];]
Out: {5.5393168, Null}
Вместо этого создайте его с помощью ArrayPad
и IdentityMatrix
:
In: AbsoluteTiming[ArrayPad[IdentityMatrix[k], {{0, 2L-k}, {0, 2L-k}}
Out: {0.0140008, Null}
Это действительно не работает для k = 0, но вы можете использовать специальный случай, если вам это нужно. Кроме того, в зависимости от размера k это может быть быстрее или медленнее. Это всегда быстрее, чем таблица [...], хотя.
Вы даже можете записать это с помощью SparseArray
:
In: AbsoluteTiming[SparseArray[{i_, i_} /; i <= k -> 1, {2 L, 2 L}];]
Out: {0.0040002, Null}
Я мог бы продолжить некоторые другие вещи, но я боюсь, если я это сделаю, я сделаю этот ответ необоснованно большим. Я накопил ряд методов для создания этих различных матриц и списков за время, потраченное на оптимизацию кода. Базовый код, с которым я работал, занял более 6 дней, чтобы выполнить один расчет, и теперь требуется всего 6 часов, чтобы сделать то же самое.
Я посмотрю, смогу ли я выбрать общие методы, которые я придумал, и просто вставлять их в блокнот для использования.
TL; DR: Кажется, что в этих случаях функциональный путь превосходит процедурный путь. Но при компиляции процедурный код превосходит функциональный код.