Mailman 3 January 2021 - SciPy-Dev

Global Optimization Benchmarks
by Andrea Gavana July 11, 2021

July 11, 2021

Dear SciPy Developers & Users, long time no see :-) . I thought to start 2021 with a bit of a bang, to try and forget how bad 2020 has been... So I am happy to present you with a revamped version of the Global Optimization Benchmarks from my previous exercise in 2013. This new set of benchmarks pretty much superseeds - and greatly expands - the previous analysis that you can find at this location: http://infinity77.net/global_optimization/ . The approach I have taken this time is to select as many benchmark test suites as possible: most of them are characterized by test function *generators*, from which we can actually create almost an unlimited number of unique test problems. Biggest news are: 1. This whole exercise is made up of *6,825* test problems divided across *16* different test suites: most of these problems are of low dimensionality (2 to 6 variables) with a few benchmarks extending to 9+ variables. With all the sensitivities performed during this exercise on those benchmarks, the overall grand total number of functions evaluations stands at *3,859,786,025* - close to *4 billion*. Not bad. 2. A couple of "new" optimization algorithms I have ported to Python: - MCS: Multilevel Coordinate Search <https://www.mat.univie.ac.at/~neum/software/mcs/>, it’s my translation to Python of the original Matlab code from A. Neumaier and W. Huyer (giving then for free also GLS and MINQ) I have added a few, minor improvements compared to the original implementation. - BiteOpt: BITmask Evolution OPTimization <https://github.com/avaneev/biteopt> , I have converted the C++ code into Python and added a few, minor modifications. Enough chatting for now. The 13 tested algorithms are described here: http://infinity77.net/go_2021/ High level description & results of the 16 benchmarks: http://infinity77.net/go_2021/thebenchmarks.html Each benchmark test suite has its own dedicated page, with more detailed results and sensitivities. List of tested algorithms: 1. *AMPGO*: Adaptive Memory Programming for Global Optimization: this is my Python implementation of the algorithm described here: http://leeds-faculty.colorado.edu/glover/fred%20pubs/416%20-%20AMP%20(TS)%2… I have added a few improvements here and there based on my Master Thesis work on the standard Tunnelling Algorithm of Levy, Montalvo and Gomez. After AMPGO was integrated in lmfit <https://github.com/lmfit/lmfit-py>, I have improved it even more - in my opinion. 2. *BasinHopping*: Basin hopping is a random algorithm which attempts to find the global minimum of a smooth scalar function of one or more variables. The algorithm was originally described by David Wales: http://www-wales.ch.cam.ac.uk/ BasinHopping is now part of the standard SciPy distribution. 3. *BiteOpt*: BITmask Evolution OPTimization, based on the algorithm presented in this GitHub link: https://github.com/avaneev/biteopt I have converted the C++ code into Python and added a few, minor modifications. 4. *CMA-ES*: Covariance Matrix Adaptation Evolution Strategy, based on the following algorithm: http://www.lri.fr/~hansen/cmaesintro.html http://www.lri.fr/~hansen/cmaes_inmatlab.html#python (Python code for the algorithm) 5. *CRS2*: Controlled Random Search with Local Mutation, as implemented in the NLOpt <http://ab-initio.mit.edu/wiki/index.php/NLopt> package: http://ab-initio.mit.edu/wiki/index.php/NLopt_Algorithms#Controlled_Random_… 6. *DE*: Differential Evolution, described in the following page: http://www1.icsi.berkeley.edu/~storn/code.html DE is now part of the standard SciPy distribution, and I have taken the implementation as it stands in SciPy: https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.differe… 7. *DIRECT*: the DIviding RECTangles procedure, described in: https://www.tol-project.org/export/2776/tolp/OfficialTolArchiveNetwork/NonL… http://ab-initio.mit.edu/wiki/index.php/NLopt_Algorithms#DIRECT_and_DIRECT-L (Python code for the algorithm) 8. *DualAnnealing*: the Dual Annealing algorithm, taken directly from the SciPy implementation: https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.dual_an… 9. *LeapFrog*: the Leap Frog procedure, which I have been recommended for use, taken from: https://github.com/flythereddflagg/lpfgopt 10. *MCS*: Multilevel Coordinate Search, it’s my translation to Python of the original Matlab code from A. Neumaier and W. Huyer (giving then for free also GLS <https://www.mat.univie.ac.at/~neum/software/ls/> and MINQ <https://www.mat.univie.ac.at/~neum/software/minq/>): https://www.mat.univie.ac.at/~neum/software/mcs/ I have added a few, minor improvements compared to the original implementation. See the MCS <http://infinity77.net/go_2021/mcs.html#mcs> section for a quick and dirty comparison between the Matlab code and my Python conversion. 11. *PSWARM*: Particle Swarm optimization algorithm, it has been described in many online papers. I have used a compiled version of the C source code from: http://www.norg.uminho.pt/aivaz/pswarm/ 12. *SCE*: Shuffled Complex Evolution, described in: Duan, Q., S. Sorooshian, and V. Gupta, Effective and efficient global optimization for conceptual rainfall-runoff models, Water Resour. Res., 28, 1015-1031, 1992. The version I used was graciously made available by Matthias Cuntz via a personal e-mail. 13. *SHGO*: Simplicial Homology Global Optimization, taken directly from the SciPy implementation: https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.shgo.ht… List of benchmark test suites: 1. SciPy Extended <http://infinity77.net/go_2021/scipy_extended.html#scipy-extended>: 235 multivariate problems (where the number of independent variables ranges from 2 to 17), again with multiple local/global minima. I have added about 40 new functions to the standard SciPy benchmarks <https://github.com/scipy/scipy/tree/master/benchmarks/benchmarks/go_benchma…> and fixed a few bugs in the existing benchmark models in the SciPy repository. 2. GKLS <http://infinity77.net/go_2021/gkls.html#gkls>: 1,500 test functions, with dimensionality varying from 2 to 6, generated with the super famous GKLS Test Functions Generator <https://arxiv.org/abs/1103.2695>. I have taken the original C code (available at http://netlib.org/toms/) and converted it to Python. 3. GlobOpt <http://infinity77.net/go_2021/globopt.html#globopt>: 288 tough problems, with dimensionality varying from 2 to 5, created with another test function generator which I arbitrarily named “GlobOpt”: https://www.researchgate.net/publication/225566516_A_new_class_of_test_func… . The original code is in C++ and I have bridged it to Python using Cython. *Many thanks* go to Professor Marco Locatelli for providing an updated copy of the C++ source code. 4. MMTFG <http://infinity77.net/go_2021/mmtfg.html#mmtfg>: sort-of an acronym for “Multi-Modal Test Function with multiple Global minima”, this test suite implements the work of Jani Ronkkonen: https://www.researchgate.net/publication/220265526_A_Generator_for_Multimod… . It contains 981 test problems with dimensionality varying from 2 to 4. The original code is in C and I have bridge it to Python using Cython. 5. GOTPY <http://infinity77.net/go_2021/gotpy.html#gotpy>: a generator of benchmark functions using the Bocharov-Feldbaum “Method-Min”, containing 400 test problems with dimensionality varying from 2 to 5. I have taken the Python implementation from https://github.com/redb0/gotpy and improved it in terms of runtime. Original paper from http://www.mathnet.ru/php/archive.phtml?wshow=paper&jrnid=at&paperid=11985&… . 6. Huygens <http://infinity77.net/go_2021/huygens.html#huygens>: this benchmark suite is very different from the rest, as it uses a “fractal” approach to generate test functions. It is based on the work of Cara MacNish on Fractal Functions <https://www.cs.bham.ac.uk/research/projects/ecb/displayDataset.php?id=150>. The original code is in Java, and at the beginning I just converted it to Python: given it was slow as a turtle, I have re-implemented it in Fortran and wrapped it using f2py <http://www.scipy.org/F2py>, then generating 600 2-dimensional test problems out of it. 7. LGMVG <http://infinity77.net/go_2021/lgmvg.html#lgmvg>: not sure about the meaning of the acronym, but the implementation follows the “Max-Set of Gaussians Landscape Generator” described in http://boyuan.global-optimization.com/LGMVG/index.htm . Source code is given in Matlab, but it’s fairly easy to convert it to Python. This test suite contains 304 problems with dimensionality varying from 2 to 5. 8. NgLi <http://infinity77.net/go_2021/ngli.html#ngli>: Stemming from the work of Chi-Kong Ng and Duan Li, this is a test problem generator for unconstrained optimization, but it’s fairly easy to assign bound constraints to it. The methodology is described in https://www.sciencedirect.com/science/article/pii/S0305054814001774 , while the Matlab source code can be found in http://www1.se.cuhk.edu.hk/~ckng/generator/ . I have used the Matlab script to generate 240 problems with dimensionality varying from 2 to 5 by outputting the generator parameters in text files, then used Python to create the objective functions based on those parameters and the benchmark methodology. 9. MPM2 <http://infinity77.net/go_2021/mpm2.html#mpm2>: Implementing the “Multiple Peaks Model 2”, there is a Python implementation at https://github.com/jakobbossek/smoof/blob/master/inst/mpm2.py . This is a test problem generator also used in the smoof <https://github.com/jakobbossek/smoof> library, I have taken the code almost as is and generated 480 benchmark functions with dimensionality varying from 2 to 5. 10. RandomFields <http://infinity77.net/go_2021/randomfields.html#randomfields>: as described in https://www.researchgate.net/publication/301940420_Global_optimization_test… , it generates benchmark functions by “smoothing” one or more multidimensional discrete random fields and composing them. No source code is given, but the implementation is fairly straightforward from the article itself. 11. NIST <http://infinity77.net/go_2021/nist.html#nist>: not exactly the realm of Global Optimization solvers, but the NIST StRD <https://www.itl.nist.gov/div898/strd/nls/nls_main.shtml> dataset can be used to generate a single objective function as “sum of squares”. I have used the NIST dataset as implemented in lmfit <https://github.com/lmfit/lmfit-py/tree/master/NIST_STRD>, thus creating 27 test problems with dimensionality ranging from 2 to 9. 12. GlobalLib <http://infinity77.net/go_2021/globallib.html#globallib>: Arnold Neumaier maintains <https://www.mat.univie.ac.at/~neum/glopt/coconut/Benchmark/Benchmark.html> a suite of test problems termed “COCONUT Benchmark” and Sahinidis has converted the GlobalLib and PricentonLib AMPL/GAMS dataset into C/Fortran code (http://archimedes.cheme.cmu.edu/?q=dfocomp ). I have used a simple C parser to convert the benchmarks from C to Python. The global minima are taken from Sahinidis <http://thales.cheme.cmu.edu/dfo/comparison/comp.html> or from Neumaier or refined using the NEOS server <https://neos-server.org/neos/> when the accuracy of the reported minima is too low. The suite contains 181 test functions with dimensionality varying between 2 and 9. 13. CVMG <http://infinity77.net/go_2021/cvmg.html#cvmg>: another “landscape generator”, I had to dig it out using the Wayback Machine at http://web.archive.org/web/20100612044104/https://www.cs.uwyo.edu/~wspears/… , the acronym stands for “Continuous Valued Multimodality Generator”. Source code is in C++ but it’s fairly easy to port it to Python. In addition to the original implementation (that uses the Sigmoid <https://en.wikipedia.org/wiki/Sigmoid_function> as a softmax/transformation function) I have added a few others to create varied landscapes. 360 test problems have been generated, with dimensionality ranging from 2 to 5. 14. NLSE <http://infinity77.net/go_2021/nlse.html#nlse>: again, not really the realm of Global optimization solvers, but Nonlinear Systems of Equations can be transformed to single objective functions to optimize. I have drawn from many different sources (Publications <https://www.sciencedirect.com/science/article/pii/S0898122111004299>, ALIAS/COPRIN <http://www-sop.inria.fr/coprin/logiciels/ALIAS/Benches/node1.html#SECTION00…> and many others) to create 44 systems of nonlinear equations with dimensionality ranging from 2 to 8. 15. Schoen <http://infinity77.net/go_2021/schoen.html#schoen>: based on the early work of Fabio Schoen and his short note <https://link.springer.com/article/10.1007%2FBF01096734> on a simple but interesting idea on a test function generator, I have taken the C code in the note and converted it into Python, thus creating 285 benchmark functions with dimensionality ranging from 2 to 6. *Many thanks* go to Professor Fabio Schoen for providing an updated copy of the source code and for the email communications. 16. Robust <http://infinity77.net/go_2021/robust.html#robust>: the last benchmark test suite for this exercise, it is actually composed of 5 different kind-of analytical test function generators, containing deceptive, multimodal, flat functions depending on the settings. Matlab source code is available at http://www.alimirjalili.com/RO.html , I simply converted it to Python and created 420 benchmark functions with dimensionality ranging from 2 to 6. Enjoy, and Happy 2021 :-) . Andrea.

