Mailman 3 March 2010 - 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

np.savetxt: apply patch in enhancement ticket 1079 to add headers?
by Skipper Seabold 08 Jun '10

08 Jun '10

Hi all, I am assuming that this is ok to request via the list... Could we discuss or could someone apply the patch in enhancement ticket 1079? http://projects.scipy.org/numpy/ticket/1079 I needed this functionality recently, and this is a quick and easy fix that may have been overlooked. There is also another enhancement request about this here: http://projects.scipy.org/numpy/ticket/1236 The only thing that I can think of that might need to be added is a test to see that the header length is the same as the number of columns, but really that might just be up to the user to supply the right headers. It might also be nice to have a header = True, that uses the field names for a structured array, but I can live without that. Cheers, Skipper

3 11

namespace recommendation in scipy tutorial and docs for np.emath
by josef.pktd＠gmail.com 15 May '10

15 May '10

I'm going through some scipy tickets and I found http://projects.scipy.org/scipy/ticket/882 I don't know much about the old usage of mixing (mixed up) namespaces. But I think it would be good if someone could check the namespace convention/recommendation in the basic scipy tutorial http://docs.scipy.org/doc/scipy/reference/tutorial/basic.html When I tried to figure out the difference between scipy.* and numpy.*, I also found that the functions in numpy.emath (numpy.scipy) have good docstrings, but they are not included in the docs (at least I don't find them), because http://docs.scipy.org/numpy/docs/numpy-docs/reference/routines.emath.rst/#n… is still a stub page without function links. http://docs.scipy.org/doc/numpy/reference/routines.emath.html Josef

2 2

requesting commit rights
by Ralf Gommers 12 Apr '10

12 Apr '10

Hi, I would like to get commit rights for both scipy and numpy. The things I want to work on and think I need them for: 1. making releases (including related build issues and tests) 2. docs (merging wiki edits, doc machinery for some scipy modules) Thanks, Ralf

2 4

few patches for your consideration
by Yaroslav Halchenko 05 Apr '10

05 Apr '10

Hi Everyone, Being lazy and git spoiled I am exposing tentative BUG: and ENH: commits upon your consideration: http://github.com/yarikoptic/scipy3/commit/0cef22c0342faf8e70d84dd4008dc5ee… BUG: use correct string formatting in raising depending on info from lapack http://github.com/yarikoptic/scipy3/commit/4feb14edcb9d6efdb307fc14c351aab4… ENH: a bit more of information upon failure in LinAlgError's from lapack P.S. If I get a RW commit to SVN I would be happy to dcommit, or if you instruct so I could file 'Pull request' into David's github clone, but since it is just git-svned, not sure if it wouldn't create more problems than would provide benefits. -- .-. =------------------------------ /v\ ----------------------------= Keep in touch // \\ (yoh(a)|www.)onerussian.com Yaroslav Halchenko /( )\ ICQ#: 60653192 Linux User ^^-^^ [175555]

2 2

Speaking of tickets...
by Warren Weckesser 01 Apr '10

01 Apr '10

Just curious: what are the plans for the "Statistic Review" tickets in scipy.stats, from April, 2006? Can any of these be closed? Warren

4 4

GSoC: Speed Up 2to3 in Pattern Matching
by gombiuda JHL 31 Mar '10

31 Mar '10

