Problem Summary: The custom Python library PYFEAST (https://github.com/empter/PYFEAST) is a library which allows users to call the Fortran FEAST libraries (https://arxiv.org/pdf/2002.04807.pdf) which efficiently solve eigenvalue and eigenvector problems on massive matrices. In the configuration available on GitHub it is assumed that FEAST was compiled using the Intel-MKL library and so in the setup.py file (https://github.com/empter/PYFEAST/blob/master/cython_feast/setup.py) of PYFEAST the extra_link_args includes various Intel-MKL libraries. This is a problem due to the cripple AMD feature which is built into Intel-MKL (https://en.wikipedia.org/wiki/Math_Kernel_Library) which can no longer be easily bypassed due to Intel removing certain debugging commands in the most recent versions. My goal has been to re-write the compile commands for FEAST and modify the PYFEAST setup.py file so as to make it processor agnostic by using LAPACK and BLAS, or OpenBLAS libraries. This has worked, but the FEAST-Python library generated by the PYFEAST setup.py file does not parallelize and so runs unacceptably slow when dealing with large arrays. If possible, I want to know what I need to change in the compile/make commands as well as the setup.py file in order to make the final FEAST-Python library run with the maximum number of threads available. I know that this is possible since the SciPy eigensolver scipy.linalg automatically parallelizes to the maximum number of threads available to the system when called by Jupyter Notebook. Also, according to the author of PYFEAST the version made with Intel-MKL also automatically parallelizes to the maximum number of threads available on the system when called from Jupyter Notebooks.
Steps Taken so Far: In order to compile FEAST without Intel-MKL (for which no example is provided within the documentation, see page 6 of https://arxiv.org/pdf/2002.04807.pdf) I needed to use the following command:
make F90=gfortran MPI=mpich MKL=no feast
This command (when run in the directory with the FEAST makefile) creates the FEAST library libfeast.a and can also create the library libpfeast.a when the feast argument is replaced with pfeast. As you can see in the above command I am using MPICH for parallelization and gfortran as the FORTRAN compiler.
Next, I compile the Python library. However, in order to do this without Intel-MKL I again need to make several modifications to the code. These modifications to ext_modules=[Extension(...)] are the following: (using OpenBLAS)
Extension("feast",
sources=["feast.pyx"],
include_dirs=[numpy.get_include(),"../../../FEAST/FEAST/4.0/include"],
library_dirs=["../../../FEAST/FEAST/4.0/lib/x64","../../../../../usr/lib/x86_64-linux-gnu"],
libraries=["feast","openblas"],
extra_link_args=['-lopenblas','-lgomp','-lpthread','-lm','-ldl','-lgfortran'],
extra_compile_args=['-lopenblas','-lgomp','-lpthread','-lm','-ldl'],
extra_objects=['../../../../../usr/lib/x86_64-linux-gnu/libopenblas.a']
)
Changes regarding directory library references are due to my version of FEAST being installed in a different way from that of the author of PYFEAST. To understand exactly what I changed see the original code here. I also added the arguments extra_compile_args=[...] and extra_objects=[...] which for some reason where needed in order to get PYFEAST to work without the Intel-MKL libraries (an explanation for why this might be would be appreciated). After doing all of this and running the command python3 setup.py build_ext install I was able to successfully generate the python compatible FEAST library. Note that I didn't make any changes to feast.pyx which is the Cython file from which the setup.py file builds the library feast.pyx file source code.
After getting the library to compile I then started testing it in a Jupyter Notebook. First, I replicated the array from the helloworld.f90 example in FEAST user guide and then ran the FEAST eigensolver using the python commands
obj=feast.HSolver(S,M0=4,Em=0.1,Ex=5.0,which=0)
obj.setfpm(cp=16,eps=14,ml=3,it=1)
This worked perfectly! However, when I attempted to find the eigenvalues of larger matrices it turned out that the FEAST functions were only running on a single thread!
It is not clear to me why this is happening. Although there are some problems with the FEAST code when dealing with sparse matrices (see top of page 7 FEAST User Guide) I don't think these should prevent parallelization for FEAST libraries compiled without Intel-MKL. According to the author the Intel-MKL version automatically parallelizes to all available threads when called from Jupyter Notebooks just like scipy.linalg! I have tried the following commands within Jupyter Notebook to set the OpenBLAS or mpich number of threads: (none of these have changed how the FEAST functions execute when called)
os.environ["OMP_NUM_THREADS"]=12os.environ["OPENBLAS_NUM_THREADS"]=12!export OMP_NUM_THREADS=12!export OPENBLAS_NUM_THREADS=12
Conclusion: I don't know why my version of PYFEAST which avoids using Intel-MKL doesn't parallelize when called. I know other eigensolver programs which use LAPACK and BLAS (like scipy.linalg) do automatically parallelize when called from Jupyter Notebooks and I want to know why and how so I can reproduce that behavior. Is there something wrong with how I compiled FEAST or PYFEAST? Please let me know your thoughts. At this point I just want to make this work just to give a metaphorical middle-finger to Intel. Thank you for your assistance!