first commit

.gitignore

@@ -0,0 +1,6 @@
+*.*~
+*.pyc
+dist
+build
+MANIFEST
+*.egg-info

README

@@ -0,0 +1,57 @@
+======================================================================
+``pyswarm``: Particle swarm optimization (PSO) with constraint support
+======================================================================
+
+The ``pyswarm`` package is a gradient-free optimization package for python
+that supports design constraints.
+
+What's New
+==========
+
+This is the initial release, but it should work out of the box. The syntax is
+similar to the other SciPy optimization routines (the ones like ``fmin_slsqp``).
+Currently, only a single objective is supported, but with any number of
+constraints (or none).
+
+Requirements
+============
+
+- NumPy
+
+Installation and download
+=========================
+
+See the `package homepage`_ for helpful hints relating to downloading
+and installing pyswarm.
+
+Source Code
+===========
+
+The latest, bleeding-edge, but working, `code
+<https://github.com/tisimst/pyDOE/tree/master/pyswarm>`_
+and `documentation source
+<https://github.com/tisimst/pyswarm/tree/master/doc/>`_ are
+available `on GitHub <https://github.com/tisimst/pyswarm/>`_.
+
+Contact
+=======
+
+Any feedback, questions, bug reports, or success stores should
+be sent to the `author`_. I'd love to hear from you!
+
+License
+=======
+
+This package is provided under two licenses:
+
+1. The *BSD License*
+2. Any other that the author approves (just ask!)
+
+References
+==========
+
+- `Particle swarm optimization`_ on Wikipedia
+
+.. _author: mailto:tisimst@gmail.com
+.. _Particle swarm optimization: http://en.wikipedia.org/wiki/Particle_swarm_optimization
+.. _package homepage: http://pythonhosted.org/pyswarm

pyswarm/__init__.py

@@ -0,0 +1,16 @@
+"""
+================================================================================
+pyswarm: Particl swarm optimization (PSO) with constraint support
+================================================================================
+
+Author: Abraham Lee
+Copyright: 2013
+
+"""
+# from __future__ import absolute_import
+
+__author__ = 'Abraham Lee'
+__version__ = '0.5'
+
+from pyswarm.pso import *
+

pyswarm/pso.py

