issue when calling pyshtools.spectralanalysis.SHReturnTapersMap #444
Comments
|
In case it helps, I can confirm that I still get the error when running within a 'clean' conda environment created as follows: conda create --name shtest -c conda-forge pyshtools which installs as follows: Channels:
Package Planenvironment location: /Users/simon/anaconda3/envs/shtest added / updated specs: The following packages will be downloaded: |
|
Could you let me know if this error occurs every time for the same value of lmax on a given machine? Or if this only occurs sometimes? |
|
It turns out that I get the same error using |
|
I ran 5 iterations on Apple M2 Max. |
|
Here are the results of some tests, for lmax=5-130
So, the error occurs on all machines, but at slightly different values of I am not really sure what this means. The good thing is that it is repeatable on a given machine, so it probably isn't a memory related problem. I suspect that the problem might be related to slightly different ways of compiling LAPACK/BLAS that cause the computations to be performed slightly differently (such as a different optimization, or some option like fast-math). To be honest, I am a little surprised that determining the eigenvectors of a 1000x1000 matrix works at all! For such large matrices, it is possible that there is a better way of doing this than the way I chose. If this is important to you, then perhaps we should look into it. |
|
Thanks for this, I can add the following:
One other thing I noticed - the execution time appears to be significantly different between cases. The above config takes close to double the time of the config with python3.11 / pyshtools from conda. For now, it is helpful to know that the issue is reproducible and that I can get better results with the different python versions and installation from source. (it was a little frustrating to upgrade from an intel Mac only to find the problem appeared worse based on the setup I was using initially on M2) As for importance - my specific use case involves Earth's lithospheric magnetic field. Since LCS-1 is defined to degree 185, it would be nice to be able to do localisations to this degree. However, for the masks that I want to use (different from the test case) I am now able to get results for LMAX > 140, a lot higher than with the conda install. (In one case LMAX = 170 works, the failure point is somewhere between 170 and 185 - impressive indeed). Therefore, I am happy to work around this for now, and don't have a sense of how important this would be for other users. If you do look into this at some point I'd be happy to help run further tests. |
Description of the problem
Following up on question originally asked on shtools gitter...
When using 'pyshtools.spectralanalysis.SHReturnTapersMap', I often run into problems when I set the value of lmax to be large. A simple example to illustrate the problem is as follows:
import numpy as np
import pyshtools
topo_coeffs = pyshtools.datasets.Earth.Earth2014.tbi(lmax=300)
topo = topo_coeffs.expand(extend=False)
mask_dh = pyshtools.SHGrid.from_array(np.array(topo.data < 0).astype(int))
lmax = 70
tapers, eigvals = pyshtools.spectralanalysis.SHReturnTapersMap(mask_dh.data, lmax)
On my current machine (An M2 Mac), my experience is that if I set lmax to be a fairly small number, then the code runs without problem (in this case, lmax=60 is fine). However, with lmax=70 (or higher), the code fails with the following message:
Error --- EigValVecSym
Problem determining eigenvalues and eigenvectors of tridiagonal matrix.
DSTEGR info = 15
SHToolsError Traceback (most recent call last)
Cell In[6], line 10
7 mask_dh = pyshtools.SHGrid.from_array(np.array(topo.data < 0).astype(int))
9 lmax = 70
---> 10 tapers, eigvals = pyshtools.spectralanalysis.SHReturnTapersMap(mask_dh.data, lmax)
File ~/anaconda3/envs/pygmt10/lib/python3.10/site-packages/pyshtools/backends/shtools.py:212, in _raise_errors..wrapped_func(*args, **kwargs)
210 returned_values = func(*args, **kwargs)
211 if returned_values[0] != 0:
--> 212 raise SHToolsError(_shtools_status_message(returned_values[0]))
213 elif len(returned_values) == 2:
214 return returned_values[1]
SHToolsError: Unhandled Fortran 95 error.
Note that I've experienced the same problem on various machines, but the precise value of lmax that triggers the error can vary from one to another.
System information
Mac OS v14.1.2, M2 Max Chip
Python 3.10.13
pyshtools version: 'v4.10.4', installed through conda-forge
typing 'conda list' gives the following entries for LAPACK
liblapack 3.9.0 19_osxarm64_openblas conda-forge
liblapacke 3.9.0 19_osxarm64_openblas conda-forge
(the full list is very long, can post if relevant)
fortran compiler is probably(?) gfortran
libgfortran 5.0.0 13_2_0_hd922786_1 conda-forge
libgfortran-devel_osx-arm64 13.2.0 h5d7a38c_1 conda-forge
libgfortran5 13.2.0 hf226fd6_1 conda-forge
The text was updated successfully, but these errors were encountered: