ASKSAGE: Sage Q&A Forum - RSS feedhttps://ask.sagemath.org/questions/Q&A Forum for SageenCopyright Sage, 2010. Some rights reserved under creative commons license.Mon, 14 May 2018 11:04:41 +0200inverse_laplace of a fraction whose denominator has real roots (or complex)https://ask.sagemath.org/question/42331/inverse_laplace-of-a-fraction-whose-denominator-has-real-roots-or-complex/ Hi, this is the simplest and very important method of Laplace inversion . It can be applied to any Laplace transform
(LT), by starting with a Pade approximation of the LT, then partial fractions and inversion. For some reason, a program I had written last year stopped working
def Ruin(Fx, rho, x, u):#assumes rational survival function Fx, i.e. a hyperexponential density
var('s')
L_F=laplace(Fx,x,s)
m1=L_F(s=0) #some algebraic manipulations known as Pollaczek-Khinchine formula
fe=L_F/m1
Fe=factor((1-fe)/s)
L_rui=rho*Fe/(1-rho*fe)
C=ComplexField(53);
dec=Frac(C['s'])(L_rui).partial_fraction_decomposition();
n=len(dec)
par=[inverse_laplace(dec[1][i],s,u) for i in [0..n]];
rui=sum(par)
return rui, L_rui, fe
It's pretty easy to repair the program when roots are rational -- see question "partial_fraction_decomposition" with possibly "complex roots", again. For nonrational roots, I proposed there the following test case
var('s,u')
R.<s> = PolynomialRing(QQbar)
F = R.fraction_field()
L=3/4*(19*s^2 + 156*s + 284)/(19*s^3 + 174*s^2 + 422*s + 228)
whole,LL=L.partial_fraction_decomposition()
show(LL[0])
inverse_laplace(LL[0],s,u)
Thanks in advance :)florinMon, 14 May 2018 11:04:41 +0200https://ask.sagemath.org/question/42331/"partial_fraction_decomposition" with possibly "complex roots", againhttps://ask.sagemath.org/question/42220/partial_fraction_decomposition-with-possibly-complex-roots-again/ Hi I come back to this question, even though it's been answered before, since I am still not able to make it work. I should mention maybe that I am just trying to teach undergraduate students to invert
rational Laplace transforms (for myself I am able to afford a Mathematica licence, but it would be nice to be able to show students that such simple things may be done nowadays for free). Following a previous answer, I tried
L = 2*(s + 3)/(3*s^2 + 13*s + 10)
Ks=FractionField(PolynomialRing(CC,names='s'))
Lr=Ks(L)
Of course, with quadratic rational roots, I could do this by hand, but the purpose here is to do it when you do not know the roots. The first two commands work, but the third has error message
TypeError: ('cannot convert {!r}/{!r} to an element of {}', 2*(s + 3)/(3*s^2 + 13*s + 10), 1.00000000000000, Fraction Field of Univariate Polynomial Ring in s over Complex Field with 53 bits of precision)
I should add that a different attempt which used to work last year
C=ComplexField(53);
dec=Frac(C['s'])(Lrs).partial_fraction_decomposition();
gets now same error message. So, this is probably due to an "improvement" of sageflorinMon, 30 Apr 2018 10:58:31 +0200https://ask.sagemath.org/question/42220/"partial_fraction_decomposition" with "complex roots"https://ask.sagemath.org/question/41893/partial_fraction_decomposition-with-complex-roots/ Hi six months ago I used a hint from Zimmerman's book
var('s')
L=2*(s + 3)/(3*s^2 + 13*s + 10)
C=ComplexField(53)
dec=Frac(C['s'])(L).partial_fraction_decomposition()
but this does not work anymore (see below). How to do this? Thanks
----> 2 dec=Frac(C['s'])(Lrs).partial_fraction_decomposition();
TypeError: ('cannot convert {!r}/{!r} to an element of {}', 2*(s + 3)/(3*s^2 + 13*s + 10), 1.00000000000000, Fraction Field of Univariate Polynomial Ring in s over Complex Field with 53 bits of precision)
/opt/sagemath-8.1/src/sage/structure/parent.pyx in sage.structure.parent.Parent.__call__ (build/cythonized/sage/structure/parent.c:9641)()
937 if mor is not None:
938 if no_extra_args:
--> 939 return mor._call_(x)
940 else:
941 return mor._call_with_args(x, args, kwds)
/opt/sagemath-8.1/src/sage/structure/coerce_maps.pyx in sage.structure.coerce_maps.DefaultConvertMap_unique._call_ (build/cythonized/sage/structure/coerce_maps.c:4928)()
152 print(type(C), C)
153 print(type(C._element_constructor), C._element_constructor)
--> 154 raise
155
156 cpdef Element _call_with_args(self, x, args=(), kwds={}):
/opt/sagemath-8.1/src/sage/structure/coerce_maps.pyx in sage.structure.coerce_maps.DefaultConvertMap_unique._call_ (build/cythonized/sage/structure/coerce_maps.c:4796)()
147 cdef Parent C = self._codomain
148 try:
--> 149 return C._element_constructor(x)
150 except Exception:
151 if print_warnings:
/opt/sagemath-8.1/local/lib/python2.7/site-packages/sage/rings/fraction_field.py in _element_constructor_(self, x, y, coerce)
616 except AttributeError:
617 raise TypeError("cannot convert {!r}/{!r} to an element of {}",
--> 618 x0, y0, self)
619 try:
620 return self._element_class(self, x, y, coerce=coerce)
florinFri, 06 Apr 2018 21:42:27 +0200https://ask.sagemath.org/question/41893/Change Precision of complex_roots()https://ask.sagemath.org/question/39987/change-precision-of-complex_roots/I am trying to find the complex roots of the polynomial
poly = x^7 - 6*x^6 + 15*x^5 - 20*x^4 + 15*x^3 - 6*x^2 + x
But when I do poly.complex_roots(), the system gives:
PariError: overflow in expo()
Apparently there are options for how much precision you want when computing roots -- one option is to use Pari, which is the high-precision option, and the other NumPy, which is the low-precision option. The default is set to use Pari, which apparently overloads when I try to compute the roots of this polynomial (and many others as well, this polynomial is just one example).
How do I change the complex_roots() function to get lower-precision roots?cshiringThu, 07 Dec 2017 06:10:52 +0100https://ask.sagemath.org/question/39987/Change Precision on complex_roots()https://ask.sagemath.org/question/39988/change-precision-on-complex_roots/ I am trying to find the complex roots of the polynomial
poly = x^7 - 6*x^6 + 15*x^5 - 20*x^4 + 15*x^3 - 6*x^2 + x
But when I do poly.complex_roots(), the system gives:
PariError: overflow in expo()
Apparently there are options for how much precision you want when computing roots -- one option is to use Pari, which is the high-precision option, and the other NumPy, which is the low-precision option. The default is set to use Pari, which apparently overloads when I try to compute the roots of this polynomial (and many others as well, this polynomial is just one example).
How do I change the complex_roots() function to get lower-precision roots?
Alternatively, how do I deal with the PariError?
cshiringThu, 07 Dec 2017 06:13:28 +0100https://ask.sagemath.org/question/39988/Finding roots of equation made of Bessel functions, some of which have complex argumentshttps://ask.sagemath.org/question/39863/finding-roots-of-equation-made-of-bessel-functions-some-of-which-have-complex-arguments/I want to find the roots (values or w) for an equation of spherical Bessel functions for L = 1 and \Omega^2 = 14.28 (capital Omega, not w):
My code is:
L = 1
ome = 14.28 #MeV
w = var('w')
Arg1 = sqrt(ome*(1-w^2))
Arg2 = I*sqrt(ome)*w
a = -I*sqrt(1-w^2)
b = spherical_bessel_J(L-1, Arg1)
c = spherical_bessel_J(L, Arg1)
d = spherical_bessel_J(L, Arg2)+I*spherical_bessel_Y(L, Arg2)
e = spherical_bessel_J(L-1, Arg2)+I*spherical_bessel_Y(L-1, Arg2)
eq = a*(b/c)*(d/e)-w == 0
eq.find_root(0, 1)
But I get an error on line 12, which specifies:
"TypeError: unable to convert 3.77888872553824*I to float; use abs() or real_part() as desired"3ndingFri, 01 Dec 2017 05:55:27 +0100https://ask.sagemath.org/question/39863/Finding complex roots numerically using sagehttps://ask.sagemath.org/question/8205/finding-complex-roots-numerically-using-sage/Can sage find complex roots numerically for equations like
cos(cos(cos(cos(x)))) == sin(sin(sin(sin(x))))
if so, how?ebsSat, 02 Jul 2011 21:04:50 +0200https://ask.sagemath.org/question/8205/