def __call__(self, x, y, id=None, best=0, k=False):
     super(VerboseLoggingMonitor,self).__call__(x, y, id, best, k=k)
     if self._vyinterval is not numpy.inf and \
        int((self._step-1) % self._vyinterval) == 0:
         if not list_or_tuple_or_ndarray(y):
             who = ''
             y = " %f" % self._ik(self._y[-1], k)
         elif self._all:
             who = ''
             y = " %s" % self._ik(self._y[-1], k)
         else:
             who = ' best'
             y = " %f" % self._ik(self._y[-1][best], k)
         msg = "Generation %d has%s Chi-Squared:%s" % (self._step-1,who,y)
         if id is not None: msg = "[id: %d] " % (id) + msg
         print(msg)
     if self._vxinterval is not numpy.inf and \
        int((self._step-1) % self._vxinterval) == 0:
         if not list_or_tuple_or_ndarray(x):
             who = ''
             x = " %f" % self._x[-1]
         elif self._all:
             who = ''
             x = "\n %s" % self._x[-1]
         else:
             who = ' best'
             x = "\n %s" % self._x[-1][best]
         msg = "Generation %d has%s fit parameters:%s" % (self._step-1,who,x)
         if id is not None: msg = "[id: %d] " % (id) + msg
         print(msg)
     return

def generate_penalty(conditions, ptype=None, **kwds):
    """Converts a penalty constraint function to a mystic.penalty function.

Inputs:
    conditions -- a penalty constraint function, or list of constraint functions
    ptype -- a mystic.penalty type, or a list of mystic.penalty types
        of the same length as the given conditions

    For example:
        >>> constraints = '''
        ...     x2 = x0/2.
        ...     x0 >= 0.'''
        >>> ineqf,eqf = generate_conditions(constraints, nvars=3)
        >>> penalty = generate_penalty((ineqf,eqf))
        >>> penalty([1.,2.,0.])
        25.0
        >>> penalty([1.,2.,0.5])
        0.0

Additional Inputs:
    k -- penalty multiplier
    h -- iterative multiplier
"""
    # allow for single condition, list of conditions, or nested list
    if not list_or_tuple_or_ndarray(conditions):
        conditions = list((conditions,))
    else: pass #XXX: should be fine...
    conditions = list(flatten(conditions))

    # allow for single ptype, list of ptypes, or nested list
    if ptype is None:
        ptype = []
        from mystic.penalty import quadratic_equality, quadratic_inequality
        for condition in conditions:
            if 'inequality' in condition.__name__: 
                ptype.append(quadratic_inequality)
            else:
                ptype.append(quadratic_equality)
    elif not list_or_tuple_or_ndarray(ptype):
        ptype = list((ptype,))*len(conditions)
    else: pass #XXX: is already a list, should be the same len as conditions
    ptype = list(flatten(ptype))

    # iterate through penalties, building a compound penalty function
    pf = lambda x:0.0
    pfdoc = ""
    for penalty, condition in zip(ptype, conditions):
        pfdoc += "%s: %s\n" % (penalty.__name__, condition.__doc__)
        apply = penalty(condition, **kwds)
        pf = apply(pf)
    pf.__doc__ = pfdoc.rstrip('\n')
    pf.__name__ = 'penalty'
    return pf

 def restore(variables, mystring):
     if list_or_tuple_or_ndarray(variables):
         vars = get_variables(mystring,'_')
         indices = [int(v.strip('_')) for v in vars]
         for i in range(len(vars)):
             mystring = mystring.replace(vars[i],variables[indices[i]])
     return mystring

def get_variables(constraints, variables='x'):
    """extract a list of the string variable names from constraints string

Inputs:
    constraints -- a string of symbolic constraints, with one constraint
        equation per line. Constraints can be equality and/or inequality
        constraints. Standard python syntax should be followed (with the
        math and numpy modules already imported).

    For example:
        >>> constraints = '''
        ...     x1 + x2 = x3*4
        ...     x3 = x2*x4'''
        >>> get_variables(constraints)
        ['x1', 'x2', 'x3', 'x4'] 

Additional Inputs:
    variables -- desired variable name. Default is 'x'. A list of variable
        name strings is also accepted for when desired variable names
        don't have the same base, and can include variables that are not
        found in the constraints equation string.

    For example:
        >>> constraints = '''              
        ...     y = min(u,v) - z*sin(x)
        ...     z = x**2 + 1.0 
        ...     u = v*z'''
        >>> get_variables(constraints, list('pqrstuvwxyz'))
        ['u', 'v', 'x', 'y', 'z']
"""
    if list_or_tuple_or_ndarray(variables):
        equations = replace_variables(constraints,variables,'_')
        vars = get_variables(equations,'_')
        indices = [int(v.strip('_')) for v in vars]
        varnamelist = []
        from numpy import sort
        for i in sort(indices):
            varnamelist.append(variables[i])
        return varnamelist

    import re
    target = variables+'[0-9]+'
    varnamelist = []
    equation_list = constraints.splitlines()
    for equation in equation_list:
        vars = re.findall(target,equation)
        for var in vars:
            if var not in varnamelist:
                varnamelist.append(var)
    return varnamelist

