Hi,
I hit a bit of a set-back today with my integration code. Specifically,
I am getting wrong results, because of improper treatment of branch cuts.
As a first example consider
sin(sqrt(z)) - sin(sqrt(w)),
where we want to substitute z=exp(I*pi)=-1 and w=exp(-I*pi)=-1. If we do
this naively, the answer comes out as zero, but what we (at least I, in
my integration code) mean is 2*sin(I).
I have for now been alleviating this problem (since I started the
integration code) by explicitely tracking the branch we are on by
(essentially) substituting early z=|z|*exp(I*arg(z)/pi*my_pi), where
my_pi is a dummy that is not recognised by any of the simplifications
code. When doing this, the above expression gets changed to
sin(exp(I*my_pi/2)) - sin(exp(-I*my_pi/2))
and in the end we can safely substitute pi for my_pi and get a correctly
computed branch.
However, this does not work for all functions. Generically, this will
fail when we replace sin(sqrt(z)) above by a non-elementary function
(i.e. one that is not recognised by the simplification machinery)
branched at the origin. A particular example I have hit is uppergamma(0,
z). As you may know, the function uppergamma(s, z) (for fixed s) is an
entire function when s is a positive integer, satisfies uppergamma(s, z)
= z**s * F(s, z) for an entire function F if s is not an integer, and
satisfies uppergamma(s, exp(2*pi*I*n)*z) = uppergamma(s, z) +
2*pi*n/factorial(-s) for non-positive integral s.
My question thus is how to put this information into code. Notice that
all functions in question come from meijer g-functions, but I don't
think there is a chance of understanding the monodromy action in this
generality (but I would love to be wrong!). Assuming this is true, the
best I could come up with would be to add
._eval_rewrite_as_{un,}branched functions to special functions, perhaps
even with a default implementation raising NotImplementedError. This
would be somewhat akin to the ._eval_rewrite_as{in,}tractable used for
gruntz' algorithm.
[Since the answer is likely messy, any code should obviously try to
avoid calling these functions at any cost, but that is something I can
take care of locally.]
In more detail, somewhere in functions/elementary/complexes.py I would
add principal_branch(z) which represents z with its argument converted
to range (-pi, pi], and argument_period(z) which then satisfies z =
principal_branch(z)*exp(I*2*pi*argument_period(z)). Then e.g.
uppergamma._rewrite_as_unbranched could return an expression involving
these functions. The integration code would use a sequence of operations
as follows to compute the branch of the result:
# TODO try to avoid this when possible
res = res.rewrite('unbranched')
res = res.subs(my_pi, pi)
res = res.rewrite('branched')
My real questions are then as follows, I guess:
1. Does anyone see a better way around the issue?
2. Any suggestions for a better interface?
Thanks,
Tom
--
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.