6 21

ENH: improve RBF interpolation
by Trever Hines Feb. 28, 2021

Feb. 28, 2021

Hello scipy developers, I would like to contribute code to scipy to address some issues regarding RBF interpolation. The code can be found on my branch here <https://github.com/treverhines/scipy/blob/improve-rbf-interpolation/scipy/i…>. My contribution consists of two new classes for scattered N-D interpolation: 1. `RBFInterpolator`: This is intended to be a replacement for `Rbf` that addresses issues mentioned in 9904 <https://github.com/scipy/scipy/issues/9904> and 4790 <https://github.com/scipy/scipy/issues/4790>. Namely, the major differences with `Rbf` are 1) the usage is similar to `NearestNDInterpolator` and `LinearNDInterpolator` making it easier to swap out different interpolation methods, 2) the sign of the smoothing parameter is correct (see page 10 of these lecture notes <http://pages.stat.wisc.edu/~wahba/stat860public/lect/lect8/lect8.pdf>), and 3) the interpolant includes polynomial terms. For some RBF choices (values of “linear”, “thin_plate”, “cubic”, “quintic”, or “multiquadratic” for `function` in `Rbf`), the additional polynomial terms are needed to ensure that the interpolation problem is well-posed (see theorem 3.2.7 in this document <http://amadeus.math.iit.edu/~fass/603_ch3.pdf>). Without the additional polynomial terms for these RBFs, I have noticed that some values for the smoothing parameter (with the corrected sign) result in an obviously erroneous interpolant. Even when the chosen RBF does not require additional polynomial terms, they still can improve the quality of the interpolant. In particular, the polynomial terms are able to accommodate shifts or linear trends in data, which the RBFs tend to struggle with by themselves. 2. `KNearestRBFInterpolator`: This class performs RBF interpolation using only the k nearest data points to each interpolation point (which was suggested in 5180 <https://github.com/scipy/scipy/issues/5180>). This class is useful when there are too many observations for `RBFInterpolator` (on the order of tens of thousands) and you want an interpolant that *looks* smoother than what you get with `NearestNDInterpolator` or `LinearNDInterpolator`. My concern with interpolation using the k nearest neighbors is that it is a bit of an ad hoc strategy to work around computational limitations. That being said, I have seen a similar strategy used in the Kriging <https://pro.arcgis.com/en/pro-app/latest/help/analysis/geostatistical-analy…> world (Kriging is a form of RBF interpolation). I would appreciate your feedback on whether you think these would be valuable contributions to scipy. If so, I will make the pull request after adding benchmarks, unit tests, and more docs. Thanks, Trever

