Sympy Solver - Limitations & Suggestions:

[Brett Smith]

Hello All,

New to the forum, but I have used Sympy on a few projects in the past.

I am currently enrolled in a machine design course, and often times we have to solve large systems of equations.

I do not really like how any of the solver we use in class work, so I thought I'd right my own.
If it works well enough Id like to add a GUI front end and incorporate commonly used engineering functions.

I have had success compiling systems of linear equations up to about 40 equations & 40 unknowns.

However adding any non-linear equations seems to break its back...

I was wondering if anybody here had a better way of going about this project.

Right now it just reads a text file full of equations, parses the strings into sympy Functions, and then runs sympy.solve, and prints the solution to a new file.

So currently operations looks like:

import solver

eqn = [] #holds all functions in file

symb = [] #holds unique symbols

sln = []  #holds all solution sets

og  = [] #holds original, unparsed functions for printing to output file

in_file = "problem_file.txt"

out_file = "solution_file.txt"

eqn , og = solver.parse_input(in_file)

symb = solver.sort_symbols(eqn)

sln = solver.solver(eqn)

print_output(out_file,sln,og)

Ill post the code below for easy reading, but I'll also attach a few demo files.

Any ideas on how to make this more powerful (i.e. more equations, or more non-linear equations)???

Running this in my sublime editor a couple times also slows my whole system WAY down from time to time (Ubuntu 14.04)

Thanks Gang!!

from sympy import *
from sympy.parsing.sympy_parser import parse_expr
from IPython.display import display_pretty
import string
import time
import os

init_printing(use_latex = True)

#Vars
sig_figs = 6


def sort_symbols(_functions): #Returns unique Sympy sybols from array of functions
_symbols  = []
temp = []
symbol_str = []
symbols_sorted = []

for i in _functions:
_symbols.append(i.atoms(Symbol))

for i in _symbols: #Puts all occurences of all symbols in 1 list
temp += (set.union(i))

temp = list(set(temp))         #Eliminate duplicates

for i in temp: #Convert to string so symbols may be sorted
   symbol_str.append(str(i))

symbol_str = sorted(symbol_str)

for i in symbol_str:   #Convert string back to sympy symbols
symbols_sorted.append(parse_expr(i))

return symbols_sorted

def parse_input(_file): #Converts file into array of Sympy Functions
_functions = []  
original_functions = []
lines = []
with open(_file,'r+') as f:
lines = f.readlines()
for line in lines:
line = line[:-1]
comment = line[0] == '#'
if not comment:
if "#" in line: #remove comments
split = string.find(line,"#")
line = line[0:split]

original_functions.append(line) #Store original syntax

line = string.replace(line," ","") #remove white space
line = string.replace(line," ","")

if "^" in line: #Convert conventional equation syntax in sympy syntax
line = string.replace(line,"^","**")

if "=" in line:
split = string.find(line,'=')
line = "Eq(" + line[0:split] + "," + line[split + 1:len(line)] + ")"
temp = parse_expr(line)
_functions.append(temp)
return _functions, original_functions


def solver(_functions):
_solution = solve(_functions[0:len(_functions)],manual=True)     #Pass functions to solver #,symbols[0:len(symbols)] #passes symbols
print _solution
return _solution

def print_output(_solutions):

print "The submitted functions:"
for i in original_functions:
pretty_print(i)

print "Have the following unique variables"
for i in sort_symbols():
pretty_print(i)

print "And have solutions set(s):"
for i in _solutions:
for u in i:
print str(u) + " = " + str(Float(i[u],sig_figs))

def print_to_file(f,_solutions,_functions):
##Print File Header
_file = open(f,'w+')
header = "Brett's Sympy Solver!\nSolution Created : " + time.strftime("%d/%m/%Y") + "\n"
_file.write(header)
header = "Equation File : " + os.getcwd() + "/" + in_file + "\n\n" 
_file.write(header)

##Some nice underscoring
max_length = 0

for i in _solutions:
for u in i:
length = len(str(Float(i[u],sig_figs))+" = "+ str(u)) #find longest equation
if length > max_length:
max_length = length

