What would you ask Mathematica to do on a big system?

Question

Our hardware supplier has given us the opportunity to test Mathematica on a big Windows system. The specs will be 1 to 4 TB of working memory and 12 sockets with 18 cores each. I got 1 full week of computer time to do my testing. The question is what should I ask Mathematica to do because it is more a test of Mathematica and not a test of the hardware.

I could load a complete database in memory using the association structure and then test search and consolidate questions.

Another option is to run an algorithm that benefits from parallel computations.

The goal is to write a test plan on what to test and how to test it.

I welcome your suggestions on what to evaluate in Mathematica to measure how useful many processors and almost unlimited working memory is?

Kind regards,

Pieter van Bijnen

note: it is one system (no grid involved)

Answer 1

I've always wondered about the scalability of MathLink (now officially "Wolfram Symbolic Transfer Protocol"). This is the protocol used by Mathematica to communicate between the front end and the kernel, and the basis of the Parallel` package. It has quite low bandwidth and high latency relative to, for example, MPI libraries. I also wonder how many MathLink connections the library can handle before it falls apart completely. This seems an ideal opportunity to test that, and so to explore the limits of Mathematica's scaling. You will need a great number of kernel licences to be able to do this, so you will probably have to ask WRI if they will grant you a temporary licence with unlimited kernels to be able to do these tests. I think it is in their interest to know the answer, so hopefully they will.