3 3

NumPy 1.20.0 released
by Charles R Harris Jan. 30, 2021

Jan. 30, 2021

Hi All, On behalf of the NumPy team I am pleased to announce the release of NumPy 1.20.0rc2. This NumPy release is the largest to date, containing some 684 merged pull requests contributed by 184 people. See the list of highlights below. The Python versions supported for this release are 3.7-3.9, support for Python 3.6 has been dropped. Wheels can be downloaded from PyPI <https://pypi.org/project/numpy/1.20.0/>; source archives, release notes, and wheel hashes are available on Github <https://github.com/numpy/numpy/releases/tag/v1.20.0>. Linux users will need pip >= 0.19.3 in order to install manylinux2010 and manylinux2014 wheels. *Highlights* - Annotations for NumPy functions. This work is ongoing and improvements can be expected pending feedback from users. - Wider use of SIMD to increase execution speed of ufuncs. Much work has been done in introducing universal functions that will ease use of modern features across different hardware platforms. This work is ongoing. - Preliminary work in changing the dtype and casting implementations in order to provide an easier path to extending dtypes. This work is ongoing but enough has been done to allow experimentation and feedback. - Extensive documentation improvements comprising some 185 PR merges. This work is ongoing and part of the larger project to improve NumPy's online presence and usefulness to new users. - Further cleanups related to removing Python 2.7. This improves code readability and removes technical debt. - Preliminary support for the upcoming Cython 3.0. *Contributors* A total of 184 people contributed to this release. People with a "+" by their names contributed a patch for the first time. * Aaron Meurer + * Abhilash Barigidad + * Abhinav Reddy + * Abhishek Singh + * Al-Baraa El-Hag + * Albert Villanova del Moral + * Alex Leontiev + * Alex Rockhill + * Alex Rogozhnikov * Alexander Belopolsky * Alexander Kuhn-Regnier + * Allen Downey + * Andras Deak * Andrea Olivo + * Andrew Eckart + * Anirudh Subramanian * Anthony Byuraev + * Antonio Larrosa + * Ashutosh Singh + * Bangcheng Yang + * Bas van Beek + * Ben Derrett + * Ben Elliston + * Ben Nathanson + * Bernie Gray + * Bharat Medasani + * Bharat Raghunathan * Bijesh Mohan + * Bradley Dice + * Brandon David + * Brandt Bucher * Brian Soto + * Brigitta Sipocz * Cameron Blocker + * Carl Leake + * Charles Harris * Chris Brown + * Chris Vavaliaris + * Christoph Gohlke * Chunlin Fang * CloseChoice + * Daniel G. A. Smith + * Daniel Hrisca * Daniel Vanzo + * David Pitchford + * Davide Dal Bosco + * Derek Homeier * Dima Kogan + * Dmitry Kutlenkov + * Douglas Fenstermacher + * Dustin Spicuzza + * E. Madison Bray + * Elia Franzella + * Enrique Matías Sánchez + * Erfan Nariman | Veneficus + * Eric Larson * Eric Moore * Eric Wieser * Erik M. Bray * EthanCJ-git + * Etienne Guesnet + * FX Coudert + * Felix Divo * Frankie Robertson + * Ganesh Kathiresan * Gengxin Xie * Gerry Manoim + * Guilherme Leobas * Hassan Kibirige * Hugo Mendes + * Hugo van Kemenade * Ian Thomas + * InessaPawson + * Isabela Presedo-Floyd + * Isuru Fernando * Jakob Jakobson + * Jakub Wilk * James Myatt + * Jesse Li + * John Hagen + * John Zwinck * Joseph Fox-Rabinovitz * Josh Wilson * Jovial Joe Jayarson + * Julia Signell + * Jun Kudo + * Karan Dhir + * Kaspar Thommen + * Kerem Hallaç * Kevin Moore + * Kevin Sheppard * Klaus Zimmermann + * LSchroefl + * Laurie + * Laurie Stephey + * Levi Stovall + * Lisa Schwetlick + * Lukas Geiger + * Madhulika Jain Chambers + * Matthias Bussonnier * Matti Picus * Melissa Weber Mendonça * Michael Hirsch * Nick R. Papior * Nikola Forró * Noman Arshad + * Paul YS Lee + * Pauli Virtanen * Paweł Redzyński + * Peter Andreas Entschev * Peter Bell * Philippe Ombredanne + * Phoenix Meadowlark + * Piotr Gaiński * Raghav Khanna + * Raghuveer Devulapalli * Rajas Rade + * Rakesh Vasudevan * Ralf Gommers * Raphael Kruse + * Rashmi K A + * Robert Kern * Rohit Sanjay + * Roman Yurchak * Ross Barnowski * Royston E Tauro + * Ryan C Cooper + * Ryan Soklaski * Safouane Chergui + * Sahil Siddiq + * Sarthak Vineet Kumar + * Sayed Adel * Sebastian Berg * Sergei Vorfolomeev + * Seth Troisi * Sidhant Bansal + * Simon Gasse * Simon Graham + * Stefan Appelhoff + * Stefan Behnel + * Stefan van der Walt * Steve Dower * Steve Joachim + * Steven Pitman + * Stuart Archibald * Sturla Molden * Susan Chang + * Takanori H + * Tapajyoti Bose + * Thomas A Caswell * Tina Oberoi * Tirth Patel * Tobias Pitters + * Tomoki, Karatsu + * Tyler Reddy * Veniamin Petrenko + * Wansoo Kim + * Warren Weckesser * Wei Yang + * Wojciech Rzadkowski * Yang Hau + * Yogesh Raisinghani + * Yu Feng * Yuya Unno + * Zac Hatfield-Dodds * Zuhair Ali-Khan + * @abhilash42 + * @danbeibei + * @dojafrat * @dpitch40 + * @forfun + * @iamsoto + * @jbrockmendel + * @leeyspaul + * @mitch + * @prateek arora + * @serge-sans-paille + * @skywalker + * @stphnlyd + * @xoviat * @谭九鼎 + * @JMFT + * @Jack + * @Neal C + Cheers, Charles Harris

