I have a problem with Laplace inversion and my function is not numerically stable for the Laplace inverse, but I do not understand the cause of this problem.
Here is my code and graph of this problem. Does anyone have any idea? I have also tried all 3 Dehoog
Stehfest
and Talbot
, they all have a different type of instabillity
import mpmath as mp
import numpy as np
import matplotlib.pyplot as plt
def Laplace_Max(s):
#mp.dps = 100
A= s**(5/4)/ (
0.00480931 +
0.0244077*s**(1/4) +
0.0129056*mp.exp(-35*s)*s**(1/4) +
0.00329167*mp.exp(0.707997*s)*s**(1/4) * mp.gammainc(0.0, 35.708*s,mp.inf, regularized=True) -
0.00530593*mp.gammainc(1.25, 35*s,mp.inf, regularized=True)
)
return A
t_R_inv = np.linspace(1,1200,12000)
Max_R=np.zeros(len(t_R_inv))
for i in range(len(t_R_inv)):
Max_R[i]=mp.invertlaplace(Laplace_Max, t_R_inv[i], method = 'stehfest', dps = 100, degree = 1000)
plt.plot(Max_R)
Update The inverse laplace of this function!
$\frac{1}{s}^{\frac{5}{4}}\cdot (0.00480931 + 0.02440766 \cdot s^{\frac{1}{4}} + 0.0129056 \cdot e^{-35 \cdot s} \cdot s^{\frac{1}{4}}+ 0.003291670 \cdot e^{0.7079967\cdot s} \cdot s^{\frac{1}{4}} \cdot \Gamma[0, 35.7079967 \cdot s] - 0.00530593 \cdot\Gamma[1.25, 35\cdot s])$