Thanks guys, this is great stuff. You might want to memorialize this on a
wiki or something for newcomers to SymPy.
So a little bit about my background with Python: I'm a pretty experienced
Python programming, I understand operator overloading, I see your
@_sympifyit decorators all over the place (I'm guessing they're to convert 2
to Integer(2), as Aaron explained). In fact, I started writing my own
Python symbolics math library years ago.
One follow-up question is:
Is there an advantage to having your own classes for Pow, Add, Mul, and
other operators? Why didn't you just absorb those in classes for
expressions, functions, etc, and have an expression + expression return an
expression?
Here's an example. Specifically, x^2, instead of being Pow(x, 2), would be
something like:
MonoTerm (self.units = {x: 1}).__pow__(2) becomes MonoTerm, self.units = {x:
2}.
For example, digging up my old code, when I was making a symbolics math
library, I started doing it this way:
class MonoTerm:
"""A MonoTerm is an instance of terms multiplied together, but without
any sums or differences. As an example, 2 could be a MonoTerm, as could 3
z, but 2 x + 3 y must be a PolyTerm. MonoTerms can only be operated on by
MonoTerms with the same dimensional units."""
def __init__(self, string=''):
# instance variables
self.number=0
self.units={}
# initialize to the string, with the become method
self.become(string)
def become(self, string=''):
"""Handles the parsing of the MonoTerm into number and units."""
if string=='':
# nothing should happen; this is the "empty term"
return
if m.search('^\s*('+patterns.number+')', string):
self.number = float(m[0])
else:
self.number = 1.0 # 'kg' == '1.0 kg'
unitSearch =
re.findall('('+patterns.unit+')\s*(?:(?:\^|\*\*)\s*('+patterns.number+'))?',
string)
for unit, exponent in unitSearch:
if unit == '':
continue
if exponent == '':
self.units[unit]=1.0
else:
self.units[unit]=float(exponent)
def __add__(self, other):
"""Adds two MonoTerms together, provided they have the same
units."""
if self.units != other.units:
# can only subtract MonoTerms with MonoTerms
raise Exception('Cannot add two MonoTerms with incompatible
units. Use PolyTerms instead.')
try:
newTerm = MonoTerm()
newTerm.number = self.number + other.number
newTerm.units = self.units
newTerm.prune()
return newTerm
except:
return other+self
(etc)
class PolyTerm:
"""A PolyTerm is a term made up of the sum of two or more MonoTerms."""
def __init__(self, string=''):
# instance variables
self.terms={}
self.become(string)
def become(self, string=''):
# assumes the string is well-formed various terms split by +
# Actually, for now, assume the string is really a MonoTerm... this
will be updated later.
newTerm = MonoTerm(string)
self.terms[newTerm.getLabel()]=newTerm
def __add__(self, other):
"""Adds two PolyTerms together."""
newTerm = copy.deepcopy(self)
for term in other:
if term.getLabel() in newTerm.terms:
newTerm.terms[term.getLabel()] =
newTerm.terms[term.getLabel()] + term
else:
newTerm.terms[term.getLabel()] = copy.deepcopy(term)
newTerm.prune()
return newTerm
(etc)
term1 = PolyTerm('1 x y')
term2 = PolyTerm('1 y z')
term3 = PolyTerm('2 x y^2 z')
term4 = PolyTerm('1 y^2')
print ((term1+term2)**2 - term3) / term4 # ((xy+yz)^2 - 2xy^2z)/y^2 == x^2 +
z^2
--
Jeff Pickhardt
(650) 530-0036
[email protected]
On Thu, Jun 23, 2011 at 10:33 PM, Aaron Meurer <[email protected]> wrote:
> Well, you picked a particularly tricky example because it uses both
> the main core and the polys.
>
> So let's start with my_expresion (the core part). You already know
> that x is defined to be a Symbol object (because you did this
> yourself). You should also know that sin and cos are Function
> objects, which you can think of as container objects that hold
> arguments and which also have various mathematical relations defined
> on them.
>
> As you may know, Python let's you override the behavior of the
> built-in operations *, +, /, -, etc. on your own objects. So all
> objects have method __mul__, __add__, etc. defined on them. When you
> call x + y, this reduces to x.__add__(y). In SymPy, a.__add__(b) is
> converted to Add(a, b). The same is true for __mul__ and Mul. So in
> sin(x)**2 + 2*sin(x) + 1, sin(x) creates a sin object (when it is
> called). This reduces to
>
> sin(x).__pow__(2) + sin(x).__rmul__(2) + 1
> Pow(x, 2).__add__(Mul(2, sin(x))).__add__(1)
> Add(Pow(x, 2), Mul(2, sin(x)), 1)
>
> which is what you will get if you call srepr(my_expression). Note
> that 2*sin(x) actually calls __rmul__. That's because 2 (type int)
> doesn't know how to multiply sin(x) (type sin), so sin(x)'s __rmul__
> method is called. This brings up an important point. All Python
> types are converted in this process to SymPy types, through the
> sympify() function. So, for example, sin(x).__pow__(2) reduces to
> sin(x).__pow__(sympify(2)), which results in
> sin(x).__pow__(Integer(2)).
>
> Now, Mul, Pow, and Add's __new__ methods contain logic to do automatic
> simplifications like x + 2*x => 3*x or x*x => x**2. In this case, no
> such simplifications needed to be applied.
>
> All the args for any SymPy expression are stored in expr.args. So you
> would have
>
> >>> my_expression.args
> (1, sin(x)**2, 2*sin(x))
> >>> my_expression.args[1].args
> (sin(x), 2)
>
> If you want to see the code for all of this, you can look at the files
> in sympy/core. The __op__ stuff is mostly in basic.py and expr.py.
> The flattening routines are in add.py, pow.py, and mul.py. The code
> for functions that sin is built off of is in function.py.
>
> Now, to the factoring part. my_expression.factor() is a shortcut to
> factor(my_expression). Because factoring is a polynomial algorithm,
> the expression has to be converted to a Poly first. Poly is able to
> represent polynomials with arbitrary symbolic "generators". In this
> case, it determines that it should use sin(x) as a generator, so it
> creates Poly(sin(x)**2 + 2*sin(x) + 1, sin(x)), which you can think of
> as a wrapper around the polynomial y**2 + 2*y + 1, where y is set to
> be sin(x) (for the purposes of Poly, it does not matter what the
> generators are, other than that coefficients cannot contain symbols
> from them, so you can think of it in this way).
>
> Then, it calls the factorization algorithm on Poly. If you are
> interested in how this works, I suggest you read the code and the
> papers referenced there. In this case, it is able to factor the
> polynomial using the squarefree factorization algorithm, which is
> actually not too difficult to understand. In the general case, it
> uses a complicated multivariate factorization algorithm that factors
> any multivariate polynomial into irreducibles.
>
> Anyway, Poly.factorlist returns something like [(Poly(sin(x) + 1),
> 2)]. factor() converts this into a normal SymPy expression (also
> called a Basic expression or Expr expression) by passing it to Mul and
> Pow (something like Mul(*[Pow(b, e) for b, e in expr]) would do it, I
> think).
>
> All the Poly stuff lives in sympy/polys. If you are interested, I can
> explain a little how they work (internal representation, etc.). The
> code for the Poly class lives in polytools.py, though the actual
> factorization algorithm lives in sqftools.py and factortools.py.
>
> And one thing that I didn't mention (and maybe you didn't even think
> of) is the printing. You do not use pretty printing, so the printing
> is rather simple (just recursively print the objects using the proper
> operators). If you are interested, you can look at the code in
> sympy/printing.
>
> Let me know if this made sense, and if there are bits that you still
> would like to know about. Also, remember that SymPy is written in a
> fairly modular way, so it's completely unnecessary to know how a
> module works unless you want to work on that module specifically
> (e.g., you don't need to how the core works to write some
> simplification algorithm, like in simplify.py).
>
> Aaron Meurer
>
> On Thu, Jun 23, 2011 at 10:42 PM, Jeff Pickhardt <[email protected]>
> wrote:
> > Hey guys,
> >
> > I'm reading through the SymPy code and, well, it's somewhat overwhelming
> if
> > you're new to the project because there's so much going on. (That's a
> good
> > thing too - it means its robust!)
> >
> > Can someone help explain how this works?
> >
> >>>> from sympy import *
> >>>> x = Symbol("x")
> >>>> my_expression = sin(x)**2 + 2*sin(x) + 1
> >>>> my_expression.factor()
> > (1 + sin(x))**2
> >>>>
> >
> > For instance, what data structures happen when I create my_expression,
> what
> > happens when I factor it, etc. A high-level walk through would help. I
> see
> > there's stuff going on at polytools.py, and I think _symbolic_factor gets
> > called. It's just confusing to keep everything in my head when I don't
> yet
> > have a high level understanding of how sympy expressions and what not
> > actually work.
> >
> > Also, the new quantum mechanics stuff looks really cool. Wish I had had
> > that a few years ago!
> >
> > Thanks,
> > Jeff
> >
> > --
> > Jeff Pickhardt
> > (650) 530-0036
> > [email protected]
> >
> > --
> > You received this message because you are subscribed to the Google Groups
> > "sympy" group.
> > To post to this group, send email to [email protected].
> > To unsubscribe from this group, send email to
> > [email protected].
> > For more options, visit this group at
> > http://groups.google.com/group/sympy?hl=en.
> >
>
> --
> You received this message because you are subscribed to the Google Groups
> "sympy" group.
> To post to this group, send email to [email protected].
> To unsubscribe from this group, send email to
> [email protected].
> For more options, visit this group at
> http://groups.google.com/group/sympy?hl=en.
>
>
--
You received this message because you are subscribed to the Google Groups
"sympy" group.
To post to this group, send email to [email protected].
To unsubscribe from this group, send email to
[email protected].
For more options, visit this group at
http://groups.google.com/group/sympy?hl=en.
