2012-05-09

Mathematica Tips – 2

1.Module[]中注意加分号, 否则在计算时会出现意想不到的错误;而且格式化代码时会乱掉.


2.提取列表中的某一项元素出现的位置
Position[{a, b, a, a, b, c, b}, b]


3.Reap 和 Sow 用法:
In[3]:= Reap[Sow[1, {x, x}]; Sow[2, y]; Sow[4, x], {x, x, y}]
Out[3]= {4, {{{1, 1, 4}}, {{1, 1, 4}}, {{2}}}}

关于Reap和Sow以及Scan的一个巧妙用法:
In[18]:= partition[l_, v_, comp_] :=
Flatten /@
  Reap[Scan[
     Which[comp[#1, v], Sow[#1, less], comp[v, #1], Sow[#1, large],
       True, Sow[#1, equ]] &, l], {large, equ, less}][[2]]

In[19]:= partition[{3, 5, 7, 9, 2, 4, 6, 8, 3, 4}, 4, Less]
Out[19]= {{5,7,9,6,8},{4,4},{3,2,3}}


4.不要轻易加脚标:
Subscript[s, k] = Array[0 &, nK];
Subscript[s, k][[1]] = 1;


5.Map(/@) 和 Apply(@@) 对比:

Map不是替换,而是把函数映射到列表的每一个元素, 即把元素分别传递给函数的变量:
In[295]:= f /@ {1, 2, 3, 4, 5}
Out[295]= {f[1], f[2], f[3], f[4], f[5]}

Apply 自动把列表中的元素传递给多变量函数(也就是它仅仅替换函数头),@@@形式自动作用于第一层
In[52]:= Mod[#1, #2] & @@@ {{10, 4}, {5, 2}}
Out[52]= {2,1}

同样的效果, 二者结合方式(借助纯函数):
In[41]:= Apply[Mod, #] & /@ {{10, 4}, {5, 2}}
Out[41]= {2,1}

6.Evaluate函数很有用:
    In[1]:= ch = ChebyshevT[5, x]
    Out[1]= 16 x^5-20 x^3+5 x

    In[4]:= Function[x, Evaluate[ch]]
    Out[4]= x\[Function]16 x^5-20 x^3+5 x

    In[5]:= %[10]
    Out[5]= 1580050
   
   
7.Compile 函数可以提高数值运算速度! (当然采用内部函数是最快的)
它不但可以处理数学表达式,还可以处理各种简单的 Mathematica 程序. 例如,Compile 可以处理条件和控制流结构.
对比:
1.使用Compile:
In[40]:= newtonIteration :=
  Compile[{x, {n, _Integer}}, Module[{t}, t = x;
    Do[t = (t + x/t)/2, {n}]; t]];
newtonIteration[2.4, 666666] // Timing

Out[41]= {0.094,1.54919}

2.不使用Compile:
In[42]:= newtonIteration2[x_, n_] := Module[{t}, t = x;
   Do[t = (t + x/t)/2, {n}]; t];
newtonIteration2[2.4, 666666] // Timing

Out[43]= {2.137,1.54919}

运行时间相差20倍.


8.编写高效代码最重要的方式之一就是要避免显式部分引用,特别是在内部循环中.

如果将要进行实数的操作,一定要确保使用实数进行初始化.

混合的符号/数值矩阵通常比操作数值矩阵要慢.

一个整数矩阵将使用符号计算技术,其速度较慢, 但可以给出精确解。


9.条件迭代:
ff = # + 1 &;
NestWhileList[ff, 2, # < 5 &]

双层迭代的一个例子:
Mp = Table[1, {1}, {4}]
Module[{i, k},
For[i = 1, i < 5, i++,
   Mpp = {};
   Do[AppendTo[Mpp, Sin[k + i]];
    , {k, 4}];
   AppendTo[Mp, Mpp]];]
Mp

 

10.迭代例子2:

intN = 10;
intK = 200;
(*注意是矩阵*)pki = RandomInteger[{1, 2 intK}, {1, intK}];
mA = RandomInteger[{1, 2 intK}, {intK, intN}];
ak = RandomInteger[{1, 2 intK}, {intK, intN}];
y = RandomInteger[{1, 2 intK}, intN];
(*mP[list_]:=DiagonalMatrix[list];*)
(*Subscript[R, y](i)更新*)
Ryi[i_] := mA\[HermitianConjugate].(DiagonalMatrix[pki[[i]]]).mA;
Ryi[1] = RandomInteger[{1, 2 intN}, {intN, intN}];
(*Subscript[p, k](i)更新*)
pkki[k_, i_] := (Abs[ak[[k]]\[Conjugate].Inverse[Ryi[i - 1]].y])^2/(ak[[k]]\[Conjugate].Inverse[Ryi[i - 1]].ak[[k]])^2;


For[i = 1; endCond = 1, endCond > 10^-4(*i<=11*), i++,
bufMp = {};
Do[AppendTo[bufMp,(*N[ pkki[k,i+1] ]*)pkki[k, i + 1]];
  , {k, intK}];
(*AppendTo使用需要注意层次*)
AppendTo[pki, bufMp];
endCond = Norm[pki[[i + 1]] - pki[[i]]]/Norm[pki[[i]]];
]
vectorPk = pki[[i]]

没有评论:

发表评论