Archive for March, 2008


Question: Laplace’s Equation confined to a circle

March 26th, 2008

Question: Solve Laplace’s equation inside a circle of radius a.

\nabla ^2 u = \frac{1}{r} \frac{\partial}{\partial r}(r \frac{\partial u}{\partial r}) + \frac{1}{r^2}\frac{\partial^2 u}{\partial \theta^2} = 0

Answer: Heat Distribution in a Rod

March 26th, 2008

Question: A 1D rod of length L has an initial heat distribution of

u(x,0) = - cos(\frac{8 \pi x}{L})

If the rod has insulated ends (\frac{\partial u}{\partial x}(0,t) = \frac{\partial u}{\partial x}(L,t) = 0) and obeys the heat equation

\frac{\partial u(x,t)}{\partial t} = k \frac{\partial^2 u(x,t)}{\partial x^2},

What is the heat distribution of the rod as a function of time?

Answer:
First assume that the solution to the PDE

\frac{\partial u(x,t)}{\partial t} = k \frac{\partial^2 u(x,t)}{\partial x^2},

has the form

u(x,t) = \phi(x) G (t)

Therefore, we can reduce the equation to the following:

\frac{\partial u(x,t)}{\partial t} = k \frac{\partial^2 u(x,t)}{\partial x^2}
\phi(x) \frac{\partial G(t)}{\partial t} = k G(t) \frac{\partial^2 \phi(x)}{\partial x^2}
\frac{1}{k G(t)} \frac{\partial G(t)}{\partial t} =\frac{1}{\phi(x)} \frac{\partial^2 \phi(x)}{\partial x^2}

Since this equation is true for all x and t, therefore both sides of the equation must equal a constant.

\frac{1}{k G(t)} \frac{\partial G(t)}{\partial t} =\frac{1}{\phi(x)} \frac{\partial^2 \phi(x)}{\partial x^2} = - \lambda

which can be written as two ODEs

\frac{d G(t)}{dt} = - \lambda k G(t)
\frac{d^2 \phi(x)}{dx^2} = - \lambda \phi(x)

Solving the first-order differential \frac{d G(t)}{dt} = - \lambda k G(t) is trivial and can easily be shown to have the solution:

G(t) = C e^{-\lambda k t}

Next, we need to solve \frac{d^2 \phi(x)}{dx^2} = - \lambda \phi(x)

The second-order differential has 3 different cases (\lambda = 0, \lambda > 0, \lambda < 0).

Case \lambda = 0
The ODE for this case would be

\frac{d^2 \phi(x)}{dx^2} =0

which has the solution

\phi(x) = A_0 x + A_1

If we apply the BC \frac{d \phi}{dx}(0)=0, we get

\frac{d \phi}{dx}(0) = A_0 = 0

which implies

\phi(x) = A_1

Case \lambda > 0
The ODE for this case would be

\frac{d^2 \phi(x)}{dx^2} = - \lambda \phi(x)

which has the solution

\phi(x) = C_1 cos(\sqrt{\lambda}x) + C_2 sin(\sqrt{\lambda}x)

Next, we apply the BC \frac{d \phi}{dx}(0)=0.

\frac{d \phi}{dx}(0) = \sqrt{\lambda}(-C_1 sin(\sqrt{\lambda}0) + C_2 cos(\sqrt{\lambda}0))
\frac{d \phi}{dx}(0) = \sqrt{\lambda}C_2

Since \lambda > 0, this implies C_2 = 0. Therefore,

\phi(x) = C_1 cos(\sqrt{\lambda}x)

If we apply the 2nd BC \frac{d \phi}{dx}(L)=0, we find the following

\frac{d \phi}{dx}(L) = -C_1 \sqrt{\lambda} sin(\sqrt{\lambda}L) = 0

Therefore, for non-trivial solutions \sqrt{\lambda}L = n \pi where n=1,2,3 \ldots.

This means,

\phi(x) = C_1 cos(\frac{n \pi x}{L})

is a solution.

Case \lambda < 0
The ODE for this case would be

\frac{d^2 \phi(x)}{dx^2} = s \phi(x) where s = -\lambda and s > 0

which has the solution

\phi(x) = B_1 e^{\sqrt{s}x} + B_2 e^{-\sqrt{s}x}

If we apply the BC \frac{d \phi}{dx}(0)=0,we find

\frac{d\phi}{dx}(0) = B_1 \sqrt{s} e^{\sqrt{s}0} - B_2 \sqrt{s} e^{-\sqrt{s}0} = 0
\sqrt{s}(B_1 - B_2) = 0

Since s > 0, we know that B_1 = B_2 which means

\phi(x) = B_1 (e^{\sqrt{s}x} + e^{-\sqrt{s}x})

is a solution.

If we apply the 2nd BC \frac{d \phi}{dx}(L)=0, we find

\frac{d\phi}{dx}(L) = B_1 \sqrt{s} (e^{\sqrt{s}L} - e^{-\sqrt{s}L}) = 0

Since s > 0 and e^{\sqrt{s}L} \neq e^{-\sqrt{s}L}, we know that B_1 =0.

This means that there is no solution where \lambda < 0.


Since the PDE will satify any linear combination of the above solutions, we find that

u(x,t) = A_0 + \sum_{n=1}^{\infty} A_n cos(\frac{n \pi x}{L}) e^{-(\frac{n\pi}{L})^2 kt}

If we apply the IC u(x,0) = - cos(\frac{8 \pi x}{L}) to the solution, we find that A_8 = -3 and every other A_n is zero.

Therefore, the solution would be

u(x,t) =-3 cos(\frac{8 \pi x}{L}) e^{-(\frac{8\pi}{L})^2 kt}