1 0

merged scipy.stats.qmc with quasi-Monte Carlo functionality
by Ralf Gommers Jan. 30, 2021

Jan. 30, 2021

Hi all, I think this is worth an announcement: we just merged https://github.com/scipy/scipy/pull/10844 which adds a new submodule with quasi-Monte Carlo functionality, `scipy.stats.qmc`. See http://scipy.github.io/devdocs/stats.qmc.html for details. This was in the works for two years. Thanks to Pamphile Roy and Max Balandat for contributing this module, and to the many reviewers and QMC domain experts (see https://github.com/scipy/scipy/pull/10844#issuecomment-770231427) who pitched in. Cheers, Ralf

1 0

PSF Scientific WG grants
by Stefan van der Walt Jan. 27, 2021

Jan. 27, 2021

Hi all, I'd like to make you aware of this call for funding from the PSF: https://pyfound.blogspot.com/2020/12/psf-scientific-working-group-announces… Funding is for up to 4000 USD. The deadline has been extended (the blog post will soon be updated accordingly); I recommend submitting requests within the next two weeks. Stéfan

1 0

Proposal to add Takagi factorisation
by Marcus P Jan. 23, 2021

Jan. 23, 2021

Hi, I have written some code to perform a Takagi factorisation of complex symmetric matrices (similar to singular value decomposition but it only exists for symmetric matrices and decomposes a matrix A as A = U @ S @ U.T where U is unitary and S is the diagonal matrix of singular values), based on an algorithm published in doi.org/10.1103/PhysRevA.94.062109. I saw that SciPy does not currently have this and I think that it would be a useful addition. If the community is interested, I would be happy to contribute the code. Question about preferred name: It is sometimes known as Autonne—Takagi or Autonne decomposition (this is the name used in Matrix Analysis by Horn and Johnson). It was first discovered by Autonne for non-singular matrices and later rediscovered by Takagi for both singular and non-singular matrices. >From what I can tell, the more common name is Takagi (see e.g. the paper above, arxiv.org/abs/physics/0607103, and doi.org/10.1134/S0965542512010034), so I think that name makes more sense. Are there any opinions about this? Regards, Marcus

