Mailman 3 May 2011 - SciPy-Dev

SLSQP Constrained Optimizer Status
by Rob Falck 24 Nov '13

24 Nov '13

I'm currently implementing the Sequential Least Squares Quadratic Programming (SLSQP) optimizer, by Dieter Kraft, for use in Scipy. The Fortran code being wrapped with F2PY is here: http://www.netlib.org/toms/733 (its use within Scipy has been cleared) SLSQP provides for bounds on the independent variables, as well as equality and inequality constraint functions, which is a capability that doesn't exist in scipy.optimize. Currently the code works, although the constraint normals are being generated by approximation. I'm working on a way to pass these in. I think the most elegant way will be a single function that returns the matrix of constraint normals. For a demonstration of what the code can do, here is an optimization of f(x,y) = 2xy + 2x - x**2 - 2y**2 Example 14.2 in Chapra & Canale gives the maximum as x=2.0, y=1.0. The unbounded optimization tests find this solution. As expected, its faster when derivatives are provided rather than approximated. Unbounded optimization. Derivatives approximated. Elapsed time: 1.45792961121 ms Results [[1.9999999515712266, 0.99999996181577444], -1.9999999999999984, 4, 0, 'Optimization terminated successfully.'] Unbounded optimization. Derivatives provided. Elapsed time: 1.03211402893 ms Results [[1.9999999515712266, 0.99999996181577444], -1.9999999999999984, 4, 0, 'Optimization terminated successfully.'] The following example uses an equality constraint to find the optimal when x=y. Bound optimization. Derivatives approximated. Elapsed time: 1.384973526 ms Results [[0.99999996845920858, 0.99999996845920858], -0.99999999999999889, 4, 0, 'Optimization terminated successfully.'] I've tried to conform to the syntax used by the other optimizers in scipy.optimize. The function definition and doc string are below. If anyone is interested in testing it out, let me know. def fmin_slsqp( func, x0 , eqcons=[], f_eqcons=None, ieqcons=[], f_ieqcons=None, bounds = [], fprime = None, fprime_cons=None,args = (), iter = 100, acc = 1.0E-6, iprint = 1, full_output = 0, epsilon = _epsilon ): """ Minimize a function using Sequential Least SQuares Programming Description: Python interface function for the SLSQP Optimization subroutine originally implemented by Dieter Kraft. Inputs: func - Objective function (in the form func(x, *args)) x0 - Initial guess for the independent variable(s). eqcons - A list of functions of length n such that eqcons[j](x0,*args) == 0.0 in a successfully optimized problem f_eqcons - A function of the form f_eqcons(x, *args) that returns an array in which each element must equal 0.0 in a successfully optimized problem. If f_eqcons is specified, eqcons is ignored. ieqcons - A list of functions of length n such that ieqcons[j](x0,*args) >= 0.0 in a successfully optimized problem f_ieqcons - A function of the form f_ieqcons(x0, *args) that returns an array in which each element must be greater or equal to 0.0 in a successfully optimized problem. If f_ieqcons is specified, ieqcons is ignored. bounds - A list of tuples specifying the lower and upper bound for each independent variable [ (xl0, xu0), (xl1, xu1), ...] fprime - A function that evaluates the partial derivatives of func fprime_cons - A function of the form f(x, *args) that returns the m by n array of constraint normals. If not provided, the normals will be approximated. Equality constraint normals precede inequality constraint normals. args - A sequence of additional arguments passed to func and fprime iter - The maximum number of iterations (int) acc - Requested accuracy (float) iprint - The verbosity of fmin_slsqp. iprint <= 0 : Silent operation iprint == 1 : Print summary upon completion (default) iprint >= 2 : Print status of each iterate and summary full_output - If 0, return only the minimizer of func (default). Otherwise, output final objective function and summary information. epsilon - The step size for finite-difference derivative estimates. Outputs: ( x, { fx, gnorm, its, imode, smode }) x - The final minimizer of func. fx - The final value of the objective function. its - The number of iterations. imode - The exit mode from the optimizer, as an integer. smode - A string describing the exit mode from the optimizer. Exit modes are defined as follows: -1 : Gradient evaluation required (g & a) 0 : Optimization terminated successfully. 1 : Function evaluation required (f & c) 2 : Number of equality constraints larger than number of independent variables 3 : More than 3*n iterations in LSQ subproblem 4 : Inequality constraints incompatible 5 : Singular matrix E in LSQ subproblem 6 : Singular matrix C in LSQ subproblem 7 : Rank-deficient equality constraint subproblem HFTI 8 : Positive directional derivative for linesearch 9 : Iteration limit exceeded -- - Rob Falck

4 16

Platform-dependent tests?
by Warren Weckesser 01 Jun '11

01 Jun '11

Jeff Armstrong has enhanced the Schur decomposition function scipy.linalg.schur to expose the ability to sort the eigenvalues (more accurately, group them according to a boolean function); see https://github.com/scipy/scipy/pull/23. schur() is a wrapper for the DGEES Fortran library function. According to its documentation, if the matrix is poorly scaled, the eigenvalues can fail to satisfy the sorting condition *after* they've been sorted, because of rounding errors. The function returns an error code when this occurs. Jeff's code checks for this condition, and he also added a couple unit tests for it. Unfortunately, the tests fail on my computer--that is, the error condition does not happen. This is not surprising, as the error condition relies on the behavior of rounding errors. So, to the question: is there a way to make these tests platform-dependent, so that they only run on a platform/architecture/whatever where the tests are known to trigger the desired condition? Or is the result likely to depend on so many other factors (compiler and optimization settings, third party math library, phase of the moon, etc) that it is hopeless to try? Warren

2 1

ANN: SfePy 2011.2
by Robert Cimrman 31 May '11

31 May '11

I am pleased to announce release 2011.2 of SfePy. Description ----------- SfePy (simple finite elements in Python) is a software for solving systems of coupled partial differential equations by the finite element method. The code is based on NumPy and SciPy packages. It is distributed under the new BSD license. Home page: http://sfepy.org Mailing lists, issue tracking: http://code.google.com/p/sfepy/ Git (source) repository: http://github.com/sfepy Documentation: http://docs.sfepy.org/doc Highlights of this release -------------------------- - experimental implementation of terms aiming at easier usage and definition of new terms - Mooney-Rivlin membrane term - update build system to use exclusively setup.py - allow switching boundary conditions on/off depending on time - support for variable time step solvers For more information on this release, see http://sfepy.googlecode.com/svn/web/releases/2011.2_RELEASE_NOTES.txt (full release notes, rather long and technical). Best regards, Robert Cimrman and Vladimír Lukeš

1 0

where and how I can submit the translation.
by 江大伟 29 May '11

29 May '11

*Hi guys* I have done some translation on the first part of the Scipy Reference Guide. But I do not know where and how I can submit.

2 1

numpy,scipy installation using mkl
by akshar bhosale 29 May '11

29 May '11

Hi, Is this right forum for doubts abt numpy/scipy installation? Please find our issue below .. we have machine having intel xeon x7350 processors(8 nos) and RHEL 5.2 x86_64 with kernel 2.6.18-92.el5. We have following configuration : /opt/intel/Compiler/11.0/069/ mkl/lib/em64t Now we want to install numpy and scipy as an user in my home directory. Following are the libraries build inside MKL. libfftw2x_cdft_DOUBLE.a libmkl_blacs_sgimpt_ilp64.a libmkl_intel_ilp64.a libmkl_pgi_thread.so libmkl_vml_mc2.so libfftw2xc_intel.a libmkl_blacs_sgimpt_lp64.a libmkl_intel_ilp64.so libmkl_scalapack.a libmkl_vml_mc3.so libfftw2xf_intel.a libmkl_blas95.a libmkl_intel_lp64.a libmkl_scalapack_ilp64.a libmkl_vml_mc.so libfftw3xc_intel.a libmkl_cdft.a libmkl_intel_lp64.so libmkl_scalapack_ilp64.so libmkl_vml_p4n.so libfftw3xf_intel.a libmkl_cdft_core.a libmkl_intel_sp2dp.a libmkl_scalapack_lp64.a locale libmkl_blacs_ilp64.a libmkl_core.a libmkl_intel_sp2dp.so libmkl_scalapack_lp64.so mkl77_blas.mod libmkl_blacs_intelmpi20_ilp64.a libmkl_core.so libmkl_intel_thread.a libmkl_sequential.a mkl77_lapack1.mod libmkl_blacs_intelmpi20_lp64.a libmkl_def.so libmkl_intel_thread.so libmkl_sequential.so mkl77_lapack.mod libmkl_blacs_intelmpi_ilp64.a libmkl_em64t.a libmkl_lapack95.a libmkl.so mkl95_blas.mod libmkl_blacs_intelmpi_ilp64.so libmkl_gf_ilp64.a libmkl_lapack.a libmkl_solver.a mkl95_lapack.mod libmkl_blacs_intelmpi_lp64.a libmkl_gf_ilp64.so libmkl_lapack.so libmkl_solver_ilp64.a mkl95_precision.mod libmkl_blacs_intelmpi_lp64.so libmkl_gf_lp64.a libmkl_mc3.so libmkl_solver_ilp64_sequential.a libmkl_blacs_lp64.a libmkl_gf_lp64.so libmkl_mc.so libmkl_solver_lp64.a libmkl_blacs_openmpi_ilp64.a libmkl_gnu_thread.a libmkl_p4n.so libmkl_solver_lp64_sequential.a libmkl_blacs_openmpi_lp64.a libmkl_gnu_thread.so libmkl_pgi_thread.a libmkl_vml_def.so Version we are trying to build of numpy and scipy are : numpy : 1.6.0b2 and scipy : 0.9.0 We have configured python as a user in my home directory with version : 2.6.6 our machine has python : 2.4.3 We want to know the exact procedure for installing it using MKL as we are facing lot of issues while installing the same. Please help us.

1 0

Entropy from empirical high-dimensional data
by Gael Varoquaux 26 May '11

26 May '11

Hi list, I am looking at estimating entropy and conditional entropy from data for which I have only access to observations, and not the underlying probabilistic laws. With low dimensional data, I would simply use an empirical estimate of the probabilities by converting each observation to its quantile, and then apply the standard formula for entropy (for instance using scipy.stats.entropy). However, I have high-dimensional data (~100 features, and 30000 observations). Not only is it harder to convert observations to probabilities in the empirical law, but I am also worried of curse of dimensionality effects: density estimation in high-dimension is a difficult problem. Does anybody has advices, or code in Python to point to, for this task? Cheers, Gaël

2 3

Used "Automatic Merge" on my github pull request...
by Warren Weckesser 26 May '11

26 May '11

I used the "Automatic merge" button on my own pull request on github, but afterwards discovered that it uses --no-ff, so my single commit also resulted in a "merge" commit: https://github.com/scipy/scipy/commit/cf04a2b8dd4cf258413687ec146883ea5ab19… Should I try to get rid of that merge? If so, how? (The git book is by my side, but I suspect an answer will show up here faster than I can find it.) Warren

3 2

Question about subpackage/submodule API
by Warren Weckesser 23 May '11

23 May '11

The recent discussion about what to call the module in Pim Schellart's nice Lomb-Scargle contribution brought up a question for me that I have not seen discussed before. In many cases, a subpackage contains many modules, and the assorted functions and classes in these modules are made available in the (sub)package-level namespace by importing them from within __init__.py. For example, scipy/stats contains distributions.py, kde.py, morestats.py, mstats.py, rv.py, stats.py and vonmises.py. Except for mstats.py, the contents of all these are imported into the __init__.py namespace and included in __all__, so users can say, for example, >>> from scipy.stats import bayes_mvs instead of >>> from scipy.stats.more_stats import bayes_mvs Is this an *intentional* definition of an API? For example, is it "wrong" for a user to refer to the submodule 'more_stats' explicitly when importing? I think the answer is yes, because I don't think any promise is being made that the submodule won't be renamed or refactored as scipy development progresses. But I'd like to be sure that that is everyone else's understanding. If that is not the case, then changes that I have made in the past have violated the deprecation policy. For example, a while back I moved the functions in signal/filter_design.py that were related to FIR filters to their own module, fir_filter_design.py. By making appropriate changes to signal/__init__.py, all the function were still importable from scipy.signal. So, from my point of view, I didn't change the public API and no deprecation was necessary. But if scipy.signal.filter_design is part of the public API, then I should have made sure that something like >>> from scipy.signal.filter_design import firwin still worked after the refactor. Note: I'm not saying that *all* submodules should be private and have their objects exposed only through __init__.py. I'm just looking for some clarification of existing policy. Warren

6 12

Comments on optimize.newton function
by Gökhan Sever 23 May '11

23 May '11

Hello, Could we add "tolf" argument into the newton function signature? I am guessing this should match the tolf argument in IDL newton ( http://star.pst.qub.ac.uk/idl/NEWTON.html) TOLF: Set the convergence criterion on the function values. The default value is 1.0 x 10-4. def newton(func, x0, fprime=None, args=(), tol=1.48e-8, maxiter=50, tolf= 1.e-4) In second method part of the newton function (from https://github.com/scipy/scipy/blob/master/scipy/optimize/zeros.py) # Secant method p0 = x0 if x0 >= 0: p1 = x0*(1 + 1e-4) + 1e-4 else: p1 = x0*(1 + 1e-4) - 1e-4 without increasing the tolf or (1e-4) in these statements I can't get a proper root solution for my function. The reason that I am experiencing with newton is because fsolve seems slower comparing to the newton for scalar root finding for a given function. Consider this example run Sage v4.6.1 notebook: %cython cpdef double myfunc(double x): return x**3 + 2*x - 1 timeit('scipy.optimize.newton(myfunc,1)') 625 loops, best of 3: 22.1 µs per loop timeit('scipy.optimize.fsolve(myfunc,1)') 625 loops, best of 3: 86.5 µs per loop I am also experimenting to Cythonize the Newton secant method, which the Cython written version shows significant speed-ups comparing to the Python version. Without going any further, I would like to know if there is any Cythonized code around for newton or any other approach to make fsolve faster? Thanks. -- Gökhan

3 20

scikit-morph
by Nathan Faggian 23 May '11

23 May '11

Greetings! I would like to start a project working on image registration algorithms in python. I am proposing a new sci-kit called scikit-morph because I believe registration is too specific for inclusion in scipy directly. My short term goals are to: 1) implement a basic image registration framework - perhaps similar to the design used by ITK. 2) Implement image resampling methods for use in the optimisation (c++). 3) Use py.test to thoroughly test all code. I have started by setting up a github project, here: http://github.com/nfaggian/scikit-morph How I can get this going in a way that meshes well with existing scikit structures? I have looked at both the image and learn scikits, is it a matter of simply using the same directory structures?? Any help on getting a skeleton project going would be great - also very open to people joining in - please email me directly if you would like to help out! Kind Regards, Nathan Faggian

5 7