1 | initial version |
The working code for n=3 is
A=[1/6,2/6,3/6];ex=[2,6,4];rho=5/8;var('x s u')
F=sum(A[i]*exp(-x*ex[i]) for i in [0..2])#sum of exponentials
R.<s> = PolynomialRing(QQbar)
FF = R.fraction_field()
L_F=sum(A[i] /(s+ex[i]) for i in [0..2]);#Laplace transform of F
#Compute Pollackek-Khinchine formula L_rui for the Laplace transform of ruin probability
m1=L_F(s=0)
fe=L_F/m1
Fe=(1-fe)/s
show(Fe) #note the s is simplified when working in QQbar
L_rui=rho*Fe/(1-rho*fe)
R.change_ring(RR)
FF = R.fraction_field()
wh,LL=L_rui.partial_fraction_decomposition()
show(LL[0])
with output
[0.5332274783298694?s+1.087668109036214?,0.07804278965136865?s+2.967947295921641?,0.01372973201876195?s+5.575963542410567?]
I need now to "remove" somehow the ? before taking
inverse_laplace(LL[0],s,u)
It would be nice if something like
inverse_laplace(LL[0].change_ring(RR),s,u)
or inverse_laplace(RR(LL[0]),s,u)
would work