[Date Index] [Thread Index] [Author Index]

RE: Fast (compiled) routine for element testing and replacement in large matrices?

• To: mathgroup at smc.vnet.net
• Subject: [mg23958] RE: [mg23929] Fast (compiled) routine for element testing and replacement in large matrices?
• From: "David Park" <djmp at earthlink.net>
• Date: Sun, 18 Jun 2000 03:00:41 -0400 (EDT)
• Sender: owner-wri-mathgroup at wolfram.com

> -----Original Message-----
> From: Gareth J. Russell [mailto:russell at cerc.columbia.edu]
To: mathgroup at smc.vnet.net

> Hello,
>
> My first post. I am moving from another software package to
> Mathematica, and
> I often work with very large matrices in the context of simulations, so
> speed is important. Here is a problem for which I would like a
> serious speed
> increase.
>
> Given a matrix of real numbers, I want to generate a new matrix where
> values from a certain range are replaced with 1, and the remainder with
> zero. Matrix size is often on the order of 1000 by 1000 (they are always
> square).
>
> I defined a function that does the comparison using If[],
> and made it listable:
>
> d2wreal[w_,l_,u_] := If[(w>l)&&(D<=u),1,0]
> SetAttributes[d2wreal,Listable]
>
> These are typical times:
>
> a = Table[Random[], {i, 1, 100}, {j, 1, 100}];
> t = Timing[b = d2wreal[a, 0.2, 0.4]][[1]]
> 0.46666666666666856 Second
>
> a = Table[Random[], {i, 1, 1000}, {j, 1, 1000}];
> t = Timing[b = d2wreal[a, 0.2, 0.4]][[1]]
> 65.55 Second
>
> I tried to generate some compiled alternatives myself, but I couldn't find
> one that was faster. I reckon a 10x speed-up should be possible, possibly
> more. My other software takes only 0.75 seconds for the 1000 by
> 1000 matrix
> function.
>
> Can anyone suggest a way to speed this up?
>
> Thanks,
>
> Gareth

Gareth,

I didn't quite understand your d2wreal routine so rewrote it as follows.

op[l_, u_][x_] := If[l <= x <= u, 1, 0]

This generates a random test matrix.

a = Table[Random[], {i, 1, 1000}, {j, 1, 1000}];

This, more or less, corresponds to your method, and on my computer gives:

Map[op[0.2, 0.4], a, {2}]; // Timing
{17.2 Second, Null}

After playing around, I came up with the following method which is about
three times as fast.

Unprotect[LessEqual];
SetAttributes[LessEqual, Listable];
Protect[LessEqual];

LessEqual[0.2, a, 0.4] /. {False -> 0, True -> 1}; // Timing
{5.82 Second, Null}

But with any luck at all, some of the speed demons in MathGroup will come up
with something even better.

David Park
djmp at earthlink.net
http://home.earthlink.net/~djmp/