def _solve_nonlinear(constraints, variables='x', target=None, **kwds):
    """Build a constraints function given a string of nonlinear constraints.
Returns a constraints function. 

Inputs:
    constraints -- a string of symbolic constraints, with one constraint
        equation per line. Constraints must be equality constraints only.
        Standard python syntax should be followed (with the math and numpy
        modules already imported).

    For example:
        >>> constraints = '''x1 = x3*3. + x0*x2'''
        >>> print _solve_nonlinear(constraints)
        x0 = (x1 - 3.0*x3)/x2
        >>> constraints = '''
        ...     spread([x0,x1]) - 1.0 = mean([x0,x1])   
        ...     mean([x0,x1,x2]) = x2'''
        >>> print _solve_nonlinear(constraints)
        x0 = -0.5 + 0.5*x2
        x1 = 0.5 + 1.5*x2

Additional Inputs:
    variables -- desired variable name. Default is 'x'. A list of variable
        name strings is also accepted for when desired variable names
        don't have the same base, and can include variables that are not
        found in the constraints equation string.
    target -- list providing the order for which the variables will be solved.
        If there are "N" constraint equations, the first "N" variables given
        will be selected as the dependent variables. By default, increasing
        order is used.

    For example:
        >>> constraints = '''
        ...     spread([x0,x1]) - 1.0 = mean([x0,x1])   
        ...     mean([x0,x1,x2]) = x2'''
        >>> print _solve_nonlinear(constraints, target=['x1'])
        x1 = -0.833333333333333 + 0.166666666666667*x2
        x0 = -0.5 + 0.5*x2

Further Inputs:
    locals -- a dictionary of additional variables used in the symbolic
        constraints equations, and their desired values.
"""
    nvars = None
    permute = False # if True, return all permutations
    warn = True  # if True, don't supress warnings
    verbose = False # if True, print details from _classify_variables
    #-----------------------undocumented-------------------------------
    permute = kwds['permute'] if 'permute' in kwds else permute
    warn = kwds['warn'] if 'warn' in kwds else warn
    verbose = kwds['verbose'] if 'verbose' in kwds else verbose
    #------------------------------------------------------------------
    if target in [None, False]:
        target = []
    elif isinstance(target, str):
        target = target.split(',')
    else:
        target = list(target) # not the best for ndarray, but should work

    from mystic.symbolic import replace_variables, get_variables
    if list_or_tuple_or_ndarray(variables):
        if nvars is not None: variables = variables[:nvars]
        constraints = replace_variables(constraints, variables, '_')
        varname = '_'
        ndim = len(variables)
    else:
        varname = variables # varname used below instead of variables
        myvar = get_variables(constraints, variables)
        if myvar: ndim = max([int(v.strip(varname)) for v in myvar]) + 1
        else: ndim = 0
    if nvars is not None: ndim = nvars

    # create function to replace "_" with original variables
    def restore(variables, mystring):
        if list_or_tuple_or_ndarray(variables):
            vars = get_variables(mystring,'_')
            indices = [int(v.strip('_')) for v in vars]
            for i in range(len(vars)):
                mystring = mystring.replace(vars[i],variables[indices[i]])
        return mystring

    locals = kwds['locals'] if 'locals' in kwds else None
    if locals is None: locals = {}

    eqns = constraints.splitlines()
    # Remove empty strings:
    actual_eqns = []
    for j in range(len(eqns)):
        if eqns[j].strip():
           actual_eqns.append(eqns[j].strip())
    orig_eqns = actual_eqns[:]

    neqns = len(actual_eqns)

    xperms = [varname+str(i) for i in range(ndim)]
    if target:
        [target.remove(i) for i in target if i not in xperms]
        [target.append(i) for i in xperms if i not in target]
        _target = []
        [_target.append(i) for i in target if i not in _target]
#.........这里部分代码省略.........

