Ask Your Question

Issue on numerical_integral with a Pari/GP function

asked 2021-03-18 22:37:39 +0200

jb gravatar image

updated 2021-03-18 22:41:37 +0200

Hello there,

I tried to evaluate the integral of a Pari/GP function but failed so far. Here is an example code on (SageMath version 8.9, Release Date: 2019-09-29 and windows 10):

sage: lchi4 = DirichletGroup(4).list()[1].lfunction(); lchi4
PARI L-function associated to Dirichlet character modulo 4 of conductor 4 mapping 3 |--> -1
sage: numerical_integral(lambda t: lchi4(t).real(), 2,3)
SystemError                               Traceback (most recent call last)
/opt/sagemath-8.9/local/lib/python2.7/site-packages/sage/all_cmdline.pyc in <module>()
----> 1 numerical_integral(lambda t: lchi4(t).real(), Integer(2),Integer(3))

/opt/sagemath-8.9/local/lib/python2.7/site-packages/sage/calculus/integration.pyx in 
sage.calculus.integration.numerical_integral (build/cythonized/sage/calculus/integration.c:4061)()
353          _b = b
354          W = <gsl_integration_workspace*> gsl_integration_workspace_alloc(n)
--> 355          sig_on()
356          gsl_integration_qag(&F,_a,_b,eps_abs,eps_rel,n,rule,W,&result,&abs_err)
357          sig_off()

SystemError: calling remove_from_pari_stack() inside sig_on()

However, for the Riemann zeta function, numerical_integral works fine as below:

sage: numerical_integral(lambda t: zeta(t), 2,3)
(1.3675256886839795, 1.518258506343328e-14)

It looks like it would work for lchi4 if it can be coerced into a symbolic expression like zeta(x) as shown below:

sage: type(lchi4(3))
<type 'sage.rings.complex_number.ComplexNumber'>
sage: type(zeta(3))
<type 'sage.symbolic.expression.Expression'>

Can you please let me know how to evaluate the integral numerically for a Pari/GP function as lchi4?

Thank you in advance.

edit retag flag offensive close merge delete

1 Answer

Sort by ยป oldest newest most voted

answered 2021-03-23 00:20:01 +0200

jb gravatar image

updated 2021-03-23 00:29:16 +0200

After searching Google for this issue, I've found myself an answer:

sage: def lchi4(s):
....:     return RR(gp('lfun(lfuncreate(-4), %f)' %s))
sage: lchi4(1)
sage: numerical_integral(lambda t: lchi4(t), 2,3)
(0.9465779691617932, 1.163268688464501e-09)

I hope this might help for those who have the same problem. Please note that this is not the solution for DirichletGroup(4).list()[1].lfunction() and it only works after changing GP element into RR as shown above.

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

1 follower


Asked: 2021-03-18 22:37:39 +0200

Seen: 25 times

Last updated: Mar 23