Round symbolic expressions

59 views
Skip to first unread message

Cristóvão Sousa

unread,
Oct 24, 2010, 6:41:54 AM10/24/10
to sage-s...@googlegroups.com
Hi,

Is it possible to round reals inside a symbolic expression?

Example,

sage: x = var('x')
sage: r = 0.1*0.1 - 0.01
sage: sr = 0.1*0.1*x - 0.01*x
sage: r
1.73472347597681e-18
sage: sr
(1.73472347597681e-18)*x

sage: round(r,10)
0.0

is OK, but

sage: round(sr,10)

doesn't output 0.0 but raises an error instead (and I understand why, as the
expression evaluation really depends on the value of x).

But, is there any way of round reals even if they are inside a symbolic
expression?

Thanks,
Cristóvão Sousa

Cristóvão Sousa

unread,
Oct 23, 2010, 6:50:44 PM10/23/10
to sage-s...@googlegroups.com

Simon King

unread,
Oct 24, 2010, 7:34:47 AM10/24/10
to sage-support
Hi Cristóvão!

On 24 Okt., 00:50, Cristóvão Sousa <cris...@gmail.com> wrote:
> But, is there any way of round reals even if they are inside a symbolic
> expression?

I guess it would be needed to define a recursive function for that
purpose, using operator and operands of a symbolic expression. Such
as:

sage: def symround(x, ndigits=0):
....: if hasattr(x,'operator') and x.operator():
....: return x.operator()(*map(lambda
y:symround(y,ndigits=ndigits), x.operands()))
....: try:
....: return round(x,ndigits=ndigits)
....: except TypeError:
....: return x
....:
sage: sr = (1.01*x^2+(0.1*0.1-0.01)*x)*2.0001*sin(3.01*x)
sage: sr
(2.02010100000000*x^2 +
(3.46962042430121e-18)*x)*sin(3.01000000000000*x)
sage: symround(sr)
2.0*x^2.0*sin(3.0*x)
sage: symround(sr,ndigits=10)
2.020101*x^2.0*sin(3.01*x)

There might be ways to do it more efficient, though. I don't know.

Cheers,
Simon

Cristóvão Sousa

unread,
Oct 24, 2010, 9:30:12 PM10/24/10
to sage-support
On Oct 24, 12:34 pm, Simon King <simon.k...@nuigalway.ie> wrote:
> I guess it would be needed to define a recursive function for that
> purpose, using operator and operands of a symbolic expression. Such
> as:
>
> sage: def symround(x, ndigits=0):
> ....:     if hasattr(x,'operator') and x.operator():
> ....:         return x.operator()(*map(lambda
> y:symround(y,ndigits=ndigits), x.operands()))
> ....:     try:
> ....:         return round(x,ndigits=ndigits)
> ....:     except TypeError:
> ....:         return x

Thank you Simon!

It resolves my problem, and the efficiency is not an issue.

It just has a minor bug when the operator has more than two operands,
like x+y+z, but I'll try to fix it as I got the picture now.

Thanks,
Cristóvão Sousa

Simon King

unread,
Oct 25, 2010, 3:41:07 AM10/25/10
to sage-support
Hi Cristóvão!

On 25 Okt., 03:30, Cristóvão Sousa <cris...@gmail.com> wrote:
> ...
> It just has a minor bug when the operator has more than two operands,
> like x+y+z, but I'll try to fix it as I got the picture now.

Yes, to my surprise, the "add" operator only accepts two arguments,
but the list of operands of a sum may be longer:
sage: s = 0.001*x^2+0.01*x+0.1*sin(1.01*x)+1
sage: s.operands()
[0.00100000000000000*x^2, 0.0100000000000000*x,
0.100000000000000*sin(1.01000000000000*x), 1]
sage: s.operator()(*s.operands())
---------------------------------------------------------------------------
TypeError Traceback (most recent call
last)

/home/king/<ipython console> in <module>()

TypeError: op_add expected 2 arguments, got 4

I think that one should have
s.operator()(*s.operands()) == s
for any symbolic expression s.

@symbolic experts (Burcin et al):
Is it really necessary that x.operator() returns None and x.operands()
returns []? What about an identity operator?

Is it really necessary to break s.operands() into smaller pieces in
order to reconstruct a sum? Why can' op_add accept an argument list of
arbitrary length, in particular since the list of operands of a sum
can be longer than two?