def _solve_linear(constraints, variables='x', target=None, **kwds):
    """Solve a system of symbolic linear constraints equations.

Inputs:
    constraints -- a string of symbolic constraints, with one constraint
        equation per line. Constraints must be equality constraints only.
        Standard python syntax should be followed (with the math and numpy
        modules already imported).

    For example:
        >>> constraints = '''
        ...     x0 - x2 = 2.
        ...     x2 = x3*2.'''
        >>> print _solve_linear(constraints)
        x2 = 2.0*x3
        x0 = 2.0 + 2.0*x3

Additional Inputs:
    variables -- desired variable name. Default is 'x'. A list of variable
        name strings is also accepted for when desired variable names
        don't have the same base, and can include variables that are not
        found in the constraints equation string.
    target -- list providing the order for which the variables will be solved.
        If there are "N" constraint equations, the first "N" variables given
        will be selected as the dependent variables. By default, increasing
        order is used.

    For example:
        >>> constraints = '''
        ...     x0 - x2 = 2.
        ...     x2 = x3*2.'''
        >>> print _solve_linear(constraints, target=['x3','x2'])
        x3 = -1.0 + 0.5*x0
        x2 = -2.0 + x0

Further Inputs:
    locals -- a dictionary of additional variables used in the symbolic
        constraints equations, and their desired values.
"""
    nvars = None
    permute = False # if True, return all permutations
    warn = True  # if True, don't supress warnings
    verbose = False  # if True, print debug info
    #-----------------------undocumented-------------------------------
    permute = kwds['permute'] if 'permute' in kwds else permute
    warn = kwds['warn'] if 'warn' in kwds else warn
    verbose = kwds['verbose'] if 'verbose' in kwds else verbose
    #------------------------------------------------------------------
    if target in [None, False]:
        target = []
    elif isinstance(target, str):
        target = target.split(',')
    else:
        target = list(target) # not the best for ndarray, but should work

    from mystic.symbolic import replace_variables, get_variables
    if list_or_tuple_or_ndarray(variables):
        if nvars is not None: variables = variables[:nvars]
        _constraints = replace_variables(constraints, variables, '_')
        varname = '_'
        ndim = len(variables)
        for i in range(len(target)):
            if variables.count(target[i]):
                target[i] = replace_variables(target[i],variables,markers='_')
    else:
        _constraints = constraints
        varname = variables # varname used below instead of variables
        myvar = get_variables(constraints, variables)
        if myvar: ndim = max([int(v.strip(varname)) for v in myvar]) + 1
        else: ndim = 0
    if nvars is not None: ndim = nvars

    # create function to replace "_" with original variables
    def restore(variables, mystring):
        if list_or_tuple_or_ndarray(variables):
            vars = get_variables(mystring,'_')
            indices = [int(v.strip('_')) for v in vars]
            for i in range(len(vars)):
                mystring = mystring.replace(vars[i],variables[indices[i]])
        return mystring

    # default is _locals with sympy imported
    _locals = {}
    locals = kwds['locals'] if 'locals' in kwds else None
    if locals is None: locals = {}
    # if sympy not installed, return original constraints
    try:
        code = """from sympy import Eq, Symbol;"""
        code += """from sympy import solve as symsol;"""
        code = compile(code, '<string>', 'exec')
        exec code in _locals
    except ImportError: # Equation will not be simplified."
        if warn: print "Warning: sympy not installed."
        return constraints

    # default is _locals with numpy and math imported
    # numpy throws an 'AttributeError', but math passes error to sympy
    code = """from numpy import *; from math import *;""" # prefer math
    code += """from numpy import mean as average;""" # use np.mean not average
    code += """from numpy import var as variance;""" # look like mystic.math
#.........这里部分代码省略.........