@@ -0,0 +1,250 @@
+import numpy as np
+
+def pso(func, lb, ub, ieqcons=[], f_ieqcons=None, args=(), kwargs={}, 
+        swarmsize=100, omega=0.5, phip=0.5, phig=0.5, maxiter=100, 
+        minstep=1e-8, minfunc=1e-8, debug=False):
+    """
+    Perform a particle swarm optimization (PSO)
+   
+    Parameters
+    ==========
+    func : function
+        The function to be minimized
+    lb : array
+        The lower bounds of the design variable(s)
+    ub : array
+        The upper bounds of the design variable(s)
+   
+    Optional
+    ========
+    ieqcons : list
+        A list of functions of length n such that ieqcons[j](x,*args) >= 0.0 in 
+        a successfully optimized problem (Default: [])
+    f_ieqcons : function
+        Returns a 1-D 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 (Default: None)
+    args : tuple
+        Additional arguments passed to objective and constraint functions
+        (Default: empty tuple)
+    kwargs : dict
+        Additional keyword arguments passed to objective and constraint 
+        functions (Default: empty dict)
+    swarmsize : int
+        The number of particles in the swarm (Default: 100)
+    omega : scalar
+        Particle velocity scaling factor (Default: 0.5)
+    phip : scalar
+        Scaling factor to search away from the particle's best known position
+        (Default: 0.5)
+    phig : scalar
+        Scaling factor to search away from the swarm's best known position
+        (Default: 0.5)
+    maxiter : int
+        The maximum number of iterations for the swarm to search (Default: 100)
+    minstep : scalar
+        The minimum stepsize of swarm's best position (Default: 1e-8)
+    minfunc : scalar
+        The minimum change of swarm's best objective value (Default: 1e-8)
+    debug : boolean
+        If True, progress statements will be displayed every iteration
+        (Default: False)
+   
+    Returns
+    =======
+    g : The swarm's best known position (optimal design)
+   
+    """
+
+    assert len(lb)==len(ub), 'Lower- and upper-bounds must be the same length'
+    assert hasattr(func, '__call__'), 'Invalid function handle'
+    lb = np.array(lb)
+    ub = np.array(ub)
+    assert np.all(ub>lb), 'All upper-bound values must be greater than lower-bound values'
+
+    vhigh = np.abs(ub - lb)
+    vlow = -vhigh
+
+    # Check for constraint function(s) #########################################
+    obj = lambda (x): func(x, *args, **kwargs)
+    if f_ieqcons is None:
+        if not len(ieqcons):
+            if debug:
+                print 'No constraints given.'
+            cons = lambda (x): np.array([0])
+        else:
+            if debug:
+                print 'Converting ieqcons to a single constraint function'
+            cons = lambda (x): np.array([y(x, *args, **kwargs) for y in ieqcons])
+    else:
+        if debug:
+            print 'Single constraint function given in f_ieqcons'
+        cons = lambda (x): np.array(f_ieqcons(x, *args, **kwargs))
+
+    def is_feasible(x):
+        check = np.all(cons(x)>=0)
+        return check
+
+    # Initialize the particle swarm ############################################
+    S = swarmsize
+    D = len(lb)  # the number of dimensions each particle has
+    x = np.random.rand(S, D)  # particle positions
+    v = np.zeros_like(x)  # particle velocities
+    p = np.zeros_like(x)  # best particle positions
+    fp = np.zeros(S)  # best particle function values
+    g = []  # best swarm position
+    fg = 1e100  # artificial best swarm position starting value
+
+    for i in xrange(S):
+        # Initialize the particle's position
+        x[i, :] = lb + x[i, :]*(ub - lb)
+
+        # Initialize the particle's best known position
+        p[i, :] = x[i, :]
+
+        # Calculate the objective's value at the current particle's
+        fp[i] = obj(p[i, :])
+
+        # If the current particle's position is better than the swarm's,
+        # update the best swarm position
+        if fp[i]<fg and is_feasible(p[i, :]):
+            fg = fp[i]
+            g = p[i, :].copy()
+
+        # Initialize the particle's velocity
+        v[i, :] = vlow + np.random.rand(D)*(vhigh - vlow)
+
+    # Iterate until termination criterion met ##################################
+    it = 1
+    while it<=maxiter:
+        rp = np.random.uniform(size=(S, D))
+        rg = np.random.uniform(size=(S, D))
+        for i in xrange(S):
+
+            # Update the particle's velocity
+            v[i, :] = omega*v[i, :] + phip*rp[i, :]*(p[i, :] - x[i, :]) + \
+                      phig*rg[i, :]*(g - x[i, :])
+
+            # Update the particle's position, correcting lower and upper bound 
+            # violations, then update the objective function value
+            x[i, :] = x[i, :] + v[i, :]
+            mark1 = x[i, :]<lb
+            mark2 = x[i, :]>ub
+            x[i, mark1] = lb[mark1]
+            x[i, mark2] = ub[mark2]
+            fx = obj(x[i, :])
+
+            # Compare particle's best position (if constraints are satisfied)
+            if fx<fp[i] and is_feasible(x[i, :]):
+                p[i, :] = x[i, :].copy()
+                fp[i] = fx
+
+                # Compare swarm's best position to current particle's position
+                # (Can only get here if constraints are satisfied)
+                if fx<fg:
+                    if debug:
+                        print 'New best for swarm at iteration %d:'%it, x[i, :], fx
+
+                    tmp = x[i, :].copy()
+                    stepsize = np.sqrt(np.sum((g-tmp)**2))
+                    if np.abs(fg - fx)<=minfunc:
+                        print 'Stopping search: Swarm best objective change less than:', minfunc
+                        return tmp, fx
+                    elif stepsize<=minstep:
+                        print 'Stopping search: Swarm best position change less than:', minstep
+                        return tmp, fx
+                    else:
+                        g = tmp.copy()
+                        fg = fx
+
+        if debug:
+            print 'Best after iteration %d'%it, g, fg
+        it += 1
+
+    print 'Stopping search: maximum iterations reached -->', maxiter
+
+    if g is []:
+        print 'No feasible point found'
+    return g, fg
+
+if __name__=='__main__':
+    print '*'*65
+    print 'Example minimization of 4th-order banana function (no constraints)'
+    def myfunc(x):
+        x1 = x[0]
+        x2 = x[1]
+        return x1**4 - 2*x2*x1**2 + x2**2 + x1**2 - 2*x1 + 5
+
+    lb = [-3, -1]
+    ub = [2, 6]
+
+    xopt1, fopt1 = pso(myfunc, lb, ub)
+
+    print 'The optimum is at:'
+    print '   ', xopt1
+    print 'Optimal function value:'
+    print '    myfunc: ', fopt1
+
+    print '*'*65
+    print 'Example minimization of 4th-order banana function (with constraint)'
+    def mycon(x):
+        x1 = x[0]
+        x2 = x[1]
+        return [-(x1 + 0.25)**2 + 0.75*x2]
+
+    xopt2, fopt2 = pso(myfunc, lb, ub, f_ieqcons=mycon)
+
+    print 'The optimum is at:'
+    print '    ', xopt2
+    print 'Optimal function value:'
+    print '    myfunc: ', fopt2
+    print '    mycon : ', mycon(xopt2)
+
+    print '*'*65
+    print 'Example minimization of twobar truss weight, subject to:'
+    print '  Yield Stress < 100 psi'
+    print '  Yield Stress < Buckling Stress'
+    print '  Deflection   < 0.25 inches'
+    def weight(x, *args):
+        H, d, t = x
+        B, rho, E, P = args
+        return rho*2*np.pi*d*t*np.sqrt((B/2)**2 + H**2)
+
+    def stress(x, *args):
+        H, d, t = x
+        B, rho, E, P = args
+        return (P*np.sqrt((B/2)**2 + H**2))/(2*t*np.pi*d*H)
+
+    def buckling_stress(x, *args):
+        H, d, t = x
+        B, rho, E, P = args
+        return (np.pi**2*E*(d**2 + t**2))/(8*((B/2)**2 + H**2))
+
+    def deflection(x, *args):
+        H, d, t = x
+        B, rho, E, P = args
+        return (P*np.sqrt((B/2)**2 + H**2)**3)/(2*t*np.pi*d*H**2*E)
+
+    def mycons(x, *args):
+        strs = stress(x, *args)
+        buck = buckling_stress(x, *args)
+        defl = deflection(x, *args)
+        return [100 - strs, buck - strs, 0.25 - defl]
+
+    B = 60
+    rho = 0.3
+    E = 30000
+    P = 66
+    args = (B, rho, E, P)
+    lb = [10, 1, 0.01]
+    ub = [30, 3, 0.25]
+    xopt4, fopt4 = pso(weight, lb, ub, f_ieqcons=mycons, args=args)
+
+    print 'The optimum is at:'
+    print '   ', xopt4
+    print 'Optimal function value:'
+    print '    weight         :', fopt4
+    print '    stress         :', stress(xopt4, *args)
+    print '    buckling stress:', buckling_stress(xopt4, *args)
+    print '    deflection     :', deflection(xopt4, *args)
+

