View source for Partial Differential Equations/The Laplacian and Laplace's Equation

==The Laplacian and Laplace's Equation==

By now, you've most likely grown sick of the one dimensional transient diffusion PDE we've been playing with:


:<math>\frac{\partial u}{\partial t} = \alpha \frac{\partial^2 u}{\partial x^2}\,</math>


Make no mistake: we're not nearly done with this stupid thing; but for the sake of variety let's introduce a fresh new equation and, even though it's not strictly a separation of variables concept, a really cool quantity called the Laplacian. You'll like this chapter; it has many pretty pictures in it.

[[Image:PDEBook0012.jpg|frame|Graph of <math>u(x) = x^3 - x</math>.]]

===The Laplacian===

The Laplacian is a linear operator in Euclidean n-space.

We'll start with the 3D Cartesian "version". Let <math>u = u(x, y, z, t)</math>. The '''Laplacian''' of <math>u</math> is defined and notated as:


:<math>\nabla^2 u = \Delta u = \operatorname{div} \ \operatorname{grad} \ u = \frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2} + \frac{\partial^2 u}{\partial z^2}\,</math>


It's the sum of the nonmixed second derivatives with respect to the ''Cartesian'' space variables. The "del squared" notation is preferred since the capital delta can be confused with increments and differences, and <math>\mbox{div grad }u</math> is too long and doesn't involve pretty math symbols. The Laplacian is also known as the Laplace operator or Laplace's operator, not to be confused with the Laplace transform.

In one dimension, recall that the second derivative measures concavity. Suppose <math>y = f(x)</math>; <B>if <math>f''(x)</math> is positive, <math>y</math> is concave up</B>, and <B>if <math>f''(x)</math> is negative, <math>y</math> is concave down</B>. The Laplacian may be thought of as a generalization of the concavity concept to multivariate functions.

This idea is demonstrated at the right, in one dimension: <math>u(x) = x^3 - x</math>. To the left of <math>x = 0</math>, the Laplacian (simply the second derivative here) is negative, and the graph is concave down. At <math>x = 0</math>, the curve inflects and the Laplacian is <math>0</math>. To the right of <math>x = 0</math>, the Laplacian is positive and the graph is concave up.



Concavity may or may not do it for you. Thankfully, there's another very important view of the Laplacian, with deep implications for any equation it shows itself in: the Laplacian compares the value of <math>u</math> at some point in space to the average of the values of <math>u</math> in the neighborhood of the same point. The three cases are:

:*If <math>u</math> is greater at some point then the average of its neighbors, <math>\scriptstyle \nabla^2 u < 0</math>.
:*If <math>u</math> is at some point equal to the average of its neighbors, <math>\scriptstyle \nabla^2 u = 0</math>.
:*If <math>u</math> is smaller at some point then the average of its neighbors, <math>\scriptstyle \nabla^2 u > 0</math>.

So the laplacian may be thought of as, at some point <math>(x_0, y_0, z_0)</math>:


:<math>\nabla^2 u \Big|_{x = x_o, y = y_o, z = z_o} \approx \mbox{average of u in the neighborhood of } (x_0, y_0, z_0) \ - \ u(x_0, y_0, z_0, t) \,</math>


[[Image:PDEBook0013.png|frame|The neighborhood of <math>(x_0, y_0, z_0)</math>.]]

The neighborhood of some point is defined as the open set that lies within some Euclidean distance δ from the point. Referring to the picture at right (a 3D example), the neighborhood of the point <math>(x_0, y_0, z_0)</math> is the shaded region which satisfies:


:<math>(x - x_0)^2 + (y - y_0)^2 + (z - z_0)^2 < \delta^2 \,</math>


Note that our one dimensional transient diffusion equation involves the Laplacian:


:<math>\frac{\partial u}{\partial t} = \alpha \frac{\partial^2 u}{\partial x^2} = \alpha \nabla^2 u\,</math>


With this mentality, let's examine the behavior of this very important PDE. On the left is the time derivative and on the right is the Laplacian. This equation is saying that:

:'''The rate of change of <math>u</math> at some point is proportional to the difference between the average value of <math>u</math> around that point and the value of <math>u</math> at that point.'''

For example, if there's at some position a "hot spot" where <math>u</math> is on average greater then its neighbors, the Laplacian will be negative and thus the time derivative will be negative, this will cause <math>u</math> to decrease at that position, "cooling" it down. This is illustrated below. The arrows reflect upon the magnitude of the Laplacian and, by grace of the time derivative, the direction the curve will move.

[[Image:PDEBook0014.png|frame|center|Visualization of transient diffusion.]]

It's worth noting that in 3D, this equation fully describes the flow of heat in a homogeneous solid that's not generating it's own heat (like too much electricity through a narrow wire would).

===Laplace's Equation===

Laplace's equation describes a steady state condition, and this is what it looks like:


:<math>\nabla^2 u = 0\,</math>


Solutions of this equation are called ''harmonic functions''. Some things to note:

:*Time is absent. This equation describes a steady state condition.
:*The absence of time implies the absence of an IC, so we'll be dealing with BVPs rather then IBVPs.
:*In one dimension, this is the ODE of a straight line passing through the boundaries at their specified values.
:*All functions that satisfy this equation in some domain are analytic (informally, an analytic function is equal to its Taylor expansion) in that domain.
:*Despite appearances, solutions of Laplace's equation are generally not minimal surfaces.
:*Laplace's equation is linear.

Laplace's equation is separable in the Cartesian (and almost any other) coordinate system. So, we shouldn't have too much problem solving it if the BCs involved aren't too convoluted.

===Laplace's Equation on a Square: Cartesian Coordinates===

[[Image:PDEBook0015.png|frame|Steady state conditions on a square.]]

Imagine a 1 x 1 square plate that's insulated top and bottom and has constant temperatures applied at its uninsulated edges, visualized to the right. Heat is flowing in and out of this thing steadily through the edges only, and since it's "thin" and "insulated", the temperature may be given as <math>u(x, y)</math>. This is the first time we venture into two spatial coordinates, note the absence of time.

Let's make up a BVP, referring to the picture:


:<math>\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2} = 0\,</math>


:<math>u(0, y) = 0 \qquad \mbox{(Edge D)}\,</math>


:<math>u(1, y) = 0 \qquad \mbox{(Edge B)}\,</math>


:<math>u(x, 0) = 0 \qquad \mbox{(Edge A)}\,</math>


:<math>u(x, 1) = 1 \qquad \mbox{(Edge C)}\,</math>


So we have one nonhomogeneous BC. Assume that <math>u(x, y) = X(x) Y(y)</math>:


:<math>\frac{\partial^2}{\partial x^2}(X(x) Y(y)) + \frac{\partial^2}{\partial y^2}(X(x) Y(y)) = 0\,</math>


:<math>X''(x)Y(y) + X(x)Y''(y) = 0\,</math>


:<math>-\frac{X''(x)}{X(x)} = \frac{Y''(y)}{Y(y)} = k^2\,</math>


::<math>\Big\Downarrow</math>


:<math>X''(x) = -k^2 X(x)\,</math>


:<math>Y''(y) = k^2 Y(y)\,</math>


As with before, calling the separation constant <math>k^2</math> in favor of just <math>k</math> (or something) happens to make the problem easier to solve. Note that the negative sign was kept for the <math>X(x)</math> equation: again, these choices happen to make things simpler. Solving each equation and combining them back into <math>u(x, y)</math>:


:<math>X(x) = C_1 \cos(k x) + C_2 \sin(k x)\,</math>


:<math>Y(y) = C_3 e^{k y} + C_4 e^{-k y}\,</math>


::<math>\Big\Downarrow</math>


:<math>u(x, y) = X(x) Y(y)\,</math>


:<math>u(x, y) = \left(C_1 \cos(k x) + C_2 \sin(k x)\right)\left(C_3 e^{k y} + C_4 e^{-k y}\right)\,</math>


At edge D:


:<math>u(0, y) = 0\,</math>


:<math>\left(C_1 \cos(k \cdot 0) + C_2 \sin(k \cdot 0)\right)\left(C_3 e^{k y} + C_4 e^{-k y}\right) = 0\,</math>


:<math>C_1 \left(C_3 e^{k y} + C_4 e^{-k y}\right) = 0 \Rightarrow C_1 = 0\,</math>


::<math>\Big\Downarrow</math>


:<math>u(x, y) = C_2 \sin(k x)\left(C_3 e^{k y} + C_4 e^{-k y}\right)\,</math>


Note that the constants can be merged, but we won't do it so that a point can be made in a moment. At edge A:


:<math>u(x, 0) = 0\,</math>


:<math>C_2 \sin(k x)\left(C_3 e^{k \cdot 0} + C_4 e^{-k \cdot 0}\right) = 0\,</math>


Taking <math>C_2</math> as <math>0</math> would satisfy this particular BC, however this would yield a plane solution of <math>u(x, y) = 0</math>, which can't satisfy the temperature at edge C. This is why the constants weren't merged a few steps ago, to make it obvious that <math>C_2</math> may not be <math>0</math>. So, we instead take <math>C_3 = -C_4</math> to satisfy the above, and then combine the three constants into one, call it <math>B</math>:


:<math>u(x, y) = C_2 \sin(k x)\left(C_3 e^{k y} + C_4 e^{-k y}\right) \qquad ; \ C_3 = -C_4\,</math>


::<math>\Big\Downarrow</math>


:<math>u(x, y) = B \sin(k x)\left(e^{k y} - e^{-k y}\right)\,</math>


Now look at edge B:


:<math>u(1, y) = 0\,</math>


:<math>B \sin(k)\left(e^{k y} - e^{-k y}\right) = 0\,</math>


It should go without saying by now that <math>B</math> can't be zero, since this would yield <math>u(x, y) = 0</math> which couldn't satisfy the nonzero BC. Instead, we can take <math>k = n \pi</math>:


:<math>u(x, y) = B \sin(n \pi x)\left(e^{n \pi y} - e^{-n \pi y}\right)\,</math>


As of now, this solution will satisfy 3 of the 4 BCs. All that is left is edge C, the nonhomogeneous BC. 


:<math>u(x, 1) = 1\,</math>


:<math>B \sin(n \pi x)\left(e^{n \pi y} - e^{-n \pi y}\right) = 1\,</math>


Neither <math>B</math> nor <math>n</math> can be contorted to fit this BC.

Since Laplace's equation is linear, a linear combination of solutions to the PDE is also a solution to the PDE. Another thing to note: since the BCs (so far) are homogeneous, we can add the solutions without worrying about nonzero boundaries adding up.

Though <math>u(x, y)</math> as shown above will not solve this problem, we can try summing (based on <math>n</math>) solutions to form a linear combination which might solve the BVP as a whole:


:<math>u(x, y) = \sum_{n=1}^{\infty} u_n(x, y)\,</math>


:<math>u_n(x, y) = B_n \sin(n \pi x)\left(e^{n \pi y} - e^{-n \pi y}\right)\,</math>


::<math>\Big\Downarrow</math>


:<math>u(x, y) = \sum_{n=1}^{\infty} B_n \sin(n \pi x)\left(e^{n \pi y} - e^{-n \pi y}\right)\,</math>


Assuming this form is correct (review [[Partial Differential Equations/Parallel Plate Flow: Realistic IC| Parallel Plate Flow: Realistic IC]] for motivation), let's again try applying the last BC:


:<math>u(x, 1) = 1\,</math>


:<math>\sum_{n=1}^{\infty} B_n \sin(n \pi x)\left(e^{n \pi \cdot 1} - e^{-n \pi \cdot 1}\right) = 1\,</math>


It looks like it needs Fourier series methodology. Finding <math>B_n</math> via orthogonality should solve this problem:


:<math>2 \sin(m \pi x) \cdot \sum_{n=1}^{\infty} B_n \sin(n \pi x)\left(e^{n \pi} - e^{-n \pi}\right) =
2 \sin(m \pi x) \cdot 1\,</math>


:<math>\sum_{n=1}^{\infty} B_n \cdot 2 \sin(m \pi x) \sin(n \pi x)\left(e^{n \pi} - e^{-n \pi}\right) = 2 \sin(m \pi x) \,</math>


:<math>\int_{0}^{1} \sum_{n=1}^{\infty} B_n \cdot 2 \sin(m \pi x) \sin(n \pi x)
\left(e^{n \pi} - e^{-n \pi}\right) dx = \int_{0}^{1} 2 \sin(m \pi x) dx\,</math>


:<math>\sum_{n=1}^{\infty} B_n \int_{0}^{1} 2 \sin(m \pi x) \sin(n \pi x) dx
\left(e^{n \pi} - e^{-n \pi}\right) dx = 2 \frac{1 - \cos(m \pi)}{m \pi}\,</math>


[[Image:PDEBook0016.png|frame|25 term partial sum of the series solution.]]
:<math>B_m \left(e^{m \pi} - e^{-m \pi}\right) = 2 \frac{1 - \cos(m \pi)}{m \pi}\,</math>


:<math>B_n = 2 \frac{1 - (-1)^n}{n \pi \left(e^{n \pi} - e^{-n \pi}\right)}\,</math>


<math>m</math> was changed to <math>n</math> in the last step. Also, for integer <math>n</math>, <math>\cos(n \pi) = (-1)^n</math>. Note that a Fourier sine expansion has been done. The solution to the BVP can finally be assembled:


:<math>u(x, y) = \sum_{n=1}^{\infty} B_n \sin(n \pi x)\left(e^{n \pi y} - e^{-n \pi y}\right)\,</math>


:<math>u(x, y) = \sum_{n=1}^{\infty} 2 \frac{1 - (-1)^n}{n \pi \left(e^{n \pi} - e^{-n \pi}\right)} \sin(n \pi x)\left(e^{n \pi y} - e^{-n \pi y}\right)\,</math>


That solves it!

It's finally time to mention that the BCs are discontinuous at the points <math>(0, 1)</math> and <math>(1, 1)</math>. As a result, the series should converge slowly at those points. This is clear from the plot at right: it's a 25 term partial sum (note that half of the terms are <math>0</math>), and it looks perfect except at <math>y = 1</math>, especially near the discontinuities at <math>x = 0</math> and <math>x = 1</math>.

===Laplace's Equation on a Circle: Polar Coordinates===

Now, we'll specify the value of <math>u</math> on a circular boundary. A circle can be represented in Cartesian coordinates without too much trouble; however, it would result in nonlinear BCs which would render the approach useless. Instead, polar coordinates <math>(r, \theta)</math> should be used, since in such a system the equation of a circle is very simple. In order for this to be realized, a polar representation of the Laplacian is necessary. Without going in to the details just yet, the Laplacian is given in (2D) polar coordinates:


:<math>\nabla^2 u = \frac{\partial^2 u}{\partial r^2} + \frac{1}{r} \frac{\partial u}{\partial r} + \frac{1}{r^2} \frac{\partial^2 u}{\partial \theta^2}\,</math>


This result may be derived using differentials and the chain rule; it's not difficult but it's a little long. In these coordinates Laplace's equation reads:


:<math>\frac{\partial^2 u}{\partial r^2} + \frac{1}{r} \frac{\partial u}{\partial r} + \frac{1}{r^2} \frac{\partial^2 u}{\partial \theta^2} = 0\,</math>


Note that in going from Cartesian to polar coordinates, a price was paid: though still linear, Laplace's equation now has ''variable coefficients''. This implies that after separation at least one of the ODEs will have variable coefficients as well.


Let's make up the following BVP, letting '''<math>u = u(r, \theta)</math>''':


:<math>\frac{\partial^2 u}{\partial r^2} + \frac{1}{r} \frac{\partial u}{\partial r} + \frac{1}{r^2} \frac{\partial^2 u}{\partial \theta^2} = 0\,</math>


:<math>u(1, \theta) = \phi(\theta)\,</math>


:<math>u(r, \theta) \ \mbox{is bounded inside the unit circle.}\,</math>


This could represent a physical problem analogous to the previous one: replace the square plate with a disc. Note the apparent absence of sufficient BC to obtain a unique solution. The funny looking statement that u is bounded inside the domain of interest turns out to be the key to getting a unique solution, and it often shows itself in polar coordinates. It "makes up" for the "lack" of BCs. To separate, we as usual incorrectly assume that <math>u(r, \theta) = R(r) \Theta(\theta)</math>:


:<math>\frac{\partial^2}{\partial r^2}(R(r)\Theta(\theta)) + \frac{1}{r} \frac{\partial}{\partial r}(R(r)\Theta(\theta)) + \frac{1}{r^2} \frac{\partial^2}{\partial \theta^2}(R(r)\Theta(\theta)) = 0\,</math>


:<math>R''(r)\Theta(\theta) + \frac{1}{r} R'(r)\Theta(\theta) + \frac{1}{r^2} R(r)\Theta''(\theta) = 0\,</math>


:<math>r^2 R''(r) + r R'(r) + R(r) \frac{\Theta''(\theta)}{\Theta(\theta)} = 0\,</math>


:<math>-\frac{\Theta''(\theta)}{\Theta(\theta)} = \frac{r^2 R''(r) + r R'(r)}{R(r)} = k^2\,</math>


::<math>\Big\Downarrow</math>


:<math>r^2 R''(r) + r R'(r) - k^2 R(r) = 0\,</math>


:<math>\Theta''(\theta) + k^2 \Theta(\theta) = 0\,</math>


Once again, the way the negative sign and the separation constant are arranged makes the solution easier later on. These decisions are made mostly by trial and error.

The <math>R</math> equation is probably one you've never seen before, it's a ''special case'' of the '''Euler differential equation''' (not to be confused with the Euler-Lagrange differential equation). There are a couple of ways to solve it, the most general method would be to change the variables so that an equation with constant coefficients is obtained. An easier way would be to note the pattern in the order of the coefficients and the order of the derivatives, and from there guess a power solution. Either way, the general solution to this simple case of Euler's ODE is given as:


:<math>R(r) = C_1 r^k + C_2 r^{-k}\,</math>


This is a very good example problem since it goes to show that PDE problems very often turn into obscure ODE problems; we got lucky this time since the solution for <math>R</math> was rather simple though its ODE looked pretty bad at first sight. The solution to the <math>\Theta(\theta)</math> equation is:


:<math>\Theta(\theta) = C_3 \cos(k \theta) + C_4 \sin(k \theta)\,</math>


Combining:


:<math>u(r, \theta) = R(r)\Theta(\theta)\,</math>


:<math>u(r, \theta) = \left(C_1 r^k + C_2 r^{-k}\right)(C_3 \cos(k \theta) + C_4 \sin(k \theta))\,</math>


Now, this is where the English sentence condition stating that u must be bounded in the domain of interest may be invoked. As <math>\scriptstyle r  \to 0</math>, the term involving <math>r^{-k}</math> is unbounded. The '''only''' way to fix this is to take <math>C_2 = 0</math>. Note that if this problem were solved between two concentric circles, this term would be nonzero and very important. With that term gone, constants can be merged:


:<math>u(r, \theta) = r^k (A \cos(k \theta) + B \sin(k \theta))\,</math>


Only one condition remains: <math>u = \phi(\theta)</math> on <math>r = 1</math>, yet there are 3 constants. Let's say for now that:


:<math>\phi(\theta) = 2 \cos(4 \theta) + 3 \sin(4 \theta)\,</math>


Then, it's a simple matter of equating coefficients to obtain:


:<math>u(r, \theta) = r^4 (2 \cos(4 \theta) + 3 \sin(4 \theta))\,</math>


Now, let's make the frequencies differ:


:<math>\phi(\theta) = 4 \cos(3 \theta) + \sin(10 \theta)\,</math>


Equating coefficients won't work. However, if the IC were broken up into individual terms, the sum of the solution to the terms just happens to solve the BVP as a whole:


:<math>\phi_3(\theta) = 4 \cos(3 \theta) \quad , \quad \phi_{10}(\theta) = \sin(10 \theta)\,</math>


:<math>u_3(r, \theta) = r^3 4 \cos(3 \theta) \quad , \quad u_{10}(r, \theta) = r^{10} \sin(10 \theta)\,</math>


:<math>u(r, \theta) = u_3(r, \theta) + u_{10}(r, \theta) = r^3 4 \cos(3 \pi \theta) + r^{10} \sin(10 \theta)\,</math>


Verify that the solution above is really equal to the BC at <math>r = 1</math>:


:<math>u(1, \theta) = 1^3 4 \cos(3 \theta) + 1^{10} \sin(10 \theta)) = 4 \cos(3 \theta) + \sin(10 \theta) = \phi(\theta)\,</math>


And, since Laplace's equation is linear, this must solve the PDE as well. What all of this implies is that, if some generic function <math>\phi(\theta)</math> may be expressed as a sum of sinusoids with angular frequencies given by <math>n</math>, all that is needed is a linear combination of the appropriate sum. Notated:


:<math>u(r, \theta) = \sum_{n = 1}^{\infty} r^{n} (A_n \cos(n \theta) + B_n \sin(n \theta))\,</math>


To identify the coefficients, substitute the BC:


:<math>u(1, \theta) = \sum_{n = 1}^{\infty} 1^{n} (A_n \cos(n \theta) + B_n \sin(n \theta))\,</math>


:<math>\phi(\theta) = \sum_{n = 1}^{\infty} (A_n \cos(n \theta) + B_n \sin(n  \theta))\,</math>


The coefficients <math>A_n</math> and <math>B_n</math> may be determined by a (full) Fourier expansion on <math>0 \leq \theta \leq 2 \pi</math>. Note that it's implied that <math>phi(\theta)</math> ''must'' have period <math>2 \pi</math> since we are solving this in a domain (a circle specifically) where <math>0 \leq \theta \leq 2 \pi</math>.

You probably don't like infinite series solutions. Well, it happens that through a variety of manipulations it's possible to express the full solution of this particular problem as:


:<math>u(r, \theta) = \frac{1}{2 \pi}\int_{0}^{2 \pi} \phi(\psi) \frac{1 - r^2}{1 - 2 r \cos(\theta - \psi) + r^2}\ d\psi\,</math>


This is called Poisson's integral formula.

===Derivation of the Laplacian in Polar Coordinates===

Though not necessarily a PDEs concept, it is very important for anyone studying this kind of math to be comfortable with going from one coordinate system to the next. What follows is a long derivation of the Laplacian in 2D polar coordinates using the multivariable chain rule and the concept of differentials. Know, however, that there are really many ways to do this.

Three definitions are all we need to begin:


:<math>\nabla^2 u = \frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2}\,</math>


:<math>x = r \cos(\theta)\,</math>


:<math>y = r \sin(\theta)\,</math>


If it's known that <math>u = u(r, \theta) = u(r(x, y), \theta(x, y))</math>, then the chain rule may be used to express derivatives in terms of <math>r</math> and <math>\theta</math> alone. Two applications will be necessary to obtain the second derivatives. Manipulating operators as if they meant something on their own:


:<math>\frac{\partial}{\partial x} = \frac{\partial}{\partial r}\frac{\partial r}{\partial x} + 
\frac{\partial}{\partial \theta}\frac{\partial \theta}{\partial x}\,</math>


Applying this to itself, treating the underlined bit as a unit dependent on <math>r</math> and <math>\theta</math>:


:<math>\frac{\partial^2}{\partial x^2} = \frac{\partial}{\partial x}\left(\frac{\partial}{\partial x}\right) = 
\frac{\partial}{\partial x}\left(
\underline{
\frac{\partial}{\partial r}\frac{\partial r}{\partial x} + 
\frac{\partial}{\partial \theta}\frac{\partial \theta}{\partial x}
}\right)
\,</math>


:<math>\frac{\partial^2}{\partial x^2} = 

\frac{\partial}{\partial r}\left(
\underline{\frac{\partial}{\partial r}\frac{\partial r}{\partial x} + \frac{\partial}{\partial \theta}\frac{\partial \theta}{\partial x}}
\right)\frac{\partial r}{\partial x}
+
\frac{\partial}{\partial \theta}\left(
\underline{\frac{\partial}{\partial r}\frac{\partial r}{\partial x} + \frac{\partial}{\partial \theta}\frac{\partial \theta}{\partial x}}
\right)\frac{\partial \theta}{\partial x}
\,</math>


:<math>\frac{\partial^2}{\partial x^2} =

\left( 
\frac{\partial^2}{\partial r^2}\frac{\partial r}{\partial x} + 
\frac{\partial}{\partial r}\frac{\partial^2 r}{\partial r \partial x} + 
\frac{\partial^2}{\partial r \partial \theta}\frac{\partial \theta}{\partial x} + 
\frac{\partial}{\partial \theta}\frac{\partial^2 \theta}{\partial r \partial x}
\right)\frac{\partial r}{\partial x}
+
\left( 
\frac{\partial^2}{\partial \theta \partial r}\frac{\partial r}{\partial x} + 
\frac{\partial}{\partial r}\frac{\partial^2 r}{\partial \theta \partial x} + 
\frac{\partial^2}{\partial \theta^2}\frac{\partial \theta}{\partial x} + 
\frac{\partial}{\partial \theta}\frac{\partial^2 \theta}{\partial \theta \partial x}
\right)\frac{\partial \theta}{\partial x}

\,</math>


The above mess may be quickly simplified a little by manipulating the funny looking derivatives:


:<math>
\frac{\partial^2 r}{\partial r \partial x} = 
\frac{\partial^2}{\partial x \partial r}(r) = 
\frac{\partial}{\partial x}\left(\frac{\partial}{\partial r}(r)\right) = 
\frac{\partial}{\partial x}(1) = 0
\,</math>


:<math>
\frac{\partial^2 \theta}{\partial \theta \partial x} = 
\frac{\partial^2}{\partial x \partial \theta}(\theta) = 
\frac{\partial}{\partial x}\left(\frac{\partial}{\partial \theta}(\theta)\right) = 
\frac{\partial}{\partial x}(1) = 0
\,</math>


::<math>\Big\Downarrow</math>


:<math>\frac{\partial^2}{\partial x^2} =

\left( 
\frac{\partial^2}{\partial r^2}\frac{\partial r}{\partial x} + 
\frac{\partial^2}{\partial r \partial \theta}\frac{\partial \theta}{\partial x} + 
\frac{\partial}{\partial \theta}\frac{\partial^2 \theta}{\partial r \partial x}
\right)\frac{\partial r}{\partial x}
+
\left( 
\frac{\partial^2}{\partial \theta \partial r}\frac{\partial r}{\partial x} + 
\frac{\partial}{\partial r}\frac{\partial^2 r}{\partial \theta \partial x} + 
\frac{\partial^2}{\partial \theta^2}\frac{\partial \theta}{\partial x}
\right)\frac{\partial \theta}{\partial x}

\,</math>


This may be made slightly easier to work with if a few changes are made to the way some of the derivatives are written. Also, the <math>y</math> variable follows analogously:


:<math>\frac{\partial^2}{\partial x^2} =
\left( 
\frac{\partial^2}{\partial r^2}\frac{\partial r}{\partial x} + 
\frac{\partial^2}{\partial r \partial \theta}\frac{\partial \theta}{\partial x} + 
\frac{\partial}{\partial \theta}\frac{\partial}{\partial r}\left(\frac{\partial \theta}{\partial x}\right)
\right)\frac{\partial r}{\partial x}
+
\left( 
\frac{\partial^2}{\partial \theta \partial r}\frac{\partial r}{\partial x} + 
\frac{\partial}{\partial r}\frac{\partial}{\partial \theta}\left(\frac{\partial r}{\partial x}\right) + 
\frac{\partial^2}{\partial \theta^2}\frac{\partial \theta}{\partial x}
\right)\frac{\partial \theta}{\partial x}
\,</math>


:<math>\frac{\partial^2}{\partial y^2} =
\left( 
\frac{\partial^2}{\partial r^2}\frac{\partial r}{\partial y} + 
\frac{\partial^2}{\partial r \partial \theta}\frac{\partial \theta}{\partial y} + 
\frac{\partial}{\partial \theta}\frac{\partial}{\partial r}\left(\frac{\partial \theta}{\partial y}\right)
\right)\frac{\partial r}{\partial y}
+
\left( 
\frac{\partial^2}{\partial \theta \partial r}\frac{\partial r}{\partial y} + 
\frac{\partial}{\partial r}\frac{\partial}{\partial \theta}\left(\frac{\partial r}{\partial y}\right) + 
\frac{\partial^2}{\partial \theta^2}\frac{\partial \theta}{\partial y}
\right)\frac{\partial \theta}{\partial y}
\,</math>


Now we need to obtain expressions for some of the derivatives appearing above. The most direct path would use the concept of differentials. If:


:<math>x = r \cos(\theta) \quad , \quad y = r \sin(\theta)\,</math>


Then:


:<math>dx = dr \cos(\theta) - r \sin(\theta) d\theta \quad , \quad dy = dr \sin(\theta) + r \cos(\theta) dr\,</math>


Solving simultaneously for <math>dr</math> and <math>d \theta</math> gives:


:<math>dr = \cos(\theta) dx + \sin(\theta) dy \quad , 
\quad d\theta = -\frac{\sin(\theta)}{r} dx + \frac{\cos(\theta)}{r} dy\,</math>


If <math>z = z(x, y)</math>, then the '''total differential''' is given as:


:<math>dz = \frac{\partial z}{\partial x}dx + \frac{\partial z}{\partial y}dy\,</math>


Note that the two previous equations are of this form (recall that <math>r = r(x, y)</math> and <math>\theta = \theta(x, y)</math>, just like <math>z</math> above), which means that:


:<math>dr = \frac{\partial r}{\partial x}dx + \frac{\partial r}{\partial y}dy \quad , 
\quad d\theta = \frac{\partial \theta}{\partial x}dx + \frac{\partial \theta}{\partial y}dy\,</math>


Equating coefficients quickly yields a bunch of derivatives:


:<math>
\frac{\partial r}{\partial x} = \cos(\theta) \quad , \quad
\frac{\partial r}{\partial y} = \sin(\theta) \quad , \quad
\frac{\partial \theta}{\partial x} = -\frac{\sin(\theta)}{r} \quad , \quad
\frac{\partial \theta}{\partial y} = \frac{\cos(\theta)}{r}
\,</math>


There's an easier but more abstract way to obtain the derivatives above that may be overkill but is worth mentioning anyway. The Jacobian of the functions <math>x(r, \theta)</math> and <math>y(r, \theta)</math> is:


:<math>

\frac{\partial (x, y)}{\partial (r, \theta)} = 

\begin{bmatrix}
\cfrac{\partial x}{\partial r} & \cfrac{\partial x}{\partial \theta} \\
\cfrac{\partial y}{\partial r} & \cfrac{\partial y}{\partial \theta} \\
\end{bmatrix} = 

\begin{bmatrix}
\cos(\theta) & -r \sin(\theta) \\
\sin(\theta) & r \cos(\theta) \\
\end{bmatrix}

</math>


Note that the Jacobian is a compact representation of the coefficients of the total derivative; using <math>\mathbf{z} = \mathbf{z}(\mathbf{x})</math> as an example (bold indicating vectors):


:<math>d\mathbf{z} = \frac{\partial (z_1, z_2, ... \ , z_n)}{\partial (x_1, x_2, ... \ , x_n)} \cdot d\mathbf{x}</math>


So, it follows then that the derivatives that we're interested in may be obtained by inverting the Jacobian matrix:


:<math>
\frac{\partial (r, \theta)}{\partial (x, y)} = \left(\frac{\partial (x, y)}{\partial (r, \theta)}\right)^{-1}
</math>


:<math>
\begin{bmatrix}
\cfrac{\partial r}{\partial x} & \cfrac{\partial r}{\partial y} \\
\cfrac{\partial \theta}{\partial x} & \cfrac{\partial \theta}{\partial y} \\
\end{bmatrix} = 

\begin{bmatrix}
\cfrac{\partial x}{\partial r} & \cfrac{\partial x}{\partial \theta} \\
\cfrac{\partial y}{\partial r} & \cfrac{\partial y}{\partial \theta} \\
\end{bmatrix}^{-1}
</math>


:<math>
\begin{bmatrix}
\cfrac{\partial r}{\partial x} & \cfrac{\partial r}{\partial y} \\
\cfrac{\partial \theta}{\partial x} & \cfrac{\partial \theta}{\partial y} \\
\end{bmatrix} = 


\begin{bmatrix}
\cos(\theta) & -r \sin(\theta) \\
\sin(\theta) & -r \cos(\theta) \\
\end{bmatrix}^{-1}
</math>


:<math>
\begin{bmatrix}
\cfrac{\partial r}{\partial x} & \cfrac{\partial r}{\partial y} \\
\cfrac{\partial \theta}{\partial x} & \cfrac{\partial \theta}{\partial y} \\
\end{bmatrix} = 

\begin{bmatrix}
\cos(\theta) & \sin(\theta) \\
-\cfrac{\sin(\theta)}{r} & \cfrac{\cos(\theta)}{r} \\
\end{bmatrix}
</math>


Though somewhat obscure, this is very convenient and it's just one of the many utilities of the Jacobian matrix. An interesting bit of insight is gained: coordinate changes are senseless unless the Jacobian is invertible everywhere except at isolated points, stated another way the determinant of the Jacobian matrix must be nonzero, otherwise the coordinate change is not one-to-one (note that the determinant will be zero at <math>r = 0</math> in this example. An isolated point such as this is not problematic.).

Either path you take, there should now be enough information to evaluate the Cartesian second derivatives. Working on <math>x</math>:


:<math>\frac{\partial^2}{\partial x^2} =
\left( 
\frac{\partial^2}{\partial r^2}\cos(\theta) - 
\frac{\partial^2}{\partial r \partial \theta}\frac{\sin(\theta)}{r} - 
\frac{\partial}{\partial \theta}\frac{\partial}{\partial r}\left(\frac{\sin(\theta)}{r}\right)
\right)\cos(\theta)
-
\left( 
\frac{\partial^2}{\partial \theta \partial r}\cos(\theta) + 
\frac{\partial}{\partial r}\frac{\partial}{\partial \theta}\left(\cos(\theta)\right) - 
\frac{\partial^2}{\partial \theta^2}\frac{\sin(\theta)}{r}
\right)\frac{\sin(\theta)}{r}
\,</math>


:<math>\frac{\partial^2}{\partial x^2} =
\cos(\theta)^2\frac{\partial^2}{\partial r^2} - 
\frac{\sin(\theta)\cos(\theta)}{r}\frac{\partial^2}{\partial r \partial \theta} + 
\frac{\sin(\theta)\cos(\theta)}{r^2}\frac{\partial}{\partial \theta} - 
\frac{\sin(\theta)\cos(\theta)}{r}\frac{\partial^2}{\partial \theta \partial r} + 
\frac{\sin(\theta)^2}{r}\frac{\partial}{\partial r} + 
\frac{\sin(\theta)^2}{r^2}\frac{\partial^2}{\partial \theta^2}
\,</math>


:<math>\frac{\partial^2}{\partial x^2} =
\cos(\theta)^2\frac{\partial^2}{\partial r^2} - 
\frac{2\sin(\theta)\cos(\theta)}{r}\frac{\partial^2}{\partial r \partial \theta} + 
\frac{\sin(\theta)\cos(\theta)}{r^2}\frac{\partial}{\partial \theta} + 
\frac{\sin(\theta)^2}{r}\frac{\partial}{\partial r} + 
\frac{\sin(\theta)^2}{r^2}\frac{\partial^2}{\partial \theta^2}
\,</math>


Proceeding similarly for <math>y</math>:


:<math>\frac{\partial^2}{\partial y^2} =
\left( 
\frac{\partial^2}{\partial r^2}\sin(\theta) + 
\frac{\partial^2}{\partial r \partial \theta}\frac{\cos(\theta)}{r} + 
\frac{\partial}{\partial \theta}\frac{\partial}{\partial r}\left(\frac{\cos(\theta)}{r}\right)
\right)\sin(\theta)
+
\left( 
\frac{\partial^2}{\partial \theta \partial r}\sin(\theta) + 
\frac{\partial}{\partial r}\frac{\partial}{\partial \theta}\left(\sin(\theta)\right) + 
\frac{\partial^2}{\partial \theta^2}\frac{\cos(\theta)}{r}
\right)\frac{\cos(\theta)}{r}
\,</math>


:<math>\frac{\partial^2}{\partial y^2} =
\sin(\theta)^2\frac{\partial^2}{\partial r^2} + 
\frac{\sin(\theta)\cos(\theta)}{r}\frac{\partial^2}{\partial r \partial \theta} - 
\frac{\sin(\theta)\cos(\theta)}{r^2}\frac{\partial}{\partial \theta} + 
\frac{\sin(\theta)\cos(\theta)}{r}\frac{\partial^2}{\partial \theta \partial r} + 
\frac{\cos(\theta)^2}{r}\frac{\partial}{\partial r} + 
\frac{\cos(\theta)^2}{r^2}\frac{\partial^2}{\partial \theta^2}
\,</math>


:<math>\frac{\partial^2}{\partial y^2} =
\sin(\theta)^2\frac{\partial^2}{\partial r^2} + 
\frac{2\sin(\theta)\cos(\theta)}{r}\frac{\partial^2}{\partial r \partial \theta} - 
\frac{\sin(\theta)\cos(\theta)}{r^2}\frac{\partial}{\partial \theta} + 
\frac{\cos(\theta)^2}{r}\frac{\partial}{\partial r} + 
\frac{\cos(\theta)^2}{r^2}\frac{\partial^2}{\partial \theta^2}
\,</math>


Now, add these tirelessly hand crafted differential operators and watch the result collapse into just 3 nontrigonometric terms:


:<math>\nabla^2 = \frac{\partial^2}{\partial x^2} + \frac{\partial^2}{\partial y^2}\,</math>


:<math>\nabla^2 = 
\cos(\theta)^2\frac{\partial^2}{\partial r^2} - 
\frac{2\sin(\theta)\cos(\theta)}{r}\frac{\partial^2}{\partial r \partial \theta} + 
\frac{\sin(\theta)\cos(\theta)}{r^2}\frac{\partial}{\partial \theta} + 
\frac{\sin(\theta)^2}{r}\frac{\partial}{\partial r} + 
\frac{\sin(\theta)^2}{r^2}\frac{\partial^2}{\partial \theta^2}
 + \,</math>

:<math>\qquad
\sin(\theta)^2\frac{\partial^2}{\partial r^2} + 
\frac{2\sin(\theta)\cos(\theta)}{r}\frac{\partial^2}{\partial r \partial \theta} - 
\frac{\sin(\theta)\cos(\theta)}{r^2}\frac{\partial}{\partial \theta} + 
\frac{\cos(\theta)^2}{r}\frac{\partial}{\partial r} + 
\frac{\cos(\theta)^2}{r^2}\frac{\partial^2}{\partial \theta^2}
\,</math>


:<math>\nabla^2 = 
(\sin(\theta)^2 + \cos(\theta)^2)\frac{\partial^2}{\partial r^2} + 
\left(-\frac{2\sin(\theta)\cos(\theta)}{r} + \frac{2\sin(\theta)\cos(\theta)}{r}\right)\frac{\partial^2}{\partial r \partial \theta} + \,</math>

:<math>
\left(\frac{\sin(\theta)\cos(\theta)}{r^2} - \frac{\sin(\theta)\cos(\theta)}{r^2}\right)\frac{\partial}{\partial \theta} + 
\left(\frac{\sin(\theta)^2}{r} + \frac{\cos(\theta)^2}{r}\right)\frac{\partial}{\partial r} + 
\left(\frac{\sin(\theta)^2}{r^2} + \frac{\cos(\theta)^2}{r^2}\right)\frac{\partial^2}{\partial \theta^2}
\,</math>


:<math>\nabla^2 = \frac{\partial^2}{\partial r^2} + \frac{1}{r}\frac{\partial}{\partial r} + \frac{1}{r^2}\frac{\partial^2}{\partial \theta^2}\,</math>


That was a lot of work. To save trouble, here is the Laplacian in other two other popular coordinate systems:


'''Cylindrical''':

:<math>\nabla^2 u = \frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial u}{\partial r}\right) + \frac{1}{r^2}\frac{\partial^2 u}{\partial \theta^2} + \frac{\partial^2 u}{\partial z^2}\,</math>


'''Spherical''':

:<math>\nabla^2 u = \frac{1}{r^2}\frac{\partial}{\partial r}\left(r^2\frac{\partial u}{\partial r}\right) + 
\frac{1}{r^2 \sin(\theta)}\frac{\partial}{\partial \theta}\left(\sin(\theta) \frac{\partial u}{\partial \theta}\right) + 
\frac{1}{r^2 \sin(\theta)^2}\frac{\partial^2 u}{\partial \phi^2}
\,</math>


Derivatives have been combined wherever possible (not done previously).

===Concluding Remarks===

This was a long, involved chapter. It should be clear that the solutions derived work only for very simple geometries, other geometries may be worked with by grace of '''conformal mappings'''.

The Laplacian (and variations of it) is a very important quantity and it's behaviour is worth knowing like the back of your hand. A sampling of important equations that involve the Laplacian:

:*The Navier Stokes equations.
:*The diffusion equation.
:*Laplace's equation.
:*Poisson's equation.
:*The Helmholtz equation.
:*The Schrödinger equation.
:*The wave equation.


There's a couple of other operators that are similar to (though less important than) the Laplacian, which deserve mention:


* '''Biharmonic operator''', in three Cartesian dimensions:


:<math>\nabla^4 = \frac{\partial^4 u}{\partial x^4} + \frac{\partial^4 u}{\partial y^4} + \frac{\partial^4 u}{\partial z^4} + 2\frac{\partial^4 u}{\partial x^2 \partial y^2} + 2\frac{\partial^4 u}{\partial y^2 \partial z^2} + 2\frac{\partial^4 u}{\partial x^2 \partial z^2}\,</math>


The biharmonic equation is useful in linear elastic theory, for example it can describe "creeping" fluid flow:


:<math>\nabla^4 \psi  = 0\,</math>


* '''d'Alembertian''':


:<math>\Box^2 = \nabla^2 - \frac{1}{c^2}\frac{\partial^2}{\partial t^2}\,</math>


The wave equation may be expressed using the d'Alembertian:


:<math>\Box^2 u = 0\,</math>


Though expressed with the Laplacian is more popular:


:<math>\nabla^2 u = \frac{1}{c^2}\frac{\partial^2 u}{\partial t^2}\,</math>