def _solve_single(constraint, variables='x', target=None, **kwds):
    """Solve a symbolic constraints equation for a single variable.

Inputs:
    constraint -- a string of symbolic constraints. Only a single constraint
        equation should be provided, and must be an equality constraint. 
        Standard python syntax should be followed (with the math and numpy
        modules already imported).

    For example:
        >>> equation = "x1 - 3. = x0*x2"
        >>> print _solve_single(equation)
        x0 = -(3.0 - x1)/x2

Additional Inputs:
    variables -- desired variable name. Default is 'x'. A list of variable
        name strings is also accepted for when desired variable names
        don't have the same base, and can include variables that are not
        found in the constraints equation string.
    target -- list providing the order for which the variables will be solved.
        By default, increasing order is used.

    For example:
        >>> equation = "x1 - 3. = x0*x2"
        >>> print _solve_single(equation, target='x1')
        x1 = 3.0 + x0*x2

Further Inputs:
    locals -- a dictionary of additional variables used in the symbolic
        constraints equations, and their desired values.
""" #XXX: an very similar version of this code is found in _solve_linear XXX#
    # for now, we abort on multi-line equations or inequalities
    if len(constraint.replace('==','=').split('=')) != 2:
        raise NotImplementedError, "requires a single valid equation" 
    if ">" in constraint or "<" in constraint:
        raise NotImplementedError, "cannot simplify inequalities" 

    nvars = None
    permute = False # if True, return all permutations
    warn = True  # if True, don't supress warnings
    verbose = False  # if True, print debug info
    #-----------------------undocumented-------------------------------
    permute = kwds['permute'] if 'permute' in kwds else permute
    warn = kwds['warn'] if 'warn' in kwds else warn
    verbose = kwds['verbose'] if 'verbose' in kwds else verbose
    #------------------------------------------------------------------
    if target in [None, False]:
        target = []
    elif isinstance(target, str):
        target = target.split(',')
    else:
        target = list(target) # not the best for ndarray, but should work

    from mystic.symbolic import replace_variables, get_variables
    if list_or_tuple_or_ndarray(variables):
        if nvars is not None: variables = variables[:nvars]
        constraints = replace_variables(constraint, variables, markers='_')
        varname = '_'
        ndim = len(variables)
        for i in range(len(target)):
            if variables.count(target[i]):
                target[i] = replace_variables(target[i],variables,markers='_')
    else:
        constraints = constraint # constraints used below
        varname = variables # varname used below instead of variables
        myvar = get_variables(constraint, variables)
        if myvar: ndim = max([int(v.strip(varname)) for v in myvar]) + 1
        else: ndim = 0
    if nvars is not None: ndim = nvars

    # create function to replace "_" with original variables
    def restore(variables, mystring):
        if list_or_tuple_or_ndarray(variables):
            vars = get_variables(mystring,'_')
            indices = [int(v.strip('_')) for v in vars]
            for i in range(len(vars)):
                mystring = mystring.replace(vars[i],variables[indices[i]])
        return mystring

    # default is _locals with sympy imported
    _locals = {}
    locals = kwds['locals'] if 'locals' in kwds else None
    if locals is None: locals = {}
    try:
        code = """from sympy import Eq, Symbol;"""
        code += """from sympy import solve as symsol;"""
        code = compile(code, '<string>', 'exec')
        exec code in _locals
    except ImportError: # Equation will not be simplified."
        if warn: print "Warning: sympy not installed."
        return constraint

    # default is _locals with numpy and math imported
    # numpy throws an 'AttributeError', but math passes error to sympy
    code = """from numpy import *; from math import *;""" # prefer math
    code += """from numpy import mean as average;""" # use np.mean not average
    code += """from numpy import var as variance;""" # look like mystic.math
    code += """from numpy import ptp as spread;"""   # look like mystic.math
    code = compile(code, '<string>', 'exec')
    exec code in _locals
#.........这里部分代码省略.........

