# How to compile Map and MapAt

I have a long list (say, million of elements, like this one:

{{{1.0, 2.0}, 3}, {{4.0, 5.0}, 6},...}


I would like to compile a function powq[list_,q_] that would replace this list by a list in which the terms 3,6,.... are replaced by 3^q, 6^q,....:

{{{1.0, 2.0}, 3^q}, {{4.0, 5.0}, 6^q},...}


For a general function f, interactively, I can do it as follows:

mf[list_]:=MapAt[f, list, 2]


Then Map[mf,list] produces the desired result:

{{{1.0, 2.0}, f}, {{4.0, 5.0}, f}}


But I am having problems with defining

powq:=Compile[{{list,_Real,2},{q_Real]}},....]

• Try using Set[] (powq = (* stuff *)) instead of SetDelayed[] (powq := (* stuff *)). – J. M.'s discontentment Sep 29 '12 at 10:49
• Welcome to Mathematica.SE, arkajad! I formatted your code for you, to make it easier for people to read. This page is helpful if you would like to know how to do this yourself. – Verbeia Sep 29 '12 at 11:12
• @arkajad Are you sure, you you have understood the usage of MapAt? Check for instance MapAt[f, {{{9, 2, 2}, 5}, {{5, 9, 4}, 1}, {{3, 9, 8}, 5}, {{4, 4, 5}, 2}, {{5, 6, 0}, 7}}, 2]. Additionally, your array is ragged (meaning it's not a rect. matrix) which does not play well with Compile. – halirutan Sep 29 '12 at 12:30
• Sorry, I skimmed over your question too fast to see your Map/MapAt construct. – halirutan Sep 29 '12 at 14:08

I know you explicitly asked for an answer using Compile but as stated in my comment, I'm not sure this is required. Additionally, I don't think it is possible as you expect it. Your list is a ragged array, which means it is a non-rectangular tensor. To my knowledge it is not possible to use it with Compile. Even the simplest example fails, which does nothing more than returning the input

testF = Compile[{{l, _Real, 2}}, l];
testF[{{{9, 2, 2}, 5}, {{5, 9, 4}, 1}, {{3, 9, 8}, 5}, {{4, 4, 5},
2}, {{5, 6, 0}, 7}}]

(*
During evaluation of In:= CompiledFunction::cfta: Argument
{{{9,2,2},5},{{5,9,4},1},{{3,9,8},5},{{4,4,5},2},{{5,6,0},7}} at position
1 should be a rank 2 tensor of machine-size real numbers. >>

{{{9, 2, 2}, 5}, {{5, 9, 4}, 1}, {{3, 9, 8}, 5}, {{4, 4, 5},
2}, {{5, 6, 0}, 7}}
*)


This would mean you have to extract/transform your array which is IMO (and as you'll see at the end) a great part of the work.

## Different approaches with timings

Here is the one million sample data and a real q

With[{n = 1000000},
data = Transpose[{RandomReal[{0, 10}, {n, 3}],
RandomReal[{0, 10}, n]}];
];
q = 1.2;


First, we measure your approach by substitution a real function f

f[x_] := x^q;
mf[list_] := MapAt[f, list, 2];
First@AbsoluteTiming[Map[mf, data]]
(* 1.289277 *)


Next, lets test a simple approach using Apply and a pure function

Function[{l, b}, {l, b^q}] @@@
{{{9, 2, 2}, 5}, {{5, 9, 4}, 1}, {{3, 9, 8}, 5}, {{4, 4, 5}, 2}, {{5, 6, 0}, 7}}
(*
{{{9,2,2},5^q},{{5,9,4},1},{{3,9,8},5^q},{{4,4,5},2^q},{{5,6,0},7^q}}
*)


Testing this on the sample data

First@AbsoluteTiming[Function[{l, b}, {l, b^q}] @@@ data]
(* 1.145594 *)


A faster one is a rule based approach

First@AbsoluteTiming[data /. {a_, b_} :> {a, b^q}]
(* 0.799231 *)


Then follows a Map approach which is kind of similar to what you did but without explicitly using MapAt. I use a Function instead

First@AbsoluteTiming[Map[{First[#], Last[#]^q} &, data]]
(* 0.654500 *)


Even faster is to transpose the array and separate the numbers you want to process from the rest. Then you use the Listable attribute of taking the power

First@AbsoluteTiming[Transpose[Transpose[data] /. {a_, b_} :> {a, b^q}]]
(* 0.210061 *)


The last attempt gains a lot of speed due to the vectorizaion of b^q. Therefore, we should investigate in this part. Defining a measure function which repeats the measure, we could time how long it takes to raise all numbers of a one million vector to the power q

SetAttributes[measure, {HoldFirst}];
measure[expr_, ntimes_] :=
Mean[Table[First[AbsoluteTiming[expr]], {ntimes}]];

vec = RandomReal[{0, 1}, 1000000];

measure[vec^q, 10]
(* 0.068889 *)


This executes the power 10 times and takes the mean of all execution times. With compile to "C" and the vectorization of this function we may gain some speed. This does of course not take into account that you need time to compile the function.

pow = Compile[{{a, _Real, 0}, {q, _Real, 0}}, a^q,
CompilationTarget -> "C", RuntimeAttributes -> {Listable},
Parallelization -> True, RuntimeOptions -> "Speed"]

measure[pow[vec,q], 10]
(* 0.029345 *)


So we won about 40 milliseconds. This gain can be seen instantly in the overall approach

measure[Transpose[Transpose[data] /. {a_, b_} :> {a, b^q}], 10]
measure[Transpose[Transpose[data] /. {a_, b_} :> {a, pow[b,q]}], 10]

(* 0.203806 *)

(* 0.162934 *)


but I doubt it is worth in a real setting. I believe (but I'm not sure of course) that you should stick with calling the power on the whole vector (whether as Mathematica or as compiled function) because it is pretty fast. This leaves the transformation (which is done here with Transpose and the rule) for further improvement.

Choice of q

It should be noted that the execution speed varies heavily with the choice of the power

measure[vec^2.1, 10]
(* 0.069307 *)

measure[vec^2.0, 10]
(* 0.016818 *)


You pointed out in the comment what you have tried (this should have be placed in the question!). Although it did not work out your basic idea was to sequentially go through the array and raise the elements to the power q. This cannot compete with a vector-operation. Lets assume we have a normal vector and we follow your approach as pointed out in the comment. I use Table to iterate through lst:

pow1 = Compile[{{lst, _Real, 1}, {q, _Real}}, Module[{i = 0},
Table[lst[[i]]^q, {i, Length[lst]}]], CompilationTarget -> "C", RuntimeOptions -> "Speed"]


Compare the execution speeds with the normal Mathematica power and my compiled first pow version

measure[#[vec, q], 10] & /@ {Power, pow, pow1}

(* {0.069090, 0.025614, 0.070855} *)


Such an approach is even slower than the built-in Power!

• Thanks a lot. Yes, I realized that my ragged list causes a problem. For instance I tried fn := Compile[{{lst, _Real, 2}, {q, _Real}}, Module[{copy = lst, i = 1, dummy}, For[i = 1, i <= Length[lst], i++, dummy = copy[[i, 2]]; copy[[i, 2]] = dummy^q]; copy], CompilationTarget -> "C", RuntimeOptions -> "Speed"] and it is not working. I will try your approaches as I will have to call this function lot of times (calculating generalized fractal dimension). – arkajad Sep 29 '12 at 16:52
• Thanks once more. I have tried the second way, that in my comment, later, after reading through the earlier Forum entries. Well, I think the lesson is that I better have two separate uniform lists than one "ragged" list. After all first I am constructing these lists from data, so I can as well construct two separate lists of equal length. Manipulating later on with two lists should not be much more troublesome and memory consuming than with one list - though it will be "less elegant". Thanks for all help. I was happy to receive such quick and serious answers. They helped. – arkajad Sep 29 '12 at 19:27