Question: Heat Distribution in a Rod

March 17th, 2008

Question: A 1D rod of length L has an initial heat distribution of

u(x,0) = - cos(\frac{8 \pi x}{L})

If the rod has insulated ends (\frac{\partial u}{\partial x}(0,t) = \frac{\partial u}{\partial x}(L,t) = 0) and obeys the heat equation

\frac{\partial u(x,t)}{\partial t} = k \frac{\partial^2 u(x,t)}{\partial x^2},

What is the heat distribution of the rod as a function of time?

Answer: Solve Laplace’s Equation

March 17th, 2008

Question: Find the solutions to Laplace’s Equation:

\frac{\partial^2 u(x,y)}{\partial x^2} + \frac{\partial^2 u(x,y)}{\partial y^2} = 0

Answer:
First assume that the solution to the PDE

\frac{\partial^2 u(x,y)}{\partial x^2} + \frac{\partial^2 u(x,y)}{\partial y^2} = 0

has the form

u(x,y) = \phi(x) \xi (y)

Therefore, we can reduce the equation to the following:

\xi (y) \frac{\partial^2 \phi(x)}{\partial x^2} + \phi(x) \frac{\partial^2 \xi (y)}{\partial y^2} = 0
\frac{1}{\phi (x)} \frac{\partial^2 \phi(x)}{\partial x^2} = \frac{-1}{\xi (y)}\frac{\partial^2 \xi (y)}{\partial y^2}

Since this equation is true for all x and y, therefore both sides of the equation must equal a constant.

\frac{1}{\phi (x)} \frac{\partial^2 \phi(x)}{\partial x^2} = \frac{-1}{\xi (y)}\frac{\partial^2 \xi (y)}{\partial y^2} = \lambda

This implies that we need to solve two ODEs.

 \frac{\partial^2 \phi(x)}{\partial x^2} =\lambda\phi (x)~~~~~~\frac{\partial^2 \xi (y)}{\partial y^2} = -\lambda\xi (y)

Since the solutions to the two ODEs will be very similar, I will solve the ODE

 \frac{\partial^2 \Psi}{\partial r^2} =\kappa \Psi(r)

and apply the results to the two ODEs.

There are 3 cases which we need to solve (\kappa > 0, \kappa = 0, and \kappa < 0).

Case \kappa = 0
The ODE for this case would be

 \frac{\partial^2 \Psi}{\partial r^2} =0

which has the solution

\Psi = C_1 r + C_2

Case \kappa < 0 and \kappa > 0
The ODE

 \frac{\partial^2 \Psi}{\partial r^2} =\kappa \Psi(r)

which will have the solution

\Psi = B_1 e^{\sqrt{\kappa}r} + B_2 e^{-\sqrt{\kappa}r}

Please note that \sqrt{\kappa} will be an imaginary number when \kappa < 0.


Therefore, if we apply the above solution, we can find the functions that solve the ODEs

 \frac{\partial^2 \phi(x)}{\partial x^2} =\lambda\phi (x)~~~~~~\frac{\partial^2 \xi (y)}{\partial y^2} = -\lambda\xi (y)

which would be

\phi (x)= A_1 e^{i\sqrt{\lambda}x} + A_2 e^{-i\sqrt{\lambda}x}
\xi  (y)= B_1 e^{\sqrt{\lambda}y} + B_2 e^{-\sqrt{\lambda}y} when \lambda < 0

\phi (x)= A_1 e^{\sqrt{\lambda}x} + A_2 e^{-\sqrt{\lambda}x}
\xi  (y)= B_1 e^{i\sqrt{\lambda}y} + B_2 e^{-i\sqrt{\lambda}y} when \lambda > 0

\phi (x)= C_1 x + C_2
\xi  (y)= D_1 y + D_2 when \lambda = 0

Therefore, the solution would have the form

u(x,y) = \begin{cases}(A_1 e^{i\sqrt{\lambda}x} + A_2 e^{-i\sqrt{\lambda}x})(B_1 e^{\sqrt{\lambda}y} + B_2 e^{-\sqrt{\lambda}y}) & for~\lambda < 0 \\ (A_1 e^{\sqrt{\lambda}x} + A_2 e^{-\sqrt{\lambda}x})(B_1 e^{i\sqrt{\lambda}y} + B_2 e^{-i\sqrt{\lambda}y}) & for~\lambda > 0 \\ (C_1 x + C_2)(D_1 y + D_2) & for~\lambda = 0\end{cases}

Any superposition of the above equation will satisfy Laplace’s equation. In order to reduce this solution more, we would need to be given Boundary and Initial Conditions.

Physics Simulations

March 11th, 2008

Recently, I started studying ODEs and PDEs over again. I have always found this subject interesting, but I haven’t done much with it since University (besides the odd Classical Mechanics ODE). You should expect a lot more questions about Differential Equations in the near future.

Besides solving Differential Equations analytically, I have been looking at numerical approaches to solving these problems. If you use Maple or Matlab, it is fairly straightforward to solve Differential Equations numerically. However, these approaches are not really optimal for larger problems (especially problems that require parallelization). In order to solve these larger problems, it usually requires using a “real” language.

Initially, I started looking for numerical libraries for Perl. The only library that looked decent was Cephes which seems a bit incomplete. The Python library called NumPy seems a lot better, but this would require me to learn Python.

Anyways, while searching the web, I came across a cool little site called “My Physics Lab”. It contains a bunch applets that have physics simulations. The coolest one, I think, is the double pendulum. It is a simple example of chaotic motion.

Also, if you know of any good numerical libraries, please leave a comment (unless it is in C/C++ or Fortran).