def _classify_variables(constraints, variables='x', nvars=None): 
    """Takes a string of constraint equations and determines which variables
are dependent, independent, and unconstrained. Assumes there are no duplicate
equations. Returns a dictionary with keys: 'dependent', 'independent', and
'unconstrained', and with values that enumerate the variables that match
each variable type.

Inputs:
    constraints -- a string of symbolic constraints, with one constraint
        equation per line. Constraints must be equality constraints only.
        Standard python syntax should be followed (with the math and numpy
        modules already imported).

    For example:
        >>> constraints = '''
        ...     x0 = x4**2
        ...     x2 = x3 + x4'''
        >>> _classify_variables(constraints, nvars=5)
        {'dependent':['x0','x2'], 'independent':['x3','x4'], 'unconstrained':['x1']}
        >>> constraints = '''
        ...     x0 = x4**2
        ...     x4 - x3 = 0.
        ...     x4 - x0 = x2'''
        >>> _classify_variables(constraints, nvars=5)
        {'dependent': ['x0','x2','x4'], 'independent': ['x3'], 'unconstrained': ['x1']}

Additional Inputs:
    nvars -- number of variables. Includes variables not explicitly
        given by the constraint equations (e.g. 'x1' in the example above).
    variables -- desired variable name. Default is 'x'. A list of variable
        name strings is also accepted for when desired variable names
        don't have the same base, and can include variables that are not
        found in the constraints equation string.
"""
    if ">" in constraints or "<" in constraints:
        raise NotImplementedError, "cannot classify inequalities" 

    from mystic.symbolic import replace_variables, get_variables
    #XXX: use solve? or first if not in form xi = ... ?
    if list_or_tuple_or_ndarray(variables):
        if nvars is not None: variables = variables[:nvars]
        constraints = replace_variables(constraints, variables)
        varname = '$'
        ndim = len(variables)
    else:
        varname = variables # varname used below instead of variables
        myvar = get_variables(constraints, variables)
        if myvar: ndim = max([int(v.strip(varname)) for v in myvar]) + 1
        else: ndim = 0
    if nvars is not None: ndim = nvars

    eqns = constraints.splitlines()
    indices = range(ndim)
    dep = []
    indep = []
    for eqn in eqns: # find which variables are used
        if eqn:
            for var in range(ndim):
                if indices.count(var) != 0:
                    if eqn.find(varname + str(var)) != -1:
                        indep.append(var)
                        indices.remove(var)
    indep.sort()
    _dep = []
    for eqn in eqns: # find which variables are on the LHS
        if eqn:
            split = eqn.split('=')
            for var in indep:
                if split[0].find(varname + str(var)) != -1:
                    _dep.append(var)
                    indep.remove(var)
                    break
    _dep.sort()
    indep = _dep + indep # prefer variables found on LHS
    for eqn in eqns: # find one dependent variable per equation
        _dep = []
        _indep = indep[:]
        if eqn:
            for var in _indep:
                if eqn.find(varname + str(var)) != -1:
                    _dep.append(var)
                    _indep.remove(var)
        if _dep:
            dep.append(_dep[0])
            indep.remove(_dep[0])
    #FIXME: 'equivalent' equations not ignored (e.g. x2=x2; or x2=1, 2*x2=2)
    """These are good:
    >>> constraints = '''
    ...     x0 = x4**2
    ...     x2 - x4 - x3 = 0.'''
    >>> _classify_variables(constraints, nvars=5)
    {'dependent': ['x0','x2'], 'independent': ['x3','x4'], 'unconstrained': ['x1']}
    >>> constraints = '''
    ...     x0 + x2 = 0.
    ...     x0 + 2*x2 = 0.'''
    >>> _classify_variables(constraints, nvars=5)
    {'dependent': ['x0','x2'], 'independent': [], 'unconstrained': ['x1','x3','x4']}

    This is a bug:
    >>> constraints = '''
#.........这里部分代码省略.........

def _prepare_sympy(constraints, variables='x', nvars=None):
    """Parse an equation string and prepare input for sympy. Returns a tuple
of sympy-specific input: (code for variable declaration, left side of equation
string, right side of equation string, list of variables, and the number of
sympy equations).

Inputs:
    constraints -- a string of symbolic constraints, with one constraint
        equation per line. Constraints must be equality constraints only.
        Standard python syntax should be followed (with the math and numpy
        modules already imported).

    For example:
        >>> constraints = '''
        ...     x0 = x4**2
        ...     x4 - x3 = 0.
        ...     x4 - x0 = x2'''
        >>> code, lhs, rhs, vars, neqn = _prepare_sympy(constraints, nvars=5)
        >>> print code
        x0=Symbol('x0')
        x1=Symbol('x1')
        x2=Symbol('x2')
        x3=Symbol('x3')
        x4=Symbol('x4')
        rand = Symbol('rand')
        >>> print lhs, rhs
        ['x0 ', 'x4 - x3 ', 'x4 - x0 '] [' x4**2', ' 0.', ' x2']
        print "%s in %s eqns" % (vars, neqn)        
        x0,x1,x2,x3,x4, in 3 eqns

Additional Inputs:
    nvars -- number of variables. Includes variables not explicitly
        given by the constraint equations (e.g. 'x1' in the example above).
    variables -- desired variable name. Default is 'x'. A list of variable
        name strings is also accepted for when desired variable names
        don't have the same base, and can include variables that are not
        found in the constraints equation string.
"""
    if ">" in constraints or "<" in constraints:
        raise NotImplementedError, "cannot simplify inequalities" 

    from mystic.symbolic import replace_variables, get_variables
    #XXX: if constraints contain x0,x1,x3 for 'x', should x2 be in code,xlist?
    if list_or_tuple_or_ndarray(variables):
        if nvars is not None: variables = variables[:nvars]
        constraints = replace_variables(constraints, variables, markers='_')
        varname = '_'
        ndim = len(variables)
    else:
        varname = variables # varname used below instead of variables
        myvar = get_variables(constraints, variables)
        if myvar: ndim = max([int(v.strip(varname)) for v in myvar]) + 1
        else: ndim = 0
    if nvars is not None: ndim = nvars

    # split constraints_str into lists of left hand sides and right hand sides
    eacheqn = constraints.splitlines()
    neqns = 0
    left = []
    right = []
    for eq in eacheqn: #XXX: Le/Ge instead of Eq; Max/Min... (NotImplemented ?)
        splitlist = eq.replace('==','=').split('=') #FIXME: no inequalities
        if len(splitlist) == 2:   #FIXME: first convert >/< to min/max ?

            # If equation is blank on one side, raise error.
            if len(splitlist[0].strip()) == 0 or len(splitlist[1].strip()) == 0:
                print eq, "is not an equation!" # Raise exception?
            else:
                left.append(splitlist[0])
                right.append(splitlist[1])
                neqns += 1

        # If equation doesn't have one equal sign, raise error.
        if len(splitlist) != 2 and len(splitlist) != 1:
            print eq, "is not an equation!" # Raise exception?

    # First create list of x variables
    xlist = ""
    for i in range(ndim):
        xn = varname + str(i)
        xlist += xn + ","

    # Start constructing the code string
    code = ""
    for i in range(ndim):
        xn = varname + str(i)
        code += xn + '=' + "Symbol('" + xn + "')\n"

    code += "rand = Symbol('rand')\n"
    return code, left, right, xlist, neqns

def replace_variables(constraints, variables=None, markers='$'):
    """Replace variables in constraints string with a marker.
