# Laplace Transform and complex contours

1. Nov 2, 2005

### saltydog In the HW section, the following equation was proposed:

$$ty(t)=\int_0^t \tau^{1-\alpha}y(t-\tau)d\tau;\quad \int_0^\infty y(x)dx=1$$

Using Laplace Transforms:

$$Y(p)=\mathcal{L}\left\{y(t)\right\}$$

and noting the RHS is the inverse transform of a convolution, we obtain the IVP:

$$-Y^{'}=Y\frac{\Gamma(\alpha)}{p^\alpha};\quad Y(0)=\int_0^\infty e^{-0t}y(t)dt=1$$

solving for $\alpha \ne 1$:

$$Y(p)=Exp\left[-\frac{\Gamma(\alpha)}{1-\alpha}s^{1-\alpha}\right]$$

I (Mathematica) can invert this only for $\alpha=1/2$:

$$y(t)=\frac{e^{-\pi/t}}{t^{3/2}}$$

I'd like to know how to invert it for other values of alpha. My only recourse is to evaluate the complex contour for the inverse transform:

$$\frac{1}{2\pi i}\int_{\gamma-i\infty}^{\gamma+i\infty} Y(p)e^{pt}dp$$

I'm not too good at this. I know I need to learn more about Complex Analysis and contours specifically, but for the time being, can someone help me through the integral for alpha=1/2:

$$\frac{1}{2\pi i}\int_{\gamma-i\infty}^{\gamma+i\infty} \frac{e^{st}}{e^{2\sqrt{\pi s}}}ds=\frac{e^{-\pi/t}}{t^{3/2}}$$

Last edited: Nov 2, 2005
2. 3. Nov 2, 2005

### Tide Briefly: Make a change of variables using $s = \sigma^2$, complete the square in the exponential, scale out the t (which gives you the scaling of $t^{-3/2}$ immediately) and note that the change of variables requires your new contour to asymptote to the lines $e^{i \pi / 4}$ and $e^{-i \pi / 4}$. It should be apparent that you're dealing with the Fresnel integral functions with the usual results.

4. Nov 2, 2005

### saltydog Thanks Tide. I'll work with it. May take a while. 5. Nov 2, 2005

### saltydog Tide, it's not happening for me. Doing as you suggest, (completing the square), I get:

$$2e^{-\pi/t^2}\int e^{t(u-\frac{\sqrt{\pi}}{t})^2}udu$$

I don't see where to go from there and I really don't understand how the limits are changed as you suggest.

But I know where my problem lies and I know how to correct it: I don't have a sufficient understanding of Complex Analysis to even ask the question, and in order to approach this problem I would first have to work on other more simple problems. That's Ok. I have the time.

I did however go to the Library today and searched several tables of inverse transforms. It turns out that only for the case of alpha=1/2 is an inverse transform published. I suspect then even when I do master the contour integral, a closed-form solution for other values of alpha will NOT be attainable and thus I will have to resort to numerical methods of computing the inverse transforms. I know absolutely nothing about such numerical methods but I'm interested in learning about them and arriving at a complete solution of this problem in terms of alpha.

I also purchased a book on Complex Analysis and will take my own advice in here: Start on page 1. 6. Nov 3, 2005

### Tide salty,

I just went through the transformations and I get

$$\frac {2 e^{-\pi / t}}{t^{3/2}} \int e^{\left( \sigma - \sqrt {\pi / t}\right)^2} \sigma d\sigma$$

The original contour for the inverse Laplace transform effectively goes from $R e^{-i\pi/2}$ to $R e^{i\pi/2}$ with $R \rightarrow \infty$. Since $\sigma = \sqrt s$ the new limits will be from $\sqrt R e^{-i\pi/4}$ to $\sqrt R e^{i\pi/4}$. Looking at it more carefully shows one more change of variables is needed (a translation) to simplify the argument of the exponential.

When you're all one with that, you should have the integral I suggested in my previous post. For now you just have to go back and recheck your transformations.

Last edited: Nov 3, 2005
7. Nov 3, 2005

### Tide P.S. I didn't include the factor $\frac {1}{2\pi i}$ in my previous post.

8. Nov 3, 2005

### saltydog Tide, I cannot arrive at the answer you propose and when I check the integral with real values, your substitution does not lead to the numerical answer:

Let's just for now consider a real integral and let t=3:

$$\int_0^1 \frac{e^{st}}{e^{2\sqrt{\pi s}}}ds$$

when I numerically integrate this I obtain:

$$\int_0^1 \frac{e^{st}}{e^{2\sqrt{\pi s}}}ds\approx 0.428599$$

When I let:

$$\sigma^2=s$$

then the limits in $\sigma$ are the same.

Completing the square as you suggest, I arrive at:

$$2e^{-\pi/t}\int_0^1 e^{t(\sigma-\sqrt{\pi}/t)^2}\sigma d\sigma$$

numerically integrating that I obtain the same result as above. However, when I numerically integrate your solution:

$$\frac {2 e^{-\pi / t}}{t^{3/2}} \int_0^1 e^{\left( \sigma - \sqrt {\pi / t}\right)^2} \sigma d\sigma$$

I obtain 0.0830375. I suspect I'm making a dumb mistake and can handle you correcting me.

Last edited: Nov 3, 2005
9. Nov 3, 2005

### Tide salty,

There are two problems. First, I didn't do the transformations quite right so the $t^(3/2)$ should just be $t$:

$$\frac {2 e^{-\pi / t}}{t} \int e^{\left( \sigma - \sqrt {\pi / t}\right)^2} \sigma d\sigma$$

Also, if you want do a numerical comparison you have to set the limits properly. Remember, there was a scaling that eliminated the factor of t from the exponential which looks like $\sigma = \sqrt {t} \sigma'$ so if you set t = 3 then the upper limit in the transformed integral will be $\sqrt {3}$.

10. Nov 3, 2005

### saltydog Tide,

As you requested, I've achieved oneness with the integral. Here it is with your substitutions in case others are following this:

Above, after letting:

$$\sigma^2=s$$

and completing the square, we obtain:

$$2e^{-\pi/t}\int e^{t(\sigma-\sqrt{\pi}/t)^2}\sigma d\sigma$$

Now, in order to remove the t in the exponent, we let:

$$v=\sqrt{t}\sigma$$

so that:

$$dv=\sqrt{t}d\sigma,\quad \sigma=\frac{v}{\sqrt{t}},\quad d\sigma=\frac{dv}{\sqrt{t}}$$

Substituting this scalling factor into the exponent:

\begin{align*} t\left[\frac{v^2}{t}-\frac{2v\sqrt{\pi}}{t\sqrt{t}}+\frac{\pi}{t^2}\right]&= v^2-2v\sqrt{\pi/t}+\frac{\pi}{t} \\ &=(v-\sqrt{\pi/t})^2 \end{align}

substituting this into the integral:

$$2e^{-\pi/t}\int e^{(v-\sqrt{\pi/t})^2}\left(\frac{v}{\sqrt{t}}\right)\frac{dv}{\sqrt{t}}$$

Simplifying:

$$\frac{2e^{-\pi/t}}{t}\int e^{(v-\sqrt{\pi/t})^2}vdv;\quad v=\sqrt{t}\sigma;\quad \sigma^2=s$$

I'll start working with the complex analysis part . . .

Last edited: Nov 3, 2005
11. Nov 4, 2005

### saltydog Tide, anyone else too you know, can you please review the analysis below and let me know if I understand it so far. Thanks!

The original integral was:

$$\int_{\gamma-i\infty}^{\gamma+i\infty} \frac{e^{st}}{e^{2\sqrt{\pi s}}}ds$$

This can be re-written as:

$$\mathop\lim\limits_{y\to\infty}\int_{\gamma-yi}^{\gamma+yi} \frac{e^{st}}{e^{2\sqrt{\pi s}}}ds$$

Gamma is taken to be ANY positive constant such that the contour is to the right of any singularities of the integrand. I believe the integrand is entire so we can just take gamma to be 1. Converting to polar coordinates, the upper and lower limits become:

$$s_1=1-yi=Re^{-i\theta}$$

$$s_2=1+yi=Re^{i\theta}$$

Now, as y goes to infinity, theta approaches pi/2 and R goes to infinity. In the first transformation we had:

$$\sigma^2=s$$

Thus in terms of sigma, the limits become:

$$\sigma_1=\sqrt{Re^{-\pi i/2}}=\sqrt{R}e^{-\pi i/4}$$

$$\sigma_2=\sqrt{Re^{\pi i/2}}=\sqrt{R}e^{\pi i/4}$$

with R approaching infinity.

In the second transformation, we have:

$$v=\sqrt{t}\sigma$$

In terms of v, then the limits become:

$$v_1=\sqrt{Rt}e^{-\pi i/4}$$

$$v_2=\sqrt{Rt}e^{\pi i/4}$$

Thus arriving at the final form of the inverse transform:

$$y(t)=\frac{1}{2\pi i}\frac{2e^{-\pi/t}}{t}\mathop\lim\limits_{R\to\infty}\int_{\sqrt{Rt}e^{-\pi i/4}}^{\sqrt{Rt}e^{\pi i/4}} e^{(v-\sqrt{\pi/t})}^2vdv$$

I've already done some preliminary calculations of the integral and the results do indeed tend to the desired result:

$$y(t)=\frac{e^{-\pi/t}}{t^{3/2}}$$

However I cannot yet provide a proof and plan to work on this next. 12. Nov 4, 2005

### saltydog For the following:

$$y(t)=\frac{1}{2\pi i}\frac{2e^{-\pi/t}}{t}\mathop\lim\limits_{R\to\infty}\int_{\sqrt{RT}e^{-\pi i/4}}^{\sqrt{RT}e^{\pi i/4}} e^{(v-\sqrt{\pi/t})}^2vdv$$

First consider:

$$\int e^{(x-k)^2} xdx$$

Letting:

$$u=x,\quad dv=e^{(x-k)^2}dx$$

so that:

$$du=dx,\quad v=\frac{\sqrt{\pi}}{2}\text{Erfi}(x-k)$$

where:

$$\text{Erfi}(z)=\frac{\text{Erf}(iz)}{i}$$

and therefore:

$$\int e^{(x-k)}^2 xdx=x\frac{\sqrt{\pi}}{2}\text{Erfi}(x-k)-\frac{\sqrt{\pi}}{2}\int \text{Erfi}(x-k)dx$$

Now:

$$\int \text{Erfi}(x-k)dx=\frac{1}{\sqrt{\pi}}e^{(x-k)^2}-k\text{Erfi}(x-k)+x\text{Erfi}(x-k)$$

and thus:

$$\int e^{(x-k)^2}xdx=\frac{\sqrt{\pi}}{2}\left[\frac{e^{(x-k)^2}}{\sqrt{\pi}}+k\text{Erfi}(x-k)\right]$$

Therefore if we let:

$$k=\sqrt{\pi/t}$$

we obtain:

$$\int e^{(v-\sqrt{\pi/t})^2}vdv=\frac{\sqrt{\pi}}{2}\left[\frac{e^{(v-\sqrt{\pi/t})^2}}{\sqrt{\pi}}+\sqrt{\pi/t}\text{Erfi}(v-\sqrt{\pi/t})\right]$$

This then reduces the inverse transform to the following expression:

$$y(t)=\frac{e^{-\pi/t}}{2\pi i t}\left[\mathop\lim\limits_{R\to\infty}\left(e^{(v-\sqrt{\pi/t})^2}+\frac{\pi}{t^{1/2}}\text{Erfi}(v-\sqrt{\pi/t})\right)_{\sqrt{Rt}e^{-\pi i/4}}^{\sqrt{Rt}e^{\pi i/4}}\right]$$

Edit: Oh yea: Equal rights for special functions.

Last edited: Nov 4, 2005
13. Nov 5, 2005

### saltydog I wish to analyze:

$$y(t)=\frac{e^{-\pi/t}}{2\pi i t}\left[\mathop\lim\limits_{R\to\infty}\left(e^{(v-\sqrt{\pi/t})^2}+\frac{\pi}{t^{1/2}}\text{Erfi}(v-\sqrt{\pi/t})\right)_{\sqrt{Rt}e^{-\pi i/4}}^{\sqrt{Rt}e^{\pi i/4}}\right]$$

in two parts.

First consider the limit of the first term:

$$\mathop\lim\limits_{R\to\infty}\left[e^{(v-\sqrt{\pi/t})^2}\right]_{\sqrt{RT}e^{-\pi i/4}}^{\sqrt{Rt}e^{\pi i/4}}$$

Converting the upper limit to polar form:

\begin{align*} \sqrt{Rt}e^{\pi i/4}&=\sqrt{Rt}\left(Cos(\pi/4)+iSin(\pi/4)\right) \\ &=\sqrt{\frac{Rt}{2}}+i\sqrt{\frac{Rt}{2}} \end{align}

so that:

\begin{align*} (v-\sqrt{\pi/t})^2&=\left[\left(\sqrt{Rt/2}-\sqrt{\pi/2}\right)+i\sqrt{Rt/2}\right]^2 \\ &=\left(-2\sqrt{Rt/2}+\pi/t\right)+i\left(Rt-2\sqrt{R\pi/2}\right) \end{align}

The upper limit then becomes:

$$\mathop\lim\limits_{R\to\infty} Exp\left[\left(-2\sqrt{R\pi/2}+\pi/t\right)+i\left(Rt-2\sqrt{R\pi/2}\right)\right]$$

Note the real part of the expression is negative and as R goes to infinity, the coefficient on the Euler expansion goes to zero and thus the limit is zero.

The same analysis can be made for the lower limit. This limit is therefore zero. We can now write the inverse transform as:

$$y(t)=\frac{e^{-\pi/t}}{2\pi i t}\mathop\lim\limits_{R\to\infty}\left\{\frac{\pi}{t^{1/2}}\text{Erfi}\left(v-\sqrt{\pi/t}\right)_{\sqrt{Rt}e^{-\pi i/4}}^{\sqrt{Rt}e^{\pi i/4}}\right\}$$

or:

$$y(t)=\frac{e^{-\pi/t}}{2 i t^{3/2}}\mathop\lim\limits_{R\to\infty}\left\{\text{Erfi}\left(v-\sqrt{\pi/t}\right)_{\sqrt{Rt}e^{-\pi i/4}}^{\sqrt{Rt}e^{\pi i/4}}\right\}$$

Last edited: Nov 5, 2005
14. Nov 5, 2005

### saltydog I wish to finally evaluate:

$$y(t)=\frac{e^{-\pi/t}}{2 i t^{3/2}}\mathop\lim\limits_{R\to\infty}\left\{\text{Erfi}\left(v-\sqrt{\pi/t}\right)_{\sqrt{Rt}e^{-\pi i/4}}^{\sqrt{Rt}e^{\pi i/4}}\right\}\tag{1}$$

First note:

$$\text{Erfi}(z)=\frac{\text{Erf(iz)}}{i}=\frac{1}{i}\frac{2}{\sqrt{\pi}}\int_0^{iz} e^{-w^2}dw$$

Using this relation, I can re-write the above limit as:

$$\frac{2}{i\sqrt{\pi}}\left\{\mathop\lim\limits_{R\to\infty}\int_0^{i(\sqrt{Rt}e^{\pi i/4})} e^{-w^2}dw-\mathop\lim\limits_{R\to\infty}\int_0^{i(\sqrt{Rt}e^{-\pi i/4})} e^{-w^2}dw\right\}$$

Using Euler's relation, I convert the upper limits to real and imaginary parts and multiply by i:

\begin{align*} i\sqrt{Rt}e^{\pi i/4}&=i\left[\left(\sqrt{Rt/2}-\sqrt{\pi/2}\right)+i\sqrt{Rt/2}\right] \\ &=-\sqrt{Rt/2}+i\left(\sqrt{Rt/2}-\sqrt{\pi/t}\right) \end{align}

Now, as R goes to infinity, this quantity goes to $-\infty+i\infty$.

And the lower limit:

\begin{align*} i\sqrt{Rt}e^{-\pi i/4}&=i\left[\sqrt{Rt/2}+\left(\sqrt{Rt/2}-\sqrt{\pi/2}\right)\right] \\ &=\sqrt{Rt/2}+i\left(\sqrt{Rt/2}-\sqrt{\pi/2}\right) \end{align}

as R goes to infinity, this limit goes to $\infty+i\infty$.

Thus the integrals can be rewritten as:

$$\frac{1}{i}\frac{2}{\sqrt{\pi}}\int_0^{-\infty+i\infty} e^{-w^2}dw- \frac{1}{i}\frac{2}{\sqrt{\pi}}\int_0^{\infty+i\infty} e^{-w^2}dw$$

I choose to evaluate these line integrals along a parametric path using the following relation from Complex Analysis:

$$\int_C f(z)dz=\int_c udx-vdy+i\int_c vdx+udy$$

Thus the first step is to express the integrand in terms of a real and imaginary part. Relying on Euler's formula:

\begin{align*} e^{-w^2}&=e^{-(x+yi)^2} \\ &=e^{-(x^2-y^2+2xyi)} \\ &=e^{-(x^2-y^2)}Cos(2xy)+ie^{-(x^2-y^2)}(-Sin(2xy))\\ &=u+iv \end{align}

For the first integral, the path goes from the origin to the point $-\infty+i\infty$

Thus if I let:

$$x=-t,\quad dx=-dt$$

$$y=t,\quad dy=dt$$

Making this substitution, note that the exponent term is 1 and we have for the first integral (with t going from 0 to infinity):

\begin{align*} \int_0^{-\infty+i\infty} e^{-w^2}dw&=\left(-\int_0^\infty Cos(-2t^2)dt+\int_0^\infty Sin(-2t^2)dt\right) \\ &+i\left(\int_0^\infty Cos(-2t^2)dt+\int_0^\infty Sin(-2t^2)dt\right) \\ &=\left(-\frac{\sqrt{\pi}}{4}-\frac{\sqrt{\pi}}{4}\right)+i\left(\frac{\sqrt{\pi}}{4}-\frac{\sqrt{\pi}}{4}\right) \\ &=-\frac{\sqrt{\pi}}{2} \end{align}

The same analysis can be made (with the path going from 0 to $\infty+i\infty$ for the second integral with the results that the line integral in that case is:

$$\frac{\sqrt{\pi}}{2}$$

Substituting these values into (1), I obtain:

$$\frac{1}{i}\frac{2}{\sqrt{\pi}}\int_0^{-\infty+i\infty} e^{-w^2}dw=\frac{1}{i}\frac{2}{\sqrt{\pi}}(-\frac{\sqrt{\pi}}{2})=-\frac{1}{i}\frac{i}{i}=i$$

and:

$$\frac{1}{i}\frac{2}{\sqrt{\pi}}\int_0^{\infty+i\infty} e^{-w^2}dw=\frac{1}{i}\frac{2}{\sqrt{\pi}}(\frac{\sqrt{\pi}}{2})=\frac{1}{i}\frac{i}{i}=-i$$

So that the limit is $2i$.

Substituting this into the last expression of the inverse transform, we obtain:

$$y(x)=\frac{e^{-\pi/t}}{2i t^{3/2}}(2i)=\frac{e^{-\pi/t}}{t^{3/2}}$$

Thanks for helping me Tide. Couldn't have done this without your help. 15. Nov 5, 2005

### Tide Way to go, Salty!

16. Nov 6, 2005

### saltydog So I wish to return to the subject of differential (integral) equations and the numerical evaluation of inverse Laplace Transforms. While going though the analysis above, I realized I can numerically evaluate an integral of a complex variable in the same (similar) way as for a real variable! Take for example the contour expression for the inverse transform:

$$y(x)=\frac{1}{2\pi i}\int_{1-i\infty}^{1+i\infty}\frac{e^{st}}{e^{2\sqrt{\pi s}}} ds$$

For the moment, let's consider the case of t=1. In that case, this integral is EXACTLY equal to:

$$y(1)=e^{-\pi}\approx 0.0432139$$

Now, for starters, suppose I just calculate a ROUGH Reimann sum for this integral:

$$\frac{1}{2\pi i}\sum_1^{1000} \frac{e^{st}}{e^{2\sqrt{\pi s}}}\Delta c$$

along the straight-line contour from the points $1-10i$ to $1+10i$.

This is my Mathematica code:

Code (Text):
sum = 0;
tval = 1;
imax = 10;
cval = 1 - i max;
nmax = 1000;
deltai = 2 imax/nmax i;
For[n = 1, n <= nmax, n++,
cval += deltai;
$$sum+=N\left[\left(\frac{e^{\text{cval tval}}}{e^{2\sqrt{\pi\text{cval}}}}\right) \text{deltai}\left]$$
];
$$sum=\frac{1}{2\pi i} sum$$
The results are:

$$0.0434844+1.31273X10^{-6} i$$

Granted there is a slight complex component but the real part agrees with the actual value to 3 decimal places and this is only a rough estimate.

My suspicion is that this will be the case ONLY for alpha=1/2 and that any other value of alpha in the original integral equation will result in a significant complex part meaning that the only REAL solution to the integral equation is the one with alpha=1/2.

Tide, you still with me?

Last edited: Nov 6, 2005
17. Nov 6, 2005

### Tide If real solutions exist to the original integral equation and if the Laplace transforms exist then the inversion must lead to a real result. A numerical approximation may lead to an imaginary part but there are two issues with a numerical approximation: (a) finite step size and (b) infinite domain. It may take that entire "trip" to infinity to eliminate the small imaginary part but it will not be there in the end! :)

18. Nov 6, 2005

### saltydog Numerical calculations

Hey Tide, I find the following results very encouraging: The first plot is the solution of the integral equation:

$$y(t)=\frac{e^{-\pi/t}}{t^{3/2}}$$

The second plot is a numerical evaluation of the inverse transform using a simple Riemann sum over the path $1-20i$ to $1+20i[/tex] sampling 1000 points. The third plot is a superposition of the two. Not bad for a rough approximation. Edit: I neglected the "small" residual imaginary component of the numerical calculations. #### Attached Files: • ###### sol1.JPG File size: 4.5 KB Views: 146 • ###### sol2.JPG File size: 4.7 KB Views: 166 • ###### sol3.JPG File size: 4.7 KB Views: 153 Last edited: Nov 6, 2005 19. Nov 6, 2005 ### Tide That's pretty good! If you want to take it a little further, you could do some asymptotic analysis on the integral. When t is large, most of the contribution to the integral (say your first eq. in #11) is from the region about a "saddle point" (e.g. [itex]\int e^{\phi(s, t)} ds$) has such a point where $d \phi(s, t) /ds = 0$ and you integrate along the path of steepest descent through the saddle point.

In principle, the asymptotic expansion of the integral will overlap your numerical evaluation over some range of large t values. I haven't looked but you should be able to find a lot of references to the "method of steepest descent" for details.

Last edited: Nov 6, 2005
20. Nov 7, 2005

### saltydog Thanks Tide. I'll check that out. Here are the plots for alpha =0.5 to alpha=0.1 with the largest one being alpha=0.1, then 0.2, 0.3, 0.4, and the smallest, 0.5.

I've yet to show that these other plots satisfy the integral equation however and plan to work on a back-substitution method to check them (within some level of accuracy).

I'm sure other numerical methods exists for evaluating integral equations, however, to pursue such would have denied me the priviledge (and learning experience) of working on this very interesting approach with Laplace Transforms. #### Attached Files:

• ###### comboplot.JPG
File size:
8.6 KB
Views:
144
21. Nov 8, 2005

### saltydog I perfected a method of back-substitution of the numeric data into the integral equation in an effort to numerically check the solution for alpha=0.1. The inverse transform is:

$$y(t)=\frac{1}{2\pi i}\mathop\lim\limits_{a\to\infty}\int_{1-ai}^{1+ai} e^{zt}\text{Exp}\left[-\frac{\Gamma[0.1]}{1-0.1}z^{1-0.1}\right]dz$$

This can be numerically integrated directly until the results reach a desired level of accuracy. The integral then is evaluated for a range of t. The results of these calculations is shown in plot 1. Mathematica has a nice feature which can represents the data points as a "Interpolation function":

$$\text{alpha01 = Interpolation[pointlist]}$$

alpha01 is now treated as a function.

The next step is to numerically evaluate the RHS of the original integral equation:

$$\int_0^t \tau^{\alpha-1}\text{alpha01}(t-\tau)d\tau$$

or:

$$rhs[t]=\int_0^t (t-z)^{\alpha-1}\text{alpha01}(z)dz$$

The LHS of the integral equation is:

$$lhs[t]=\text{alpha01[t]t}$$

Subtracting:

$$\text{dif[t]=lhs[t]-rhs[t]}$$

The results of this difference is shown in the second plot showing that the difference is no larger than 10^-6.

This then gives me some confidence that the numerical results for alpha=0.1 do indeed satisfy the integral equation.

#### Attached Files:

File size:
5.1 KB
Views:
145
• ###### alpha01dif.JPG
File size:
5.1 KB
Views:
142
Last edited: Nov 8, 2005