We carry out a series a basic experiments to compare Python related packages (Python, NumPy) and compilers (GNU Fortran, Intel Fortran). We also add Matlab and Java in our study. The primary objective of this exercise is to determine how NumPy performs with respect to the other packages and compilers. All the calculations were carried out in dali.

Problem 1

This example shows the importance of avoiding loops (as far as possible) when manipulating arrays in NumPy. We have a 5000x5000x3 matrix A and we want to perform the following operations on A:

             A(i,j,1) = A(i,j,2)

             A(i,j,3) = A(i,j,1)

             A(i,j,2) = A(i,j,3)

For instance, in Python the code looks like:

for i in range(5000):
    for j in range(5000):
        a[i,j,0] = a[i,j,1]
        a[i,j,2] = a[i,j,0]
        a[i,j,1] = a[i,j,2]

We record the elapsed time needed to do the array assignments. The results are summarized on the table below.

In Fortran, NumPy and Matlab we did not use the loop to carry out the desired operations. We instead use array indexing. Note how slow was Python and how efficient was NumPy.

Problem 2

Here, we want to multiply two randomly generated nxn matrices A and B:

C=AxB

For NumPy and Matlab, we use the predefined matrix multiplication functions whereas in Fortran, we wrote the code to perform the multiplication. The elapsed times presented here only measure the times spent on the multiplication (as the size of the matrix varies).

(*) We did not use any compilation option and only link with "-lmkl_intel_lp64 -lmkl_intel_thread -lmkl_core -liomp5 -lpthread".

The above table suggests that built-in functions are more appropriate to perform matrix multiplication. DGEMM is far more efficient. It is important to note that DGEMM is more suitable for large size matrices. If for instance n=100, the function matmul out performs DGEMM. An interesting discussion on the performance of DGEMM and matmul using the Intel Fortran compiler can be read at:

• How to calculate a multiplication of two matrices efficiently?

NumPy is very slow for matrix multiplication. It was suggested that it should be built against ATLAS BLAS for better performance.

Problem 3

We find the numerical solution of the 2D Laplace equation:

We use the Jacobi iterative solver. We are interested in two finite difference schemes:

The second-order scheme:

and the fourth-order compact scheme:

The Jacobi iterative solver stops when the difference of two consecutive approximations falls below 10^{-6}.

Second-Order Scheme

Fourth-Order Compact Scheme

Once again, we observe that Python is not suitable to manipulate arrays. NumPy and Matlab have comparable results whereas the Intel Fortran compiler displays the best performance. Java did not use array indexing like NumPy, Matlab and Fortran, but did better than NumPy and Matlab.

The results presented above are consistent with the ones done by other groups:

• numerical computing: matlab vs python+numpy+weave

• PerformancePython