Returns a modified constraints string.

Inputs:
    constraints -- a string of symbolic constraints, with one constraint
        equation per line. Constraints can be equality and/or inequality
        constraints. Standard python syntax should be followed (with the
        math and numpy modules already imported).
    variables -- list of variable name strings. The variable names will
        be replaced in the order that they are provided, where if the
        default marker "$i" is used, the first variable will be replaced
        with "$0", the second with "$1", and so on.

    For example:
        >>> variables = ['spam', 'eggs']
        >>> constraints = '''spam + eggs - 42'''
        >>> print replace_variables(constraints, variables, 'x')
        'x0 + x1 - 42'

Additional Inputs:
    markers -- desired variable name. Default is '$'. A list of variable
        name strings is also accepted for when desired variable names
        don't have the same base.

    For example:
        >>> variables = ['x1','x2','x3']
        >>> constraints = "min(x1*x2) - sin(x3)"
        >>> print replace_variables(constraints, variables, ['x','y','z'])
        'min(x*y) - sin(z)'
"""
    if variables is None: variables = []
    elif isinstance(variables, str): variables = list((variables,))

    # substitite one list of strings for another
    if list_or_tuple_or_ndarray(markers):
        equations = replace_variables(constraints,variables,'_')
        vars = get_variables(equations,'_')
        indices = [int(v.strip('_')) for v in vars]
        for i in range(len(vars)):
            equations = equations.replace(vars[i],markers[indices[i]])
        return equations

    # Sort by decreasing length of variable name, so that if one variable name 
    # is a substring of another, that won't be a problem. 
    variablescopy = variables[:]
    def comparator(x, y):
        return len(y) - len(x)
    variablescopy.sort(comparator)

    # Figure out which index goes with which variable.
    indices = []
    for item in variablescopy:
        indices.append(variables.index(item))

    # Default is markers='$', as '$' is not a special symbol in Python,
    # and it is unlikely a user will choose it for a variable name.
    if markers in variables:
        marker = '_$$$$$$$$$$' # even less likely...
    else:
        marker = markers

    '''Bug demonstrated here:
    >>> equation = """x3 = max(y,x) + x"""
    >>> vars = ['x','y','z','x3']
    >>> print replace_variables(equation,vars)
    $4 = ma$1($2,$1) + $1
    ''' #FIXME: don't parse if __name__ in __builtins__, globals, or locals?
    for i in indices: #FIXME: or better, use 're' pattern matching
        constraints = constraints.replace(variables[i], marker + str(i))
    return constraints.replace(marker, markers)

def generate_constraint(conditions, ctype=None, **kwds):
    """Converts a constraint solver to a mystic.constraints function.

Inputs:
    conditions -- a constraint solver, or list of constraint solvers
    ctype -- a mystic.constraints type, or a list of mystic.constraints types
        of the same length as the given conditions

NOTES:
    This simple constraint generator doesn't check for conflicts in conditions,
    but simply applies conditions in the given order. This constraint generator
    assumes that a single variable has been isolated on the left-hand side
    of each constraints equation, thus all constraints are of the form
    "x_i = f(x)". This solver picks speed over robustness, and thus relies on
    the user to formulate the constraints so that they do not conflict.

    For example:
        >>> constraints = '''
        ...     x0 = cos(x1) + 2.
        ...     x1 = x2*2.'''
        >>> solv = generate_solvers(constraints)
        >>> constraint = generate_constraint(solv)
        >>> constraint([1.0, 0.0, 1.0])
        [1.5838531634528576, 2.0, 1.0]

    Standard python math conventions are used. For example, if an 'int'
    is used in a constraint equation, one or more variable may be evaluate
    to an 'int' -- this can affect solved values for the variables.

    For example:
        >>> constraints = '''
        ...     x2 = x0/2.
        ...     x0 >= 0.'''
        >>> solv = generate_solvers(constraints, nvars=3)
        >>> print solv[0].__doc__
        'x[2] = x[0]/2.'
        >>> print solv[1].__doc__
        'x[0] = max(0., x[0])'
        >>> constraint = generate_constraint(solv)
        >>> constraint([1,2,3])
        [1, 2, 0.5]
        >>> constraint([-1,2,-3])
        [0.0, 2, 0.0]
