wolfram-mathematica - 数学 : using simplify to do common sub-expression elimination and reduction in strength

Question

所以最近我一直在玩弄 Mathematica 的模式匹配和术语重写如何在编译器优化中得到很好的利用……试图高度优化作为循环内部部分的短代码块。减少计算表达式所需的工作量的两种常见方法是识别多次出现的子表达式并存储结果，然后在后续点使用存储的结果来节省工作。另一种方法是尽可能使用更便宜的操作。例如，我的理解是求平方根比加法和乘法需要更多的时钟周期。需要明确的是，我对评估表达式所需的浮点运算成本感兴趣，而不是 Mathematica 评估它需要多长时间。

我的第一个想法是我会使用 Mathematica 的 simplify 解决开发问题。功能。可以指定一个复杂度函数来比较两个表达式的相对简单性。我打算为相关的算术运算使用权重创建一个，并将 LeafCount 添加到表达式中，以说明所需的分配操作。这解决了强度方面的降低问题，但让我绊倒的是消除常见的子表达式。

我正在考虑将公共(public)子表达式消除添加到简化使用的可能转换函数中。但是对于一个大表达式，可能有许多可能的子表达式可以被替换，并且在你看到表达式之前不可能知道它们是什么。我编写了一个提供可能替换的函数，但您指定的转换函数似乎只需要返回一个可能的转换，至少来自文档中的示例。关于如何绕过这个限制的任何想法？有没有人更好地了解简化如何使用可能暗示前进方向的转换函数？

我想，在幕后 Simplify 正在做一些动态编程，尝试对表达式的不同部分进行不同的简化，并返回复杂度分数最低的那个。我是否会更好地尝试使用常见的代数简化(例如因子和收集)自己进行这种动态规划？

编辑:我添加了生成可能要删除的子表达式的代码

(*traverses entire expression tree storing each node*)
AllSubExpressions[x_, accum_] := Module[{result, i, len},
  len = Length[x];
  result = Append[accum, x];
  If[LeafCount[x] > 1,
   For[i = 1, i <= len, i++,
     result = ToSubExpressions2[x[[i]], result];
     ];
   ];
  Return[Sort[result, LeafCount[#1] > LeafCount[#2] &]]
  ]

CommonSubExpressions[statements_] := Module[{common, subexpressions},
subexpressions = AllSubExpressions[statements, {}];
  (*get the unique set of sub expressions*)
  common = DeleteDuplicates[subexpressions];
  (*remove constants from the list*)
  common = Select[common, LeafCount[#] > 1 &];
  (*only keep subexpressions that occur more than once*)
  common = Select[common, Count[subexpressions, #] > 1 &];
  (*output the list of possible subexpressions to replace with the \
number of occurrences*)
  Return[common];
  ]

一旦从 CommonSubExpressions 返回的列表中选择了一个公共(public)子表达式，执行替换的函数如下所示。

eliminateCSE[statements_, expr_] := Module[{temp},
temp = Unique["r"];
Prepend[ReplaceAll[statements, expr -> temp], temp[expr]]
]

冒着这个问题变长的风险，我会放一些示例代码。我认为尝试优化的体面表达式将是经典的 Runge-Kutta求解微分方程的方法。

Input:
nextY=statements[y + 1/6 h (f[t, n] + 2 f[0.5 h + t, y + 0.5 h f[t, n]] + 
    2 f[0.5 h + t, y + 0.5 h f[0.5 h + t, y + 0.5 h f[t, n]]] + 
    f[h + t, 
     y + h f[0.5 h + t, y + 0.5 h f[0.5 h + t, y + 0.5 h f[t, n]]]])];
possibleTransformations=CommonSubExpressions[nextY]
transformed=eliminateCSE[nextY, First[possibleTransformations]]

Output:
{f[0.5 h + t, y + 0.5 h f[0.5 h + t, y + 0.5 h f[t, n]]], 
y + 0.5 h f[0.5 h + t, y + 0.5 h f[t, n]], 
0.5 h f[0.5 h + t, y + 0.5 h f[t, n]], 
f[0.5 h + t, y + 0.5 h f[t, n]], y + 0.5 h f[t, n], 0.5 h f[t, n], 
0.5 h + t, f[t, n], 0.5 h}

statements[r1[f[0.5 h + t, y + 0.5 h f[0.5 h + t, y + 0.5 h f[t, n]]]], 
y + 1/6 h (2 r1 + f[t, n] + 2 f[0.5 h + t, y + 0.5 h f[t, n]] + 
 f[h + t, h r1 + y])]

最后，判断不同表达式的相对成本的代码如下。在这一点上，权重是概念性的，因为这仍然是我正在研究的一个领域。

Input:
cost[e_] := 
 Total[MapThread[
   Count[e, #1, Infinity, Heads -> True]*#2 &, {{Plus, Times, Sqrt, 
     f}, {1, 2, 5, 10}}]]
cost[transformed]

Output:
100

Answer 1

这里也有一些程序是作者在这里实现的:http://stoney.sb.org/wordpress/2009/06/converting-symbolic-mathematica-expressions-to-c-code/

我将它打包成一个 *.M 文件并修复了一个错误(如果表达式没有重复的子表达式，它就会死掉)，我正在尝试查找作者的联系信息，看看我是否可以将他修改后的代码上传到 pastebin 或任何地方.

编辑:我已获得作者的上传许可并将其粘贴在此处:http://pastebin.com/fjYiR0B3