Help solving this heat differential equation -u'' = (x-1.5) ; delta

75 Views Asked by At

1.4 Numerical Methods

Solve the boundary value problem (BVP) by modeling the temperature in a well conducting metal rod of length 3 in the presence of a point heat source in the middle of the rod in the stationary state if both ends of the rod are in contact with ice, which keeps the temperature of the two ends at zero. Choose natural units, so that the heat equation is

$$ −u′′ = δ(x − 1.5) $$

After solving the problem analytically also solve it numerically with step size h = 1. Calculate the error at the mesh points x = 1 and x = 2.


Where I'm at:

So the goal is to find.

$$ u(x) = f(x) $$

I'm using:

$$ u'(x = 0) = u(x = 0) = 0 $$

I apply the Laplace transform:

$$ L\{-u''\} = L\{δ(x − 1.5)\} $$

$$ = -s^2U(s) = e^{-1.5s} $$

$$ = U(s) = \frac{e^{-1.5s}}{-s^2} $$

finally Laplace inverse:

$$ L^{-1}\{U(s)\} = u(x) = $$

$$ L^{-1}\{ \frac{e^{-1.5s}}{-s^2} \} = ? $$

Which I'm not sure how to evaluate. I understand $$ e^{-1.5s} $$ is a time-shift. But I'm lost at how to evaluate this.

Any help is SUPER appreciated because I've been working on this problem for weeks!

Please include all steps and/or any reference material. (including how you are applying which laplace properties).

I'm open to solving this problem anyway, it's just this is the way I understood it.

Below I attached an image of a solution I found, but I didn't understand it fully.

Solution I found

1

There are 1 best solutions below

0
On

Let's write $F(s)=-\exp(-a t)$ (with $a=1.5$) and $G(s)=s^{-2}$. We have $f(t)=\mathcal{L}^{-1}\{F(s)\}=-\delta(t-a)$ and $g(t)=\mathcal{L}^{-1}\{G(s)\}=t \cdot \Theta(t)$ (with $\Theta(t)$ being the Heaviside step function) see wikipedia "Laplace transform" where these standard inversions are tabulated.

An identity (again, check wikipedia) for the inverse Laplace transformation of the product of two functions says \begin{equation} \mathcal{L}^{-1}\{F(s)G(s)\}=\int_0^t f(\tau)g(t-\tau)\,\rm{d}\tau\,. \end{equation} In your case, this means \begin{align} \mathcal{L}^{-1}\{-\exp(-a s)/s^2\}&=-\int_0^x \delta(x'-a)\cdot(x-x') \cdot \Theta (x-x')\,{\rm d}x'\\ &=-(x-a) \cdot \Theta(x-a) \end{align} The figure below compares a numerical Laplace inversion of $-\exp(-a s)/s^2$ to the above analytical result:

enter image description here

I created it with the code below.

import numpy as np
import matplotlib.pyplot as plt
from mpmath import *

def f(x,a):
    if x<a:
        f=0
    else:
        f=-(x-a)
    return f

def F(s,a):
    return -mp.exp(-a*s)/s**2

mp.dps = 15
def fnum(x,a):
    fx = lambda s: F(s,a)
    return invertlaplace(fx,x,method='talbot')

a=1.5
x  = np.linspace(0,4,100)
x2 = np.linspace(a,4,10)

fig, ax = plt.subplots()
ax.plot(x,[f(i,a) for i in x],  'r',lw=2,label='analytical')
ax.plot(x2,[fnum(i,a) for i in x2],'r',lw=0,marker='o',label='numerical')
ax.set(xlabel=r'$x$', ylabel=r'$u(x)$') 
ax.legend()