"""
    # allow for single condition, list of conditions, or nested list
    if not list_or_tuple_or_ndarray(conditions):
        conditions = list((conditions,))
    else: pass #XXX: should be fine...
    conditions = list(flatten(conditions))

    # allow for single ctype, list of ctypes, or nested list
    if ctype is None:
        from mystic.coupler import inner #XXX: outer ?
        ctype = list((inner,))*len(conditions)
    elif not list_or_tuple_or_ndarray(ctype):
        ctype = list((ctype,))*len(conditions)
    else: pass #XXX: is already a list, should be the same len as conditions
    ctype = list(flatten(ctype))

    # iterate through solvers, building a compound constraints solver
    cf = lambda x:x
    cfdoc = ""
    for wrapper, condition in zip(ctype, conditions):
        cfdoc += "%s: %s\n" % (wrapper.__name__, condition.__doc__)
        apply = wrapper(condition, **kwds)
        cf = apply(cf)
    cf.__doc__ = cfdoc.rstrip('\n')
    cf.__name__ = 'constraint'
    return cf

def constraints_parser(constraints, variables='x', nvars=None):
    """parse symbolic constraints into a tuple of constraints solver equations.
The left-hand side of each constraint must be simplified to support assignment.

Inputs:
    constraints -- a string of symbolic constraints, with one constraint
        equation per line. Constraints can be equality and/or inequality
        constraints. Standard python syntax should be followed (with the
        math and numpy modules already imported).

    For example:
        >>> constraints = '''
        ...     x2 = x0/2.
        ...     x0 >= 0.'''
        >>> constraints_parser(constraints, nvars=3)
        ('x[2] = x[0]/2.', 'x[0] = max(0., x[0])')

Additional Inputs:
    nvars -- number of variables. Includes variables not explicitly
        given by the constraint equations (e.g. 'x2' in the example above).
    variables -- desired variable name. Default is 'x'. A list of variable
        name strings is also accepted for when desired variable names
        don't have the same base, and can include variables that are not
        found in the constraints equation string.
    """
   #from mystic.tools import src
   #ndim = len(get_variables(src(func), variables))
    if list_or_tuple_or_ndarray(variables):
        if nvars != None: variables = variables[:nvars]
        constraints = replace_variables(constraints, variables)
        varname = '$'
        ndim = len(variables)
    else:
        varname = variables # varname used below instead of variables
        myvar = get_variables(constraints, variables)
        if myvar: ndim = max([int(v.strip(varname)) for v in myvar]) + 1
        else: ndim = 0
    if nvars != None: ndim = nvars

    # Parse the constraints string
    lines = constraints.splitlines()
    parsed = [] #XXX: in penalty_parser is eqconstraints, ineqconstraints
    for line in lines:
        if line.strip():
            fixed = line
            # Iterate in reverse in case ndim > 9.
            indices = list(range(ndim))
            indices.reverse()
            for i in indices:
                fixed = fixed.replace(varname + str(i), 'x[' + str(i) + ']') 
            constraint = fixed.strip()

            # Replace 'ptp', 'average', and 'var' (uses mystic, not numpy)
            if constraint.find('ptp(') != -1:
                constraint = constraint.replace('ptp(', 'spread(')
            if constraint.find('average(') != -1:
                constraint = constraint.replace('average(', 'mean(')
            if constraint.find('var(') != -1:
                constraint = constraint.replace('var(', 'variance(')
            if constraint.find('prod(') != -1:
                constraint = constraint.replace('prod(', 'product(')

            #XXX: below this line the code is different than penalty_parser
            # convert "<" to min(LHS, RHS) and ">" to max(LHS,RHS)
            split = constraint.split('>')
            expression = '%(lhs)s = max(%(rhs)s, %(lhs)s)'
            if len(split) == 1: # didn't contain '>'
                split = constraint.split('<')
                expression = '%(lhs)s = min(%(rhs)s, %(lhs)s)'
            if len(split) == 1: # didn't contain '>' or '<'
                split = constraint.split('=')
                expression = '%(lhs)s = %(rhs)s'
            if len(split) == 1: # didn't contain '>', '<', or '='
                print "Invalid constraint: ", constraint
            eqn = {'lhs':split[0].rstrip('=').strip(), \
                   'rhs':split[-1].lstrip('=').strip()}
            expression = expression % eqn

            # allow mystic.math.measures impose_* on LHS
            lhs,rhs = expression.split('=')
            if lhs.find('spread(') != -1:
              lhs = lhs.split('spread')[-1]
              rhs = ' impose_spread( (' + rhs.lstrip() + '),' + lhs + ')'
            if lhs.find('mean(') != -1:
              lhs = lhs.split('mean')[-1]
              rhs = ' impose_mean( (' + rhs.lstrip() + '),' + lhs + ')'
            if lhs.find('variance(') != -1:
              lhs = lhs.split('variance')[-1]
              rhs = ' impose_variance( (' + rhs.lstrip() + '),' + lhs + ')'
            if lhs.find('sum(') != -1:
              lhs = lhs.split('sum')[-1]
              rhs = ' impose_sum( (' + rhs.lstrip() + '),' + lhs + ')'
            if lhs.find('product(') != -1:
              lhs = lhs.split('product')[-1]
              rhs = ' impose_product( (' + rhs.lstrip() + '),' + lhs + ')'
            expression = "=".join([lhs,rhs])

            parsed.append(expression)

    return tuple(parsed)

def penalty_parser(constraints, variables='x', nvars=None):
    """parse symbolic constraints into penalty constraints.
