10
$\begingroup$

I need to generate a Hadamard matrix with dimensions of $4096$, but it seems Mathematica is much slower than MATLAB when generating such matrices.

This program

(MatrixPlot @ HadamardMatrix[4096]) // AbsoluteTiming

cost about 36 seconds on my computer to finish. While for MATLAB,

 H = hadamard(4096); imshow(H)

cost about 0.7 seconds.

I tested the pure process of generating the H matrix. The code is listed below.

Mathematica cost 20 seconds:

HadamardMatrix[4096, WorkingPrecision -> MachinePrecision] // AbsoluteTiming

MATLAB cost 0.12 second:

hadamard(4096)

I was wondering if there is any method to increase the speed or efficiency of Mathematica on this problem?

Jim has provided a way of changing the Method:

HadamardMatrix[4096, Method -> "BitComplement"] // AbsoluteTiming
HadamardMatrix[4096, Method -> "Sequency"] // AbsoluteTiming

It improves the speed a little bit, but still much slower than MATLAB enter image description here

$\endgroup$
  • 2
    $\begingroup$ While I don’t think that this is the cause of the slow down, you still shouldn’t include the plotting function in your timing. Secondarily, consider also that the Hadamard matrix produced by Mathematica is normalized by a factor of $1/\sqrt{n}$, whereas Matlab produces a matrix of $\pm 1$. Having said that, if the timing in both systems is accurate, the difference is very large indeed. $\endgroup$ – MarcoB Jun 28 '18 at 3:09
  • $\begingroup$ Does using HadamardMatrix[4096, WorkingPrecision -> MachinePrecision] alter the timing significantly? (I’m on mobile and can’t test it myself). $\endgroup$ – MarcoB Jun 28 '18 at 3:19
  • 1
    $\begingroup$ Also, Mathematica produces two flavors: Method->"BitComplement" and Method->"Sequency". The former is much faster. Which does Matlab produce? $\endgroup$ – JimB Jun 28 '18 at 3:23
  • 1
    $\begingroup$ @MarcoB,@JimB,It should be right that it should be blamed on the Sqrt[n] here. Method->"BitComplement" and Method->"Sequency" can improve this problem a little bit. However,it seems WorkingPrecision -> MachinePrecision cann't improve the behavoir too much. I' wondering whether there is a solution to force the internal calculation to use machine precision numbers, just as Matlab does. $\endgroup$ – yulinlinyu Jun 28 '18 at 4:55
  • $\begingroup$ @MarcoB change the working precision is not a very effective way. It still uses an awfully long time to produce the result. $\endgroup$ – BNHSX Jun 28 '18 at 5:09
13
$\begingroup$

The MATLAB matrix is similar to the BitComplement type. Avoiding HadamardMatrix makes the construction slightly faster.

r1 = Nest[ArrayFlatten[{{#, #}, {#, -#}}] &, {{1}}, Log2[4096]]; // AbsoluteTiming
r2 = Sign[HadamardMatrix[4096, Method -> "BitComplement",
                         WorkingPrecision -> MachinePrecision]]; // AbsoluteTiming
r1 === r2

{0.14940275, Null}
{0.36388055, Null}
True

Regarding the creation of an image, you should not use MatrixPlot. It constructs a Graphics which is awfully slow. Rather you can create an Image

AbsoluteTiming[
 r1 = Nest[ArrayFlatten[{{#, #}, {#, -#}}] &, {{1}}, Log2[4096]];
 Image[r1]
]

Mathematica graphics

$\endgroup$
  • $\begingroup$ That's much faster! Thanks $\endgroup$ – BNHSX Jun 28 '18 at 10:13
  • 2
    $\begingroup$ This is a little bit faster: Nest[KroneckerProduct[{{1, 1}, {1, -1}}, #] &, {{1}}, BitLength[4096] - 1] $\endgroup$ – J. M. is away Sep 27 '18 at 5:30
6
$\begingroup$

Some spelunking of the code for HadamardMatrix[] reveals why it is so slow for the setting Method -> "Sequency". (Coolwater's answer already shows the method corresponding to Method -> "BitComplement").

Altho the function is internally using a CompiledFunction[] for the purpose, it is slow because of redundant computations to IntegerDigits[]. I have the reorganized the algorithm being used in the code that follows:

hadmat[n_Integer?Positive, prec_: ∞] /; BitAnd[n, n - 1] == 0 := 
   With[{bits = IntegerDigits[Range[0, n - 1], 2, BitLength[n] - 1]},
        ((-1)^(ListCorrelate[{{1, 1}}, Reverse[bits, 2], 1, {0}].Transpose[bits]))/
        Sqrt[N[n, prec]]]

where I have used the classical bit-twiddling test for a power of 2.

Even if it is uncompiled, hadmat[] handily outperforms the built-in:

AbsoluteTiming[t1 = hadmat[4096];]
   {4.54179, Null}

AbsoluteTiming[t2 = HadamardMatrix[4096];]
   {37.9819, Null}

t1 === t2
   True

For giggles, I wrote a compiled routine for generating the Hadamard matrix:

hmc = Compile[{{n, _Integer}, {l, _Integer}},
              With[{bits = Table[IntegerDigits[k - 1, 2, l], {k, n}]},
                   (-1)^(bits.Table[Switch[l - j - k + 1, 0, 1, 1, 1, _, 0],
                                           {j, l}, {k, l}].Transpose[bits])]];

hadm2[n_Integer?Positive, prec_: ∞] /; BitAnd[n, n - 1] == 0 := 
    hmc[n, BitLength[n] - 1]/Sqrt[N[n, prec]]

where convolution is replaced with multiplication by an antitriangular matrix. As expected, this is slower than the method based on ListCorrelate[], but this method is still faster than the built-in.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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