Cheers,
Simon

Burcin Erocal

unread,
Oct 25, 2010, 5:04:32 AM10/25/10
to sage-s...@googlegroups.com
Hi Simon,

On Mon, 25 Oct 2010 00:41:07 -0700 (PDT)
Simon King <simon...@nuigalway.ie> wrote:

> @symbolic experts (Burcin et al):
> Is it really necessary that x.operator() returns None and x.operands()
> returns []? What about an identity operator?

The operator and operands don't have a meaning if you have a single
variable. I agree that the current behavior is confusing:

sage: x.operands()
[]
sage: x.operator() is None
True

I suggest we raise a ValueError when there is no operator or operands.
This is already done for iterators of symbolic expressions in #7537:

http://trac.sagemath.org/sage_trac/ticket/7537

Can you open a ticket to do the same for operands() and operator()?

> Is it really necessary to break s.operands() into smaller pieces in
> order to reconstruct a sum? Why can' op_add accept an argument list of
> arbitrary length, in particular since the list of operands of a sum
> can be longer than two?

You are right, for add and mul we should return a function that can
handle multiple arguments. I am not sure if the top level sum() and
prod() functions would be suitable here though.

Can you open a ticket for this as well?


Thank you.

Burcin

Simon King

unread,
Oct 25, 2010, 8:09:06 AM10/25/10
to sage-support
Hi Burcin,

On 25 Okt., 11:04, Burcin Erocal <bur...@erocal.org> wrote:
> I suggest we raise a ValueError when there is no operator or operands.
> This is already done for iterators of symbolic expressions in #7537:
>
> http://trac.sagemath.org/sage_trac/ticket/7537
>
> Can you open a ticket to do the same for operands() and operator()?
>
> > Is it really necessary to break s.operands() into smaller pieces in
> > order to reconstruct a sum? Why can' op_add accept an argument list of
> > arbitrary length, in particular since the list of operands of a sum
> > can be longer than two?
>
> You are right, for add and mul we should return a function that can
> handle multiple arguments. I am not sure if the top level sum() and
> prod() functions would be suitable here though.
>
> Can you open a ticket for this as well?

I only opened one ticket, namely #10169, aiming at making
s==s.operator()(*s.operands()) work uniformely.

Cheers,
Simon

Burcin Erocal

unread,
Oct 25, 2010, 8:39:52 AM10/25/10
to sage-s...@googlegroups.com
Hi Simon,

On Mon, 25 Oct 2010 05:09:06 -0700 (PDT)
Simon King <simon...@uni-jena.de> wrote:

> I only opened one ticket, namely #10169, aiming at making
> s==s.operator()(*s.operands()) work uniformely.

I don't think this makes sense for variables, numeric objects, or
constants, in other words, objects for which .operator() returns None
at the moment.

If we return an identity operator for these cases, how do you plan to
test for it in your code:

sage: def symround(x, ndigits=0):
....: if hasattr(x,'operator') and x.operator():
....: return x.operator()(*map(lambda
y:symround(y,ndigits=ndigits), x.operands()))
....: try:
....: return round(x,ndigits=ndigits)
....: except TypeError:
....: return x

....:


Isn't it simpler (and faster) to check for None?


Cheers,
Burcin

Simon King

unread,
Oct 25, 2010, 9:09:16 AM10/25/10
to sage-support
Hi Burcin!

On 25 Okt., 14:39, Burcin Erocal <bur...@erocal.org> wrote:
> If we return an identity operator for these cases, how do you plan to
> test for it in your code:

Something like this:

L = x.operands()
if len(L)>1:
return x.operator()(*map(lambda ..., L))
else:
try:
return x.operator()(round(L[0],...))
except TypeError:
return x

Cheers,
Simon

Simon King

unread,
Oct 25, 2010, 9:15:21 AM10/25/10
to sage-support
PS:

I know that testing "is None" is faster than "len(L)>1" and wouldn't
insist that there *has* to be an identity operator. One has to
consider two different cases anyway.

However, if there *is* an operator, s==s.operator()(*s.operands())
should hold.

Cheers,
Simon

Burcin Erocal

unread,
Oct 25, 2010, 9:26:37 AM10/25/10
to sage-s...@googlegroups.com
Hi Simon,

