Ask Your Question
1

possible bug in residue function

asked 2020-01-21 10:57:16 +0200

rue82 gravatar image

I compute the residue of a function in two ways, using the fact that it has simple poles. It seems residue method gives the wrong anser (psi) to be compared to correct answer (chi).

   br(x) = 1-x
   q1,q2,q3,m = var('q1,q2,q3,m')
   assume(q1>0,q2>0,q3>0,m>0)
   q4 = (q1*q2*q3)^-1
   k = 2
   rho = [1,q4]
   X = [var("x%d" % i) for i in range(k)]
   chi1 = prod([ br(X[j]/m)/(br(X[j])*X[j]) for j in range(k)])
   chi2 = prod([ prod([ br(q1*q2*X[i]/X[j])*br(q1*q3*X[i]/X[j])*br(q2*q3*X[i]/X[j]) for i in range(k) if i != j]) for j in range(k)])
   chi3 = prod([ prod([ br(q1*X[i]/X[j])*br(q2*X[i]/X[j])*br(q3*X[i]/X[j])*br(q4*X[i]/X[j]) for i in range(k) if i != j]) for j in range (k)])
   chi = (chi1*chi2/chi3).factor()
   psi = chi
   for xi,rhoi in zip(X,rho):
        psi = psi.residue(xi==rhoi).combine().factor()
        chi = (chi*(xi-rhoi)).factor().subs({xi: rhoi})

Now we can print chi and psi and see that they're different. Is this a bug of residue or series?

edit retag flag offensive close merge delete

1 Answer

Sort by ยป oldest newest most voted
1

answered 2020-01-21 18:21:30 +0200

Emmanuel Charpentier gravatar image

updated 2020-01-22 21:03:46 +0200

WorksForMe(TM) in Sage 9.1?beta0 ... when done manually:

sage: list(zip(X,rho))
[(x0, 1), (x1, 1/(q1*q2*q3))]
sage: psi.residue(x0==1).combine().factor()
(q1*q2*x1 - 1)*(q1*q3*x1 - 1)*(q2*q3*x1 - 1)*(q1*q2 - x1)*(q1*q3 - x1)*(q2*q3 - x1)*(m - x1)*(m - 1)*q1^2*q2^2*q3^2/((q1*q2*q3*x1 - 1)*(q1*q2*q3 - x1)*(q1*x1 - 1)*(q2*x1 - 1)*(q3*x1 - 1)*m^2*(q1 - x1)*(q2 - x1)*(q3 - x1)*(x1 - 1))
sage: (chi*(x0-1)).factor().subs(x0==1)
(q1*q2*x1 - 1)*(q1*q3*x1 - 1)*(q2*q3*x1 - 1)*(q1*q2 - x1)*(q1*q3 - x1)*(q2*q3 - x1)*(m - x1)*(m - 1)*q1^2*q2^2*q3^2/((q1*q2*q3*x1 - 1)*(q1*q2*q3 - x1)*(q1*x1 - 1)*(q2*x1 - 1)*(q3*x1 - 1)*m^2*(q1 - x1)*(q2 - x1)*(q3 - x1)*(x1 - 1))
sage: bool(psi.residue(x0==1).combine().factor()==(chi*(x0-1)).factor().subs(x0==
....: 1))
True
sage: psi.residue(x1==1/(q1*q2*q3)).combine().factor()
0
sage: (chi*(x1-1/(q1*q2*q3))).factor().subs(x1==1/(q1*q2*q3))
0
sage: bool(psi.residue(x1==1/(q1*q2*q3)).combine().factor()==(chi*(x1-1/(q1*q2*q3
....: ))).factor().subs(x1==1/(q1*q2*q3)))
True

I think that the problem resides in your two last lines of code:

    psi = psi.residue(xi==rhoi).combine().factor()
    chi = (chi*(xi-rhoi)).factor().subs({xi: rhoi})

which scratch the previous values of psi and chi...

And, BTW, psi has a native definition in sage: from psi?:

Signature:      psi(x, *args, **kwds)
Docstring:     
   The digamma function, psi(x), is the logarithmic derivative of the
   gamma function.

      psi(x) = frac{d}{dx} log(Gamma(x)) =
      frac{Gamma'(x)}{Gamma(x)}

   We represent the n-th derivative of the digamma function with
   psi(n, x) or psi(n, x).

Unless you have a pressing necessity of redefining it, it might be wiser to leave psi's definition untouched...

EDIT :

However, my last two lines are exactly as they shuold be: I want to iterate over the residue operation

Okay. Let's replace your last five lines by:

chi = [(chi1*chi2/chi3).factor()]
psi = chi.copy()
for xi,rhoi in zip(X,rho):
    psi.insert(0,psi[0].residue(xi==rhoi).combine().factor())
    chi.insert(0,(chi[0]*(xi-rhoi)).factor().subs({xi: rhoi}))

Then check each step:

sage: list(map(lambda u,v:bool(u==v), psi, chi))
[False, True, True]

The first iteration gives the expected result, the second cannot be proved True.

edit flag offensive delete link more

Comments

The equalities can be checked with:

list(map(lambda xi,rhoi:bool(psi.residue(xi==rhoi).combine().factor() == (chi*(xi-rhoi)).factor().subs({xi: rhoi})), X, rho))

HTH,

Emmanuel Charpentier gravatar imageEmmanuel Charpentier ( 2020-01-22 06:35:27 +0200 )edit

Thanks for the tip about not redefining psi. However, my last two lines are exactly as they shuold be: I want to iterate over the residue operation. The problem I wanted to point out is that residue/series fails to give the correct result. Do you agree? Btw, do email notifications work for you?

rue82 gravatar imagerue82 ( 2020-01-22 10:55:21 +0200 )edit

In particular, if you substitute e.g. q1=3, q2=5, q3=7, m=11 you should see that chi[0] and psi[0] are different.

rue82 gravatar imagerue82 ( 2020-01-24 11:49:02 +0200 )edit

Your Answer

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

Add Answer

Question Tools

1 follower

Stats

Asked: 2020-01-21 10:57:16 +0200

Seen: 193 times

Last updated: Jan 22 '20