数学建模社区-数学中国
标题: 10个编写快速运行的Mathematica代码的小诀窍 [打印本页]
作者: Asdmath2020 时间: 2020-3-16 16:28
标题: 10个编写快速运行的Mathematica代码的小诀窍
本帖最后由 Asdmath2020 于 2020-3-16 16:32 编辑
2 K& ^- m$ q% l" B! T9 g( r8 w
}6 U" j- `% x! |3 ~1 `; B当我听到人们说Mathematica不够快的时候，我通常会提出想要看一下这段令他们烦恼的代码，然后会发现，其实并不是Mathematica本身的表现不够好，而是Mathematica没有被最优使用。我觉得我应该和大家分享一下我在优化Mathematica代码时首先会看的一些内容。* k$ K/ |% U+ P& |
% ?6 E6 ~# a/ W# Z% }- W8 q; K
01 如果可以的话尽量尽早使用浮点数" M ~" Y) J; \ F {
% C3 m, G* a2 {* Z我最常看到的导致代码变慢的问题是，程序员会不经意地让Mathematica做超出需要的细致的事情。没必要的代数精确是其中最常见的问题。
在多数与数字相关的软件中，是不需要这么精确的代数的。1/3和0.33333333333333是一样的。当你碰到特别严重的在数字上不稳定的问题时这个差异可能会被放大的特别明显，但是，在大多数情况中，浮点数已经足够使用了，而且最重要的是，浮点数运算更快。Mathematica中，任何小于16位的小数都被看作是机器浮点数，所以如果更想要速度而可以舍弃一些精确性的时候，记得用小数（比如，三分之一输入为1./3.）。以下是一个例子，可以看到使用浮点数是精确数运行速度的50.6倍。在这个例子中，两个数字的使用得到的是同一个结果。
# |! o8 h0 U3 \
5 ^! _9 d' R3 j3 ]" U
在符号运算中也是这样。如果你不是很在意符号式的结果，并且计算的稳定性也不是问题的话，那么尽快使用数值作为替代。比如，求解下面的二项式符号计算时，在使用数值作为替代之前，这个代码可能会让Mathematica生成长达五页的中间符号表达式。
+ P5 s+ U& W2 V! p; Q
* X% K8 V% g( p9 u6 E: Z
6 |% G3 P. P/ i0 E/ [, p9 @9 ~# F
但是如果先用数值替代，那么Solve会使用更快的数值方法。3 s/ |% E, Z. c( ?! A4 L( P8 y& f
% a0 {+ I8 M5 L* R7 k* K
: f2 V }0 L# g! `# R5 j( y当用数据列表工作时，使用实数的方法必须保持一致。只要一个精确的数值就可以让整个数据组处于一个更灵活但是缺乏效率的形式中。/ \7 N5 V' U0 D5 i6 J$ p
) [3 c8 d% P* d! r7 \! B; M; h
02学会Compile( r, p4 g& g% E& S. d2 j0 b5 k D' S
% _$ _- m# y3 ]
Compile函数接受Mathematica的代码，并让你预先声明输入参数的类型（比如实数、复数等）和结构（如数值、列表、矩阵等）。这虽然失去了Mathematica语言灵活性的优势，但是可以免于担心类似于“如果参数是符号怎么办？”的问题，Mathematica也可以最优化程序并创建一个字节码在虚拟器上运行。并不是所有东西都可以被编译，且简单的代码可能不会有太大效果，但是那种复杂的低阶数字代码速度可以得到大大的提升。
下面是一个例子：
Z! z8 r- K1 z( Q) h# u) {/ k6 w
& g( [+ [! [% s# ^1 M$ @
使用Compile可以比Function的运行速度提高80倍。
: @: k" H& o, Z
3 s8 s3 y! H) F* v9 B. J/ q
但是我们可以在Compile函数中加入一些代码的可并行性质，这样可以生成更好的结果。
$ J# a$ j' V2 \: R7 g- R2 C
5 v/ a) S# d( l在我的双核处理器电脑上，我的运行结果比原本快150倍，如果是多核处理器那么效果会更加明显。- t. {) [: m- `
但是要注意，很多Mathematica函数比如Table、Plot、NIntegrate等会自动编译它们的参数，这样的话你使用上述方法可能不会看到任何速度上的提升。
; o! H9 K1 k* }! Y# N9 a
02.5使用Compile生成C代码
M) y# S, ?$ Y3 T& p$ I& U6 n( m3 ]
另外，如果你的代码可编译，你还可以使用选项CompilationTarget->“C”来生成C代码，调用你的C编码器并将其汇编成一个DLL，并把这个DLL链接回Mathematica，都是自动操作的。在编译阶段，DLL直接在CPU上运行而非Mathematica的虚拟器，所以会更快得到结果。
9 j8 ^% E! I7 s( J
0 {$ U8 Z8 Z( N5 [( ]
03使用内置函数
+ N8 e% C' u6 K1 r, E! z+ X
$ a" j1 b! t3 N! B6 w) BMathematica有很多函数。起码半数以上的人可能不会坐下来学习所有函数。所以当我看见有些人会写一些代码而没有意识到其实Mathematica知道怎么做这些操作的时候，我一点也不意外。这种重复操作不仅是浪费时间，而且公司是花钱请程序员来开发研究运行这些操作的最有效方法，所以内置的函数一般是非常快的。
如果你发现有些结果很接近了但是不完全对的时候，此时可以检查选项和参数，通常它们会概括可以覆盖很多特殊用法或者专有应用的函数。
下面举一个这样的例子。如果我有一个一百万2x2矩阵的列表，我想把该列表转换成一百万个包含四个元素的列表，概念上来说最简单的方法是用Map把已经用Flatten扁平化过的数据进行映射即可。
- J4 p9 i1 X: d' _; W m1 O
5 n6 T# U; E1 T$ t
但是Flatten本身知道怎么把整个步骤完成，你只要说明数据结构的第二层和第三层应该被合并而第一层不动就可以了。说明这种细节的内容可能相对来说是比较细致的工作，但是只需要使用Flatten就能完成整个扁平化工作可以让整个进程比你自己手动做这些程序要快将近4倍。- S8 {) ~- d. G i* ?1 P
2 w+ O/ E) n. P" {4 j- b1 k; A
所以记住：在运行代码之前在帮助菜单里先搜索一遍。
, i% G+ t4 y5 t6 I, i, @# b7 ~
, Y* P* c3 Q( m: d04使用Wolfram Workbench
, A. G0 N' P8 m% _5 i2 D
9 m( d* @0 b* H, @4 q/ A. hMathematica对于某些种类的编程错误容忍度很高——如果你忘记在正确的时候初始化一个变量，Mathematica会以符号的模式顺利运行，而并不会有循环计算或者预料之外的数据类型出现。如果你只想要一个答案的话这个功能是很棒的，但是这也会让你没有得到最优的解答。
Workbench会在几个方面帮助你。首先它会帮你排除程序问题，并把大型的代码项目组织得更好，整齐易读的代码会让程序员更好地写优秀的代码。但是最关键的功能在于分析器会告诉你是哪一行代码用光了时间，而且会告诉你调用这些代码用了多少时间。
看下这个例子，一个很可怕的执行斐波那契数的方法。如果你没有考虑到数列的双重递归，你可能会惊讶计算fib[35]怎么会需要22秒钟（大约和内置函数计算Fibonacci[1000000000]所有208,987,639位数字需要的时间一样）（请看诀窍3）。
( U5 f2 t! d8 M8 }. H在分析器中运行这个代码可以解释这个现象的原因。主要规则被援引9,227,464次，fib[1]的值被请求18,454,929次。& I& a6 t3 h7 ]0 ?; P) `
学习代码能做什么，而不是想当然，会让你眼界大开。
. `! K0 J! e O; g3 l
05记住你将来会需要用到的值
) X C2 h7 C2 x, A. F' Q. S+ V: _$ X4 D& G, Z9 c" K
这个编程诀窍对任何语言都管用。Mathematica认为你想知道的是这个：
; N5 z4 B/ e/ ?
5 U3 |: ^& `) [! }8 b8 J. I/ }
& N s/ b# {$ D" _, q6 p+ C. N这省去了用任何值调用 f 的结果，这样的话如果再用相同数值调用 f，Mathematica不需要再算一遍。这里你就是用内存换取计算速度，所以如果你的函数要用大量不同数值调用而不太重复的时候这个方法可能不合适。但是如果输入的范围有限，那么这个方法就很有用了。以下就是如何拯救我刚才提到的来解释诀窍3的例子的方法。可以把第一条规则改成这样：
. h4 t. l/ Z; U! U) i2 i
( @$ P2 h5 R1 l2 q m
然后速度立刻就可以提升，因为fib[35]现在只需要用主要的规则运算33次。查询之前的结果可以防止循环递归fib[1]的问题。
1 Y! q( N+ H, ~6 j. V8 d
) g2 z+ u0 c& w6 S9 e6 u06并行
* U) g* O% M0 L: d9 X4 O) q/ J: V% S1 L0 F7 c# |; t* `
有很多Mathematica的操作都会自动在本地核中并行运行（大部分是代数、图像处理和统计），如果需要手动的话，Compile也可以。但是对于其他操作来说，或者如果你想在远程硬件上并行操作，你可以试用内置的并行编程架构来完成。
有一个这样工具的集合，但是都是为非常独立的任务服务的，比如ParallelTable,ParallelMap，ParallelTry，还有很多。每个这样的小工具都可以自动进行通信、工作管理和收集结果。发送任务和回收结果需要一点时间，所以在减少时间和增加时间上会有需要一个取舍。你的Mathematica有四个计算内核，如果你有额外的CPU可使用的话，还可以通过gridMathematica在此基础上提高这一性能。这里由于我用的是双核电脑，ParallelTable实际将我的运算时间缩少了一半。如果有更多CPU则会得到更好的结果。
" w; G, A! Z. I8 b
, y5 W2 c/ c) N5 y- G$ X3 E任何Mathematica可以做的事情都可以以并行方法运行。比如，你可以给远程硬件发送一个并任务集合，每个任务都在CPU或GPU中编译和运行。' r. n/ a7 T8 {* j
1 ^6 s! Z) V+ ^3 n8 T
06.5想想CUDALink和OPENCLLink
" |+ d I9 @0 k/ {, g3 T& x, X% ?" F Y$ B
如果你有GPU硬件，有一些用批量并行运行方法可以做的非常快的事情。除非这些最优化CUDA函数恰好可以做你想要它们做的事情，否则你还要做一点额外的工作，但是CUDALink和OpenCLLink 工具可以为你自动化很多繁琐的细节。
# w! K8 T% w' y" W0 \07使用Sow和Reap累积大量数据（不是AppendTo）
4 P' A+ ^) ]4 U7 Q4 x& e7 d7 s: u( h
因为Mathematica数据结构的灵活性，AppendTo不会假设你要追加的是一个数字，因为你要追加的可能是一个文件、音频或者图像等。所以AppendTo必须为所有数据创建一个新的副本，并重新调整架构以适应新追加的信息。当数据累积的时候这个过程会变得越来越慢。（而且构建data=Append[data,value]与AppendTo一样。）
尝试使用Sow和Reap。Sow会舍弃你想要累积的值，而Reap收集它们并一次性在末尾建立一个数据对象。下列范例是等价的：
. M0 y$ \" P+ u4 f! f% t
08使用Block或With而非Module4 P% x, M1 O6 s+ J3 n, H
% q7 A# i" b. N( f
Block，With和Module都是本地化构建的工具，但是属性上有些小区别。根据我的经验，95%以上的几率在我写的代码中Block和Module是可以互相替换的，但是Block通常快一点，而在另一些例子中（Block的变量在只读状态的情况下）With会快一些。& B$ e3 T/ {% D% F0 V
6 T9 H4 U1 g* |
5 R9 G. Y( W6 O# d* U8 ~09少用模式匹配5 r; k! G ^7 Q& N8 l+ F1 M
6 C" f: N \! R5 @
模式匹配很好，可以让项目中复杂的任务变得简单一点。但是它有时会很慢，尤其是像BlankNullSequence这种比较复杂的模式（通常写作“___”）中，可能会花很长时间仔细在你的数据中搜索一些——你作为一个程序员可能已经可以判断的——不存在的模式。如果想要速度的话，那么选择范围更窄的模式，或者不用模式会更好。
比如，下面范例使用了模式，在一行代码中简洁地执行了冒泡排序：
- g. E" \$ U/ Y- s
上例概念上很简单，但是比起这个我最开始学习编程的时候就学过的列出步骤的方法来说还是要慢很多：% x; ?5 N% n. q# T- f$ \
6 Q$ J8 B6 P4 s' j8 C) m
& P, _9 s8 p; y: l1 T$ p
当然在这个例子中你可以用内置函数（参见诀窍3），这个内置函数会使用比冒泡排序更好的排序算法。
& y4 y. L b9 N- g3 u$ C0 F! u: P* [ q3 o1 h# R" Q
010尝试不同的方法6 |) [' m! \( i0 q: p/ V
* x8 L! P9 P2 D+ _
Mathematica的一个很重要的优点是，它可以用不同的方式处理同一个问题。它允许你按照你自己的想法编程，而不是为了编程语言的风格重构你的问题。但是，概念上简单和计算效率不是一件事。有时候容易懂的想法可能会需要更多的工作才能实现。
但是另一个问题是，因为Mathematica中最优化和一些绝妙算法都是自动应用的，所以很难预测什么时候Mathematica又会做出另一个绝妙的操作。比如，下例是两种计算阶乘的方法，第二种比第一种快10倍。
- c" }& w- T5 C8 e9 r# @8 Y3 B
& q1 |( x$ M2 Z4 f6 |3 N& `7 p为什么？你可能会猜可能Do的循环很慢，或者所有这些任务缓存都需要时间，或者可能第一次执行的时候有什么东西出了问题，但是实际的原因很难预料到。Time有一个很聪明的二元分离的小技巧，可以在当你有大量整数参数的情况下使用，即将循环将参数分成两个更小的乘积(1*2*…*32767)*(32768*…*65536)，而不是把这个参数从第一个用到最后一个。当然要做的乘法数量还是一样，但是不会再包括数值非常大的整数，所以平均来说，运算的速度会更快。在Mathematica中有很多这样隐藏的小魔法，而且每次新版本发布都会有更多的小技巧加入。9 M: T% G8 M( f0 K
当然最好的方法还是使用内置函数（又说到诀窍3了）：
$ j" ^+ L. q1 _. `" v6 R- A
1 B! Q* p! C: `# qMathematica可以做非常高级的计算，而且有强大的功能和极高的精确性，但是这两者并不总能兼得。我希望这些诀窍可以在快速编程、快速执行和精确结果的冲突诉求中对你有些许帮助。1 j3 K9 u" y& m K5 h- w& O
3 ~7 D( r1 y5 I2 y: P6 d) O/ i7 \+ F6 C% {6 a
光学软件供应 软件定制开发光学软件培训 光学解决方案光学仪器设备 光学镜头设计科学计算软件 机械设计软件高校辅助教学方案 BSDF测量衍射光学元件设计开发-------------------------------: m6 |, F3 b9 [% Y6 B$ ~# Z5 \
电话：19926678573邮箱：ada@asdoptics.com网址：www.asdoptics.com: ~# S0 A" s. z* l* Z- B; X
8 \- K5 v% l; g
0 }7 a! }! }6 q. D; y2 @6 |; ]
-
008.png
(19.59 KB, 下载次数: 1414)
-
07.jpg
(15.98 KB, 下载次数: 414)
-
07.jpg
(15.98 KB, 下载次数: 396)
-
007.png
(116.25 KB, 下载次数: 466)
-
05.jpg
(15.61 KB, 下载次数: 398)
作者: 去嗨嗨嗨嗨嗨 时间: 2020-3-16 17:57
大神好厉害啊！1 q) M8 y0 x: B j
| 欢迎光临 数学建模社区-数学中国 (http://www.madio.net/) |
Powered by Discuz! X2.5 |