Returns a tuple of inequality constraints and a tuple of equality constraints.

Inputs:
    constraints -- a string of symbolic constraints, with one constraint
        equation per line. Constraints can be equality and/or inequality
        constraints. Standard python syntax should be followed (with the
        math and numpy modules already imported).

    For example:
        >>> constraints = '''
        ...     x2 = x0/2.
        ...     x0 >= 0.'''
        >>> penalty_parser(constraints, nvars=3)
        (('-(x[0] - (0.))',), ('x[2] - (x[0]/2.)',))

Additional Inputs:
    nvars -- number of variables. Includes variables not explicitly
        given by the constraint equations (e.g. 'x2' in the example above).
    variables -- desired variable name. Default is 'x'. A list of variable
        name strings is also accepted for when desired variable names
        don't have the same base, and can include variables that are not
        found in the constraints equation string.
    """
   #from mystic.tools import src
   #ndim = len(get_variables(src(func), variables))
    if list_or_tuple_or_ndarray(variables):
        if nvars != None: variables = variables[:nvars]
        constraints = replace_variables(constraints, variables)
        varname = '$'
        ndim = len(variables)
    else:
        varname = variables # varname used below instead of variables
        myvar = get_variables(constraints, variables)
        if myvar: ndim = max([int(v.strip(varname)) for v in myvar]) + 1
        else: ndim = 0
    if nvars != None: ndim = nvars

    # Parse the constraints string
    lines = constraints.splitlines()
    eqconstraints = []
    ineqconstraints = []
    for line in lines:
        if line.strip():
            fixed = line
            # Iterate in reverse in case ndim > 9.
            indices = list(range(ndim))
            indices.reverse()
            for i in indices:
                fixed = fixed.replace(varname + str(i), 'x[' + str(i) + ']') 
            constraint = fixed.strip()

            # Replace 'spread', 'mean', and 'variance' (uses numpy, not mystic)
            if constraint.find('spread(') != -1:
                constraint = constraint.replace('spread(', 'ptp(')
            if constraint.find('mean(') != -1:
                constraint = constraint.replace('mean(', 'average(')
            if constraint.find('variance(') != -1:
                constraint = constraint.replace('variance(', 'var(')

            # Sorting into equality and inequality constraints, and making all
            # inequality constraints in the form expression <= 0. and all 
            # equality constraints of the form expression = 0.
            split = constraint.split('>')
            direction = '>'
            if len(split) == 1:
                split = constraint.split('<')
                direction = '<'
            if len(split) == 1:
                split = constraint.split('=')
                direction = '='
            if len(split) == 1:
                print "Invalid constraint: ", constraint
            eqn = {'lhs':split[0].rstrip('=').strip(), \
                   'rhs':split[-1].lstrip('=').strip()}
            expression = '%(lhs)s - (%(rhs)s)' % eqn
            if direction == '=':
                eqconstraints.append(expression)
            elif direction == '<':
                ineqconstraints.append(expression)
            else:
                ineqconstraints.append('-(' + expression + ')')

    return tuple(ineqconstraints), tuple(eqconstraints)