Laplace inverse in Python with mpmath

2019-07-17 02:22发布

问题:

I want to use "DE HOOG" algorithm for numerical Laplace inverse transform. I want to use the "mpmath" package and I installed it from the link:

https://github.com/klkuhlm/mpmath

Lets say I need to find the inverse Laplace transform of the below function at t=1:

f = 1/(s-1)

The inverse Laplace transform of f is : e^(t)

At t=1 the result is expected to be = e

import mpmath as mp
import numpy as np

def f(s):
    return 1 / (s-1)

t = np.linspace(0.01,0.5,10)

G = []

for i in range(0,4):
    G.append(mp.invlapdehoog(f, t[i]))

print G 

The thing is it works perfectly only when I set the range of the "i" less than 4. For example once I substitute:

for i in range(0,5): #or for i in range(0,more than 5):

I get this error:

enter image description here

Could you please help me fix this problem?

Thanks!

回答1:

The object InverseLaplaceTransform has an attribute degrees that dictates the levels of approximation necessary to achieve a given level of precision. Your copy of InverseLaplaceTransform updates degrees each time you call it with a smaller and smaller value. Eventually, degrees is so small the parameter fp has only one value, which is not enough to continue with further calculations.

Solution: edit your call to invlapdehoog to reset degrees each time. I suggest however calling invertlaplace directly rather than invlapdehoog.

for i in xrange(0,10):
    G.append(mp.invertlaplace(f, t[i], method = 'dehoog', degree = 18))

Edit: The original poster asked a related question in the comments to this solution. They asked why the computation time increases (quite drastically) with consecuitve calls to mp.invertlaplace. In short, the mp.invertlaplace is updating its attribute precision which dictates how many decimal places it should compute in calculating the inverse laplace. As with the above solution, we can pass in precision to each call to make sure we obtain the precision we want (eg - 10 decimal places):

for i in xrange(0,10):
    G.append(mp.invertlaplace(f, t[i], method = 'dehoog', dps = 10, degree = 18))

PS - you can apply the inverse laplace to all of t at once with the following snippet:

G = map( lambda x: mp.invertlaplace(f, x, method = 'dehoog', dps = 10, degree = 18), t)