for i in _functions:
length = len(str(i))
if length > max_length:
max_length = length

underscore = "_"
for i in range (0,(max_length-20)/2):
underscore += "_"

_file.write(underscore +"SOLUTIONS"+ underscore +"\n")
for i in _solutions:
for u in i:
_file.write(str(u) + " = " + str(Float(i[u],sig_figs))+"\n")

_file.write("\n" + underscore +"EQUATIONS"+ underscore +"\n")
#_file.write("Compiler Assumes All Equations = 0 Unless specified\n")
print len(_functions)
for i in _functions:
_file.write(str(i)+"\n")

_file.write("\n\n\n\n")

[Ondřej Čertík]

Hi Brett,

On Thu, Jun 11, 2015 at 12:26 PM, Brett Smith <waym...@gmail.com> wrote:
> Hello All,
> New to the forum, but I have used Sympy on a few projects in the past.
> I am currently enrolled in a machine design course, and often times we have
> to solve large systems of equations.
> I do not really like how any of the solver we use in class work, so I
> thought I'd right my own.
> If it works well enough Id like to add a GUI front end and incorporate
> commonly used engineering functions.
> I have had success compiling systems of linear equations up to about 40
> equations & 40 unknowns.
> However adding any non-linear equations seems to break its back...

Can you post the file P5-9-FOS.txt, so that we can run solver.py?

Is anything else needed to run it?

What kind of equations are those? It seems these are numerical
equations, but they are non linear, is that right? What solvers do you
use in class?

Ondrej

> --
> You received this message because you are subscribed to the Google Groups
> "sympy" group.
> To unsubscribe from this group and stop receiving emails from it, send an
> email to sympy+un...@googlegroups.com.
> To post to this group, send email to sy...@googlegroups.com.
> Visit this group at http://groups.google.com/group/sympy.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/sympy/89416b65-b1b1-40e4-b608-1d4a93be9a1f%40googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

[Brett Smith]

Thank you for getting back to me so quickly!

I am actually excited about this program, but it seems like I will need to learn a thing or two about numerical analysis.

Our class is supposed to use TK Solver, and it seems to work OK. However, I really didn't understand how it worked, or why I was getting the bugs I had. 

In addition, TK is only a numerical solver. As far as I can tell it will NOT present an answer as a function of other variables.

I feel like in the long the run, developing and building my own suite (and hopefully sharing) of functions would save me a lot of time in the future, especially if I plan on doing engineering design work.

I attached a revised version of the solver.py file, and the input file.

Everything works perfectly. I have 32 linear functions, 4 non-linear functions and 36 unknown variables. It takes my program about 1.1s to run and solve.

I would like to add the following equations to the file "P5-9-Reactions.txt":

Na = Sy/Sa

Nb = Sy/Sb

Nc = Sy/Sc

Nd = Sy/Sd

Sy = 400*10^6

Sa^2 = 3*T14^2

Sb^2 = 3*T12^2

Sc^2 = 3*T43^2

Sd^2 = 3*T32^2

T14 = F14/A

T12 = F12/A

T43 = F43/A

T32 = F32/A

A = pi*d^2/4

d = .008

It seems to me that system is still determinant, with 51 unknowns and 51 equations.

However, as far as I can tell, the system hangs/runs indefinitely. 

Because sympy is capable of delivering symbolic answers, I don't mind a long execution time, as long as it is actually converging onto a solution, not just running off to infinity...

That being said, I do not understand the underpinnings of a numerical solver, much less a symbolic run. In a class I am taking online we have gone over linear gradient descent and thats about all I know  about the subject - besides just doing symbolic analysis by hand.

Any ideas on how I can get this to solve? I have grown attached to sympy, numpy, and matplotlib. I think they are very powerful, and I would like to see more tools become free and open source :P

[Ondřej Čertík]

Hi Brett,

On Thu, Jun 11, 2015 at 5:37 PM, Brett Smith <waym...@gmail.com> wrote:
> Thank you for getting back to me so quickly!
> I am actually excited about this program, but it seems like I will need to
> learn a thing or two about numerical analysis.
> Our class is supposed to use TK Solver, and it seems to work OK. However, I

