Ask Your Question
4

Can sage find the poles of a complex function, and the residues at those poles?

asked 2011-04-11 12:49:16 -0600

Mark G gravatar image

Hi all,

I searched the SAGE reference and didn't find the word residue in realtion to functions of a complex variable. The function I'm working with is not the ratio of two polynomials in Z, which can be factored out to find the poles and zeros. It's actually like this:

f(z) = A* sinh(z)/ (B* z^2* sinh(z)+ C* z* cosh(z)-D* sin(z))

got this from laplace transforming a partial differential equation I want to solve. I'm following the technique used in a paper I'm trying to apply to a different geometry/coordinate system, and if I can get the residues at the poles, I can get a series solution to the PDE. However, I can't fill in the omitted steps in the paper due to my never having taken a course in complex variables. Hence my interest in finding the residues at the poles programmatically.

I would greatly appreciate a clue or two. Thank you.

edit retag flag offensive close merge delete

Comments

Although not directly Sage-related, perhaps you could use the trick described here [http://en.wikipedia.org/wiki/Residue_%28complex_analysis%29#Series_methods]. I know absolutely nothing about complex analysis. :-) On another note, the Sage reference manual says that there is a way to compute Laurent series (which would give you your residue) [http://www.sagemath.org/doc/reference/sage/symbolic/expression.html#sage.symbolic.expression.Expression.taylor]. But when I fiddled with it, I could only get taylor series. You could also take a look at this [http://www.sagemath.org/doc/reference/sage/rings/laurent_series_ring_element.html], but I couldn't figure anything useful out of that, either.

Kelvin Li gravatar imageKelvin Li ( 2011-04-11 16:51:27 -0600 )edit

You can find the zeroes of the denominator: "solve(f.denominator() == 0)". I don't know if there is a built-in way to compute residues from there.

John Palmieri gravatar imageJohn Palmieri ( 2011-04-13 09:12:56 -0600 )edit

Just as a data point, Mma apparently does this: http://reference.wolfram.com/mathematica/ref/Residue.html. Unsurprisingly, rjf has also written on this topic.

kcrisman gravatar imagekcrisman ( 2011-04-16 17:25:26 -0600 )edit

@John Palmieri - unfortunately, solve(f.denominator() == 0,x) doesn't help in this case, and using to_poly_solve=True gives a weird indexing error on the symbolic expression that comes out of the original non-solution that it's too late in this time zone for me to track down tonight. But I couldn't get W|A to do it, either; this probably doesn't have a closed-form solution, n'est pas?

kcrisman gravatar imagekcrisman ( 2011-04-16 18:00:34 -0600 )edit

2 answers

Sort by ยป oldest newest most voted
0

answered 2011-04-16 17:50:55 -0600

kcrisman gravatar image

Apparently Maxima includes this, and we didn't even know it! See for instance here and the official docs.

(%i1) residue (s/(s**2+a**2), s, a*%i);
(%o1)                           1/2

In Sage, that would be doable.

sage: f(x)=1/x
sage: f.maxima_methods().residue(x,0)
1
f(x)=1/(x^2+1)
f.maxima_methods().residue(x,i)
-1/2*I

Well, we should wrap this!

This does not solve the problem of finding poles. Maple apparently has this, though perhaps not for such an example as the original one.

edit flag offensive delete link more

Comments

kcrisman gravatar imagekcrisman ( 2011-04-16 17:58:52 -0600 )edit
0

answered 2011-04-16 13:41:08 -0600

parzan gravatar image

updated 2011-04-16 13:44:58 -0600

This won't work with a complicated example such as yours, but perhaps it will give you some ideas:

sage: f = sinh(x)/(sin(x)-1/2)
sage: (1/f == 0).solve(x)
[x == 1/6*pi]
sage: pole = _[0].rhs()
sage: (f*(x-pole)).limit(x=pole)
1/3*(e^(1/3*pi) - 1)*sqrt(3)*e^(-1/6*pi)
sage: f.taylor(x,pole,-1)*(x-pole)
2/3*sqrt(3)*sinh(1/6*pi)

the second line solves 1/f(x)=0, which for meromorphic functions is like solving f(x)=infinity. We then store the pole found in "pole". The last two lines calculates the residue, assuming the pole is simple (the two answers are identical).

Here is an example with a non-simple pole:

sage: f = sinh(x)/(sin(x)-1/2)-1/(x-pi/6)^2
sage: (1/f == 0).solve(x)
[x == 1/6*pi]
sage: pole = _[0].rhs()
sage: f(x=x+pole).taylor(x,0,-1).coeff(x,-1)
2/3*sqrt(3)*sinh(1/6*pi)

The substitution "x=x+pole" is needed - f.taylor(x,pole,-1).coeff(x-pole,-1), which should have given the same answer returns zero, since some automatic simplification takes place and messes things up:

sage: 1/(x-pole)                                  
-6/(pi - 6*x)

hope this helps!

edit flag offensive delete link more

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

Stats

Asked: 2011-04-11 12:49:16 -0600

Seen: 2,046 times

Last updated: Apr 16 '11