Hi all, firstly sorry for my length of the proposal. My name is Gombiuda Jiang. Here is my proposal. I am eager to listen to all kind of comments and advices. Thank you for reading my proposal. GSoC: Speed Up 2to3 in Pattern Matching Firstly, thank your for your time to have a glance of my proposal. Idea Background ============================ As PSF pay more and more attention on porting projects using python 2.x to python 3.x, developers in the community start to pay their efforts to port their projects. Luckily, there is a great tool `2to3` helping them in porting job automatically to some extent. However, `2to3` has its limitations. Within these limitations, speed of the conversion can be improved with another implementation to replace the using one. To be more detail, `2to3` is slow because of its using back-tracking of its pattern matching job. Back-tracking is easy to understand and implement, so firstly it's put into used. In addition its implementation is using recursion instead of loop. So as the payload of `2to3` is rocketing up, it eventually shows its weakness of creating more and more frames to allow function calls. As a result, it's necessary to use another algorithm to speed up the pattern matching process. Improvement ============================ As a traditional problem, pattern matching has some efficient algorithms. One of them is named `McNaughton-Yamada-Thompson` algorithm, which bases on NFA(Nondeterministic Finite Automata) theory. This theory can efficiently judge if a string(token stream is also acceptable) can be matched a given regular expression. So it reaches its basic requirements of the pattern matching in `2to3` project. After finishing implementing the algorithm and passing test cases to certificate it works well, we can, if necessary, continue to improve more using DFA theory, which may cost more memory and pattern-compiling time but work faster in matching the patterns. The details is described as below: In original implementation, a pattern string of a fixer is firstly tokenized and then compiled to a pattern syntax tree. Then it starts matching it with the syntax tree of the python 2.x source code with back-tracking. Here we use the original tokenizing result, then put it into a `NFAGenerator` to generate the NFA. At last, we put a node of the syntax tree which is generated from the python 2.x source code into the NFA, and output the result of `True` or `False`. The left process(replacement process) hasn't change. So the changed things are a new `NFAGenerator` class, a new `NFAPattern` class to represent the result of the generation result of the `NFAGenerator` class. Here is the model of the two classes: class NFAPattern(object): """ Attributes: -`content`: The original pattern to be matched -`graph`: A 2D array to store the transition diagram Methods: -`optimize(self)`: Maybe implement later using DFA or some graph reduction algorithms -`match(self, node, result)`: Travels the tree to find out all matched node and put them into result to be passed out """ def __init__(self, content): self.graph = None self.content = content def optimize(self): # Maybe implement later using DFA or some graph reduction algorithms pass def match(self, node, result=None): # Find out all matched node(as root) and put them into result to be passed out pass class NFAGenerator(object): """ Basically the same as the original one `PatternCompiler` except of the implementation of the compile_node() method. """ def __init__(self): # implement like the original one in `PatternCompiler` pass def compile_pattern(self, input, debug=False): tokens = tokenize_wrapper(input) root = self.driver.parse_tokens(tokens, debug=debug) return self.compile_node(root) def compile_node(self, node): # implement base on `McNaughton-Yamada-Thompson` algorithm to compile a pattern string into a NFA pass With this improvement, the process will be more faster than before. But there still contains protential to improve. Right, that's pararrelization. Simplely, we can imagine that we have a bunch of files to be converted. Multithread can be easily used because there is no sequence-relationship between each file. As a result, the changes would be: # Inside refactor.py import threading class RefactoringTool(object): # .... THREADS = 10 def refactor_dir(self, dir_name, write=False, doctests_only=False): """ """ cl = threading.Condition(threading.Lock()) def generate_filenames(dir_name): for dirpath, dirnames, filenames in os.walk(dir_name): for name in filenames: if not name.startswith(".") and name.endswith("py"): fullname = os.path.join(dirpath, name) yield fullname def refacotr_one_file(iter): while True: cl.wait() cl.acquire() fullname = iter.next() cl.release() if fullname: self.refactor_file(fullname, write, doctests_only) else: break iter = generate_filenames(dir_name) threads = [] for i in range(THREADS): threads.append(threading.Thread(target=refacotr_one_file, args=(iter,))) for thread in threads: thread.start() # .... Time Schedule of GSoC ============================ April 27~May 6: Communicate with the mentor to get more requirement information May 7~May 16: Read more documentations to make sure the implementation. May 17~May 26: Write the structure of the new implementation May 27~June 6: Write test cases for testing and write benchmark tool for checking out the efficiency. June 7~July 6: Write the new implementation. July 7~July 11: Testing and bug fixing. July 12~July 16: Mid-term evaluation July 17~July 26: Improvement from NFA to DFA or just reduce the graph of the NFA or add multithread support July 27~August 6: Scrub code and test August 7~August 16: Write documentations August 17: Pencil down About Me ============================ Name: Gombiuda Jiang University: Sun-Yat sen university of China Major: Software Engineering Email/Gtalk: gombiuda(a)gmail.com Knowledge Background: - 9 years programming experience with 6 years focusing on algorithm and data structure learning and 3 years software development experience. - Some experience of researching face detection and recognition using Hidden-Markov-Model algorithm. - Some experience of researching text clustering with artifical intelligence. - Fundermental courses at software school, like `Software Engineering`, `Operating System`, ... Prizes: - The first prize NOIP of primary group district and thus roll in high school with bonus point added - The first prize NOIP of high school group district and thus recommended to the Sun-Yat sen university - The third prize `Software Inovation Competition` at Sun-Yat sen university - The second prize `Commercial Software Inovation Competition` which is hold by Citizen Bank Development Experience: - A RSS news reader using C++ in 2007 - A commercial website in Citizen Bank competition using Java SSH - A educational online-judge website using python and django which is put into use in our software school. - Some small tools like file format convertor, auto voter via Internet, cheater of the Travian game, and so on. Why Me? ============================ I have been using Linux as my default system for years and I like open source very much. I do want to contribute my energy to the community continuously. I have been using python for nearly 2 years and I like it very much for is efficiency. So I read part of the source code of python and got to realize its internal structure. As a result, I'd like to take part in developing the python projects, especially the core and the built-in libraries. And with a strong algorithm ability, I think I can manage to do it! Future View ============================ I would like to join the group of developing python with the help of this project's success. If you see my protential, please let me in, then I will try my best to learn and work hard, not making you disappoint. Finally, thanks again for reading my proposal. -- Sincerely, Gombiuda Jiang