Interesting, it looks like TK Solver
(http://en.wikipedia.org/wiki/TK_Solver) was invented by Miloš
Konopásek (http://en.wikipedia.org/wiki/Milos_Konopasek), born in the
same country as I was.

I think this would be a welcome addition. We need to learn about the
algorithms and state of the art to solve these kinds of problems. If
it is a polynomial system of equations, then Gröbner basis is one such
method and SymPy has such capability. But I am not an expert in that.
It looks like you are after a numerical solution, so that's what I
would recommend at first. It might be much easier to implement. You
should open source it. Depending on which language you'll use to
implement it, it can go to SciPy or a separate library. I like that
you are trying to find or implement an open source solution. I am the
same way, thus SymPy was born.

Can you describe the symbolic structure of those equations --- what
kind of non-linearity is allowed? Because perhaps you can tell SymPy
to first eliminate the linear system, which SymPy should be able to
do, and then we just need to handle the nonlinear part. One way is to
take one equation and substitute into the rest, and so on. This only
works for some non-linear systems. If it is more non-linear, it might
not be easy to eliminate a variable from an equation, then you need
something more advanced.

E.g. these are a bit tricky:

> Sa^2 = 3*T14^2
> Sb^2 = 3*T12^2
> Sc^2 = 3*T43^2
> Sd^2 = 3*T32^2

But you can help it by solving two cases, first:

Sa = + sqrt(3) * T14
...

and the second:

Sa = - sqrt(3) * T14

And you'll have two solutions. It might be that the whole system
doesn't have a solution for one of the cases, so then you end up with
just one solution overall.

What is the application of such systems?

Ondrej

> https://groups.google.com/d/msgid/sympy/0ca6f806-1316-4903-8cc2-99732b6b9e07%40googlegroups.com.

[Ondřej Čertík]

On Thu, Jun 11, 2015 at 9:35 PM, Ondřej Čertík <ondrej...@gmail.com> wrote:
> Hi Brett,
>

Actually, you have 2 cases per each of the 4 equations, thus 2^4 = 16
cases total, you need to do all possible combinations.

[Brett Smith]

This system of equations calculates the reaction forces internal to a machine structure (in this case, vise grip pliers). Then uses the reaction forces to solve for stress,  and eventually, the factor of safety on each pin connection.

I have made the changes you have suggested and I still don't get a solution.
I was thinking that by assigning each variable twice might over define the system:

Sa = - sqrt(3) * T14 AND Sa = + sqrt(3) * T14 

Will sympy know to solve this as a pair?

Even if I comment out all the negative assignments, I still do not get a result...

A numerical result is necessary, but in addition I'd like to play with and vary the parameters, and perhaps solve large systems for variables I choose to isolate.

Looks Grobner basis might be a good place to start looking. I wanted to get all of this functionality working this weekend so I could use it for a class project....... buuuuut I might have been overly optimistic with my timeline haha. Looks like I could learn a lot from studying these kind of numerical analysis methods. 

For my sympy.solve() call I set the flat "Manual=true", does this change the method Sympy uses to solve the system? I like the idea of the computer solving equations the same way I do, but perhaps eventually it is simply not practical or even possible.

It seems like for most solvers with non-linear systems, you must give the solver some kind of starting value for it to begin its work.

I have never actually had any of the code I've wrote make it into an official distributable  package, so that would be really awesome!

Looks like you have done a lot of work for this module! Open source projects like these are primarily responsible for any success I've had both professionally and academically :P

[Ondřej Čertík]

On Fri, Jun 12, 2015 at 2:53 AM, Brett Smith <waym...@gmail.com> wrote:
> This system of equations calculates the reaction forces internal to a
> machine structure (in this case, vise grip pliers). Then uses the reaction
> forces to solve for stress, and eventually, the factor of safety on each
> pin connection.
>
> I have made the changes you have suggested and I still don't get a solution.
> I was thinking that by assigning each variable twice might over define the
> system:
>
> Sa = - sqrt(3) * T14 AND Sa = + sqrt(3) * T14
>
> Will sympy know to solve this as a pair?

No, you have to plug them there separately. As I said, there is 16 such options.

There seem to be other non-linearities, like:

F12x + F32x
F32y + F12y + P
R32x*F32y - R32y*F32x + R12x*F12y - R12y*F12x + Rp*P
F43x + F23x
F43y + F23y + Fh
R43x*F43y - R43y*F43x + R23x*F23y - R23y*F23x + Rh*Fh
F14x + F34x
F14y + F34y
R14x*F14y - R14y*F14x + R34x*F34y - R34y*F34x

aren't these missing an equal sign (=)? Or are these equal to 0? I
don't understand yet.

> Even if I comment out all the negative assignments, I still do not get a
> result...
>
> A numerical result is necessary, but in addition I'd like to play with and
> vary the parameters, and perhaps solve large systems for variables I choose
> to isolate.

You should start by numerical solution of the system. Then if you find
one, we can try to find a symbolic solution. At least we'll know that
there is one.

Ondrej

> https://groups.google.com/d/msgid/sympy/db23d376-edff-4bc7-8cab-dce9d2538ed7%40googlegroups.com.

[Brett Smith]

Yes, sympy automatically equates these equations to zero.

You really are not supposed to write the equations this way when programming with Sympy, but for this application I wanted to be able to type equations in a little easier.

I prefer:

x^2 + y = 15 

to:

Eq(x**2 + y, 15)

I have a large amount of equations to get through in a large amount of problem sets :P

I believe the syntax above is a little more familiar. The following function in the solver.py makes it understandable to sympy:

def parse_input(_file): #Converts file into array of Sympy Functions

_functions = []  

original_functions = []

lines = []

with open(_file,'r+') as f:

lines = f.readlines()

for line in lines:

line = line[:-1]

comment = line[0] == '#'

if not comment:

if "#" in line: #remove comments

split = string.find(line,"#")

line = line[0:split]

original_functions.append(line) #Store original syntax

line = string.replace(line," ","") #remove white space

line = string.replace(line," ","")

if "^" in line: #Convert conventional equation syntax in sympy syntax

line = string.replace(line,"^","**")

if "=" in line:

split = string.find(line,'=')

line = "Eq(" + line[0:split] + "," + line[split + 1:len(line)] + ")"

temp = parse_expr(line)

_functions.append(temp)

return _functions, original_functions

If its two cases for each of the four, is that not 8? I am afraid I don't understand.

Do i need to solve for each term? i.e:

> Sa = + sqrt(3)*T14

> Sa = - sqrt(3)*T14

> T14 = Sa/sqrt(3)

> T14 = -Sa/sqrt(3)

[Brett Smith]

Just thought I'd update:

I've had a decent amount of success using sympy to solve large systems of equations!

This is making school much easier in my opinion!

I have added a function to my solver that iterates through and replace any symbols that have an obvious numeric answer, and another that iterates through and substitutes symbols from equations where there are only two variables.

This seemed to tame the large system we had put together for my machine design class.

I will re-upload if anyone is interested. This is one of my first python modules; I want to learn to code them correctly. I think this could have a lot of utility.

Still learning how to manipulate data...

The "Wheel_Chair_Solver.py" takes the data from P-Lift.txt and generates solutions for S-Lift.txt

solver.py is my imported module

Anyway, I had one more question. 

I wanted to create a symbolic solution set for one variable in a system of equations.

If you solved each equation for each variable in that system it seems like you'd be poised to generate a large number of solutions for that variable - via substitution.

It seems to me that a symbolic system of equations has a finite solution set as long you don't start substituting an equation back into into themselves. Is this seem correct?

I realize this would take a lot of processing power, so i started looking in PyOpenCL as means of generate large amounts of parallel processes.

However, it seems that GPU programming is really only capable of numerical analysis, and that it would be impossible to pass sympy functions/variables to a GPU kernel?

^^I tried doing this with multiprocessing at the bottom of my solver.py file. The functions have been commented out.

Anybody know if this is possible, or where to start looking??

 

Best regards,
Brett