# Efficient replacement of conditional nested Do over elements of a matrix

I have two matrices like these:

m1=Table[Random[],{i,10},{j,10}]

m2=Table[Random[],{i,10},{j,10}]


I want to create a new matrix in which its main diagonal elements are zero and off-diagonals are determined based on

Min(m1[[i,j]],m2[[i,j]])


I came up with this Do loop:

newm=ConstantArray[0,{10,10}];
Do[
Do[
If[m1[[i,j]]<m2[[i,j]],newm[[i,j]]=m1[[i,j]],newm[[i,j]]=m2[[i,j]]],
{j,10}],
{i,10}]


What is the better way for doing this and avoiding Do loops?Is it possible to use Map and Which?

-
You can do ReplacePart[ MapThread[Min[{#1, #2}] &, {m1, m2}, 2], {i_, i_} -> 0]. – b.gatessucks Dec 28 '13 at 19:45
a possibility : Table[If[i === j, 0, Min[m1[[i, j]], m2[[i, j]]]], {i, 10}, {j, 10}] – andre Dec 28 '13 at 19:47

If I am understanding you:

SeedRandom[1];

m1 = RandomReal[1, {5, 5}];
m2 = RandomReal[1, {5, 5}];

arrayMin = Function[Null, Min[##], Listable];
arrayMin[m1, m2] (1 - IdentityMatrix[5]) // MatrixForm


$\left( \begin{array}{ccccc} 0 & 0.11142 & 0.59326 & 0.187803 & 0.169013 \\ 0.0657388 & 0 & 0.0118355 & 0.316876 & 0.700474 \\ 0.011978 & 0.391276 & 0 & 0.247495 & 0.727517 \\ 0.790566 & 0.263269 & 0.297514 & 0 & 0.0553108 \\ 0.481571 & 0.128821 & 0.203011 & 0.544772 & 0 \end{array} \right)$

-
Nice, but I note that your solution takes 30% more time than mine ;-P Funny thing, if I combine half of your solution with half of mine the result performs best: (1 - IdentityMatrix[10]) MapThread[Min, {m1, m2}, 2] – Sjoerd C. de Vries Dec 28 '13 at 20:16
..although for larger matrices the SparseArray version is again better. – Sjoerd C. de Vries Dec 28 '13 at 20:31
@Sjoerd I made no attempt to optimize this. Instead I tried to provide a simple yet practical method. +1 on yours too. – Mr.Wizard Dec 29 '13 at 1:12
SparseArray[{{i_, i_} -> 0, {_, _} -> 1}, {10, 10}] MapThread[Min, {m1, m2}, 2]

-
A shorter and slightly faster way to define the SparseArray: SparseArray[Band[{1, 1}] -> 0, {10, 10}, 1] – Mr.Wizard Dec 29 '13 at 3:19
Sjoerd, the question is now undeleted. – R. M. Jan 6 '14 at 0:06

Firstly, you might want to write

m1 = RandomReal[1, {10, 10}]
m2 = RandomReal[1, {10, 10}]


Which should be a lot faster and does not make use of the "outdated" Random[].

Procedural code is easily compilable! You could do this (suboptimal)

cfu =
Compile[
{{m1, _Real, 2}, {m2, _Real, 2}, {nnnn, _Integer}},
Block[
{newm}
,
newm = ConstantArray[0., {nnnn, nnnn}];
Do[
Do[
If[
i != j
,
If[
m1[[i, j]] < m2[[i, j]]
,
newm[[i, j]] = m1[[i, j]]
,
newm[[i, j]] = m2[[i, j]]
]

]
,
{j, nnnn}
]
,
{i, nnnn}

];
newm
]
,
CompilationTarget -> "C"
];


One little issue with this is that it you only need to use one Do. In your code, you could write.

Do[
If[m1[[i, j]] < m2[[i, j]], newm[[i, j]] = m1[[i, j]],
newm[[i, j]] = m2[[i, j]]]
,
{i, 10},
{j, 10}
]


The main improvement we can make however is to use Table instead of Do (thanks to Andre). This is both better in the uncompiled as in the compiled case. The compiled function Andre suggested is

cfu3 =
Compile[
{{m1, _Real, 2}, {m2, _Real, 2}, {nnnn, _Integer}},
Table[
If[i === j, 0.,
Min[m1[[i, j]], m2[[i, j]]]
]
,
{i, nnnn}, {j, nnnn}
],
CompilationTarget -> "C"
]


Comparison with other answers

mmmm = 1*^3;
arrayMin = Function[Null, Min[##], Listable];

(res1 = cfu[m1, m2, mmmm]) // Timing // First
(res5 = cfu3[m1, m2, mmmm]) // Timing // First
(res3 = arrayMin[m1, m2] (1 - IdentityMatrix[mmmm])) //
Timing // First
(res4 = SparseArray[{{i_, i_} -> 0, {_, _} -> 1}, {mmmm,
mmmm}] MapThread[Min, {m1, m2}, 2]) // Timing // First


0.066221 (Me)
0.033331 (Andre)

1.257168 (Mr.Wizard)
0.882195 (Sjoerd)

res1 === res2 === res3 === res4


True

-
At least as fast as cfu2, there is :cfu3 = Compile[{{m1, _Real, 2}, {m2, _Real, 2}, {nnnn, _Integer}}, Table[If[i === j, 0., Min[m1[[i, j]], m2[[i, j]]]], {i, nnnn}, {j, nnnn}], CompilationTarget -> "C"]; – andre Dec 28 '13 at 21:26
@andre ah if that is true then I am really glad. I thought MMA would not be able to figure out that table would be a tensor of integers. If you tested the timings, please tell me. Otherwise I will probably test it myself. – Jacob Akkerboom Dec 28 '13 at 23:44
"tensor of integers" ? I don't understand. I have quickly tested the timing. Cfu3 is between 1 and 3 times faster than cfu2. It depends on the run. – andre 13 mins ago – andre Dec 29 '13 at 9:06
@andre ah I guess I meant reals, huh? :P. Thanks for testing the timings! – Jacob Akkerboom Dec 29 '13 at 10:51
@andre I absorbed your comment in the answer. – Jacob Akkerboom Dec 29 '13 at 12:07

A fast method without Compile:

n = 10;
m1 = RandomReal[1, {n, n}];
m2 = RandomReal[1, {n, n}];

SparseArray[Band[{1, 1}] -> 0, {n, n}, 1] (# m1 + (1 - #) m2) &@UnitStep[m2 - m1]


Its speed is comparable to compiled solutions because it uses packed arrays.

-
I tried this same UnitStep method myself but in version 7 is it not faster than MapThread[Min, {m1, m2}, 2]. How are the timings in v9? – Mr.Wizard Dec 30 '13 at 13:58
@Mr.Wizard For n=1000 it takes about 0.16 sec on my old laptop with v9 and 0.9 sec on university server with v6 (latter unpack arrays). I don't have an access to v7, but I think it unpack arrays too. – ybeltukov Dec 30 '13 at 20:26