This initializes a list with a single element for objects which return
None for operator() now. IMHO, this approach is inefficient. In this
case, you should act on the object directly.


In any case, we should wrap the following ginac interfaces to provide a
better way of doing this:

* A way to apply a function to the operands of an expression:

http://www.ginac.de/tutorial/Applying-a-function-on-subexpressions.html

* or tree traversal:

http://www.ginac.de/tutorial/Visitors-and-tree-traversal.html


Cheers,
Burcin

Simon King

unread,
Oct 25, 2010, 9:45:01 AM10/25/10
to sage-support
Hi Burcin,

On 25 Okt., 15:26, Burcin Erocal <bur...@erocal.org> wrote:
> This initializes a list with a single element for objects which return
> None for operator() now. IMHO, this approach is inefficient. In this
> case, you should act on the object directly.

I didnt claim that it was efficient - my code was intended to be a
work-around to achieve what the OP was asking.
Directly acting on the Ginac interface, as you suggest below, is
certainly a better way.

And as I stated in my PS, I wouldn't insist on introducing an identity
operator. After all, a case distinction between leaves and non-leaves
in an expression tree can hardly be avoided.

> In any case, we should wrap the following ginac interfaces to provide a
> better way of doing this:
>
>  * A way to apply a function to the operands of an expression:
>
> http://www.ginac.de/tutorial/Applying-a-function-on-subexpressions.html
>
>  * or tree traversal:
>
> http://www.ginac.de/tutorial/Visitors-and-tree-traversal.html

+1

Cheers,
Simon

Cristóvão Sousa

unread,
Oct 25, 2010, 1:26:21 PM10/25/10
to sage-support
I made a modification so it can take 1, 3, 4, ... number of operands
nicely.

def symround(x, ndigits=0):
if hasattr(x,'operator') and x.operator():
o = map( lambda y : symround(y,ndigits=ndigits) , x.operands() )
r = o[0]
for i in xrange(1,x.nops()):
r = x.operator()(r,o[i])
return r
try:
r = round(x,ndigits=ndigits)
if r == x: return x
else: return r
except TypeError:
return x

Then extra 'if r == x' check was introduced to avoid transformation of
unrounded ints into floats (avoiding 2x to 2.0x).

I may be ugly in terms of efficiency, but where I need it, it is a
perfect workaround.

Thanks,
Cristóvão Sousa

Francois Maltey

unread,
Oct 26, 2010, 6:04:34 AM10/26/10
to sage-s...@googlegroups.com, cri...@gmail.com, Francois Maltey
Crist�v�o wrote :

I propose a generic function mapexpression in order operate in the tree
of the expression.
You can use it over leaf (ie the numbers or the variables) or over
functions (by the fctfct parameter)

I like this function because we can use this same function for numerical
transforms (for Christavo)
and rewrite rules (with trigonometric rules).

There are 2 fuzzy cases for operator add and mul because the
fct(*operand) doesn't work for add and mul.

import operator

def mapexpression (expr, fctleaf, fctfct, param, level, addDepth=0,
mulDepth=0) :
def mapex (expr, depth) : # a very local function
opor = expr.operator()
opands = expr.operands()
if (opor == None) :
return fctleaf(expr,param) # a leaf in the expression tree
if (opor == operator.add) : # recursive call thru sum
opands = map (lambda ex : mapex (ex, depth+addDepth), opands)
return add (opands)
if (opor == operator.mul) : # recursive call thru mul
opands = map (lambda ex : mapex (ex, depth+mulDepth), opands)
return mul (opands)
if (level == -1) or (level[-1] >= depth) : # recursive call over
operands
opands = map (lambda ex : mapex (ex, depth+1), opands)
if level == -1 or depth in level : # root of the subtree must be
changed
return fctfct (opor, opands, param)
return opor (*opands) # opands may or maynot be changed by a
recursive call
return mapex (expr, 0)

def nn (ex, p) :
if ex._is_numeric() or ex==pi: return ex.n(prec=p)
else : return ex

sage: mapexpression (sin(5*pi/7)+exp(i*pi/13),nn, lambda x,y,z:x(*y),20,-1)
1.7528 + 0.23932*I

Method __.n() doesn't work over partial numerical formula, but this
mapexpression allows this.

Francois

Reply all
Reply to author
Forward
0 new messages