1 0

Evolution of the generic_filter function
by Martin DRUON 31 Mar '10

31 Mar '10

Hi, first, like it is my first post here, I would like to congratulate you for your great job... thank you very much... Now, my question. It seems that the function "scipy.ndimage.generic_filter" can only operate onto scalar objets like float or int. Moreover, it returns only float objects. Do you think that, in the near future, this function will be able to treat generic objects like tuples or user classes ? In that case, we will be able to treat directly, for example, vector fields or color images. We will also be able to apply a function onto multiple arrays at the same time, with the same windowing. Then, to my mind, this method will become very powerful for image processing tools. If this evolution is not scheduled, how I can substitute this method ? For example, I would like to call a function fct(win_a, win_b) with "win_a" a window extracted from an array "a" and "win_b" the same window extracted from "b". Example : # The original fields : two 2D vector fields with 2 components for each vector field_size = (5,5) a = numpy.random.random(field_size + (2,)) b = numpy.random.random(field_size + (2,)) # A simple function : the result is the sum of the element-wise product def fct(win1,win2): return (win1 * win2).sum() # And now, the new generic_filter() which operates simultaneously onto sub-windows of multiple arrays # "c" is the result field : a 2D scalar field with the same size "field_size" c = scipy.ndimage.generic_filter( {"win1":a, "win2":b}, fct, size=(3,3), ...) Thank you. Martin

2 1

Candidates for closing
by Warren Weckesser 31 Mar '10

31 Mar '10

I think these tickets can now be closed: http://projects.scipy.org/scipy/ticket/667 http://projects.scipy.org/scipy/ticket/675 http://projects.scipy.org/scipy/ticket/996 http://projects.scipy.org/scipy/ticket/1002 http://projects.scipy.org/scipy/ticket/1128 Warren

2 2

GSoC: scipy porting to py3k
by gombiuda JHL 30 Mar '10

30 Mar '10

Hi all, I would like to take part in porting the scipy to python3.0 in GSoC. I have been a programmer for 9 years and the first 6 years I focus myself in algorithm and data structure. I use python for almost 2 years, and have read part of the python source code. I think I can manage to do this job. And what I want to ask is, how many people it needs to do the porting? I don't know how to arrange my time in the proposal. Is there anyone give me a hand. Thanks a lot. If possible, please give me some guidelines to follow. Thanks again. -- Sincerely, Gombiuda Jiang

1 0