# How to use a matrix variable for a compiled function

I want to compile the following function which takes a matrix variable.

ClusterFind =
Compile[{{Inte, _Real}, {Matrix}},
Map[If[Last@# == Inte, Take[#, 4]] &, Matrix], CompilationTarget -> "C",
RuntimeAttributes -> {Listable}, Parallelization -> True];


Here Matrix is a $m\times n$ matrix of real,integer or complex numbers. I also would like to know if it is possible to compile without knowing the dimension of the variable matrix beforehand. Can the Matrix in the above example be of PackedArray data type?

UPDATE

Thanks for the suggestions. I found the following compiled version most efficient as far as speed is concerned.

ClusterFind =
Compile[{{Inte, _Real}, {Matrix, _Complex, 2}},
Take[Select[Matrix, (Last@# == Inte) &], All, 4],
CompilationTarget -> "C"];


BR

-
On my machine, @Oleksandr's version is faster than the one you added in your edit. Is that not so on yours? –  acl Feb 3 '12 at 17:43
@acl I gave sixteen runs for clustering with sixteen integers from 2 to 16 and I took a complex matrix of dimension {750000, 5}. I tried with my function and that of Oleksandr's. Mine took 3.54 sec. and his took 4.98 sec. Funny that u get opposite results. –  PlatoManiac Feb 3 '12 at 19:09
Sorry! I meant 15 runs.. –  PlatoManiac Feb 3 '12 at 19:16
Sorry; I didn't understand from the wording you used in the original question that the matrix had only five columns. I had assumed it would have more, so I extracted the five needed columns into a temporary matrix to avoid copying the larger matrix. As it turns out this step was unnecessary, so yes, your version without the extra operation should indeed be faster. (@acl, if you found mine faster, I guess you also assumed the matrices had more than 5 columns?) –  Oleksandr R. Feb 6 '12 at 23:04
@OleksandrR yes i tried for randomly chosen sizes –  acl Feb 6 '12 at 23:19

Not only can Matrix be a PackedArray, it must be a PackedArray. However, it will be packed for you if necessary before the compiled code is called.

The following code is substantially faster than that given by acl above and does not require any post-processing of the output, but is still sub-optimal in terms of requiring CopyTensor calls and using a rather larger working set than one would think necessary. Perhaps these limitations can be lifted, but after a brief survey of possible implementations I didn't find a way better than this (though note that I didn't try anything with InternalBag).

clusterFind = Compile[{{inte, _Real, 0}, {matrix, _Complex, 2}},
Module[{tmp = matrix[[All, {1, 2, 3, 4, -1}]]},
Select[tmp, Last[#] == inte &][[All, ;; -2]]
], RuntimeAttributes -> Listable, Parallelization -> True
];


An example of the improved timings:

range = {0, 10};
data = RandomInteger[range, {1*^5, 100}];


acl's version:

Timing[
cf4[RandomInteger[range], data];
]


producing: {2.297, Null}

My version:

Timing[
clusterFind[RandomInteger[range], data];
]


which gives {0.047, Null}.

Note that the above timings are not for C-compiled versions of the two functions; compilation to C does not help very much as there is not much computational work to be done in this process anyway and most of the timing consists of copying or extracting parts of tensors. Also, I should mention that RuntimeAttributes -> Listable and Parallelization -> True do not really buy you anything here unless you are operating on a list of matrices.

-
But you should compare to my cf2 (with Select), not cf4. Your version is then still faster, but by 50%. –  acl Feb 3 '12 at 10:06
Yes, you're right. Sorry, I only compared against the compiled one, so I didn't notice that this wasn't the best performing of your suggestions. –  Oleksandr R. Feb 6 '12 at 23:06

This

ClusterFind = Compile[{{Inte, _Real}, {Matrix, _Complex, 2}},
Map[
If[Last@# == Inte, Take[#, 4]] &, Matrix],
CompilationTarget -> "C",
RuntimeAttributes -> {Listable},
Parallelization -> True
];


should take care of the matrix. However, If[Last@# == Inte, Take[#, 4]] does not always return the same type (sometimes it returns nothing, right?)

EDIT1:

something like

cf2[inte_, matrix_] :=
Select[matrix, Last[#] == inte &][[All, 1 ;; 4]]


or

cf3[inte_, matrix_] :=
matrix[[#, 1 ;; 4]] & /@ Position[matrix[[All, -1]], inte]


might do what you want (although I don't know what form your output should have). I don't know if it is fast enough, though. If not perhaps maybe something compilable to C might be set up.

EDIT2:

This works and compiles, but is suboptimal:

cf4 = Compile[{{Inte, _Real}, {Matrix, _Complex, 2}},
Module[
{lst = Most@{1. + I}},
Scan[
If[
Last@# == Inte,
lst = Join[lst, Take[#, 4]]
] &,
Matrix
];
lst
],
CompilationTarget \[Rule] "C"
];


(I stole the Most[{1.+I}] thing from halirutan). You'll need to use Partition[#,4] on the output (if I am guessing your desired output form correctly).

-
Take is compilable. –  rcollyer Feb 2 '12 at 14:49
@rcollyer ah yes! hm... –  acl Feb 2 '12 at 14:51
That list is very nice. :) –  rcollyer Feb 2 '12 at 14:52
@rcollyer I even have a printout of the sorted list in grid format on my office wall. And, now that I think of it, I routinely use Take in my own compiled code. So, I really have no excuse :) –  acl Feb 2 '12 at 14:53
@acl You are right it often will return Null` –  PlatoManiac Feb 2 '12 at 15:02