Код, показанный ниже, доступен здесь: 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, это может быть быстрее или медленнее. Это всегда быстрее, чем версия [...] Table.
Вы могли бы даже написать это, используя SparseArray:
In: AbsoluteTiming[SparseArray[{i_, i_} /; i <= k -> 1, {2 L, 2 L}];]
Out: {0.0040002, Null}
Я мог бы поговорить о некоторых других вещах, но боюсь, если я сделаю это, я сделаю этот ответ неоправданно большим Я накопил ряд методов для формирования этих различных матриц и списков за время, потраченное на оптимизацию некоторого кода. Базовый код, с которым я работал, занял более 6 дней, чтобы выполнить одно вычисление, и теперь требуется всего 6 часов, чтобы сделать то же самое.
Я посмотрю, смогу ли я выбрать общие методы, которые я придумаю, и просто сунуть их в тетрадь для использования.
TL; DR: Похоже, что в этих случаях функциональный путь превосходит процедурный. Но при компиляции процедурный код превосходит функциональный код.