Page 8 - Laplace transform numerical inversion

Bruno Josso & Leif Larsen: Laplace transform numerical inversion - June 2012 - p 6/18
3.2 Testing function
Back to the exemple in section(1) , both the asymptotic rate qD (tD ) decline solution and its Laplace
transform L [q D (tD )] = ˆqD (p) can be expressed as follows:
qD (tD ) = Γe−µt L [.]
←−−−−−−−−→
ˆ q D (p) = Γ 1
p + µ
(10)
with:
Γ =
1
1
2 ln( 4AD
eγCA ) + S
(11)
AD =
A
r2
w
(12)
µ =
2π
AD Γ
(13)
AD denotes a dimensionless area, CA is the Dietz shape factor and γ denotes Euler's constant.
One will only focus here on this unique simple but illustrative asymptotic rate decline inversion
problem. Let ˆ
F (p) = Γ
(p+µ) , with p ∈ C, be a function to be inverted. Ones assumes that (µ, Γ) ∈ R∗2
+ .
As can be found in a Laplace transform dictionary, and as reminded in equation(10), inversion
can be performed analytically leading to function f (t):
f (t) = L−1 Γ
(p + µ)
= Γe− µt with: t ∈ R (14)
When applying the "final value theorem" to ˆ
F (p), one can demonstrate that, at late time t, the
function f(t) tends towards zero:
lim
p→0
p ˆ
F (p) = lim
p→0
Γ p
(p + µ)
= 0 (15)
Furthermore, if Γ > 0 the function f (t) always remains strictly positive. Therefore, one can write
this known results under another form :
lim
t→∞
Γ e− µt = 0 (16)
∀t ∈ R , Γ e− µt > 0 (17)
For a given couple (µ, Γ) ∈ R∗2
+ , one can see in Figure 3 the function f (t) = Γe−µt , which is the
analytical result of the inversion of ˆ
F (p) = Γ
(p+µ) , as well as the result obtained by a usual numeri-
cal inversion algorithm. In Figure 3,the analytical solution is called f(t) and it appears in green,
whereas the numerical inversion result is called Fp_inv(t) and it appears in blue; both curves
nicely overlap.
When zooming on the curves at late time, when they get close to zero, one could see as before on
a similar exemple that the reconstructed function resulting from Gaver-Stehfest inversion method
becomes negative whereas it should follow an always positive exponential decline.
In the coming sections, one will test the different inversion algorithms on this exponential func-
tion f (t) given in( 14) and perform such a close late time inspection. One will examine how some
Laplace numerical inversion methods converge towards the solution to see whether this inversion
artefact can be avoided.
Beyond the numerical inversion accuracy, one will particularly focus on the fact that for all t ∈ R,
the exponential function f (t) should remain strictly positive. The target here is to exhibit a
numerical inversion method that accurately inverts ˆ
F (p) = Γ
(p+µ) and always gives strictly positive
results.