Mathematicaで以前からこれが遅くなる原因ではないか…と思っていた部分を検証したのでメモ．
概要
Mathematicaの代入
expr /. {replace1 -> value1, replace2 -> value2, ...}
は，代入されるもの・代入するものが多いと遅くなる (vectorizationできるならすべき)．
これの類推で，expr (ここでは行列) の一部を書き換える関数ReplacePart
ReplacePart[expr, {{i1,j1} -> value1, {i2,j2} -> value2, ...}]
をTableを用いて
ReplacePart[expr, Table[position[[i]] -> value[[i]], {i, Length@position}]]
のように書くのは遅くなる原因なのだろうか (7年ぐらい前に書いたコードはこれ)．
MapThreadを使うとベクトルをベクトル表記のまま同じことができる:
ReplacePart[expr, MapThread[Rule, {position, value}]]
ので，Tableを使う場合とMapThreadを使う場合で差が出るか見てみた．
結果
MapThreadのほうが1割ぐらい速かった．実行環境はやや古いiMacのMathematica14.0.0.
コード
以下コード．
pattern01[nrow_, ncol_, nreplace_, nDo_] :=
Module[{timebegin, mat, replacePosition, newValue},
(* record computational time *)
timebegin = SessionTime[] ;
(*"set seed"*)
ParallelEvaluate[SeedRandom[1]] (* enable to be replicated *);
Do[
(* a matrix that I want replace with newValue *)
mat = RandomReal[1, {nrow, ncol}];
(* positions of mat to be replaced.
the number of replace is at most nreplace. *)
replacePosition =
DeleteDuplicates@
Transpose@{RandomInteger[{1, nrow}, nreplace],
RandomInteger[{1, ncol}, nreplace]};
(* new values (taking Exp is optional) *)
newValue = Exp@RandomReal[1, Length@replacePosition];(*
replacepart with Table *)
mat =
ReplacePart[mat,
Table[replacePosition[[i]] -> newValue[[i]], {i,
Length@replacePosition}] ];,
nDo] (* end of Do *);
Print[ToString@(SessionTime[] - timebegin) <> " sec"];
Return@DateString[]
]
pattern02[nrow_, ncol_, nreplace_, nDo_] :=
Module[{timebegin, mat, replacePosition, newValue},
(* record computational time *)
timebegin = SessionTime[] ;
(*"set seed"*)
ParallelEvaluate[SeedRandom[1]] (* enable to be replicated *);
Do[
(* a matrix that I want replace with newValue *)
mat = RandomReal[1, {nrow, ncol}];
(* positions of mat to be replaced.
the number of replace is at most nreplace. *)
replacePosition =
DeleteDuplicates@
Transpose@{RandomInteger[{1, nrow}, nreplace],
RandomInteger[{1, ncol}, nreplace]};
(* new values (taking Exp is optional) *)
newValue = Exp@RandomReal[1, Length@replacePosition];(*
replacepart with MapThread *)
mat =
ReplacePart[mat, MapThread[Rule, {replacePosition, newValue}] ];,
nDo] (* end of Do *);
Print[ToString@(SessionTime[] - timebegin) <> " sec"];
Return@DateString[]
]
"Parallel Karnels" (* 使ってないけど *)
LaunchKernels[];
"run (execution requires about 20min)"
(* 16384x128行列のうち最大1024箇所を書き換える，を65536回実行 *)
pattern01[2^14, 2^7, 2^10, 2^16]
(* 1420 sec *)
pattern02[2^14, 2^7, 2^10, 2^16]
(* 1318 sec *)
その他
Parallelize@ReplacePartのようにして並列化はできなさそう．
