Рассмотрим следующий простой, иллюстрирующий пример
cf = Block[{a, x, degree = 3},
With[{expr = Product[x - a[[i]], {i, degree}]},
Compile[{{x, _Real, 0}, {a, _Real, 1}}, expr]
]
]
Это один из возможных способов передачи кода в теле оператора Compile. Он выдает ошибку Part :: partd, так как [[i]] на момент оценки не является списком.
Самое простое решение - просто проигнорировать это сообщение или выключить его. Конечно, есть и другие способы. Например, можно обойти оценку [[i]], заменив его внутри тела компиляции до его компиляции
cf = ReleaseHold[Block[{a, x, degree = 3},
With[{expr = Product[x - a[i], {i, degree}]},
Hold[Compile[{{x, _Real, 0}, {a, _Real, 1}}, expr]] /.
a[i_] :> a[[i]]]
]
]
Если скомпилированная функция содержит большой кусок кода, Hold, Release и замена в конце идут вразрез с моим представлением о прекрасном коде. Есть ли короткое и приятное решение, которое я еще не рассмотрел?
Ответ на пост сабольца
Не могли бы вы рассказать мне, почему вы используете здесь?
Да, и это связано с причиной, по которой я не могу использовать: = здесь. Я использую With, чтобы иметь что-то вроде #define
в C, что означает замену кода в том месте, где оно мне нужно. Использование: = in With задерживает оценку, и то, что видит тело Compile, не является последним фрагментом кода , который предполагается компилировать. Таким образом,
<< CompiledFunctionTools`
cf = Block[{a, x, degree = 3},
With[{expr := Product[x - a[[i]], {i, degree}]},
Compile[{{x, _Real, 0}, {a, _Real, 1}}, expr]]];
CompilePrint[cf]
показывает, что в скомпилированной функции вызывается ядро Mathematica
I4 = MainEvaluate[ Function[{x, a}, degree][ R0, T(R1)0]]
Это плохо, потому что Compile должен использовать только локальные переменные для вычисления результата.
Обновление
В этом случае решение Szabolcs работает, но оно оставляет все выражение без оценки. Позвольте мне объяснить, почему важно, чтобы выражение было расширено до его компиляции. Должен признать, мой игрушечный пример был не лучшим. Итак, давайте попробуем лучше, используя With и SetDelayed, как в решении Szabolcs
Block[{a, x}, With[
{expr := D[Product[x - a[[i]], {i, 3}], x]},
Compile[{{x, _Real, 0}, {a, _Real, 1}}, expr]
]
]
Скажем, у меня есть многочлен степени 3, и мне нужна его производная внутри компиляции. В приведенном выше коде я хочу, чтобы Mathematica вычисляла производную для неназначенных корней a [[i]], чтобы я могла использовать формулу очень часто в скомпилированном коде. Глядя на скомпилированный код выше
каждый видит, что D [..] не может быть скомпилирован так же хорошо, как Продукт, и остается неоцененным
11 R1 = MainEvaluate[ Hold[D][ R5, R0]]
Поэтому мой обновленный вопрос: возможно ли оценить кусок кода без оценки Part [] - доступ к нему лучше / лучше, чем с использованием
Block[{a, x}, With[
{expr = D[Quiet@Product[x - a[[i]], {i, 3}], x]},
Compile[{{x, _Real, 0}, {a, _Real, 1}}, expr]
]
]
Редактировать: Я положил Тихий на место, которому он принадлежит. У меня было это перед блоком кода, чтобы всем было видно, что я использовал Quiet здесь, чтобы подавить предупреждение. Как уже указывал Рубенко, в реальном коде он всегда должен быть как можно ближе к месту, к которому он принадлежит. При таком подходе вы, вероятно, не пропустите другие важные предупреждения / ошибки.
Обновление 2
Поскольку мы отходим от первоначального вопроса, нам следует перенести это обсуждение, возможно, в новую ветку. Я не знаю, кому мне следует дать лучший ответ на мой вопрос, поскольку мы обсудили Mathematica и Scoping больше, чем способы подавления проблемы [[i]].
Обновление 3
Чтобы дать окончательное решение: я просто подавляю (к сожалению, как делал все время) предупреждение [[i]] с помощью Quiet. В настоящем примере ниже я должен использовать Quiet за пределами полного блока, чтобы подавить предупреждение.
Чтобы вставить требуемый код в тело компиляции, я использую чистую функцию и даю код в качестве аргумента. Это тот же подход, который использует Майкл Тротт, например его книга чисел. Это немного похоже на предложение where
в Haskell, где вы определяете то, что использовали потом.
newtonC = Function[{degree, f, df, colors},
Compile[{{x0, _Complex, 0}, {a, _Complex, 1}},
Block[{x = x0, xn = 0.0 + 0.0 I, i = 0, maxiter = 256,
eps = 10^(-6.), zeroId = 1, j = 1},
For[i = 0, i < maxiter, ++i,
xn = x - f/(df + eps);
If[Abs[xn - x] < eps,
Break[]
];
x = xn;
];
For[j = 1, j <= degree, ++j,
If[Abs[xn - a[[j]]] < eps*10^2,
zeroId = j + 1;
Break[];
];
];
colors[[zeroId]]*(1 - (i/maxiter)^0.3)*1.5
],
CompilationTarget -> "C", RuntimeAttributes -> {Listable},
RuntimeOptions -> "Speed", Parallelization -> True]]@@
(Quiet@Block[{degree = 3, polynomial, a, x},
polynomial = HornerForm[Product[x - a[[i]], {i, degree}]];
{degree, polynomial, HornerForm[D[polynomial, x]],
List @@@ (ColorData[52, #] & /@ Range[degree + 1])}])
И эта функция теперь достаточно быстра для вычисления фрактала Ньютона полинома, в котором положение корней не фиксировано. Поэтому мы можем корректировать корни динамически.
Не стесняйтесь настроить n. Здесь он бегло до n = 756
(* ImageSize n*n, Complex plange from -b-I*b to b+I*b *)
With[{n = 256, b = 2.0},
DynamicModule[{
roots = RandomReal[{-b, b}, {3, 2}],
raster = Table[x + I y, {y, -b, b, 2 b/n}, {x, -b, b, 2 b/n}]},
LocatorPane[Dynamic[roots],
Dynamic[
Graphics[{Inset[
Image[Reverse@newtonC[raster, Complex @@@ roots], "Real"],
{-b, -b}, {1, 1}, 2 {b, b}]}, PlotRange -> {{-b, b}, {-
b, b}}, ImageSize -> {n, n}]], {{-b, -b}, {b, b}},
Appearance -> Style["\[Times]", Red, 20]
]
]
]
Teaser: