3
$\begingroup$

Suppose that there are 2000 stones in a row.
You should grind 1000 of them to make beads.
The time loads(making beads) for each stone are similar.
But which 1000 stones?

case1) the first 1000 stones.

case2) Let B = RandomSample[Range[2000], 1000].
Then, all b-th stones, where b denotes element of B.

If this is a problem of real world, then case1) is faster then case2), because
for case2),

  1. You have to read the B every time after grinding a stone.
  2. If the next stone to be processed is far away, it will take a long time to go there.
  3. B is even not sorted by size.

But, if it was a problem of mathematica,

2000 stones in a row = a list of 2000 elements
grinding a stone = applying function to an element

The time loads for case1) and time loads for case2) should be similar in my opinion.

I experienced that "just because some elements of a list are positionally close, there is no particular advantage to processing(=applying function) them at once(or consecutively)"

Also I believe that "there is a fast mathematica code for case2) always, that takes silmilar time compared to case1)"

Am I thinking correctly?

If so, can you make a fast code for the following problem ?

For L, square(=apply #^2&) each b-th element of L.(leave other element unchanged)

L = Range[40000, 59999];
B = RandomSample[Range[20000], 10000];

My first trial was

MapAt[#^2&, L, {#}&/@B] //Timing
2.85938

which is never successful.

After some other trials, I realized that MapAt is very slow, in all cases.
The performance was so bad that I felt MapAt should not be used in any case.

Finally I succeeded to shortening the time(using ReplacePart) but it is not neat and I don't think my code is as fast as professional programmer's one.

What skill do you use, when applying a function for a part of list? (faster!)

Can you help me?

$\endgroup$
1
  • 4
    $\begingroup$ Much faster: cres = lx; cres[[bx]] = lx[[bx]]*lx[[bx]]. $\endgroup$
    – Alan
    Commented Mar 5, 2022 at 18:55

4 Answers 4

7
$\begingroup$

I am sure this is not the fastest solution but it is 25 times faster.

lx = Range[40000, 59999];
bx = RandomSample[Range[20000], 10000];
ares = MapAt[#^2 &, lx, {#} & /@ bx] // Timing
bres = SubsetMap[#^2 &, lx, {#} & /@ bx] // Timing

ares[[2]] == bres[[2]]
(*True*)
{First@ares, First@bres}
{1.17001, 0.0468003}
$\endgroup$
4
  • $\begingroup$ Wow! SubsetMap is a new feature in V12! Thank you! $\endgroup$
    – imida k
    Commented Mar 5, 2022 at 17:48
  • $\begingroup$ For this particular variation of the SubsetMap syntax, you are executing a MapAt at the subset of the positions that you have chosen somehow, so there is no mental gymnastics involved. $\endgroup$
    – Syed
    Commented Mar 5, 2022 at 17:51
  • $\begingroup$ May I ask what if you were in mathematica V11 ? $\endgroup$
    – imida k
    Commented Mar 5, 2022 at 17:59
  • $\begingroup$ Please edit this requirement in to your post. I don't have M11 so I cannot say for sure, but if I find another way I will update the answer. $\endgroup$
    – Syed
    Commented Mar 5, 2022 at 18:17
8
$\begingroup$

Here are a couple more to add to Syed's & Alan's. The compiled version is much faster, if everything stays machine-size. If you're doing this repeatedly, you could compile ahead of time (to C even). You have to recompile if you change the function from #^2; and if the function to be applied is not compilable, then it obviously won't work well.

lx = Range[40000, 59999];
bx = RandomSample[Range[20000], 10000];
ares = MapAt[#^2 &, lx, {#} & /@ bx]; // AbsoluteTiming
bres = SubsetMap[#^2 &, lx, {#} & /@ bx]; // AbsoluteTiming
(cres = lx; Do[cres[[k]] = cres[[k]]^2, {k, bx}]); // AbsoluteTiming
dres = Compile[{{lx, _Integer, 1}, {bx, _Integer, 1}},
     Block[{dres = lx},
      Do[dres[[k]] = dres[[k]]^2, {k, bx}];
      dres
      ]
     ][lx, bx]; // AbsoluteTiming
(eres = lx; eres[[bx]] = eres[[bx]]^2); // AbsoluteTiming

ares == bres == cres == dres == eres
(*
  {0.633558, Null}
  {0.013454, Null}
  {0.006371, Null}
  {0.000657, Null}
  {0.000221, Null}  <-- Alan's is fastest
  True
*)
$\endgroup$
2
  • $\begingroup$ Note @Alan's depends on B or bx not having duplicate indices, which it doesn't in the problem given in the OP. $\endgroup$
    – Michael E2
    Commented Mar 5, 2022 at 19:34
  • 1
    $\begingroup$ FunctionCompile beats them all, but there is an initial long compilation delay. I have posted the code as an extension of this answer. $\endgroup$ Commented Mar 6, 2022 at 4:55
7
$\begingroup$

I just extended @Michael-E2 answer by including FunctionCompile:

lx = Range[40000, 59999];
bx = RandomSample[Range[20000], 10000];
ares = MapAt[#^2 &, lx, {#} & /@ bx]; // AbsoluteTiming
bres = SubsetMap[#^2 &, lx, {#} & /@ bx]; // AbsoluteTiming
(cres = lx; Do[cres[[k]] = cres[[k]]^2, {k, bx}]); // AbsoluteTiming
dres = Compile[{{lx, _Integer, 1}, {bx, _Integer, 1}}, 
     Block[{dres = lx}, Do[dres[[k]] = dres[[k]]^2, {k, bx}];
      dres]][lx, bx]; // AbsoluteTiming
(eres = lx; eres[[bx]] = eres[[bx]]^2); // AbsoluteTiming
ffun = FunctionCompile[
   Function[
    {Typed[lx, TypeSpecifier["PackedArray"]["MachineInteger", 1]], 
     Typed[bx, TypeSpecifier["PackedArray"]["MachineInteger", 1]]}, 
    Block[{fres = lx}, Do[fres[[k]] = fres[[k]]^2, {k, bx}]; fres]]];
fres = ffun[lx, bx]; // AbsoluteTiming

ares == bres == cres == dres == eres == fres

(*
  {1.29098,Null}
  {0.0542107,Null}
  {0.0229712,Null}
  {0.0005062,Null}
  {0.0001768,Null}
  {0.0000713,Null}
  True
*)
$\endgroup$
1
  • $\begingroup$ +1 because updating fres is so fast. Nice. I wonder how much the compile time should be counted. (The compile time is included in dres, but it's fairly negligible.) $\endgroup$
    – Michael E2
    Commented Mar 6, 2022 at 5:30
4
$\begingroup$

Another approach that is supported by the previous versions (11.3 at least) is to create a list of functions which are either (Identity which doesn't touch the data or fn which apply your function) and then MapThread this list with Construct:

Block[{temp},

(* create a list of Identity function *)
temp = ConstantArray[Identity, Length[L]];

(* replace the positions with your function *)
temp[[B]] = #^2 &;

(* apply your list of function to your data (element-wise) *)
result = MapThread[Construct, {temp, L}];

] // MaxMemoryUsed // AbsoluteTiming

(* Out: {0.0161858, 800568} *)

Result:

result == SubsetMap[#^2 &, L, {#} & /@ B]
(* Out: True *)

Another benefit is the memory footprint (could vary):

MapAt[#^2 &, L, {#} & /@ B]; // MaxMemoryUsed // AbsoluteTiming
(* Out: {1.2336, 1912832} *)

(* @Syed answer *)
SubsetMap[#^2 &, L, {#} & /@ B]; // MaxMemoryUsed // AbsoluteTiming
(* Out: {0.0305638, 6541432} *)

which is around 2 and 8 times smaller respectively.

Notes:

  • Minimum version 11.3 is because of Construct, if you happen to have an alternative, the minimum requirement would be lowered (although performance would be affected).
  • This method could easily be extended to apply multiple functions.
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.