Changeset 622 for branches


Ignore:
Timestamp:
12/26/12 17:27:23 (3 years ago)
Author:
mmckerns
Message:

merge simplify_symbolic, _build_linear, _build_nonlinear, and _build_constraint;
now is single function 'solve', which returns a solve system of equations

Location:
branches/decorate
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • branches/decorate/_symbolic.py

    r620 r622  
     1# The following code attempts to construct something like: 
     2# >>> from sympy import Eq, Symbol 
     3# >>> from sympy import solve as symsol 
     4# >>> x1 = Symbol('x1') 
     5# >>> x2 = Symbol('x2') 
     6# >>> x3 = Symbol('x3') 
     7# >>> eq1 = Eq(x2, x1 - 2.) 
     8# >>> eq2 = Eq(x2, x3*2.) 
     9# >>> soln = symsol([eq2, eq1], [x1, x2, x3]) 
     10 
    111from symbolic import * 
    212from mystic.tools import permutations 
    … …  
    4151        raise NotImplementedError, "cannot classify inequalities"  
    4252 
    43     #XXX: use simplify_symbolic? or first if not in form xi = ... ? 
     53    #XXX: use solve? or first if not in form xi = ... ? 
    4454    if list_or_tuple_or_ndarray(variables): 
    4555        if nvars: variables = variables[:nvars] 
    … …  
    122132def _prepare_sympy(constraints, variables='x', nvars=None): 
    123133    """Parse an equation string and prepare input for sympy. Returns a tuple 
    124 of sympy-specific input (code for variable declaration, left side of equation 
     134of sympy-specific input: (code for variable declaration, left side of equation 
    125135string, right side of equation string, list of variables, and the number of 
    126136sympy equations). 
    … …  
    211221 
    212222 
    213 def simplify_symbolic(constraint, variables='x', target=None, **kwds): 
    214     """Solve a single equation for each variable found within the constraint. 
    215 Returns permutations of the constraint equation, solved for each variable. 
    216 If the equation fails to simplify, the original equation is returned. 
     223def _solve_single(constraint, variables='x', target=None, **kwds): 
     224    """Solve a symbolic constraints equation for a single variable. 
    217225 
    218226Inputs: 
    … …  
    224232    For example: 
    225233        >>> equation = "x2 - 3. = x1*x3" 
    226         >>> print simplify_symbolic(equation) 
    227         x3 = -(3.0 - x2)/x1 
    228         x2 = 3.0 + x1*x3 
     234        >>> print _solve_single(equation) 
    229235        x1 = -(3.0 - x2)/x3 
    230236 
    … …  
    234240        don't have the same base, and can include variables that are not 
    235241        found in the constraints equation string. 
    236     target -- specify the variable to isolate on the left side upon return. 
     242    target -- list providing the order for which the variables will be solved. 
     243        By default, increasing order is used. 
    237244 
    238245    For example: 
    239246        >>> equation = "x2 - 3. = x1*x3" 
    240         >>> print simplify_symbolic(equation, target='x2') 
     247        >>> print _solve_single(equation, target='x2') 
    241248        x2 = 3.0 + x1*x3 
    242249 
    … …  
    244251    locals -- a dictionary of additional variables used in the symbolic 
    245252        constraints equations, and their desired values. 
    246 """ #XXX: an expanded version of this code is found in build_constraints XXX# 
     253""" #XXX: an very similar version of this code is found in _solve_linear XXX# 
    247254    # for now, we abort on multi-line equations or inequalities 
    248255    if len(constraint.replace('==','=').split('=')) != 2: 
    … …  
    251258        raise NotImplementedError, "cannot simplify inequalities"  
    252259 
     260    nvars = None 
     261    permute = False # if True, return all permutations 
     262    warn = True  # if True, don't supress warnings 
     263    verbose = False  # if True, print debug info 
     264    #-----------------------undocumented------------------------------- 
     265    if kwds.has_key('permute'): permute = kwds['permute'] 
     266    if kwds.has_key('warn'): warn = kwds['warn'] 
     267    if kwds.has_key('verbose'): verbose = kwds['verbose'] 
     268    #------------------------------------------------------------------ 
     269    if target in [None, False]: 
     270        target = [] 
     271    elif isinstance(target, str): 
     272        target = target.split(',') 
     273    else: 
     274        target = list(target) # not the best for ndarray, but should work 
     275 
    253276    if list_or_tuple_or_ndarray(variables): 
    254277        if nvars: variables = variables[:nvars] 
    … …  
    256279        varname = '_' 
    257280        ndim = len(variables) 
    258         if variables.count(target): 
    259             target = replace_variables(target, variables, markers='_') 
     281        for i in range(len(target)): 
     282            if variables.count(target[i]): 
     283                target[i] = replace_variables(target[i],variables,markers='_') 
    260284    else: 
    261285        constraints = constraint # constraints used below 
    … …  
    264288        if myvar: ndim = max([int(v.strip(varname)) for v in myvar]) 
    265289        else: ndim = 0 
    266  
    267     warn = True  # if True, don't supress warnings 
    268     verbose = False  # if True, print debug info 
    269     if kwds.has_key('warn'): warn = kwds['warn'] 
    270     if kwds.has_key('verbose'): verbose = kwds['verbose'] 
     290    if nvars: ndim = nvars 
     291 
     292    # create function to replace "_" with original variables 
     293    def restore(variables, mystring): 
     294        if list_or_tuple_or_ndarray(variables): 
     295            vars = get_variables(mystring,'_') 
     296            indices = [int(v.strip('_'))-1 for v in vars] 
     297            for i in range(len(vars)): 
     298                mystring = mystring.replace(vars[i],variables[indices[i]]) 
     299        return mystring 
    271300 
    272301    # default is _locals with sympy imported 
    … …  
    302331    eqlist = eqlist.rstrip(',') 
    303332 
    304     if not target: #XXX: the goal is solving *only one* equation 
    305        #code += 'soln = symsol([' + eqlist + '], [' + xlist + '])\n' 
    306         code += '_xlist = %s\n' % xlist 
     333    # get full list of variables in 'targeted' order 
     334    xperms = xlist.split(',')[:-1] 
     335    targeted = target[:] 
     336    [targeted.remove(i) for i in targeted if i not in xperms] 
     337    [targeted.append(i) for i in xperms if i not in targeted] 
     338    _target = [] 
     339    [_target.append(i) for i in targeted if i not in _target] 
     340    targeted = _target 
     341    targeted = tuple(targeted) 
     342 
     343    ######################################################################## 
     344    # solve each xi: symsol(single_equation, [x1,x2,...,xi,...,xn]) 
     345    # returns: {x1: f(xn,...), x2: f(xn,...), ..., xn: f(...,x1)} 
     346    if permute or not target: #XXX: the goal is solving *only one* equation 
     347        code += '_xlist = %s\n' % ','.join(targeted) 
     348        code += '_xlist = [symsol(['+eqlist+'], [i]) for i in _xlist]\n' 
    307349        code += 'soln = {}\n' 
    308         code += '[soln.update({i:symsol('+eqlist+', i)}) for i in _xlist]\n' 
    309     else: 
    310         code += 'soln = symsol(' + eqlist + ', [' + target + '])\n' 
     350        code += '[soln.update(i) for i in _xlist if i]\n' 
     351    else: 
     352        code += 'soln = symsol([' + eqlist + '], [' + target[0] + '])\n' 
     353       #code += 'soln = symsol([' + eqlist + '], [' + targeted[0] + '])\n' 
     354    ######################################################################## 
     355 
    311356    if verbose: print code 
    312357    code = compile(code, '<string>', 'exec') 
    … …  
    314359        exec code in globals(), _locals 
    315360        soln = _locals['soln'] 
     361        if not soln: 
     362            if warn: print "Warning: target variable is not valid" 
     363            soln = {} 
    316364    except NotImplementedError: # catch 'multivariate' error for older sympy 
    317         if warn: print "Warning: sympy could not simplify equation." 
    318         return constraint 
     365        if warn: print "Warning: could not simplify equation." 
     366        return constraint      #FIXME: resolve diff with _solve_linear 
    319367    except NameError, error: # catch when variable is not defined 
    320368        if warn: print "Warning:", error 
    321         return None 
     369        soln = {} 
    322370    if verbose: print soln 
    323371 
    324     #XXX Not the best way to handle multiple solutions? 
    325     if not target: 
    326         solvedstring = "" 
    327         for key, value in soln.iteritems(): 
    328             if value: 
    329                 for v in value: 
    330                     solvedstring += str(key) + ' = ' + str(v) + '\n' 
    331         if solvedstring: solvedstring = solvedstring[:-1] 
    332     #XXX: consistent to return string (above) and None (below)? 
    333     else: 
    334         if not soln: 
    335             if warn: print "Warning: target variable is not valid" 
    336             return None 
    337         solvedstring = target + ' = ' + str(soln[0]) 
    338     # replace '_' with the original variable names 
    339     if list_or_tuple_or_ndarray(variables): 
    340         vars = get_variables(solvedstring,'_') 
    341         indices = [int(v.strip('_'))-1 for v in vars] 
    342         for i in range(len(vars)): 
    343             solvedstring = solvedstring.replace(vars[i],variables[indices[i]]) 
    344     return solvedstring 
    345  
    346  
    347 def _build_linear(constraints, variables='x', nvars=None, targets=None, **kwds): 
    348     """Build a constraints function given a string of linear constraints. 
    349 Returns a constraints function.  
     372    #XXX handles multiple solutions? 
     373    soln = dict([(str(key),str(value)) for key, value in soln.iteritems()]) 
     374    soln = [(i,soln[i]) for i in targeted if i in soln] #XXX: order as targeted? 
     375     
     376    solns = []; solved = "" 
     377    for key,value in soln: 
     378        solved = str(key) + ' = ' + str(value) + '\n' 
     379        if solved: solns.append( restore(variables, solved.rstrip()) ) 
     380 
     381    if not permute: 
     382        return None if not solns[:1] else solns[0] 
     383    return tuple(solns) 
     384 
     385 
     386def _solve_linear(constraints, variables='x', target=None, **kwds): 
     387    """Solve a system of symbolic linear constraints equations. 
    350388 
    351389Inputs: 
    … …  
    356394 
    357395    For example: 
    358         >>> constraints = '''x1 = x3 + 2. 
     396        >>> constraints = ''' 
     397        ...     x1 - x3 = 2. 
    359398        ...     x3 = x4*2.''' 
    360         >>> f = _build_linear(constraints, nvars=4) 
    361         >>> f([1.0, 0.0, 0.0, 1.0]) 
    362         [4.0, 0.0, 2.0, 1.0] 
     399        >>> print _solve_linear(constraints) 
     400        x3 = 2.0*x4 
     401        x1 = 2.0 + 2.0*x4 
    363402 
    364403Additional Inputs: 
    365     nvars -- number of variables. Includes variables not explicitly 
    366         given by the constraint equations (e.g. 'x2' in the example above). 
    367404    variables -- desired variable name. Default is 'x'. A list of variable 
    368405        name strings is also accepted for when desired variable names 
    369406        don't have the same base, and can include variables that are not 
    370407        found in the constraints equation string. 
    371     targets -- list providing the order for which the variables will be solved. 
     408    target -- list providing the order for which the variables will be solved. 
    372409        If there are "N" constraint equations, the first "N" variables given 
    373410        will be selected as the dependent variables. By default, increasing 
    … …  
    375412 
    376413    For example: 
    377         >>> constraints = '''x1 = x3 + 2. 
     414        >>> constraints = ''' 
     415        ...     x1 - x3 = 2. 
    378416        ...     x3 = x4*2.''' 
    379         >>> f = _build_linear(constraints, nvars=4, targets=['x4','x3']) 
    380         >>> f([1.0, 0.0, 0.0, 1.0]) 
    381         [1.0, 0.0, -1.0, -0.5] 
     417        >>> print _solve_linear(constraints, target=['x4','x3']) 
     418        x4 = -1.0 + 0.5*x1 
     419        x3 = -2.0 + x1 
    382420 
    383421Further Inputs: 
    … …  
    385423        constraints equations, and their desired values. 
    386424""" 
    387     """ 
    388     guess -- list of parameter values proposed to solve the constraints. 
    389     lower_bounds -- list of lower bounds on solution values. 
    390     upper_bounds -- list of upper bounds on solution values. 
    391  
    392     NOTE: For optimization problems with initial values and either lower or 
    393         upper bounds, the constraints function must be formulated in a manner 
    394         such that it does not immediately violate the given bounds. 
    395  
    396     Returns a constraints function with the constraints simplified nicely. 
    397     For example, the constraints function is: 
    398  
    399     def f(params): 
    400         params[1] = 2.0*params[2] 
    401         params[0] = 2.0 + 2.0*params[2] 
    402         return params 
    403  
    404     This is good for systems of equations, which must be linear for Sympy to 
    405     simplify them. However, using the nonlinear all-purpose equation  
    406     inverter/simplifier will probably also work! It just won't get rid of 
    407     redundancies, but that should be ok. 
    408 """ 
    409     permute = False # if True, force to use 'permutations' code 
     425    nvars = None 
     426    permute = False # if True, return all permutations 
    410427    warn = True  # if True, don't supress warnings 
    411     verbose = False # if True, print tempcode and _classify_variables 
    412 #   guess = None 
    413 #   upper_bounds = None 
    414 #   lower_bounds = None 
     428    verbose = False  # if True, print debug info 
    415429    #-----------------------undocumented------------------------------- 
    416430    if kwds.has_key('permute'): permute = kwds['permute'] 
    417431    if kwds.has_key('warn'): warn = kwds['warn'] 
    418432    if kwds.has_key('verbose'): verbose = kwds['verbose'] 
    419 #   if kwds.has_key('guess'): guess = kwds['guess'] 
    420 #   if kwds.has_key('upper_bounds'): upper_bounds = kwds['upper_bounds'] 
    421 #   if kwds.has_key('lower_bounds'): lower_bounds = kwds['lower_bounds'] 
    422433    #------------------------------------------------------------------ 
    423     if targets in [None, False]: 
    424         targets = [] 
    425     elif isinstance(targets, str): 
    426         targets = targets.split(',') 
    427     else: 
    428         targets = list(targets) # not the best for ndarray, but should work 
     434    if target in [None, False]: 
     435        target = [] 
     436    elif isinstance(target, str): 
     437        target = target.split(',') 
     438    else: 
     439        target = list(target) # not the best for ndarray, but should work 
    429440 
    430441    if list_or_tuple_or_ndarray(variables): 
    431442        if nvars: variables = variables[:nvars] 
    432         constraints = replace_variables(constraints, variables, '_') 
     443        _constraints = replace_variables(constraints, variables, '_') 
    433444        varname = '_' 
    434445        ndim = len(variables) 
    435     else: 
     446        for i in range(len(target)): 
     447            if variables.count(target[i]): 
     448                target[i] = replace_variables(target[i],variables,markers='_') 
     449    else: 
     450        _constraints = constraints 
    436451        varname = variables # varname used below instead of variables 
    437452        myvar = get_variables(constraints, variables) 
    … …  
    439454        else: ndim = 0 
    440455    if nvars: ndim = nvars 
    441 #   elif guess: ndim = len(guess) 
    442 #   elif lower_bounds: ndim = len(lower_bounds) 
    443 #   elif upper_bounds: ndim = len(upper_bounds) 
    444  
    445     # The following code attempts to construct something like: 
    446     # >>> from sympy import Eq, Symbol 
    447     # >>> from sympy import solve as symsol 
    448     # >>> x1 = Symbol('x1') 
    449     # >>> x2 = Symbol('x2') 
    450     # >>> x3 = Symbol('x3') 
    451     # >>> eq1 = Eq(x2, x1 - 2.) 
    452     # >>> eq2 = Eq(x2, x3*2.) 
    453     # >>> soln = symsol([eq2, eq1], [x1, x2, x3]) 
     456 
     457    # create function to replace "_" with original variables 
     458    def restore(variables, mystring): 
     459        if list_or_tuple_or_ndarray(variables): 
     460            vars = get_variables(mystring,'_') 
     461            indices = [int(v.strip('_'))-1 for v in vars] 
     462            for i in range(len(vars)): 
     463                mystring = mystring.replace(vars[i],variables[indices[i]]) 
     464        return mystring 
    454465 
    455466    # default is _locals with sympy imported 
    … …  
    465476    except ImportError: # Equation will not be simplified." 
    466477        if warn: print "Warning: sympy not installed." 
    467         solv = generate_solvers(constraints, variables=varname, locals=locals) 
    468         return generate_constraint(solv) 
     478        return constraints 
    469479 
    470480    # default is _locals with numpy and math imported 
    … …  
    478488    _locals.update(locals) #XXX: allow this? 
    479489 
    480     code,left,right,xlist,neqns = _prepare_sympy(constraints, varname, ndim) 
     490    code,left,right,xlist,neqns = _prepare_sympy(_constraints, varname, ndim) 
    481491 
    482492    eqlist = "" 
    … …  
    487497    eqlist = eqlist.rstrip(',') 
    488498 
    489     # Figure out if trying various permutations is necessary 
    490 #   if (guess and lower_bounds) or (guess and upper_bounds): 
    491 #       permute = True 
    492  
     499    # get full list of variables in 'targeted' order 
    493500    xperms = xlist.split(',')[:-1] 
    494     if targets: 
    495         [targets.remove(i) for i in targets if i not in xperms] 
    496         [targets.append(i) for i in xperms if i not in targets] 
    497     targets = tuple(targets) 
     501    targeted = target[:] 
     502    [targeted.remove(i) for i in targeted if i not in xperms] 
     503    [targeted.append(i) for i in xperms if i not in targeted] 
     504    _target = [] 
     505    [_target.append(i) for i in targeted if i not in _target] 
     506    targeted = _target 
     507    targeted = tuple(targeted) 
    498508 
    499509    if permute: 
    500         # If there are strict bounds and initial x value, form a constraints  
    501         # function for each permutation. 
    502         # For sympy, change the order of the x variables passed to symsol() 
     510        # Form constraints equations for each permutation. 
     511        # This will change the order of the x variables passed to symsol() 
    503512        # to get different variables solved for. 
    504         solns = [] 
    505513        xperms = list(permutations(xperms)) #XXX: takes a while if nvars is ~10 
    506         if targets: 
    507             xperms.remove(targets) 
    508             xperms.insert(0, targets) 
    509         for perm in xperms:  
    510             tempcode = "" 
    511             tempcode += code 
    512             xstring = "" 
    513             for item in perm: 
    514                 xstring += item + "," 
    515             xstring = xstring.rstrip(',') 
    516             tempcode += 'soln = symsol([' + eqlist + '], [' + xstring + '])' 
    517             if verbose: print tempcode 
    518  
    519             tempcode = compile(tempcode, '<string>', 'exec') 
    520             try:  
    521                 exec tempcode in globals(), _locals 
    522                 soln = _locals['soln'] 
    523                 if soln == None and warn: 
    524                     print "Warning: sympy could not simplify equation." 
    525             except NotImplementedError: # catch 'multivariate' error 
    526                 if warn: print "Warning: sympy could not simplify equation." 
    527                 soln = None 
    528             except NameError, error: # catch when variable is not defined 
    529                 if warn: print "Warning:", error 
    530                 soln = None 
    531             if verbose: print soln 
    532  
    533             solvedstring = "" 
    534             for key, value in soln.iteritems(): 
    535                 solvedstring += str(key) + ' = ' + str(value) + '\n' 
    536             solns.append(solvedstring) 
    537  
    538         # Create strings of all permutations of the solved equations. 
    539         # Remove duplicates, then take permutations of the lines of equations 
    540         # to create equations in different orders. 
    541         noduplicates = list(set(solns))  
    542         stringperms = [] 
    543         for item in noduplicates: 
    544             spl = item.splitlines() 
    545             for perm in permutations(spl): 
    546                 permstring = "" 
    547                 for line in perm: 
    548                     permstring += line + '\n' 
    549                 stringperms.append(permstring) 
    550  
    551         # Feed each solved set of equations into generate_constraint to get  
    552         # constraints functions. Check if constraints(guess) is in the bounds. 
    553         for string in stringperms: 
    554             solv = generate_solvers(string, variables=varname, locals=locals) 
    555             cf = generate_constraint(solv) 
    556 #           if guess: compatible = isbounded(cf, guess, min=lower_bounds, \ 
    557 #                                                       max=upper_bounds) 
    558 #           else: compatible = True 
    559             compatible = True #(due to commented out 'guess') 
    560             if compatible: 
    561                 #print string # Debugging 
    562                 if verbose: 
    563                     print _classify_variables(string, varname, ndim) 
    564                 return cf # Return the first compatible constraints function 
    565             else: 
    566                 _constraints = cf # Save for later and try other permutations 
    567         # Perhaps raising an Exception is better? But sometimes it's ok... 
    568        #warning='Warning: guess incompatible with constraints and bounds.' 
    569         warning='Warning: an error occurred in building the constraints.' 
    570         if warn: print warning 
    571         if verbose: 
    572             print _classify_variables(string, varname, ndim) 
    573         return cf 
    574  
    575     else: 
    576         # If no bounds and guess, no need for permutations. 
    577         # Just form code and get whatever solution sympy comes up with. 
    578         if targets: 
    579             xlist = ','.join(targets).rstrip(',') 
    580         code += 'soln = symsol([' + eqlist + '], [' + xlist + '])' 
    581         if verbose: print code 
    582         code = compile(code, '<string>', 'exec') 
     514        if target: # put the tuple with the 'targeted' order first 
     515            xperms.remove(targeted) 
     516            xperms.insert(0, targeted) 
     517    else: 
     518        xperms = [tuple(targeted)] 
     519 
     520    solns = [] 
     521    for perm in xperms:  
     522        _code = code 
     523        xlist = ','.join(perm).rstrip(',') #XXX: if not all, use target ? 
     524        # solve dependent xi: symsol([linear_system], [x1,x2,...,xi,...,xn]) 
     525        # returns: {x1: f(xn,...), x2: f(xn,...), ...} 
     526        _code += 'soln = symsol([' + eqlist + '], [' + xlist + '])' 
     527        if verbose: print _code 
     528        _code = compile(_code, '<string>', 'exec') 
    583529        try:  
    584             exec code in globals(), _locals 
     530            exec _code in globals(), _locals 
    585531            soln = _locals['soln'] 
    586             if soln == None and warn: 
    587                 print "Warning: sympy could not simplify equation." 
     532            if not soln: 
     533                if warn: print "Warning: could not simplify equation." 
     534                soln = {} 
    588535        except NotImplementedError: # catch 'multivariate' error 
    589             if warn: print "Warning: sympy could not simplify equation." 
    590             soln = None 
     536            if warn: print "Warning: could not simplify equation." 
     537            soln = {} 
    591538        except NameError, error: # catch when variable is not defined 
    592539            if warn: print "Warning:", error 
    593             soln = None 
     540            soln = {} 
    594541        if verbose: print soln 
    595542 
    596         solvedstring = "" 
     543        solved = "" 
    597544        for key, value in soln.iteritems(): 
    598             solvedstring += str(key) + ' = ' + str(value) + '\n' 
    599  
    600         solv = generate_solvers(solvedstring, variables=varname, locals=locals) 
    601         cf = generate_constraint(solv) 
    602         if verbose: 
    603             print _classify_variables(solvedstring, varname, ndim) 
    604         return cf 
    605  
    606  
    607 def build_constraint(constraints,variables='x',nvars=None,targets=None,**kwds): 
    608     """Build a constraints function given a constraints string.  
    609 Returns a constraints function. 
     545            solved += str(key) + ' = ' + str(value) + '\n' 
     546        if solved: solns.append( restore(variables, solved.rstrip()) ) 
     547 
     548    if not permute: 
     549        return None if not solns[:1] else solns[0] 
     550 
     551    # Remove duplicates 
     552    filter = []; results = [] 
     553    for i in solns: 
     554        _eqs = '\n'.join(sorted(i.split('\n'))) 
     555        if _eqs not in filter: 
     556            filter.append(_eqs) 
     557            results.append(i) 
     558    return tuple(results) 
     559 
     560    # Create strings of all permutations of the solved equations. 
     561    # Remove duplicates, then take permutations of the lines of equations 
     562    # to create equations in different orders. 
     563#   noduplicates = [] 
     564#   [noduplicates.append(i) for i in solns if i not in noduplicates] 
     565#   stringperms = [] 
     566#   for item in noduplicates: 
     567#       spl = item.splitlines() 
     568#       for perm in permutations(spl): 
     569#           permstring = "" 
     570#           for line in perm: 
     571#               permstring += line + '\n' 
     572#           stringperms.append(permstring.rstrip()) 
     573#   return tuple(stringperms) 
     574 
     575 
     576def solve(constraints, variables='x', target=None, **kwds): 
     577    """Solve a system of symbolic constraints equations. 
    610578 
    611579Inputs: 
    … …  
    616584 
    617585    For example: 
    618         >>> constraints = '''x1 = x3 + 2. 
     586        >>> constraints = ''' 
     587        ...     x1 - x3 = 2. 
    619588        ...     x3 = x4*2.''' 
    620         >>> f = build_constraint(constraints, nvars=4) 
    621         >>> f([1.0, 0.0, 0.0, 1.0]) 
    622         [4.0, 0.0, 2.0, 1.0] 
     589        >>> print solve(constraints) 
     590        x3 = 2.0*x4 
     591        x1 = 2.0 + 2.0*x4 
    623592        >>> constraints = ''' 
    624593        ...     spread([x1,x2]) - 1.0 = mean([x1,x2])    
    625594        ...     mean([x1,x2,x3]) = x3''' 
    626         >>> f = build_constraint(constraints) 
    627         >>> f([1.0, 2.0, 3.0]) 
    628         [1.0, 5.0, 3.0] 
     595        >>> print solve(constraints) 
     596        x1 = -0.5 + 0.5*x3 
     597        x2 = 0.5 + 1.5*x3 
    629598 
    630599Additional Inputs: 
    631     nvars -- number of variables. Includes variables not explicitly 
    632         given by the constraint equations (e.g. 'x2' in the example above). 
    633600    variables -- desired variable name. Default is 'x'. A list of variable 
    634601        name strings is also accepted for when desired variable names 
    635602        don't have the same base, and can include variables that are not 
    636603        found in the constraints equation string. 
    637     targets -- list providing the order for which the variables will be solved. 
     604    target -- list providing the order for which the variables will be solved. 
    638605        If there are "N" constraint equations, the first "N" variables given 
    639606        will be selected as the dependent variables. By default, increasing 
    … …  
    641608 
    642609    For example: 
    643         >>> constraints = '''x1 = x3 + 2. 
     610        >>> constraints = ''' 
     611        ...     x1 - x3 = 2. 
    644612        ...     x3 = x4*2.''' 
    645         >>> f = build_constraint(constraints, nvars=4, targets=['x4','x3']) 
    646         >>> f([1.0, 0.0, 0.0, 1.0]) 
    647         [1.0, 0.0, -1.0, -0.5] 
     613        >>> print _solve_linear(constraints, target=['x4','x3']) 
     614        x4 = -1.0 + 0.5*x1 
     615        x3 = -2.0 + x1 
    648616 
    649617Further Inputs: 
    … …  
    651619        constraints equations, and their desired values. 
    652620""" 
    653     """ 
    654     guess -- list of parameter values proposed to solve the constraints. 
    655     lower_bounds -- list of lower bounds on solution values. 
    656     upper_bounds -- list of upper bounds on solution values. 
    657  
    658     NOTE: For optimization problems with initial values and either lower or 
    659         upper bounds, the constraints function must be formulated in a manner 
    660         such that it does not immediately violate the given bounds. 
    661 """ 
     621    #-----------------------undocumented------------------------------- 
     622   #kwds['permute'] = False # if True, return all permutations 
     623    kwds['warn'] = False  # if True, don't supress warnings 
     624    kwds['verbose'] = False  # if True, print debug info 
     625    #------------------------------------------------------------------ 
    662626    try: 
    663         return _build_linear(constraints, variables=variables, \ 
    664                              nvars=nvars, targets=targets, **kwds) 
     627        if len(constraints.replace('==','=').split('=')) <= 2: 
     628            soln = _solve_single(constraints, variables=variables, \ 
     629                                 target=target, **kwds) 
     630        else: 
     631            soln = _solve_linear(constraints, variables=variables, \ 
     632                                 target=target, **kwds) 
     633        if not soln: raise ValueError 
    665634    except: 
    666         return _build_nonlinear(constraints, variables=variables, \ 
    667                                 nvars=nvars, targets=targets, **kwds) 
    668  
    669  
    670 def _build_nonlinear(constraints,variables='x',nvars=None,targets=None,**kwds): 
     635        soln = _solve_nonlinear(constraints, variables=variables, \ 
     636                                target=target, **kwds) 
     637    return soln 
     638 
     639 
     640def _solve_nonlinear(constraints, variables='x', target=None, **kwds): 
    671641    """Build a constraints function given a string of nonlinear constraints. 
    672642Returns a constraints function.  
    … …  
    679649 
    680650    For example: 
    681         >>> constraints = '''x2 = x4*3. + (x1*x3)**x1''' 
    682         >>> f = _build_nonlinear(constraints, nvars=5) 
    683         >>> f([1.0, 1.0, 1.0, 1.0, 1.0]) 
    684         [1.0, 4.0, 1.0, 1.0, 1.0] 
     651        >>> constraints = '''x2 = x4*3. + x1*x3''' 
     652        >>> print _solve_nonlinear(constraints) 
     653        x1 = (x2 - 3.0*x4)/x3 
    685654        >>> constraints = ''' 
    686655        ...     spread([x1,x2]) - 1.0 = mean([x1,x2])    
    687656        ...     mean([x1,x2,x3]) = x3''' 
    688         >>> f = _build_nonlinear(constraints) 
    689         >>> f([1.0, 2.0, 3.0]) 
    690         [1.0, 5.0, 3.0] 
     657        >>> print _solve_nonlinear(constraints) 
     658        x1 = -0.5 + 0.5*x3 
     659        x2 = 0.5 + 1.5*x3 
    691660 
    692661Additional Inputs: 
    693     nvars -- number of variables. Includes variables not explicitly 
    694         given by the constraint equations (e.g. 'x5' in the example above). 
    695662    variables -- desired variable name. Default is 'x'. A list of variable 
    696663        name strings is also accepted for when desired variable names 
    697664        don't have the same base, and can include variables that are not 
    698665        found in the constraints equation string. 
    699     targets -- list providing the order for which the variables will be solved. 
     666    target -- list providing the order for which the variables will be solved. 
    700667        If there are "N" constraint equations, the first "N" variables given 
    701668        will be selected as the dependent variables. By default, increasing 
    … …  
    706673        ...     spread([x1,x2]) - 1.0 = mean([x1,x2])    
    707674        ...     mean([x1,x2,x3]) = x3''' 
    708         >>> f = _build_nonlinear(constraints, targets=['x2']) 
    709         >>> f([1.0, 2.0, 3.0]) 
    710         [1.0, -0.33333333333333204, 3.0] 
     675        >>> print _solve_nonlinear(constraints, target=['x2']) 
     676        x2 = -0.833333333333333 + 0.166666666666667*x3 
     677        x1 = -0.5 + 0.5*x3 
    711678 
    712679Further Inputs: 
    … …  
    714681        constraints equations, and their desired values. 
    715682""" 
    716     """ 
    717     guess -- list of parameter values proposed to solve the constraints. 
    718     lower_bounds -- list of lower bounds on solution values. 
    719     upper_bounds -- list of upper bounds on solution values. 
    720  
    721     NOTE: For optimization problems with initial values and either lower or 
    722         upper bounds, the constraints function must be formulated in a manner 
    723         such that it does not immediately violate the given bounds. 
    724 """ 
    725     permute = False # if True, force to use 'permutations' code 
    726     _all = False # if True, return all permutations 
     683    nvars = None 
     684    permute = False # if True, return all permutations 
    727685    warn = True  # if True, don't supress warnings 
    728686    verbose = False # if True, print details from _classify_variables 
    729 #   guess = None 
    730 #   upper_bounds = None 
    731 #   lower_bounds = None 
    732687    #-----------------------undocumented------------------------------- 
    733     if kwds.has_key('all'): _all = kwds['all'] 
    734688    if kwds.has_key('permute'): permute = kwds['permute'] 
    735689    if kwds.has_key('warn'): warn = kwds['warn'] 
    736690    if kwds.has_key('verbose'): verbose = kwds['verbose'] 
    737 #   if kwds.has_key('guess'): guess = kwds['guess'] 
    738 #   if kwds.has_key('upper_bounds'): upper_bounds = kwds['upper_bounds'] 
    739 #   if kwds.has_key('lower_bounds'): lower_bounds = kwds['lower_bounds'] 
    740691    #------------------------------------------------------------------ 
    741     if targets in [None, False]: 
    742         targets = [] 
    743     elif isinstance(targets, str): 
    744         targets = targets.split(',') 
    745     else: 
    746         targets = list(targets) # not the best for ndarray, but should work 
     692    if target in [None, False]: 
     693        target = [] 
     694    elif isinstance(target, str): 
     695        target = target.split(',') 
     696    else: 
     697        target = list(target) # not the best for ndarray, but should work 
    747698 
    748699    if list_or_tuple_or_ndarray(variables): 
    … …  
    757708        else: ndim = 0 
    758709    if nvars: ndim = nvars 
    759 #   elif guess: ndim = len(guess) 
    760 #   elif lower_bounds: ndim = len(lower_bounds) 
    761 #   elif upper_bounds: ndim = len(upper_bounds) 
     710 
     711    # create function to replace "_" with original variables 
     712    def restore(variables, mystring): 
     713        if list_or_tuple_or_ndarray(variables): 
     714            vars = get_variables(mystring,'_') 
     715            indices = [int(v.strip('_'))-1 for v in vars] 
     716            for i in range(len(vars)): 
     717                mystring = mystring.replace(vars[i],variables[indices[i]]) 
     718        return mystring 
    762719 
    763720    locals = kwds.get('locals', None) 
    764721    if locals is None: locals = {} 
    765  
    766     # Figure out if trying various permutations is necessary 
    767 #   if (guess and lower_bounds) or (guess and upper_bounds): 
    768 #       permute = True 
    769722 
    770723    eqns = constraints.splitlines() 
    … …  
    779732 
    780733    xperms = [varname+str(i+1) for i in range(ndim)] 
    781     if targets: 
    782         [targets.remove(i) for i in targets if i not in xperms] 
    783         [targets.append(i) for i in xperms if i not in targets] 
    784     targets = tuple(targets) 
     734    if target: 
     735        [target.remove(i) for i in target if i not in xperms] 
     736        [target.append(i) for i in xperms if i not in target] 
     737        _target = [] 
     738        [_target.append(i) for i in target if i not in _target] 
     739        target = _target 
     740    target = tuple(target) 
    785741 
    786742    xperms = list(permutations(xperms)) #XXX: takes a while if nvars is ~10 
    787     if targets: # Try the suggested order first. 
    788         xperms.remove(targets) 
    789         xperms.insert(0, targets) 
     743    if target: # Try the suggested order first. 
     744        xperms.remove(target) 
     745        xperms.insert(0, target) 
    790746 
    791747    complete_list = [] 
    … …  
    821777            # Trying to use xi as a pivot. Loop through the equations 
    822778            # looking for one containing xi. 
    823             target = usedvars[i] 
     779            _target = usedvars[i] 
    824780            for eqn in actual_eqns[i:]: 
    825                 invertedstring = simplify_symbolic(eqn, variables=varname, target=target, warn=warn) 
     781                invertedstring = _solve_single(eqn, variables=varname, target=_target, warn=warn) 
    826782                if invertedstring: 
    827783                    warn = False 
    … …  
    835791            for othereqn in othereqns: 
    836792                expression = invertedstring.split("=")[1] 
    837                 fixed = othereqn.replace(target, '(' + expression + ')') 
     793                fixed = othereqn.replace(_target, '(' + expression + ')') 
    838794                k = actual_eqns.index(othereqn) 
    839795                newsystem[k] = fixed 
    … …  
    843799        simplified = [] 
    844800        for eqn in actual_eqns: 
    845             target = usedvars[actual_eqns.index(eqn)] 
    846             simplified.append(simplify_symbolic(eqn, variables=varname, target=target, warn=warn)) 
    847  
    848         solv = generate_solvers('\n'.join(simplified), variables=varname, locals=locals) 
    849         cf = generate_constraint(solv) 
    850  
    851         if _all: 
    852             complete_list.append(cf) 
     801            _target = usedvars[actual_eqns.index(eqn)] 
     802            mysoln = _solve_single(eqn, variables=varname, target=_target, warn=warn) 
     803            if mysoln: simplified.append(mysoln) 
     804        simplified = restore(variables, '\n'.join(simplified).rstrip()) 
     805 
     806        if permute: 
     807            complete_list.append(simplified) 
    853808            continue 
    854809 
    855 #       if permute and guess: 
    856 #           try: # catches trying to square root a negative number, for example. 
    857 #               cf(guess) 
    858 #           except ValueError: 
    859 #               continue 
    860 #           compatible = isbounded(cf, guess, min=lower_bounds, \ 
    861 #                                             max=upper_bounds) 
    862 #           satisfied = issolution(cf, cf(guess)) 
    863 #           if compatible and satisfied: 
    864 #               if verbose: 
    865 #                   print _classify_variables('\n'.join(simplified), varname, ndim) 
    866 #               return cf 
    867 #       else: #not permute or not guess 
    868         if True: #(due to commented out 'guess') 
    869             if verbose: 
    870                 print _classify_variables('\n'.join(simplified), varname, ndim) 
    871             return cf 
    872  
    873 #   warning='Warning: guess incompatible with constraints and bounds.' 
     810        if verbose: 
     811            print _classify_variables(simplified, variables, ndim) 
     812        return simplified 
     813 
    874814    warning='Warning: an error occurred in building the constraints.' 
    875815    if warn: print warning 
    876816    if verbose: 
    877         print _classify_variables('\n'.join(simplified), varname, ndim) 
    878     if _all: 
    879         return complete_list 
    880     return cf 
     817        print _classify_variables(simplified, variables, ndim) 
     818    if permute: #FIXME: target='x4,x2' may order correct, while 'x2,x4' doesn't 
     819        filter = []; results = [] 
     820        for i in complete_list: 
     821            _eqs = '\n'.join(sorted(i.split('\n'))) 
     822            if _eqs and (_eqs not in filter): 
     823                filter.append(_eqs) 
     824                results.append(i) 
     825        return tuple(results) #FIXME: somehow 'rhs = xi' can be in results 
     826    return simplified 
    881827 
    882828 

  • branches/decorate/test_symbolic.py

    r616 r622  
    5353  assert almostEqual(constraint([1,2,3]), [0.0,5.0,10.0], 1e-10) 
    5454 
    55 def test_build_constraint(): 
     55def test_solve_constraint(): 
    5656 
    5757  constraints = """ 
    … …  
    6060 
    6161  from mystic.math.measures import mean, spread 
    62   constraint = build_constraint(constraints) 
     62  _constraints = solve(constraints) 
     63  solv = generate_solvers(_constraints) 
     64  constraint = generate_constraint(solv) 
    6365  x = constraint([1.0, 2.0, 3.0]) 
    6466  assert all(x) == all([1.0, 5.0, 3.0]) 
    … …  
    7173  test_numpy_penalty() 
    7274  test_generate_constraint() 
    73   test_build_constraint() 
     75  test_solve_constraint() 
    7476 
Note: See TracChangeset for help on using the changeset viewer.