1 0

Name for Page's L test
by Warren Weckesser Jan. 21, 2021

Jan. 21, 2021

Hey all, Matt Haberland has implemented Page's L test in the pull request https://github.com/scipy/scipy/pull/12531. I'd like to merge the PR, but Ralf has suggested that the name, `pagel`, is "a terrible name", and has suggested `page_test` or `page_l_test`. (See the comments at https://github.com/scipy/scipy/pull/12531#pullrequestreview-447123727.) Matt's reasoning for `pagel` is that it is consistent with the style used for many other tests in stats. I don't have a strong preference, and when that happens in such a case, I tend to go with the original author's preference. I think Ralf prefers to have "test" in the name. Some existing tests do, but many others don't. Anyone else have an opinion? Either about this specific case, or about the general question of objective criteria for a "good" name? Warren

5 8

Cannot generate very large very sparse random matrix
by Emanuele Olivetti Jan. 19, 2021

Jan. 19, 2021

Hi, I've just noticed that it is not possible to generate very large very sparse random matrices with scipy.sparse.random(). For example: scipy.sparse.random(1_000_000, 1_000_000, density = 1e-11) should create a sparse matrix with only 10 non-zero values... but instead triggers a MemoryError: ---- MemoryError Traceback (most recent call last) <ipython-input-8-eb81d3aec480> in <module> ----> 1 scipy.sparse.random(1_000_000, 1_000_000, density = 1e-11) ~/miniconda3/envs/lap/lib/python3.8/site-packages/scipy/sparse/construct.py in random(m, n, density, format, dtype, random_state, data_rvs) 787 data_rvs = partial(random_state.uniform, 0., 1.) 788 --> 789 ind = random_state.choice(mn, size=k, replace=False) 790 791 j = np.floor(ind * 1. / m).astype(tp, copy=False) mtrand.pyx in numpy.random.mtrand.RandomState.choice() mtrand.pyx in numpy.random.mtrand.RandomState.permutation() MemoryError: Unable to allocate 7.28 TiB for an array with shape (1000000000000,) and data type int64 ---- Here is the problematic line in current master branch of SciPy: https://github.com/scipy/scipy/blob/master/scipy/sparse/construct.py#L806 In short, the issue is due to random_state.choice(... replace=False) which needs to allocate the humongous array in order to pick the ten random numbers... I understand the technical difficulty of generating random numbers without replacement, but it is quite counterintuitive that in order to generate a sparse random matrix it is necessary to create an equally large but *dense* vector first. Is there a solution to this problem? Thanks in advance, Emanuele -- -- Le informazioni contenute nella presente comunicazione sono di natura privata e come tali sono da considerarsi riservate ed indirizzate esclusivamente ai destinatari indicati e per le finalità strettamente legate al relativo contenuto. Se avete ricevuto questo messaggio per errore, vi preghiamo di eliminarlo e di inviare una comunicazione all’indirizzo e-mail del mittente. -- The information transmitted is intended only for the person or entity to which it is addressed and may contain confidential and/or privileged material. If you received this in error, please contact the sender and delete the material.

2 4

Time efficient code to calculate range of Bessel Functions
by Vlad Blazhko Jan. 16, 2021

Jan. 16, 2021

Hello, I would like to contribute time efficient and numerically stable cython code for calculating a range of Bessel's functions of first and second order (jv and yv) to the subpackage special in scipy. It is useful since often you have a sum that uses jv of the one argument, but of several orders that increase by one. They can be efficiently computed through the recurrence relation. While recurrence relation works fine for the yv function, it quickly diverges for the jv function, so for the jv function I have implemented Miller's recurrence algorithm, which is stable. In my experience it is usually faster 5-6 times then using just jv and yv, and such code usually lies in the heart of simulations (in hot places). So I propose the following: Functions names: jv_range, yv_range Example of code: 1) jv_range(v_from=0.3, n=5, z=np.ndarray) -> computes values for orders 0.3, 1.3, 2.3, 3.3, 4.3 and puts them in last axis, i.e. the shape of result is z.shape + (n,) 2) jv_range(v_from=-2.3, n=6, z=np.ndarray) -> computes values for orders -2.3, -1.3, -0.3, 0.7, 1.7, 2.7 If the functions don't have enough accuracy for given orders, then they will fallback to jv and yv functions. Do you think such functions will be a good contribution to the special subpackage? Feel free to give a suggestion on the functions interface And please let me know if you would be so kind to review my PR. Best regards, Vlad Blazhko

1 0

Willing to contribute to SciPy
by 赵刚 Jan. 15, 2021

Jan. 15, 2021

Dear SciPy team, I have two questions when I read SciPy documents and source codes. 1. In the directory "scipy/sparse/linalg/isolve/", I see the implementations of some sparse iterative solver and parameter interfaces with preconditioner, such as CG/PCG, BiCGSTAB/PBiCGSTAB, GMRES/PGMRES, FGMRES/PFGMRES, etc., it's great. But it seems that the detailed preconditioner has not been implemented currently (If multigrid preconditioners had been implemented, please let me know, thanks). So I'd like to ask if I can help you implement some Krylov methods by Multigrid preconditioned (the fastest preconditioner for SPD) to contribute to the SciPy community (From the perspective of necessity and value). 2. On the other hand, I would also like to know if the other parts of SciPy need me to do some improvements and enhancements. I am very willing to make a contribution to SciPy. Best Wishes, Gang Zhao

2 1