In the Parallel` package, connections are made from a master kernel to each slave kernel (one per slave). But, for message passing implementations, it would be necessary to have one connection between each pair of kernels that could potentially wish to communicate with each other. I wouldn't be surprised if you can start 216 slave kernels (although, I also wouldn't be too shocked if this is already too many). But can you open another $216^2/2 = 23328$ MathLink connections? And if so, how quickly can you transfer data from each kernel to all the others?

I wrote some code that sets up all of the connections in this answer. The notebook is here.

Assuming you can open the connections, you might like to try some all-to-all messaging. In MPI, typically benchmarks are done for message sizes increasing as powers of two. Having run the code in the other answer so that the connections are established and the required definitions are in place, you can write

alltoall[n_Integer] := ParallelEvaluate[
   LinkWrite[Rank[$AllRanks], ConstantArray[$MaxMachineNumber, n]];
   LinkFlush@Rank[$AllRanks];
   While[! And @@ LinkReadyQ@Rank[$AllRanks], Pause[0.001]];
   LinkRead@Rank[$AllRanks];,
   Rank[$AllRanks]
  ];

Table[{2^n, First@AbsoluteTiming@alltoall[2^n]}, {n, 0, 12}]

to get a set of timings.

---

Bug?

Unfortunately, on my computer (64-bit Windows), the LinkWrite call hangs when n goes above 12 (message size more than 32KB but less than 64KB), no matter how many kernels/connections are involved. This is curious, because LinkWrite is supposed to be non-blocking, and MathLink is meant to be thread-safe, as far as I understand it. The deadlock is definitely in MathLink rather than my code, because even with random delays inserted anywhere in the sequence, it always completes correctly for messages smaller than the critical size. What is even stranger is that this size depends on whether the links have LinkProtocol -> "SharedMemory" or LinkProtocol -> "TCPIP". In the latter case it sometimes works for n up to 15, but also sometimes hangs using 100% CPU (a livelock).

I am not sure if this also happens on other platforms, but I have reported it to WRI under [CASE:3521773]. A self-contained notebook to reproduce it with just two kernels (not the arbitrary number proposed in this post) is available here.

---

Not a bug; just poor documentation.

Despite many claims in the various sources of documentation about MathLink that LinkWrite is non-blocking, I have received a response today from WRI stating that this behavior is actually expected. The solution, as I figured out during a conversation with halirutan, is to know that LinkWrite will block until LinkRead is posted at the other end of the link, if the message is big enough. This seems not to be documented anywhere, as far as I could discern, and neither halirutan nor Szabolcs had heard of it either.

But the good thing is, when LinkRead and LinkWrite are used simultaneously, an arbitrarily large(?) message can be sent. With this knowledge it is at least possible to design a working all-to-all collective communication. I will do so and update this answer in due course, so that the 17 upvotes will not have been in vain.

Answer 2

I am sure you can easily install also Linux on it and then you could contact Vladyslav Shtabovenko, the current maintainer of FeynCalc (https://github.com/vsht) and ask him about hard problems in High Energy Physics he would like to benchmark on such a King-Kong machine. Either him or somebody else could also provide you with more complicated examples of FIRE (http://science.sander.su/FIRE-Examples.htm).

Answer 3

I would choose a problem that exploits the unique power of Mathematica, in particular the natural functions involving graph theory, symbolic math, graphics, and the high compute power you have available. So I would choose some image recognition and clustering problem such as: Take some large number of images ($\sim\!\!10^8$) and perform deep learning on them but form clusters and find relationships that can be described by graphs. You can break the problem into parallel tasks, then merge appropriately.

Answer 4

4 TB memory? N-body of course;

dir = 2;

T = 10;
k = 1;
l = 10; n = l^3;
v = 6;

q[i_] := (-1)^i
m[i_] := RandomReal[5]
Rem[o_] := l (o - IntegerPart[o])

eqns = Table[
   {D[Subscript[r, i][t], {t, 2}] m[i] == Sum[
      Normalize[Subscript[r, i][t] - Subscript[r, j][t]] k q[
        i] q[j]/(Subscript[r, i][t] - 
           Subscript[r, j][t]).(Subscript[r, i][t] - 
           Subscript[r, j][t]), {j, Delete[Table[k, {k, 1, n}], i]}
      ],
    Subscript[r, i][0] == {Rem[IntegerPart[(i - 1)/l]/l], 
      Rem[(i - 1)/l], IntegerPart[(i - 1)/l^2]},
    Subscript[r, i]'[0] == 
     v Normalize[
       Subscript[r, If[i == n - dir, n, Mod[i + dir, n]]][0] - 
        Subscript[r, i][0]]},
   {i, 1, n}
   ];
sol = NDSolve[eqns, Table[Subscript[r, i][t], {i, 1, n}], {t, 0, T}, 
    AccuracyGoal -> 5, 
    Method -> {"EquationSimplification" -> "Residual"}]; // Timing
Show[
 Table[
  ParametricPlot3D[Subscript[r, i][t] /. sol, {t, 0, T}, 
   PlotStyle -> Hue[1 - i/(4 n)]], {i, 1, n}
  ],
 ListPointPlot3D[
  Flatten[Table[Subscript[r, i][t] /. sol /. t -> 0, {i, 1, n}], 1], 
  PlotStyle -> {Gray, PointSize[0.02]}],
 PlotRange -> All, Background -> Black, BoxRatios -> Automatic, 
 Boxed -> False, Axes -> False
 ]
Beep[]

Try to up l, play with the definitions of forces and strain your machine. Here is the result with l=4:

Answer 5

What I am most interested in is described in Oleksandr's answer. Here's something else I would also try:

There are several built-in functions that take advantage of parallelization without any special settings (and without using the Parallel tools framework). I would like to know how well these scale to a high number of cores.

Examples:

Matrix operations, such as addition and matrix multiplication typically take advantage of multiple cores. I would try how Mathematica performs on tasks such as eigenvalue computations for dense and sparse matrices for very large matrix sizes. Large matrices take a lot of storage space, so you'll also make use of the memory. For large matrices it'll be realistic to compute only the few largest eigenvalues.

I noticed that NIntegrate can also use multiple cores. I would try how much speedup can be achieved on slow-to-compute integrals (e.g. multidimensional ones).

The point of this test would be to see how much speedup more cores can bring. Thus you would need to adjust the number of threads Mathematica uses for a computation and see how many times faster operations complete on 200 cores than 10 cores (i.e. how much less than 20 times).

The following system options are interesting:

In[130]:= SystemOptions["ParallelOptions"]

Out[130]= {"ParallelOptions" -> {"AbortPause" -> 2., 
   "BusyWait" -> 0.01, "MachineFunctionParallelThreshold0" -> 32768, 
   "MachineFunctionParallelThreshold1" -> 32768, 
   "MachineFunctionParallelThreshold2" -> 16384, 
   "MachineFunctionParallelThreshold3" -> 16384, 
   "MathLinkTimeout" -> 15., "MKLThreadNumber" -> 4, 
   "ParallelThreadNumber" -> 4, "RecoveryMode" -> "Retry", 
   "RelaunchFailedKernels" -> True}}

I believe that both "MKLThreadNumber" and "ParallelThreadNumber" are relevant. See SetSystemOptions on how to change them. I am not sure about their exact effect, but if you decide to do the test we can figure it out together.

---

Example benchmark:

We'll use this function to set the number of cores:

setThreadCount[
  threads_] := (SetSystemOptions[
   "ParallelOptions" -> "MKLThreadNumber" -> threads];
  SetSystemOptions[
   "ParallelOptions" -> "ParallelThreadNumber" -> threads];)

Eigenvalues for a dense symmetric matrix.

In[151]:= a = RandomReal[1, {5000, 5000}];
a = a + Transpose[a];

I have a laptop with a 4-core hyperthreading CPU. It appears as an 8-core one to the OS, but it only has 4 physical cores.

In[163]:= Table[
 setThreadCount[t];
 {t, First@AbsoluteTiming[Eigenvalues[a, 20];]},
 {t, 1, 8}
 ]

Out[163]= {{1, 12.2941}, {2, 8.18082}, {3, 6.32726}, {4, 6.22842}, {5,
   6.42646}, {6, 6.43579}, {7, 6.53777}, {8, 6.60184}}

Benchmarking is hard. Part of what we see here is due to the fact that as the CPU heats up it automatically throttles its clock rate. During the last few calculations it got hot so it has to slow down. If I repeat the 4-core test after waiting for it to cool down, I see better performance:

In[164]:= setThreadCount[4]

In[165]:= First@AbsoluteTiming[Eigenvalues[a, 20];]
Out[165]= 5.7838

AFAIK this won't be a problem with servers like the one you are testing.

But do pay attention to using only AbsoluteTiming and not Timing to measure parallel calculations. Timing measures the time used by each core separately then adds them up ...

Answer 6

To take advantage of that kind of memory, you really want to do some parallel processing. Mathematica's parallel processing focuses on data-parallelism, or more simply, embarrassingly parallel problems. So you might try out various Monte Carlo simulations.

I don't recommend trying to reproduce full MPI functionality with LinkCreate et al. In the case mentioned, it sounds like all the links would be made to be SharedMemory. SharedMemory links handle indexing information through files in /tmp/MathLink (on Linux and OSX at least) and I suspect there might be some overhead there though this was improved greatly in 10.3.1. Ironically, you won't run into this issue with TCPIP connections though those are slower for more obvious reasons. Also, depending on the physical memory organization, the shared memory links between processes on different sockets may be significantly slower than the shared memory links between processes on the same socket, so obfuscating the architecture in an MPI-like fashion here may lead to oddities in performance.

As a last point, Mathematica's generic parallelization relies on the underlying operating system to appropriately give processes to cores on processors. On most modern laptops and desktops, this is straight forward. You have a chip with 8 logical (not including hyperthreading) cores and you launch eight subprocesses. The operating system sees those processes as active and distributes them to the 8 cores and they can run in parallel. I've never managed these high-socket count server systems before so I am not familiar with their process management. It might be more complex.

EDIT: Sorry to put this into my post as I don't have the reputation yet to leave a comment on someone else's post. Pertaining to Szalbolcs post below, if I am remembering correctly, MKLThreadNumber is the number of threads launched by the MKL library and ParallelThreadNumber is an override for when you are using subkernels.

Answer 7

The more I've thought about this question, the more my answer (above) has changed. Now my answer is this: Assuming there is a high-bandwidth connection to the Wolfram server, choose a problem that relies on Mathematica's superior handling of curated data. Create an enormous problem that relies on curated financial data, geographic data, biological data, historical data, etc. Some problem that relies on search through combinatorically difficult conditions. (Other languages will never be able to deal with such a problem as well as Mathematica.)

A silly example, off the top of my head:

• Find two major cities whose mutual distance in meters is equal to their difference in population, and which traveling from the larger city to the smaller city is shortest by heading eastward, but where their average temperature difference in July 1950 is the same as in June 2015.

Silly, yes... but you get the idea.

Imagine trying to solve such a problem in Matlab or C++, or ...

Answer 8

Pieter, here is one suggestion: - try a FEM calculation on a huge and complicatd structure, eg. stress around clusters of cracks :)

Answer 9

You could just try to find the next largest prime number

http://www.iflscience.com/editors-blog/largest-ever-prime-number-found-gimps

https://www.youtube.com/watch?v=tlpYjrbujG0