setup.py

@@ -0,0 +1,42 @@
+import os
+from setuptools import setup
+
+def read(fname):
+    return open(os.path.join(os.path.dirname(__file__), fname)).read()
+
+setup(
+    name='pywarm',
+    version="0.5",
+    author='Abraham Lee',
+    author_email='tisimst@gmail.com',
+    description='Particle swarm optimization (PSO) with constraint support',
+    url='https://github.com/tisimst/pyswarm',
+    license='BSD License',
+    long_description=read('README'),
+    packages=['pyswarm'],
+    install_requires=['numpy'],
+    keywords=[
+        'PSO',
+        'particle swarm optimization',
+        'optimization',
+        'python'
+        ],
+    classifiers=[
+        'Development Status :: 5 - Production/Stable',
+        'Intended Audience :: Education',
+        'Intended Audience :: Science/Research',
+        'License :: OSI Approved :: BSD License',
+        'Operating System :: OS Independent',
+        'Programming Language :: Python',
+        'Programming Language :: Python :: 2.6',
+        'Programming Language :: Python :: 2.7',
+        'Topic :: Education',
+        'Topic :: Scientific/Engineering',
+        'Topic :: Scientific/Engineering :: Mathematics',
+        'Topic :: Scientific/Engineering :: Physics',
+        'Topic :: Software Development',
+        'Topic :: Software Development :: Libraries',
+        'Topic :: Software Development :: Libraries :: Python Modules',
+        'Topic :: Utilities'
+        ]
+    )