A place to keep quick notes about Math that I keep forgetting. This is meant to be a scratch notes and cheat sheet for me to write math notes before I forget them or move them somewhere else. Can and will contain errors and/or not complete description in number of places. Use at your own risk.
\(\blacksquare \) Methods to ﬁnd Green function are
reference Wikipedia I need to make one example and apply each of the above methods on it.
\(\blacksquare \) In solving an ODE with constant coeﬃcient just use the characteristic equation to solve the solution.
\(\blacksquare \) In solving an ODE with coeﬃcients that are functions that depends on the independent variable, as in \(y^{\prime \prime }\left ( x\right ) +q\left ( x\right ) y^{\prime }\left ( x\right ) +p\left ( x\right ) y\left ( x\right ) =0\), ﬁrst classify the point \(x_{0}\) type. This means to check how \(p\left ( x\right ) \) and \(q\left ( x\right ) \) behaves at \(x_{0}\). We are talking about the ODE here, not the solution yet.
There are 3 kinds of points. \(x_{0}\) can be normal, or regular singular point, or irregular singular
point. Normal point \(x_{0}\) means \(p\left ( x\right ) \) and \(q\left ( x\right ) \) have Taylor series expansion \(y\left ( x\right ) =\sum _{n=0}^{\infty }a_{n}\left ( x-x_{0}\right ) ^{n}\) that converges to \(y\left ( x\right ) \) at
\(x_{0}\).
Regular singular point \(x_{0}\) means that the above test fails, but \(\lim _{x\rightarrow x_{0}}\left ( x-x_{0}\right ) q\left ( x\right ) \) has a convergent Taylor
series, and also that \(\lim _{x\rightarrow x_{0}}\left ( x-x_{0}\right ) ^{2}p\left ( x\right ) \) now has a convergent Taylor series at \(x_{0}\). This also means the limit
exist.
All this just means we can get rid of the singularity. i.e. \(x_{0}\) is a removable singularity. If this is the case, then the solution at \(x_{0}\) can be assumed to have a Frobenius series \(y\left ( x\right ) =\sum _{n=0}^{\infty }a_{n}\left ( x-x_{0}\right ) ^{n+\alpha }\) where \(a_{0}\neq 0\) and \(\alpha \) must be integer values.
The third type of point, is the hard one. Called irregular singular point. We can’t get rid of it using the above. So we also say the ODE has an essential singularity at \(x_{0}\) (another fancy name for irregular singular point). What this means is that we can’t approximate the solution at \(x_{0}\) using either Taylor nor Frobenius series.
If the point is an irregular singular point, then use the methods of asymptotic. See advanced mathematical methods for scientists and engineers chapter 3. For normal point, use \(y\left ( x\right ) =\sum _{n=0}^{\infty }a_{n}x^{n}\), for regular singular point use \(y\left ( x\right ) =\sum _{n=0}^{\infty }a_{n}x^{n+r}\). Remember, to solve for \(r\) ﬁrst. This should give two values. If you get one root, then use reduction of order to ﬁnd second solution.
\(\blacksquare \) Asymptotic series \(S\left ( z\right ) =c_{0}+\frac{c_{1}}{z}+\frac{c_{2}}{z^{2}}+\cdots \) is series expansion of \(f\left ( z\right ) \) which gives good and rapid approximation for large \(z\) as long as we know when to truncate \(S\left ( z\right ) \) before it becomes divergent. This is the main diﬀerence Asymptotic series expansion and Taylor series expansion.
\(S\left ( z\right ) \) is used to approximate a function for large \(z\) while Taylor (or power series) is used for local approximation or for small distance away from the point of expansion. \(S\left ( z\right ) \) will become divergent, hence it needs to be truncated at some \(n\) to use, where \(n\) is the number of terms in \(S_{n}\left ( z\right ) \). It is optimally truncated when \(n\approx \left \vert z\right \vert ^{2}\).
\(S\left ( x\right ) \) has the following two important properties
We write \(S\left ( z\right ) \sim f\left ( z\right ) \) when \(S\left ( z\right ) \) is the asymptotic series expansion of \(f\left ( z\right ) \) for large \(z\). Most common method to ﬁnd \(S\left ( z\right ) \) is by integration by parts. At least this is what we did in the class I took.
\(\blacksquare \) For Taylor series, leading behavior is \(a_{0}\) no controlling factor? For Frobenius series, leading behavior term is \(a_{0}x^{\alpha }\) and controlling factor is \(x^{\alpha }\). For asymptotic series, controlling factor is assumed to be \(e^{S\left ( x\right ) }\) always. proposed by Carlini (1817)
\(\blacksquare \) Method to ﬁnd the leading behavior of the solution \(y\left ( x\right ) \) near irregular singular point using
asymptotic is called the dominant balance method.
\(\blacksquare \) When solving \(\epsilon y^{\prime \prime }+p\left ( x\right ) y^{\prime }+q\left ( x\right ) y=0\) for very small \(\epsilon \) then use WKB method, if there is no boundary layer between the boundary conditions. If the ODE non-linear, can’t use WKB, has to use boundary layer (B.L.). Example \(\epsilon y^{\prime \prime }+yy^{\prime }-y=0\) with \(y\left ( 0\right ) =0,y\left ( 1\right ) =-2\) then use BL.
\(\blacksquare \) good exercise is to solve say \(\epsilon y^{\prime \prime }+(1+x)y^{\prime }+y=0\) with \(y\left ( 0\right ) =y\left ( 1\right ) \) using both B.L. and WKB and compare the solutions, they should come out the same. \(y\sim \frac{2}{1+x}-\exp \left ( \frac{-x}{\epsilon }-\frac{x^{2}}{2\epsilon }\right ) +O\left ( \epsilon \right ) .\) with BL had to do the matching between the outer and the inner solutions. WKB is easier. But can’t use it for non-linear ODE.
\(\blacksquare \) When there is rapid oscillation over the entire domain, WKB is better. Use WKB to solve Schrodinger equation where \(\epsilon \) becomes function of \(\hslash \) (Planck’s constant, \(6.62606957\times 10^{-34}\) m\(^{2}\)kg/s)
\(\blacksquare \) In second order ODE with non constant coeﬃcient, \(y^{\prime \prime }\left ( x\right ) +p\left ( x\right ) y^{\prime }\left ( x\right ) +q\left ( x\right ) y\left ( x\right ) =0\), if we know one solution \(y_{1}\left ( x\right ) \), then a method called the reduction of order can be used to ﬁnd the second solution \(y_{2}\left ( x\right ) \). Write \(y_{2}\left ( x\right ) =u\left ( x\right ) y_{1}\left ( x\right ) \), plug this in the ODE, and solve for \(u\left ( x\right ) \). The ﬁnal solution will be \(y\left ( x\right ) =c_{1}y_{1}\left ( x\right ) +c_{2}y_{2}\left ( x\right ) \). Now apply I.C.’s to ﬁnd \(c_{1},c_{2}\).
\(\blacksquare \) To ﬁnd particular solution to \(y^{\prime \prime }\left ( x\right ) +p\left ( x\right ) y^{\prime }\left ( x\right ) +q\left ( x\right ) y\left ( x\right ) =f\left ( x\right ) \), we can use a method called undetermined coeﬃcients. But a better method is called variation of parameters, In this method, assume \(y_{p}\left ( x\right ) =u_{1}\left ( x\right ) y_{1}\left ( x\right ) +u_{2}\left ( x\right ) y_{2}\left ( x\right ) \) where \(y_{1}\left ( x\right ) ,y_{2}\left ( x\right ) \) are the two linearly independent solutions of the homogeneous ODE and \(u_{1}\left ( x\right ) ,u_{2}\left ( x\right ) \) are to be determined. This ends up with \(u_{1}\left ( x\right ) =-\int \frac{y_{2}\left ( x\right ) f\left ( x\right ) }{W}dx\) and \(u_{2}\left ( x\right ) =\int \frac{y_{1}\left ( x\right ) f\left ( x\right ) }{W}dx\). Remember to put the ODE in standard form ﬁrst, so \(a=1\), i.e. \(ay^{\prime \prime }\left ( x\right ) +\cdots \). In here, \(W\) is the Wronskian \(W=\begin{vmatrix} y_{1}\left ( x\right ) & y_{2}\left ( x\right ) \\ y_{1}^{\prime }\left ( x\right ) & y_{2}^{\prime }\left ( x\right ) \end{vmatrix} \)
\(\blacksquare \) Two solutions of \(y^{\prime \prime }\left ( x\right ) +p\left ( x\right ) y^{\prime }\left ( x\right ) +q\left ( x\right ) y\left ( x\right ) =0\) are linearly independent if \(W\left ( x\right ) \neq 0\), where \(W\) is the Wronskian.
\(\blacksquare \) For second order linear ODE deﬁned over the whole real line, the Wronskian is either always zero, or not zero. This comes from Abel formula for Wronskian, which is \(W\left ( x\right ) =k\exp \left ( -\int \frac{B\left ( x\right ) }{A\left ( x\right ) }dx\right ) \) for ODE of form \(A\left ( x\right ) y^{\prime \prime }+B\left ( x\right ) y^{\prime }+C\left ( x\right ) y=0\). Since \(\exp \left ( -\int \frac{B\left ( x\right ) }{A\left ( x\right ) }dx\right ) >0\), then it is decided by \(k\). The constant of integration. If \(k=0\), then \(W\left ( x\right ) =0\) everywhere, else it is not zero everywhere.
\(\blacksquare \) For linear PDE, if boundary condition are time dependent, can not use separation of variables. Try Transform method (Laplace or Fourier) to solve the PDE.
\(\blacksquare \) If unable to invert Laplace analytically, try numerical inversion or asymptotic methods. Need to ﬁnd example of this.
\(\blacksquare \) Green function takes the homogeneous solution and the forcing function and constructs a particular solution. For PDE’s, we always want a symmetric Green’s function.
\(\blacksquare \) To get a symmetric Green’s function given an ODE, start by converting the ODE to a Sturm-Liouville form ﬁrst. This way the Green’s function comes out symmetric.
\(\blacksquare \) For numerical solutions of ﬁeld problems, there are basically two diﬀerent problems: Those with closed boundaries and those with open boundaries but with initial conditions. Closed boundaries are elliptical problems which can be cast in the form \(Au=f\), and the other are either hyperbolic or parabolic.
\(\blacksquare \) For numerical solution of elliptical problems, the basic layout is something like this:
Always start with trial solution \(u(x)\) such that \(u_{trial}(x)=\sum _{i=0}^{i=N}C_{i}\phi _{i}(x)\) where the \(C_{i}\) are the unknowns to be determined and the \(\phi _{i}\) are set of linearly independent functions (polynomials) in \(x\).
How to determine those \(C_{i}\) comes next. Use either residual method (Galerkin) or variational
methods (Ritz). For residual, we make a function based on the error \(R=A-u_{trial}f\). It all comes down to
solving \(\int f(R)=0\) over the domain. This is a picture
\(\blacksquare \) Geometric probability distribution. Use when you want an answer to the question: What is the probability you have to do the experiment \(N\) times to ﬁnally get the output you are looking for, given that a probability of \(p\) showing up from doing one experiment.
For example: What is the probability one has to ﬂip a fair coin \(N\) times to get a head? The answer is \(P(X=N)=(1-p)^{k-1}p\). So for a fair coin, \(p=\frac{1}{2}\) that a head will show up from one ﬂip. So the probability we have to ﬂip a coin \(10\) times to get a head is \(P(X=10)=(1-0.5)^{9}(0.5)=0.00097\) which is very low as expected.
\(\blacksquare \) To generate random variable drawn from some distribution diﬀerent from uniform distribution, by only using uniform distribution \(U(0,1)\) do this: Lets say we want to generate random number from exponential distribution with mean \(\mu \).
This distribution has \(pdf(X)=\frac{1}{\mu }e^{\frac{-x}{\mu }}\), the ﬁrst step is to ﬁnd the cdf of exponential distribution, which is known to be \(F(x)=P(X<=x)=1-e^{\frac{-x}{\mu }}\).
Now ﬁnd the inverse of this, which is \(F^{-1}(x)=-\mu \ln (1-x)\). Then generate a random number from the uniform distribution \(U(0,1)\). Let this value be called \(z\).
Now plug this value into \(F^{-1}(z)\), this gives a random number from exponential distribution, which will be \(-\mu \ \ln (1-z)\) (take the natural log of both side of \(F(x)\)).
This method can be used to generate random variables from any other distribution by knowing on \(U(0,1)\). But it requires knowing the CDF and the inverse of the CDF for the other distribution. This is called the inverse CDF method. Another method is called the rejection method
\(\blacksquare \) Given \(u\), a r.v. from uniform distribution over [0,1], then to obtain \(v\), a r.v. from uniform distribution over [A,B], then the relation is \(v=A+(B-A)u\).
\(\blacksquare \) When solving using F.E.M. is best to do everything using isoparametric element (natural coordinates), then ﬁnd the Jacobian of transformation between the natural and physical coordinates to evaluate the integrals needed. For the force function, using Gaussian quadrature method.
\(\blacksquare \) A solution to diﬀerential equation is a function that can be expressed as a convergent series. (Cauchy. Briot and Bouquet, Picard)
\(\blacksquare \) To solve a ﬁrst order ODE using integrating factor. \[ x^{\prime }(t)+p(t)x(t)=f(t) \] then as long as it is linear and \(p(t),f(t)\) are integrable functions in \(t\), then follow these steps
Integrating both sides gives \begin{align*} \ln (I) & =\int{p(t)dt}\\ I(t) & =e^{\int{p(t)dt}} \end{align*}
Where \(I(t)\) is given by (2). Hence \[ x(t)=\frac{\int{e^{\int{p(t)dt}}f(t)\,dt}+C}{e^{\int{p(t)dt}}}\] \(\blacksquare \) A polynomial is called ill-conditioned if we make small change to one of its coeﬃcients and this causes large change to one of its roots.
\(\blacksquare \) To ﬁnd rank of matrix \(A\) by hand, ﬁnd the row echelon form, then count how many zero rows there are. subtract that from number of rows, i.e. \(n\).
\(\blacksquare \) To ﬁnd the basis of the column space of \(A\), ﬁnd the row echelon form and pick the columns with the pivots, there are the basis (the linearly independent columns of \(A\)).
\(\blacksquare \) For symmetric matrix \(A\), its second norm is its spectral radius \(\rho (A)\) which is the largest eigenvalue of \(A\) (in absolute terms).
\(\blacksquare \) The eigenvalues of the inverse of matrix \(A\) is the inverse of the eigenvalues of \(A\).
\(\blacksquare \) If matrix \(A\) of order \(n\times n\), and it has \(n\) distinct eigenvalues, then it can be diagonalized \(A=V\Lambda V^{-1}\), where \[ \Lambda =\begin{pmatrix} e^{\lambda _{1}} & 0 & 0\\ 0 & \ddots & 0\\ 0 & 0 & e^{\lambda n}\end{pmatrix} \] and \(V\) is matrix that has the \(n\) eigenvectors as its columns.
\(\blacksquare \) \(\lim _{k\rightarrow \infty }\int _{x_{1}}^{x_{2}}f_{k}\left ( x\right ) dx=\int _{x_{1}}^{x_{2}}\lim _{k\rightarrow \infty }f_{k}\left ( x\right ) dx\) only if \(f_{k}\left ( x\right ) \) converges uniformly over \(\left [ x_{1},x_{2}\right ] \).
\(\blacksquare \) \(A^{3}=I\), has inﬁnite number of \(A\) solutions. Think of \(A^{3}\) as 3 rotations, each of \(120^{0}\), going back to where we started. Each rotation around a straight line. Hence inﬁnite number of solutions.
\(\blacksquare \) How to integrate \(I=\int \frac{\sqrt{x^{3}-1}}{x}\,dx\).
Let \(u=x^{3}+1\), then \(du=3x^{2}dx\) and the above becomes\[ I=\int \frac{\sqrt{u}}{3x^{3}}\,du=\frac{1}{3}\int \frac{\sqrt{u}}{u-1}\,du \] Now let \(u=\tan ^{2}v\) or \(\sqrt{u}=\tan v\), hence \(\frac{1}{2}\frac{1}{\sqrt{u}}du=\sec ^{2}v\,dv\) and the above becomes\begin{align*} I & =\frac{1}{3}\int \frac{\sqrt{u}}{\tan ^{2}v-1}\left ( 2\sqrt{u}\sec ^{2}v\right ) \,dv\\ & =\frac{2}{3}\int \frac{u}{\tan ^{2}v-1}\sec ^{2}v\,dv\\ & =\frac{2}{3}\int \frac{\tan ^{2}v}{\tan ^{2}v-1}\sec ^{2}v\,dv \end{align*}
But \(\tan ^{2}v-1=\sec ^{2}v\) hence\begin{align*} I & =\frac{2}{3}\int \tan ^{2}v\,dv\\ & =\frac{2}{3}\left ( \tan v-v\right ) \end{align*}
Substituting back\[ I=\frac{2}{3}\left ( \sqrt{u}-\arctan \left ( \sqrt{u}\right ) \right ) \] Substituting back\[ I=\frac{2}{3}\left ( \sqrt{x^{3}+1}-\arctan \left ( \sqrt{x^{3}+1}\right ) \right ) \]
\(\blacksquare \) (added Nov. 4, 2015) Made small diagram to help me remember long division terms used.
\(\blacksquare \) If a linear ODE is equidimensional, as in \(a_{n}x^{n}y^{(n)}+a_{n-1}x^{n-1}y^{(n01)}+\dots \) for example \(x^{2}y^{\prime \prime }-2y=0\) then use ansatz \(y=x^{r}\) this will give equation in \(r\) only. Solve for \(r\) and obtain \(y_{1}=x^{r_{1}},y_{2}=x^{r_{2}}\) and the solution will be \[ y=c_{1}y_{1}+c_{2}y_{2}\] For example, for the above ode, the solution is \(c_{1}x^{2}+\frac{c_{2}}{x}\). This ansatz works only if ODE is equidimensional. So can’t use it on \(xy^{\prime \prime }+y=0\) for example.
If \(r\) is multiple root, use \(x^{r},x^{r}\log (x),x^{r}(\log (x))^{2}\dots \) as solutions.
\(\blacksquare \) for \(x^{i}\), where \(i=\sqrt{-1}\), write it as \(x=e^{\log{x}}\) hence \(x^{i}=e^{i\,\log{x}}=\cos (\log{x})+i\,\sin (\log{x})\)
\(\blacksquare \) For Sturm-Liouville, write it as \(py^{\prime \prime }+p^{\prime }y^{\prime }+\left ( q+\lambda \sigma \right ) y=0\). Now it is easy to compare to. For example, \(x^{\prime \prime }+\lambda x=0\), then we see that \(p=1,q=0,\sigma =1\). And for \(xy^{\prime \prime }+x^{\prime }+\lambda xy=0\) then \(p=1,q=0,\sigma =x\). Do not divide by \(x\) ﬁrst to re-write \(xy^{\prime \prime }+x^{\prime }+\lambda xy=0\) as \(y^{\prime \prime }+\frac{1}{x}x^{\prime }+\lambda y=0\) and then compare, as this makes it hard to compare to standard form and will make error. Books write Sturm-Liouville as \(\frac{d}{dx}\left ( py^{\prime }\right ) +qy+\lambda \sigma y=0\) but I ﬁnd this confusing to compare to this way. Better to expand it to \(py^{\prime \prime }+p^{\prime }y^{\prime }+\left ( q+\lambda \sigma \right ) y=0\) ﬁrst so now it is in a more familiar form and can read it out more directly.
\(\blacksquare \) Some integral tricks: \(\int \sqrt{a^{2}-x^{2}}dx\) use \(x=a\sin \theta \). For \(\int \sqrt{a^{2}+x^{2}}dx\) use \(x=a\tan \theta \) and for \(\int \sqrt{x^{2}-a^{2}}dx\) use \(x=a\sec \theta \).
\(\blacksquare \) \(y^{\prime \prime }+x^{n}y=0\) is called Emden-Fowler form.
\(\blacksquare \) For second order ODE, boundary value problem, with eigenvalue (Sturm-Liouville), remember that having two boundary conditions is not enough to fully solve it.
One boundary condition is used to ﬁnd the ﬁrst constant of integration, and the second boundary condition is used to ﬁnd the eigenvalues.
We still need another input to ﬁnd the second constant of integration. This is normally done by giving the initial value. This problem happens as part of initial value, boundary value problem. The point is, with boundary value and eigenvalue also present, we need 3 inputs to fully solve it. Two boundary conditions is not enough.
\(\blacksquare \) If given ODE \(y^{\prime \prime }\left ( x\right ) +p\left ( x\right ) y^{\prime }\left ( x\right ) +q\left ( x\right ) y\left ( x\right ) =0\) and we are asked to classify if it is singular at \(x=\infty \), then let \(x=\frac{1}{t}\) and check what happens at \(t=0\). The \(\frac{d^{2}}{dx^{2}}\) operator becomes \(\left ( 2t^{3}\frac{d}{dt}+t^{4}\frac{d^{2}}{dt^{2}}\right ) \) and \(\frac{d}{dx}\) operator becomes \(-t^{2}\frac{d}{dt}\). And write the ode now where \(t\) is the independent variable, and follow standard operating procedures. i.e. look at \(\lim _{t\rightarrow 0}xp\left ( t\right ) \) and \(\lim _{t\rightarrow 0}x^{2}q\left ( t\right ) \) and see if these are ﬁnite or not. To see how the operator are mapped, always start with \(x=\frac{1}{t}\) then write \(\frac{d}{dx}=\frac{d}{dt}\frac{dt}{dx}\) and write \(\frac{d^{2}}{dx^{2}}=\left ( \frac{d}{dx}\right ) \left ( \frac{d}{dx}\right ) \). For example, \(\frac{d}{dx}=-t^{2}\frac{d}{dt}\) and \begin{align*} \frac{d^{2}}{dx^{2}} & =\left ( -t^{2}\frac{d}{dt}\right ) \left ( -t^{2}\frac{d}{dt}\right ) \\ & =-t^{2}\left ( -2t\frac{d}{dt}-t^{2}\frac{d^{2}}{dt^{2}}\right ) \\ & =\left ( 2t^{3}\frac{d}{dt}+t^{4}\frac{d^{2}}{dt^{2}}\right ) \end{align*}
Then the new ODE becomes \begin{align*} \left ( 2t^{3}\frac{d}{dt}+t^{4}\frac{d^{2}}{dt^{2}}\right ) y\left ( t\right ) +p\left ( t\right ) \left ( -t^{2}\frac{d}{dt}y\left ( t\right ) \right ) +q\left ( t\right ) y\left ( t\right ) & =0\\ t^{4}\frac{d^{2}}{dt^{2}}y+\left ( -t^{2}p\left ( t\right ) +2t^{3}\right ) \frac{d}{dt}y+q\left ( t\right ) y & =0\\ \frac{d^{2}}{dt^{2}}y+\frac{\left ( -p\left ( t\right ) +2t\right ) }{t^{2}}\frac{d}{dt}y+\frac{q\left ( t\right ) }{t^{4}}y & =0 \end{align*}
The above is how the ODE will always become after the transformation. Remember to change \(p\left ( x\right ) \) to \(p\left ( t\right ) \) using \(x=\frac{1}{t}\) and same for \(q\left ( x\right ) \). Now the new \(p\) is \(\frac{\left ( -p\left ( t\right ) +2t\right ) }{t^{2}}\) and the new \(q\) is \(\frac{q\left ( t\right ) }{t^{4}}\). Then do \(\lim _{t\rightarrow 0}t\frac{\left ( p\left ( t\right ) +2t^{3}\right ) }{t^{4}}\) and \(\lim _{t\rightarrow 0}t^{2}\frac{q\left ( t\right ) }{t^{4}}\) as before.
\(\blacksquare \) If the ODE \(a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+c\left ( x\right ) y=0\), and say \(0\leq x\leq 1\), and there is essential singularity at either end, then use boundary layer or WKB. But Boundary layer method works on non-linear ODE’s (and also on linear ODE) and only if the boundary layer is at end of the domain, i.e. at \(x=0\) or \(x=1\).
WKB method on the other hand, works only on linear ODE, but the singularity can be any where (i.e. inside the domain). As rule of thumb, if the ODE is linear, use WKB. If the ODE is non-linear, we must use boundary layer.
Another diﬀerence, is that with boundary layer, we need to do matching phase at the interface between the boundary layer and the outer layer in order to ﬁnd the constants of integrations. This can be tricky and is the hardest part of solving using boundary layer.
Using WKB, no matching phase is needed. We apply the boundary conditions to the whole solution obtained. See my HWs for NE 548 for problems solved from Bender and Orszag text book.
\(\blacksquare \) In numerical, to ﬁnd if a scheme will converge, check that it is stable and also check that if it is consistent.
It could also be conditionally stable, or unconditionally stable, or unstable.
To check it is consistent, this is the same as ﬁnding the LTE (local truncation error) and
checking that as the time step and the space step both go to zero, the LTE goes to zero. What
is the LTE? You take the scheme and plug in the actual solution in it. An example is better to
explain this part. Lets solve \(u_{t}=u_{xx}\). Using forward in time and centered diﬀerence in space, the
numerical scheme (explicit) is
\[ U_{j}^{n+1}=U_{j}^{n}+\frac{k}{h^{2}}\left ( U_{j-1}^{n}-2U_{j}^{n}+U_{j+1}^{n}\right ) \] The LTE is the diﬀerence between these two (error)\[ LTE=U_{j}^{n+1}-\left ( U_{j}^{n}+\frac{k}{h^{2}}\left ( U_{j-1}^{n}-2U_{j}^{n}+U_{j+1}^{n}\right ) \right ) \] Now plug-in \(u\left ( t^{n},x_{j}\right ) \) in place of \(U_{j}^{n}\) and \(u\left ( t^{n}+k,x_{j}\right ) \) in place of \(U_{j}^{n+1}\)
and plug-in \(u\left ( t^{n},x+h\right ) \) in place of \(U_{j+1}^{n}\) and plug-in \(u\left ( t^{n},x-h\right ) \) in place of \(U_{j-1}^{n}\) in the above. It becomes\begin{equation} LTE=u\left ( t+k,x_{j}\right ) -\left ( u\left ( t^{n},x_{j}\right ) +\frac{k}{h^{2}}\left ( u\left ( t,x-h\right ) -2u\left ( t^{n},x_{j}\right ) +u\left ( t,x+h\right ) \right ) \right ) \tag{1} \end{equation}
Where in the above \(k\) is the time step (also written as \(\Delta t\)) and \(h\) is the space step size. Now comes
the main trick. Expanding the term \(u\left ( t^{n}+k,x_{j}\right ) \) in Taylor, \begin{equation} u\left ( t^{n}+k,x_{j}\right ) =u\left ( t^{n},x_{j}\right ) +k\left . \frac{\partial u}{\partial t}\right \vert _{t^{n}}+\frac{k^{2}}{2}\left . \frac{\partial ^{2}u}{\partial t^{2}}\right \vert _{t^{n}}+O\left ( k^{3}\right ) \tag{2} \end{equation}
And expanding\begin{equation} u\left ( t^{n},x_{j}+h\right ) =u\left ( t^{n},x_{j}\right ) +h\left . \frac{\partial u}{\partial x}\right \vert _{x_{j}}+\frac{h^{2}}{2}\left . \frac{\partial ^{2}u}{\partial x^{2}}\right \vert _{x_{j}}+O\left ( h^{3}\right ) \tag{3} \end{equation}
And expanding\begin{equation} u\left ( t^{n},x_{j}-h\right ) =u\left ( t^{n},x_{j}\right ) -h\left . \frac{\partial u}{\partial x}\right \vert _{x_{j}}+\frac{h^{2}}{2}\left . \frac{\partial ^{2}u}{\partial x^{2}}\right \vert _{x_{j}}-O\left ( h^{3}\right ) \tag{4} \end{equation}
Now plug-in (2,3,4) back into (1). Simplifying, many things drop out, and we should obtain that
\[ LTE=O(k)+O\left ( h^{2}\right ) \] Which says that \(LTE\rightarrow 0\) as \(h\rightarrow 0,k\rightarrow 0\). Hence it is consistent.
To check it is stable, use Von Neumann method for stability. This check if the solution at next time step does not become larger than the solution at the current time step. There can be condition for this. Such as it is stable if \(k\leq \frac{h^{2}}{2}\). This says that using this scheme, it will be stable as long as time step is smaller than \(\frac{h^{2}}{2}\). This makes the time step much smaller than space step.
\(\blacksquare \) do not replace \(\sqrt{x^{2}}\) by \(x\), but by \(|x|\), since \(x=\sqrt{x^{2}}\) only for non negative \(x\).
\(\blacksquare \) For \(ax^{2}+bx+c=0\), with roots \(\alpha ,\beta \) then the relation between roots and coeﬃcients is \begin{align*} \alpha +\beta & =-\frac{b}{a}\\ \alpha \beta & =\frac{c}{a} \end{align*}
\(\blacksquare \) Leibniz rules for integration \begin{align*} \frac{d}{dx}\int _{a\left ( x\right ) }^{b\left ( x\right ) }f\left ( t\right ) dt & =f\left ( b\left ( x\right ) \right ) b^{\prime }\left ( x\right ) -f\left ( a\left ( x\right ) \right ) a^{\prime }\left ( x\right ) \\ \frac{d}{dx}\int _{a\left ( x\right ) }^{b\left ( x\right ) }f\left ( t,x\right ) dt & =f\left ( b\left ( x\right ) \right ) b^{\prime }\left ( x\right ) -f\left ( a\left ( x\right ) \right ) a^{\prime }\left ( x\right ) +\int _{a\left ( x\right ) }^{b\left ( x\right ) }\frac{\partial }{\partial x}f\left ( t,x\right ) dt \end{align*}
\(\blacksquare \) \(\int _{a}^{b}f\left ( x\right ) dx=\int _{a}^{b}f\left ( a+b-x\right ) dx\)
\(\blacksquare \) Diﬀerentiable function implies continuous. But continuous does not imply diﬀerentiable. Example is \(\left \vert x\right \vert \) function.
(Added July, 2017).
If the ODE \(M\left ( x,y\right ) +N\left ( x,y\right ) \frac{dy}{dx}=0\) has both \(M\) and \(N\) homogenous functions of same power, then this ODE can be converted to separable. Here is an example. We want to solve\begin{equation} \left ( x^{3}+8x^{2}y\right ) +\left ( 4xy^{2}-y^{3}\right ) y^{\prime }=0 \tag{1} \end{equation} The above is homogenous in \(M,N\), since the total powers of each term in them is \(3\).\[ \left ( \overset{3}{\overbrace{x^{3}}}+8\overset{3}{\overbrace{x^{2}y}}\right ) +\left ( 4\overset{3}{\overbrace{xy^{2}}}-\overset{3}{\overbrace{y^{3}}}\right ) y^{\prime }=0 \] So we look at each term in \(N\) and \(M\) and add all the powers on each \(x,y\) in them. All powers should add to same value, which is \(3\) in this case. Of course \(N,M\) should be polynomials for this to work. So one should check that they are polynomials in \(x,y\) before starting this process. Once we check \(M,N\) are homogeneous, then we let \[ y=xv \] Therefore now\begin{align} M & =x^{3}+8x^{2}\left ( xv\right ) \nonumber \\ & =x^{3}+8x^{3}v \tag{2} \end{align}
And\begin{align} N & =4x\left ( xv\right ) ^{2}-\left ( xv\right ) ^{3}\nonumber \\ & =4x^{3}v^{2}-x^{3}v^{3} \tag{3} \end{align}
And \begin{equation} y^{\prime }=v+xv^{\prime } \tag{4} \end{equation} Substituting (3,4,5) into (1) gives\begin{align*} \left ( x^{3}+8x^{3}v\right ) +\left ( 4x^{3}v^{2}-x^{3}v^{3}\right ) \left ( v+xv^{\prime }\right ) & =0\\ \left ( x^{3}+8x^{3}v\right ) +\left ( 4x^{3}v^{3}-x^{3}v^{4}\right ) +\left ( 4x^{4}v^{2}-x^{4}v^{3}\right ) v^{\prime } & =0 \end{align*}
Dividing by \(x^{3}\neq 0\) it simpliﬁes to\[ \left ( 1+8v\right ) +\left ( 4v^{3}-v^{4}\right ) +x\left ( 4v^{2}-v^{3}\right ) v^{\prime }=0 \] Which can be written as\begin{align*} x\left ( 4v^{2}-v^{3}\right ) v^{\prime } & =-\left ( \left ( 1+8v\right ) +\left ( 4v^{3}-v^{4}\right ) \right ) \\ v^{\prime } & =\frac{-\left ( \left ( 1+8v\right ) +\left ( 4v^{3}-v^{4}\right ) \right ) }{\left ( 4v^{2}-v^{3}\right ) }\left ( \frac{1}{x}\right ) \end{align*}
We see that it is now separable. We now solve this for \(v\left ( x\right ) \) by direct integration of both sides And then using \(y=xv\) ﬁnd \(y\left ( x\right ) \).
Some simple PDE’s can be solved by direct integration, here are few examples.
Example 1
\[ \frac{\partial z\left ( x,y\right ) }{\partial x}=0 \] Integrating w.r.t. \(x\)., and remembering that now constant of integration will be function of \(y\), hence\[ z\left ( x,y\right ) =f\left ( y\right ) \] Example 2\[ \frac{\partial ^{2}z\left ( x,y\right ) }{\partial x^{2}}=x \] Integrating once w.r.t. \(x\) gives\[ \frac{\partial z\left ( x,y\right ) }{\partial x}=\frac{x^{2}}{2}+f\left ( y\right ) \] Integrating again gives\[ z\left ( x,y\right ) =\frac{x^{3}}{6}+xf\left ( y\right ) +g\left ( y\right ) \] Example 3\[ \frac{\partial ^{2}z\left ( x,y\right ) }{\partial y^{2}}=y \] Integrating once w.r.t. \(y\) gives\[ \frac{\partial z\left ( x,y\right ) }{\partial y}=\frac{y^{2}}{2}+f\left ( x\right ) \] Integrating again gives\[ z\left ( x,y\right ) =\frac{y^{3}}{6}+yf\left ( x\right ) +g\left ( x\right ) \] Example 4\[ \frac{\partial ^{2}z\left ( x,y\right ) }{\partial x\partial y}=0 \] Integrating once w.r.t \(x\) gives\[ \frac{\partial z\left ( x,y\right ) }{\partial y}=f\left ( y\right ) \] Integrating again w.r.t. \(y\) gives\[ z\left ( x,y\right ) =\int f\left ( y\right ) dy+g\left ( x\right ) \] Example 5
Solve \(u_{t}+u_{x}=0\) with \(u\left ( x,1\right ) =\frac{x}{1+x^{2}}\). Let \(u\equiv u\left ( x\left ( t\right ) ,t\right ) \), therefore\[ \frac{du}{dt}=\frac{\partial u}{\partial t}+\frac{\partial u}{\partial x}\frac{dx}{dt}\] Comparing the above with the given PDE, we see that if \(\frac{dx}{dt}=1\) then \(\frac{du}{dt}=0\) or \(u\left ( x\left ( t\right ) ,t\right ) \) is constant. At \(t=1\) we are given that\begin{equation} u=\frac{x\left ( 1\right ) }{1+x\left ( 1\right ) ^{2}} \tag{1} \end{equation} To ﬁnd \(x\left ( 1\right ) \), from \(\frac{dx}{dt}=1\) we obtain that \(x\left ( t\right ) =t+c\). At \(t=1\), \(c=x\left ( 1\right ) -1\). Hence \(x\left ( t\right ) =t+x\left ( 1\right ) -1\) or \[ x\left ( 1\right ) =x\left ( t\right ) +1-t \] Hence solution from (1) becomes\[ u=\frac{x-t+1}{1+\left ( x-t+1\right ) ^{2}}\] Example 6
Solve \(u_{t}+u_{x}+u^{2}=0\).
Let \(u\equiv u\left ( x\left ( t\right ) ,t\right ) \), therefore\[ \frac{du}{dt}=\frac{\partial u}{\partial t}+\frac{\partial u}{\partial x}\frac{dx}{dt}\] Comparing the above with the given PDE, we see that if \(\frac{dx}{dt}=1\) then \(\frac{du}{dt}=-u^{2}\) or \(\frac{-1}{u}=-t+c.\) Hence\[ u=\frac{1}{t+c}\] At \(t=0\), \(c=\frac{1}{u\left ( x\left ( 0\right ) ,0\right ) }\). Let \(u\left ( x\left ( 0\right ) ,0\right ) =f\left ( x\left ( 0\right ) \right ) \). Therefore\[ u=\frac{1}{t+\frac{1}{f\left ( x\left ( 0\right ) \right ) }}\] Now we need to ﬁnd \(x\left ( 0\right ) \). From \(\frac{dx}{dt}=1\), then \(x=t+c\) or \(c=x\left ( 0\right ) \), hence \(x\left ( 0\right ) =x-t\) and the above becomes\[ u\left ( x,t\right ) =\frac{1}{t+\frac{1}{f\left ( x-t\right ) }}=\frac{f\left ( x-t\right ) }{tf\left ( x-t\right ) +1}\]
(added Oct. 20, 2016)
(added Jan. 10, 2019)
If \(y_{1},y_{2}\) are two solutions to \(ay^{\prime \prime }+by^{\prime }+cy=0\) then to show that \(c_{1}y_{1}+c_{2}y_{2}\) is also solution:\begin{align*} ay_{1}^{\prime \prime }+by_{1}^{\prime }+cy_{1} & =0\\ ay_{2}^{\prime \prime }+by_{2}^{\prime }+cy_{2} & =0 \end{align*}
Multiply the ﬁrst ODE by \(c_{1}\) and second ODE by \(c_{2}\)\begin{align*} a\left ( c_{1}y_{1}\right ) ^{\prime \prime }+b\left ( c_{1}y_{1}\right ) ^{\prime }+c\left ( c_{1}y_{1}\right ) & =0\\ a\left ( c_{2}y_{2}\right ) ^{\prime \prime }+b\left ( c_{2}y_{2}\right ) ^{\prime }+c\left ( c_{2}y_{2}\right ) & =0 \end{align*}
Add the above two equations, using linearity of diﬀerentials\[ a\left ( c_{1}y_{1}+c_{2}y_{2}\right ) ^{\prime \prime }+b\left ( c_{1}y_{1}+c_{2}y_{2}\right ) ^{\prime }+c\left ( c_{1}y_{1}+c_{2}y_{2}\right ) =0 \] Therefore \(c_{1}y_{1}+c_{2}y_{2}\) satisﬁes the original ODE. Hence solution.
Since \[ W\left ( x\right ) =\begin{vmatrix} y_{1} & y_{2}\\ y_{1}^{\prime } & y_{2}^{\prime }\end{vmatrix} =y_{1}y_{2}^{\prime }-y_{2}y_{1}^{\prime }\] Where \(y_{1},y_{2}\) are two solutions to \(ay^{\prime \prime }+by^{\prime }+cy=0.\) Write\begin{align*} ay_{1}^{\prime \prime }+py_{1}^{\prime }+cy_{1} & =0\\ ay_{2}^{\prime \prime }+py_{2}^{\prime }+cy_{2} & =0 \end{align*}
Multiply the ﬁrst ODE above by \(y_{2}\) and the second by \(y_{1}\)\begin{align*} ay_{2}y_{1}^{\prime \prime }+py_{2}y_{1}^{\prime }+cy_{2}y_{1} & =0\\ ay_{1}y_{2}^{\prime \prime }+py_{1}y_{2}^{\prime }+cy_{1}y_{2} & =0 \end{align*}
Subtract the second from the ﬁrst\begin{equation} a\left ( y_{2}y_{1}^{\prime \prime }-y_{1}y_{2}^{\prime \prime }\right ) +p\left ( y_{2}y_{1}^{\prime }-y_{1}y_{2}^{\prime }\right ) =0 \tag{1} \end{equation} But \begin{equation} p\left ( y_{2}y_{1}^{\prime }-y_{1}y_{2}^{\prime }\right ) =-pW \tag{2} \end{equation} And\begin{align} \frac{dW}{dx} & =\frac{d}{dx}\left ( y_{1}y_{2}^{\prime }-y_{2}y_{1}^{\prime }\right ) \nonumber \\ & =y_{1}^{\prime }y_{2}^{\prime }+y_{1}y_{2}^{\prime \prime }-y_{2}^{\prime }y_{1}^{\prime }-y_{2}y_{1}^{\prime \prime }\nonumber \\ & =y_{1}y_{2}^{\prime \prime }-y_{2}y_{1}^{\prime \prime } \tag{3} \end{align}
Substituting (2,3) into (1) gives the Wronskian diﬀerential equation\begin{align*} -a\left ( \frac{dW}{dx}\right ) -pW & =0\\ aW^{\prime }+pW & =0 \end{align*}
Whose solution is
\[ W\left ( x\right ) =Ce^{-\int \frac{p}{a}dx}\] Where \(C\) is constant of integration.
Remember: \(W\left ( x_{0}\right ) =0\) does not mean the two functions are linearly dependent. The functions can still be Linearly independent on other interval, It just means \(x_{0}\) can’t be in the domain of the solution for two functions to be solutions. However, if the two functions are linearly dependent, then this implies \(W=0\) everywhere. So to check if two functions are L.D., need to show that \(W=0\) everywhere.
Regular Sturm-Liouville ODE is an eigenvalue boundary value ODE. This means the ODE has an eigenvalue in it \(\lambda \), where solutions exists only for speciﬁc values of eigenvalues. The ODE is\[ \frac{d}{dx}\left ( p\left ( x\right ) \frac{dy}{dx}\right ) +q\left ( x\right ) y\left ( x\right ) +\lambda \sigma \left ( x\right ) y\left ( x\right ) =0\qquad a<x<b \]
Some books use \(-q\left ( x\right ) y\left ( x\right ) \) in the above instead of \(+q\left ( x\right ) y\left ( x\right ) \). The above also can be written as
\[ py^{\prime \prime }+p^{\prime }y^{\prime }+\left ( q+\lambda \sigma \right ) y=0 \] With the restrictions that \(p\left ( x\right ) ,q\left ( x\right ) ,\sigma \left ( x\right ) \) are real functions that are continuous everywhere over \(a\leq x\leq b\) and also we need \(p\left ( x\right ) >0,\sigma \left ( x\right ) >0\). The \(\sigma \left ( x\right ) \) is called the weight function. But this is not all. The boundary conditions must be linear homogeneous, of this form\begin{align*} \beta _{1}y\left ( a\right ) +\beta _{2}y^{\prime }\left ( a\right ) & =0\\ \beta _{3}y\left ( b\right ) +\beta _{4}y^{\prime }\left ( b\right ) & =0 \end{align*}
Where \(\beta _{i}\) are just real constants. Some of them can be zero but not all. For example, \(\beta _{1}y\left ( a\right ) =0,\beta _{4}y^{\prime }\left ( b\right ) =0\) is OK.
Boundary conditions do not have to mixed, but they can be in general. But they must be homogeneous.
Notice that periodic boundary conditions are not allowed. Well, they are allowed, but then the problem is no longer called Sturm-Liouville. The above is just the deﬁnition of the equation and its boundary conditions. Below is list of the important properties of this ODE. Each one of these properties have a proof.
\(\blacksquare \) All eigenvalues \(\lambda _{n}\) are real. No complex \(\lambda _{n}\).
\(\blacksquare \) Each eigenvalue \(\lambda _{n}\), will have one and only one real eigenfunction \(\phi _{n}\left ( x\right ) \) associated with it. (in higher dimensions, the eigenvalue problem can have more than one eigenfunction associated with one eigenvalue. For example, heat PDE in rectangle). But we can use Gram-schmidt to make these eigenfunctions linearly independent if we need to). But for exam, just worry about 1D for now.
\(\blacksquare \) There is smallest \(\lambda \), called \(\lambda _{1}\) and there are inﬁnite number of eigenvalues. Hence eigenvalues are ordered \(\lambda _{1}<\lambda _{2}<\lambda _{3}\cdots \)
\(\blacksquare \) Each eigenfunction \(\phi _{n}\left ( x\right ) \) will have \(n-1\) zeros in \(a<x<b\). Note that this does not include the end points. This means, \(\phi _{3}\left ( x\right ) \) will have two zeros inside the domain. i.e. there are two \(x\) locations where \(\phi _{3}\left ( x\right ) =0\) in \(a<x<b\).
\(\blacksquare \) Eigenfunctions \(\phi _{n}\left ( x\right ) \) make a complete set of basis. This means any piecewise continuous function \(f\left ( x\right ) \) can be represented as \(f\left ( x\right ) \sim \sum _{n=1}^{\infty }a_{n}\phi _{n}\left ( x\right ) \). This is called the generalized Fourier series.
\(\blacksquare \) This follows from (5). Each eigenfunction is orthogonal to another eigenfunction, but with the weight in there. This means \(\int _{a}^{b}\phi _{n}\left ( x\right ) \phi _{m}\left ( x\right ) \sigma \left ( x\right ) dx=0\) if \(n\neq m\).
\(\blacksquare \) Rayleigh quotient relates an eigenvalue to its eigenfunction. Starting with the SL ODE \(\left ( p\phi ^{\prime }\right ) ^{\prime }+q\phi =-\lambda \sigma \phi \), then multiplying by\(\ \phi \) and integrating, we obtain \(\int _{a}^{b}\phi \left ( p\phi ^{\prime }\right ) ^{\prime }+q\phi ^{2}dx=\) \(-\int _{a}^{b}\lambda \sigma \phi ^{2}dx\) and solving for \(\lambda \), gives \[ \lambda =\frac{\int _{a}^{b}\phi \left ( p\phi ^{\prime }\right ) ^{\prime }+q\phi ^{2}dx}{\int _{a}^{b}\phi ^{2}\sigma \left ( x\right ) dx}\] Carrying integration by parts on ﬁrst part of the integral in numerator, it becomes \[ \lambda =\frac{\left ( -p\phi \phi ^{\prime }\right ) _{b}^{a}+\int _{a}^{b}p\left ( \phi ^{\prime }\right ) ^{2}-q\phi ^{2}dx}{\int _{a}^{b}\phi ^{2}\sigma \left ( x\right ) dx}\] But this becomes much simpler when we plug-in the boundary conditions that we must use, making the above \[ \lambda =\frac{\int _{a}^{b}p\left ( \phi ^{\prime }\right ) ^{2}-q\phi ^{2}dx}{\int _{a}^{b}\phi ^{2}\sigma \left ( x\right ) dx}\] And if \(q=0,p=1\) and \(\sigma \left ( x\right ) =1\), then it becomes \(\lambda =\frac{\int _{a}^{b}\left ( \phi ^{\prime }\right ) ^{2}dx}{\int _{a}^{b}\phi ^{2}dx}.\) Rayleigh quotient is useful to show that \(\lambda \) can be positive without solving for \(\phi \) and also used to estimate the value of the minimum \(\lambda \) by replacing \(\phi \) by a trial function and actually solving for \(\lambda \) using \(\lambda =\frac{\int _{a}^{b}\phi \left ( p\phi ^{\prime }\right ) ^{\prime }+q\phi ^{2}dx}{\int _{a}^{b}\phi ^{2}\sigma \left ( x\right ) dx}\) to obtain a numerical estimate of \(\lambda _{1}\).
\(\blacksquare \) There is symmetry relation. In operator form, let \(L\equiv \frac{d}{dx}\left ( p\left ( x\right ) \frac{d}{dx}\right ) +q\left ( x\right ) \), then we have \(\int _{a}^{b}uL\left [ v\right ] -vL\left [ u\right ] dx=0\) where \(u,v\) are any two diﬀerent eigenfunctions.
\(\blacksquare \) There is also what is called Lagrange identity, which says \(uL\left [ v\right ] -vL\left [ u\right ] =\frac{d}{dx}\left ( p\left ( uv^{\prime }-vu^{\prime }\right ) \right ) \). This comes from simply expanding \(L[v],L\left [ u\right ] \) and simplifying things. Just calculus.
\(\blacksquare \) Green formula, follows from Lagrange identity, which just gives the integral form \(\int _{a}^{b}uL\left [ v\right ] -vL\left [ u\right ] dx=\left [ p\left ( uv^{\prime }-vu^{\prime }\right ) \right ] _{a}^{b}\). But form Sturm-Liouville, we know that \(\int _{a}^{b}uL\left [ v\right ] -vL\left [ u\right ] dx=0\) (from 8). So this really just says that \(\left [ p\left ( uv^{\prime }-vu^{\prime }\right ) \right ] _{a}^{b}=0\). But we know this already from boundary conditions. So I am not sure why this is useful for Sturm-Liouville now. Since it is just saying the same thing again.
\(\blacksquare \) If \(p\left ( x\right ) =0\) at the left or right end (boundaries), then the problem is now called singular Sturm-Liouville. This is actually the important case. In regular S.L., \(p\left ( x\right ) \) must be positive everywhere. We only consider \(p\left ( x\right ) =0\) at the ends, not in the middle or any other place. When this happens, then solutions to the ODE generate special functions, such as Bessel, Legendre, Chebyshev. These special functions are solution to the Singular S.L. ODE, not the regular S.L. ODE. In addition, at the end where \(p\left ( x\right ) =0\), the boundary conditions must be bounded. Not the regular boundary conditions above. For example, if \(p\left ( a\right ) =0\) where \(x=a\) is say the left end, the the boundary conditions at the left end must be \(y\left ( a\right ) <\infty ,y^{\prime }\left ( a\right ) <\infty \). Notice that singular can be at left end, right end or at both at same time. At the end where \(p\left ( x\right ) =0\), the boundary conditions must be bounded type.
\(\blacksquare \) Symmetry relation \[ \int _{a}^{b}uL\left [ v\right ] -vL\left [ u\right ] dx=0 \] where \(u,v\) are are two eigenfunction. But remember, this does NOT mean the integrand is identically zero, or \(uL\left [ v\right ] -vL\left [ u\right ] =0\). Only when \(u,v\) happened to have same eigenvalue, only then we can say that \(uL\left [ v\right ] -vL\left [ u\right ] =0\). But for unique eigenfunctions only the integral form is zero \(\int _{a}^{b}uL\left [ v\right ] -vL\left [ u\right ] dx=0.\) This is important diﬀerence, used in the proof below, that diﬀerent eigenfunctions have diﬀerent eigenvalues (for the scalar SL case). For higher dimensions, we can have more than one eigenfunction for same eigenvalue.
Proofs for these properties follows
Given regular Sturm-Liouville (RSL) ODE\begin{align} \left ( py^{\prime }\right ) ^{\prime }+qy & =-\lambda \sigma y\tag{1}\\ py^{\prime \prime }+p^{\prime }y^{\prime }+qy & =-\lambda \sigma y\nonumber \end{align}
In operator form \[ L\equiv \frac{d}{dx}\left ( p\left ( x\right ) \frac{d}{dx}\right ) +q\left ( x\right ) \] And (1) becomes\[ L[y]=-\lambda \sigma y \] When solving RSL ode, since it is an eigenvalue ODE with associated boundary conditions, we will get inﬁnite number of non-negative eigenvalues.
For each eigenvalue, there is associated with it one eigenfunction (for 1D case). Looking at any two diﬀerent eigenfunctions, say \(u,v\), then the symmetry relation says the following\[ \int _{a}^{b}uL\left [ v\right ] dx=\int _{a}^{b}vL\left [ u\right ] dx \] Now we will show the above is true. This requires integration by parts two times. We start from the LHS expression and at the end we should end up with the integral on the RHS. Let \(I=\int _{a}^{b}uL\left [ v\right ] dx\), then \begin{align} I & =\int _{a}^{b}uL\left [ v\right ] dx\tag{1}\\ & =\int _{a}^{b}u\left ( \frac{d}{dx}\left ( p\frac{d}{dx}v\right ) +qv\right ) dx\nonumber \end{align}
Where we used \(L\left [ v\right ] =\left ( \frac{d}{dx}\left ( p\frac{d}{dx}\right ) +q\right ) v=\left ( \frac{d}{dx}\left ( p\frac{d}{dx}v\right ) +qv\right ) \) in the above. Hence\begin{equation} I=\int _{a}^{b}u\overset{dV}{\overbrace{\frac{d}{dx}\left ( p\frac{dv}{dx}\right ) }}dx+\int _{a}^{b}qvudx \tag{2} \end{equation} Now we will do integration by parts on the ﬁrst integral. Let \(I_{2}=\int _{a}^{b}u\frac{d}{dx}\left ( p\frac{dv}{dx}\right ) dx\). Using \(\int UdV=UV-\int VdU\), and if we let \(U=u,dV=\frac{d}{dx}\left ( p\frac{dv}{dx}\right ) \), then \(dU=\frac{du}{dx},V=p\frac{dv}{dx}\). Hence \[ I_{2}=\left ( up\frac{dv}{dx}\right ) _{a}^{b}-\int _{a}^{b}\overset{dV}{\overbrace{\frac{dv}{dx}}}\overset{U}{\overbrace{p\frac{du}{dx}}}dx \] We now apply integration by parts again to the second integral above. But now let \(U=p\frac{du}{dx}\) and \(dV=\frac{dv}{dx}\), hence \(dU=\frac{d}{dx}\left ( p\frac{du}{dx}\right ) \) and \(V=v\), therefore the above becomes\[ I_{2}=\left ( up\frac{dv}{dx}\right ) _{a}^{b}-\left [ \left ( p\frac{du}{dx}v\right ) _{a}^{b}-\int _{a}^{b}v\frac{d}{dx}\left ( p\frac{du}{dx}\right ) dx\right ] \] Substituting the above into (2) gives\begin{equation} I=\left ( up\frac{dv}{dx}\right ) _{a}^{b}-\left [ \left ( p\frac{du}{dx}v\right ) _{a}^{b}-\int _{a}^{b}v\frac{d}{dx}\left ( p\frac{du}{dx}\right ) dx\right ] +\int _{a}^{b}qvudx \tag{3} \end{equation} Now comes the part where the boundary conditions are important. In RSL, the boundary conditions are such that all the terms \(\left ( up\frac{dv}{dx}\right ) _{a}^{b},\left ( p\frac{du}{dx}v\right ) _{a}^{b}\) vanish. This is because the boundary conditions are\begin{align*} \beta _{1}u\left ( a\right ) +\beta _{2}u^{\prime }\left ( a\right ) & =0\\ \beta _{3}u\left ( b\right ) +\beta _{4}u^{\prime }\left ( b\right ) & =0 \end{align*}
And\begin{align*} \beta _{1}v\left ( a\right ) +\beta _{2}v^{\prime }\left ( a\right ) & =0\\ \beta _{3}v\left ( b\right ) +\beta _{4}v^{\prime }\left ( b\right ) & =0 \end{align*}
So now (3) becomes\begin{align*} I & =-\left [ -\int _{a}^{b}v\frac{d}{dx}\left ( p\frac{du}{dx}\right ) dx\right ] +\int _{a}^{b}qvudx\\ & =\int _{a}^{b}v\frac{d}{dx}\left ( p\frac{du}{dx}\right ) dx+\int _{a}^{b}qvudx\\ & =\int _{a}^{b}v\left ( \frac{d}{dx}\left ( p\frac{du}{dx}\right ) +qu\right ) dx\\ & =\int _{a}^{b}vL\left [ u\right ] dx \end{align*}
Therefore, we showed that \(\int _{a}^{b}uL\left [ v\right ] dx=\int _{a}^{b}vL\left [ u\right ] dx\). The only thing to watch for here, is which term to make \(U\) and which to make \(dV\) when making integration by parts. To remember, always start with \(uL[v]\) and make \(dV\) the term with \(\frac{dv}{dx}\) in them during integration parts.
Assume \(\lambda \) is complex. The corresponding eigenfunctions are complex also. Let \(\phi \) be corresponding eigenfunction\begin{equation} L\left [ \phi \right ] =-\lambda \sigma \phi \tag{1} \end{equation} Taking complex conjugate of both sides, and since \(\sigma \) is real, we obtain\[ \overline{L\left [ \phi \right ] }=-\overline{\lambda }\sigma \overline{\phi }\] But \(\overline{L\left [ \phi \right ] }=L\left [ \overline{\phi }\right ] \) since all the coeﬃcients of the ODE are real. The above becomes\begin{equation} L\left [ \overline{\phi }\right ] =-\overline{\lambda }\sigma \overline{\phi } \tag{2} \end{equation} But by symmetry, we know that\[ \int _{a}^{b}\phi L\left [ \overline{\phi }\right ] -\overline{\phi }L\left [ \phi \right ] dx=0 \] Substituting (1),(2) into the above gives\begin{align*} \int _{a}^{b}\phi \left ( -\overline{\lambda }\sigma \overline{\phi }\right ) -\overline{\phi }\left ( -\lambda \sigma \phi \right ) dx & =0\\ \int _{a}^{b}-\overline{\lambda }\sigma \phi \overline{\phi }+\lambda \sigma \overline{\phi }\phi dx & =0\\ \int _{a}^{b}\left ( \lambda -\overline{\lambda }\right ) \left ( \sigma \overline{\phi }\phi \right ) dx & =0\\ \left ( \lambda -\overline{\lambda }\right ) \int _{a}^{b}\sigma \overline{\phi }\phi dx & =0 \end{align*}
But \(\overline{\phi }\phi =\left \vert \phi \right \vert \) which is positive. Also the weight \(\sigma \left ( x\right ) \) is positive by deﬁnition. Hence for the above to zero, it must be that \(\lambda =\overline{\lambda }\). Which means \(\lambda \) is real. QED.
So the main tools to use in this proof: Deﬁnition of \(L\left [ \phi \right ] =-\lambda \sigma \phi \,,\) and symmetry relation and that \(\overline{L\left [ \phi \right ] }=L\left [ \overline{\phi }\right ] \). This might come up in the exam.
Now will show that there is one eigenfunction associated with each eigenvalue (Again, this is for 1D, it is possible to get more than one eigenfunction for same eigenvalue for 2D, as mentioned earlier). By contradiction, assume that \(\lambda \) has two eigenfunctions \(\phi _{1},\phi _{2}\) associated with it. Hence\begin{align*} L\left [ \phi _{1}\right ] & =-\lambda \sigma \phi _{1}\\ L\left [ \phi _{2}\right ] & =-\lambda \sigma \phi _{2} \end{align*}
From the ﬁrst equation, \(\lambda =\frac{-L\left [ \phi _{1}\right ] }{\sigma \phi _{1}}\), substituting this into the second equation gives\begin{align*} L\left [ \phi _{2}\right ] & =\frac{L\left [ \phi _{1}\right ] }{\phi _{1}}\phi _{2}\\ \phi _{1}L\left [ \phi _{2}\right ] -\phi _{2}L\left [ \phi _{1}\right ] & =0 \end{align*}
By Lagrange identity, \(\phi _{1}L\left [ \phi _{2}\right ] -\phi _{2}L\left [ \phi _{1}\right ] =\frac{d}{dx}\left ( p\left ( \phi _{1}\phi _{2}^{\prime }-\phi _{2}\phi _{1}^{\prime }\right ) \right ) \), hence this means that \begin{align*} \frac{d}{dx}\left ( p\left ( \phi _{1}\phi _{2}^{\prime }-\phi _{2}\phi _{1}^{\prime }\right ) \right ) & =0\\ p\left ( \phi _{1}\phi _{2}^{\prime }-\phi _{2}\phi _{1}^{\prime }\right ) & =c_{1} \end{align*}
Where \(c_{1}\) is some constant. This is the main diﬀerence between the above argument, and between Lagrange identity. This can be confusing. So let me talk more about this. In Lagrange identity, we write \[ \phi _{1}L\left [ \phi _{2}\right ] -\phi _{2}L\left [ \phi _{1}\right ] =\frac{d}{dx}\left ( p\left ( \phi _{1}\phi _{2}^{\prime }-\phi _{2}\phi _{1}^{\prime }\right ) \right ) \] And when \(\phi _{2},\phi _{1}\) also satisfy the SL boundary condition, only then we say that\begin{align*} \int _{a}^{b}\phi _{1}L\left [ \phi _{2}\right ] -\phi _{2}L\left [ \phi _{1}\right ] dx & =0\\ \left . p\left ( \phi _{1}\phi _{2}^{\prime }-\phi _{2}\phi _{1}^{\prime }\right ) \right \vert _{a}^{b} & =0 \end{align*}
But the above is not the same as saying \(\phi _{1}L\left [ \phi _{2}\right ] -\phi _{2}L\left [ \phi _{1}\right ] =0\). This is important to keep in mind. Only the integral form is zero for any two functions with the SL B.C..Now we continue. We showed that \(p\left ( \phi _{1}\phi _{2}^{\prime }-\phi _{2}\phi _{1}^{\prime }\right ) =c\). In SL, this constant is zero due to B.C. Hence\[ p\left ( \phi _{1}\phi _{2}^{\prime }-\phi _{2}\phi _{1}^{\prime }\right ) =0 \] But \(p>0\) by deﬁnition. Hence \(\phi _{1}\phi _{2}^{\prime }-\phi _{2}\phi _{1}^{\prime }=0\) or\begin{align*} \frac{d}{dx}\left ( \frac{\phi _{2}}{\phi _{1}}\right ) & =0\\ \frac{\phi _{2}}{\phi _{1}} & =c_{2} \end{align*}
or \(\phi _{2}=c_{2}\phi _{1}\). So the eigenfunctions are linearly dependent. One is just scaled version of the other. But eigenfunction must be linearly independent.
Hence assumption is not valid, and there can not be two linearly independent eigenfunctions for same eigenvalue. Notice also that \(\phi _{1}\phi _{2}^{\prime }-\phi _{2}\phi _{1}^{\prime }=0\) is the just the Wronskian. When it is zero, we know the functions are linearly dependent. The important part in the above proof, is that \(\phi _{1}\phi _{2}^{\prime }-\phi _{2}\phi _{1}^{\prime }=0\) only when \(\phi _{1},\phi _{2}\) happened to have same eigenvalue.
The idea of this proof is to assume the eigenfunction is complex, then show that its real part and its complex part both satisfy the ODE and the boundary conditions. But since they are both use the same eigenvalue, then the real part and the complex part must be linearly dependent. This implies the eigenfunction must be real. (think of Argand diagram)
Assume that \(\phi =U+iV\) is complex eigenfunction with real part \(U\) and complex part \(V\). Then since\[ L\left [ \phi \right ] =-\lambda \sigma \phi \] The above is just writing the Sturm-Liouville ODE in operator form, where \(L\) is the operator as above. Now we have\[ L\left [ U+iV\right ] =-\lambda \sigma \left ( U+iV\right ) \] By linearity of operator \(L\)\[ L\left [ U\right ] +iL\left [ V\right ] =-\lambda \sigma U-i\lambda \sigma V \] Which implies\begin{align*} L\left [ U\right ] & =-\lambda \sigma U\\ L\left [ V\right ] & =-\lambda \sigma V \end{align*}
So we showed that the real and complex part satisfy S.L. Now we need to show they are satisfy S.L. boundary conditions. Since\begin{align*} \beta _{1}\phi \left ( a\right ) +\beta _{2}\phi ^{\prime }\left ( a\right ) & =0\\ \beta _{3}\phi \left ( b\right ) +\beta _{4}\phi ^{\prime }\left ( b\right ) & =0 \end{align*}
Where \(x=a,x=b\) are the left and right ends of the domain. Then\begin{align*} \beta _{1}\left ( U\left ( a\right ) +iV\left ( a\right ) \right ) +\beta _{2}\left ( U^{\prime }\left ( a\right ) +iV^{\prime }\left ( a\right ) \right ) & =0\\ \beta _{3}\left ( U\left ( b\right ) +iV\left ( b\right ) \right ) +\beta _{4}\left ( U^{\prime }\left ( b\right ) +iV^{\prime }\left ( b\right ) \right ) & =0 \end{align*}
Hence\begin{align*} \beta _{1}U\left ( a\right ) +\beta _{2}U^{\prime }\left ( a\right ) & =0\\ \beta _{1}V\left ( a\right ) +\beta _{2}V^{\prime }\left ( a\right ) & =0\\ \beta _{3}U\left ( b\right ) +\beta _{4}U^{\prime }\left ( b\right ) & =0\\ \beta _{3}V\left ( b\right ) +\beta _{4}V^{\prime }\left ( b\right ) & =0 \end{align*}
So the above means both \(U\) and \(V\) satisfy the boundary conditions of S.L. But since both \(U,V\) have the same eigenvalue, then they must be linearly dependent, since we know with S.L. each eigenfunction (or one linearly dependent to it) have only one eigenvalue. This means\[ V=cU \] Where \(c\) is some constant. In other words, \begin{align*} \phi & =U+iV\\ & =U+icU\\ & =U\left ( 1+ic\right ) \\ & =c_{0}U \end{align*}
Where \(c_{0}\) is new constant. (OK it happens to be complex constant, but it is OK to do so, we always do this trick in other places, if it will make me feel better, I could take the magnitude of the constant). So all what the above says, is that we assumed \(\phi \) to be complex, and found that it is real. So it can’t be complex.
Given two diﬀerent eigenfunctions \(\phi _{1},\phi _{2}\). Hence\begin{align} L\left [ \phi _{1}\right ] & =-\lambda _{1}\sigma \phi _{1}\tag{1}\\ L\left [ \phi _{2}\right ] & =-\lambda _{2}\sigma \phi _{2} \tag{2} \end{align}
From symmetry integral relation, since these eigenfunctions also satisfy S.L. boundary conditions, we can write\[ \int _{a}^{b}\phi _{1}L\left [ \phi _{2}\right ] -\phi _{2}L\left [ \phi _{1}\right ] dx=0 \] Replacing (1,2) into the above\begin{align*} \int _{a}^{b}\phi _{1}\left ( -\lambda _{2}\sigma \phi _{2}\right ) -\phi _{2}\left ( -\lambda _{1}\sigma \phi _{1}\right ) dx & =0\\ \int _{a}^{b}-\lambda _{2}\sigma \phi _{1}\phi _{2}+\lambda _{1}\sigma \phi _{2}\phi _{1}dx & =0\\ \int _{a}^{b}\left ( \lambda _{1}-\lambda _{2}\right ) \left ( \sigma \phi _{1}\phi _{2}\right ) dx & =0\\ \left ( \lambda _{1}-\lambda _{2}\right ) \int _{a}^{b}\sigma \phi _{1}\phi _{2}dx & =0 \end{align*}
But \(\left ( \lambda _{1}-\lambda _{2}\right ) \neq 0\) since there are diﬀerent eigenvalues for diﬀerent eigenfunctions. Hence \[ \int _{a}^{b}\sigma \phi _{1}\left ( x\right ) \phi _{2}\left ( x\right ) dx=0 \] Which means \(\phi _{1,}\phi _{2}\) are orthogonal to each others with weight \(\sigma \left ( x\right ) \)
When S.L. is singular, meaning \(p=0\) at one or both ends, we end with important class of ODE’s, whose solutions are special functions (not \(\sin \) or \(\cos \)) as the case with the regular S.L. Recall that S.L. is\begin{align} \frac{d}{dx}\left ( p\left ( x\right ) \frac{dy}{dx}\right ) +q\left ( x\right ) y\left ( x\right ) +\lambda \sigma \left ( x\right ) y\left ( x\right ) & =0\tag{1}\\ py^{\prime \prime }+p^{\prime }y^{\prime }+qy & =-\lambda \sigma y\nonumber \end{align}
With the regular S.L., we say that \(p\left ( x\right ) >0\) over the whole domain, including end points. But with singular, this is not the case. Here are three important S.L. ODE’s that are singular.
Bessel equation: \begin{align*} x^{2}y^{\prime \prime }+xy+\lambda x^{2}y & =0\qquad 0<x<1\\ xy^{\prime \prime }+y+\lambda xy & =0 \end{align*}
Or in standard form\[ \left ( xy^{\prime }\right ) ^{\prime }=-\lambda xy \] Comparing to SL form, then \(p\left ( x\right ) =x,q\left ( x\right ) =0,\sigma \left ( x\right ) =x\). At \(x=0\), then \(p\left ( x\right ) =0\), which is what makes it singular. (also the weight happens to be zero), but we only care about \(p\) being zero or not, at one of the ends. So to check if S.L. is singular or not, just need to check if \(p\left ( x\right ) \) is zero or not at one of the ends.
As mentioned before, when \(p=0\) at one of the ends, we can’t use the standard B.C. for the regular S.L. instead, at the end where \(p=0\), we must use what is called bounded boundary conditions, which is in this case \(y\left ( 0\right ) <\infty ,y^{\prime }\left ( 0\right ) <\infty \). The solution to this ODE will be in terms of Bessel functions. Notice that this happened to be singular, due to the domain starting at \(x=0\). If the domain happened to be from \(x=0.5\) to say \(x=1\), then it is no longer singular S.L. but regular one.
Legendre equation \[ \left ( \left ( 1-x^{2}\right ) y^{\prime }\right ) ^{\prime }=-\lambda y\qquad 0<x<1 \] Or it can be written as\[ \left ( 1-x^{2}\right ) y^{\prime \prime }-2xy^{\prime }=-\lambda y \] We see that \(p=1-x^{2}\). And now it happens to be that \(p=0\) at \(x=1\). So it is singular at the other end compared to Bessel. In this case \(q=0,\sigma =1\). Again, the boundary conditions at \(x=1\) must now be bounded. i.e. \(y\left ( 1\right ) <\infty ,y^{\prime }\left ( 1\right ) <\infty \). On the other end, where \(p\) is not zero, we still use the standard boundary conditions \(\beta _{1}y\left ( 0\right ) +\beta _{2}y^{\prime }\left ( 0\right ) =0\).
Chebyshev equation\[ \left ( \sqrt{1-x^{2}}y^{\prime }\right ) ^{\prime }=-\lambda \sqrt{1-x^{2}}y\qquad -1<x<1 \] Or\[ \sqrt{1-x^{2}}y^{\prime \prime }-\frac{1}{2}\frac{2x}{\sqrt{1-x^{2}}}y^{\prime }=-\lambda \sqrt{1-x^{2}}y \] So we see that \(p=\sqrt{1-x^{2}},q=0,\sigma =\sqrt{1-x^{2}}\). So where is \(p=0\) here? At \(x=-1\) then \(p=\sqrt{1-1}=0\) and at \(x=-1\) then \(p=0\) also. So this is singular at both ends. So we need to use bounded boundary conditions at both ends now to solve this.\begin{align*} y\left ( -1\right ) & <\infty \\ y^{\prime }\left ( -1\right ) & <\infty \\ y\left ( +1\right ) & <\infty \\ y^{\prime }\left ( +1\right ) & <\infty \end{align*}
The solution to this singular S.L. is given in terms of special function Chebyshev.
\(\blacksquare \) Solution to \(xy^{\prime \prime }+y^{\prime }=-\lambda xy\). Set in S-L form, we get \(\left ( py^{\prime }\right ) ^{\prime }+qy=-\lambda \sigma y\), then \(p=x,q=0,\sigma =x\). For \(a<x<b\) with \(y\left ( b\right ) =0,y\left ( a\right ) <\infty \). This is singular S-L. Using Frobenius series, the solution comes out to be \(y\left ( x\right ) =\sum _{n=1}^{\infty }c_{n}J_{0}\left ( \sqrt{\lambda _{n}}x\right ) \), where \(J_{0}\) is Bessel function of ﬁrst kind. \(\sqrt{\lambda _{n}}\) are roots of \(0=J_{0}\left ( \sqrt{\lambda }b\right ) \). (check if this should be \(0=J_{0}\left ( \sqrt{\lambda }\left ( b-a\right ) \right ) \). \(c_{n}\) is constant, which still needs to be found.
\(\blacksquare \) We can ﬁnd approximate solution to S-L ODE for large eigenvalue, without solving the ODE using series method, but by using WKB. Given \(\left ( py^{\prime }\right ) ^{\prime }+qy=-\lambda \sigma y\) which is standard form S-L. Using physical optics correction, we obtain the solution as\[ y\left ( x\right ) \sim \left ( \sigma \left ( x\right ) p\left ( x\right ) \right ) ^{\frac{-1}{4}}\left ( c_{1}\cos \left ( \sqrt{\lambda }\int _{0}^{x}\sqrt{\frac{\sigma \left ( t\right ) }{p\left ( t\right ) }}dt\right ) +c_{1}\sin \left ( \sqrt{\lambda }\int _{0}^{x}\sqrt{\frac{\sigma \left ( t\right ) }{p\left ( t\right ) }}dt\right ) \right ) \] Where \(c_{1},c_{2}\) are found from boundary conditions now. The above is valid for ”large” \(\lambda \), and is found by ﬁrst letting \(\varepsilon ^{2}=\frac{1}{\lambda }\) and then assuming \(y\left ( x\right ) =\exp \left ( \frac{1}{\varepsilon }\sum _{n=0}^{\infty }\varepsilon ^{n}S_{n}\left ( x\right ) \right ) \). And working though the WKB method. Remember, WKB only works for linear homogeneous ODE and is used to estimate solution for large \(\lambda \) (or small \(\varepsilon \)).
\(\blacksquare \) in 1D, S-L is \(\left ( py^{\prime }\right ) ^{\prime }+qy=-\lambda \sigma y\) with \(0<x<L\), with B.C. given by \(\beta _{1}y\left ( 0\right ) +\beta _{2}y^{\prime }\left ( 0\right ) =0\) and \(\beta _{1}y\left ( L\right ) +\beta _{2}y^{\prime }\left ( L\right ) =0,\) and in higher dimensions, this problem becomes \(\nabla \cdot \left ( p\nabla \phi \right ) +q\phi =-\lambda \sigma \phi \) with boundary conditions now written as \(\beta _{1}\phi +\beta _{2}\left ( \nabla \phi \cdot \hat{n}\right ) =0\) on boundary \(\Gamma \). In higher dimensions, \(\phi \equiv \phi \left ( \bar{x}\right ) \).
\(\blacksquare \) Lagrange’s identity in 1D is \[ uL\left ( v\right ) -vL\left ( u\right ) =\frac{d}{dx}\left ( p\left ( u\frac{dv}{dx}-v\frac{du}{dx}\right ) \right ) \] where \(L\) above is the S-L operator, as in \(L\equiv \frac{d}{dx}\left ( p\frac{d}{dx}\right ) +q\) so that we write \(L\left ( u\right ) =-\lambda \sigma y\,\).
\(\blacksquare \) When it says eigenfunctions are normalized, it should mean that \(\int _{a}^{b}\phi _{n}^{2}\left ( x\right ) \sigma \left ( x\right ) \ dx=1\) where \(\sigma \left ( x\right ) \) is the weight function (comes from SL equation) and \(a<x<b\) is the domain. This is for 1D SL. Remember for example, that \(\int _{0}^{L}\cos ^{2}\left ( \sqrt{\lambda _{n}}x\right ) dx=\frac{L}{2}\). where \(\sqrt{\lambda _{n}}=\frac{n\pi }{L}\). Here, it is not normalized since the result is not one. The weight here is just one and here \(\phi _{n}\left ( x\right ) =\cos \left ( \frac{n\pi }{L}x\right ) \)
\(\blacksquare \) Heat PDE in cylinder is \(\frac{\partial T\left ( r,z,\theta ,t\right ) }{\partial t}=k\left ( \frac{\partial ^{2}T}{\partial r^{2}}+\frac{1}{r}\frac{\partial T}{\partial r}+\frac{\partial ^{2}T}{\partial z^{2}}+\frac{1}{r^{2}}\frac{\partial ^{2}T}{\partial \theta ^{2}}\right ) .\) For steady state, it becomes, say cylinder has length \(L\) \begin{align*} 0 & =\nabla ^{2}T\left ( r,z,\theta \right ) =\frac{\partial ^{2}T}{\partial r^{2}}+\frac{1}{r}\frac{\partial T}{\partial r}+\frac{\partial ^{2}T}{\partial z^{2}}+\frac{1}{r^{2}}\frac{\partial ^{2}T}{\partial \theta ^{2}}\\ \left \vert T\left ( 0,z,\theta \right ) \right \vert & <\infty \\ \left \vert T^{\prime }\left ( 0,z,\theta \right ) \right \vert & <\infty \\ T\left ( r,0\right ) & =f\left ( r\right ) \\ R\left ( r,L\right ) & =g\left ( r\right ) \end{align*}
And for axis symmetric (no \(\theta \) dependency), it becomes \[ 0=\nabla ^{2}T\left ( r,z\right ) =\frac{\partial ^{2}T}{\partial r^{2}}+\frac{1}{r}\frac{\partial T}{\partial r}+\frac{\partial ^{2}T}{\partial z^{2}}\] The solution is found by separation of variables. For the above is in terms of Bessel function order zero\[ T\left ( r,z\right ) =\sum _{n=1}^{\infty }J_{0}\left ( \sqrt{\lambda _{0}}r\right ) \left ( A_{n}e^{\sqrt{\lambda _{n}}z}+B_{n}e^{\sqrt{\lambda _{n}}z}\right ) \] which is found using series method \(T=\sum _{n=0}^{\infty }a_{n}r^{n+\alpha }\). The eigenvalues \(\lambda _{n}\) are zeros of \(0=J_{0}\left ( \sqrt{\lambda }r_{0}\right ) \) where \(r_{0}\) is disk radius.
\(\blacksquare \) Heat PDE in 2D polar is \(\frac{\partial T\left ( r,\theta ,t\right ) }{\partial t}=k\left ( \frac{\partial ^{2}T}{\partial r^{2}}+\frac{1}{r}\frac{\partial T}{\partial r}+\frac{1}{r^{2}}\frac{\partial ^{2}T}{\partial \theta ^{2}}\right ) \) and for steady state it becomes \(0=\frac{\partial ^{2}T}{\partial r^{2}}+\frac{1}{r}\frac{\partial T}{\partial r}+\frac{1}{r^{2}}\frac{\partial ^{2}T}{\partial \theta ^{2}}\). After separating, the \(r\) ODE becomes Euler ODE. Use guess \(R\left ( r\right ) =r^{p}\) for solution. The solution will be\[ T\left ( r,z\right ) =A_{0}+\sum _{n=1}^{\infty }A_{n}\cos \left ( n\theta \right ) r^{n}+\sum _{n=1}^{\infty }B_{n}\sin \left ( n\theta \right ) r^{n}\] \(\blacksquare \) To change \(p\left ( x\right ) y^{\prime \prime }+q\left ( x\right ) y^{\prime }+r\left ( x\right ) y\left ( x\right ) =0\) to S.L. form, multiply the ODE by \(\mu =\frac{1}{p\left ( x\right ) }\exp \left ( \int \frac{q\left ( x\right ) }{p\left ( x\right ) }dx\right ) \)
\(\blacksquare \) In using method of eigenfunction expansion to solve nonhomogeneous B.C. PDE, the eigenfunction used \(\phi _{n}\left ( x\right ) \) are the ones from solving the PDE with homogeneous B.C. For example, if given \(\frac{\partial u}{\partial t}=k\frac{\partial ^{2}u}{\partial x^{2}}\) with nonhomogeneous B.C. , say \(u\left ( 0,t\right ) =A,u\left ( L,t\right ) =B\), and we want to use \(u\left ( x,t\right ) \sim \sum _{n=1}^{\infty }b_{n}\left ( t\right ) \phi _{n}\left ( x\right ) \) on the nonhomogeneous B.C. PDE, then those \(\phi _{n}\left ( x\right ) \) are from the corresponding homogeneous B.C. PDE. They should be something like \(\sin \left ( \frac{n\pi }{L}x\right ) \) and so on.
See the table at the top. That is why we write \(\sim \) above and not \(=\) and remember not to do term-by-term diﬀerentiation in \(x\) when this is the case. We can do term-by-term diﬀerentiation only of the PDE itself also has homogeneous B.C. , even if it had a source term there also. The point is, \(\phi _{n}\left ( x\right ) \) always come from solving the homogeneous B.C. PDE. (This normally means solving a Sturm-Liouville ODE).
\(\blacksquare \) I found new relation for eigenvalue \(\lambda \), but it is probably not useful as the one the book has. Here
how to derive it
Given S-L \begin{equation} \frac{d}{dx}\left ( p\frac{dy}{dx}\right ) +qy=-\lambda \sigma y \tag{1} \end{equation}
Let \(\phi _{n}\left ( x\right ) \) be the eigenfunction associated with eigenvalue \(\lambda _{n}\). Since eigenfunctions satisfy the
ODE itself, then we can write, for any arbitrary eigenfunction (subscript removed for clarity in
what follows)\begin{align*} \frac{d}{dx}\left ( p\phi ^{\prime }\right ) +q\phi & =-\lambda \sigma \phi \\ p^{\prime }\phi ^{\prime }+p\phi ^{\prime \prime }+q\phi & =-\lambda \sigma \phi \end{align*}
Integrating both sides\begin{equation} \int _{a}^{b}p^{\prime }\phi ^{\prime }dx+\int _{a}^{b}p\phi ^{\prime \prime }dx+\int _{a}^{b}q\phi dx=-\lambda \int _{a}^{b}\sigma \phi dx \tag{2} \end{equation} Looking at \(\int _{a}^{b}p\phi ^{\prime \prime }dx\). Integrating by part. Let \(u=p,dv=\phi ^{\prime \prime }\rightarrow du=p^{\prime },v=\phi ^{\prime }\), hence\begin{align} \int _{a}^{b}p\phi ^{\prime \prime }dx & =\left [ uv\right ] _{a}^{b}-\int _{a}^{b}vdu\nonumber \\ & =\left [ p\phi ^{\prime }\right ] _{a}^{b}-\int _{a}^{b}p^{\prime }\phi ^{\prime }dx \tag{3} \end{align}
Substituting (3) into (2) gives\begin{align*} \int _{a}^{b}p^{\prime }\phi ^{\prime }dx+\left [ p\phi ^{\prime }\right ] _{0}^{L}-\int _{a}^{b}p^{\prime }\phi ^{\prime }dx+\int _{a}^{b}q\phi dx & =-\lambda \int _{a}^{b}\sigma \phi dx\\ \left [ p\phi ^{\prime }\right ] _{0}^{L}+\int _{a}^{b}q\phi dx & =-\lambda \int _{a}^{b}\sigma \phi dx \end{align*}
Hence\[ \lambda =-\frac{\left [ p\phi ^{\prime }\right ] _{a}^{b}+\int _{a}^{b}q\phi dx}{\int _{a}^{b}\sigma \phi dx}\] compare to the one in the book \(\lambda =\frac{\left ( -p\phi \phi ^{\prime }\right ) _{b}^{a}+\int _{a}^{b}p\left ( \phi ^{\prime }\right ) ^{2}-q\phi ^{2}dx}{\int _{a}^{b}\phi ^{2}\sigma \left ( x\right ) dx}\).
\(\blacksquare \) This is about using eigenfunction expansion to approximate \(f\left ( x\right ) \). Given the eigenfunctions \(\Phi _{n}\left ( x\right ) \) found by solving SL \(L\left [ y\right ] =\lambda ry\) problem, then these eigenfunctions are called complete in the mean sense, if \(\sum _{n=1}^{\infty }b_{n}\Phi _{n}\left ( x\right ) \) converges to the mean of \(f\left ( x\right ) \) at each point in \(a\leq x\leq b\).
This requires only that \(f\left ( x\right ) \) be square integrable over the domain. We can also say that the eigenfunctions are complete, in pointwise sense, if \(\sum _{n=1}^{\infty }b_{n}\Phi _{n}\left ( x\right ) \) converges to \(f\left ( x\right ) \) at each point in \(a\leq x\leq b\). This however requires that \(f\left ( x\right ) \) be continuous in \(a\leq x\leq b\) and \(f^{\prime }\left ( x\right ) \) be piecewise continuous.
So clearly the ﬁrst case is less strict than the second. This means if \(f\left ( x\right ) \) is only square integrable (and not necessarily continuos), then convergence to the mean of \(f\left ( x\right ) \) at each point can be obtained. To obtain the more strict point wise convergence, more restrictions are needed on \(f\left ( x\right ) \) as mentioned above. Clearly point wise convergence implies convergence in the mean. But not the other way around.
Problem 1 If the problem gives S-L equations, and asks to ﬁnd estimate on the smallest eigenvalue, then use Rayleigh quotient for \(\lambda \). And write \(\lambda _{\min }=\lambda _{1}\leq \frac{\left ( -p\phi \phi ^{\prime }\right ) _{b}^{a}+\int _{a}^{b}p\left ( \phi ^{\prime }\right ) ^{2}-q\phi ^{2}dx}{\int _{a}^{b}\phi ^{2}\sigma \left ( x\right ) dx}\). Now we do not need to solve the SL to ﬁnd solution \(\phi \), this is the whole point. Everything in RHS is given, except for, of course the solution \(\phi \).
Here comes the main idea: Come up with any trial \(\phi _{trial}\) and use it in place of \(\phi \). This trial function just needs to satisfy the boundary conditions, which is also given. Then all what we need to do is just evaluate the integral. Pick the simplest \(\phi _{trial}\) function of \(x\) which satisﬁes boundary conditions. All other terms \(p,q,\sigma \) we can read from the given problem.
At the end, we should get a numerical value for the integral. This is the upper limit of the lowest eigenvalue \(\lambda _{1}\)
Problem 2 We are given SL problem with boundary conditions, and asked to show that \(\lambda \geq 0\) without solving the ODE to ﬁnd \(\phi \):
Use Rayleigh quotient and argue that the denominator can’t be zero (else eigenfunction is zero, which is not possible) and it also can’t be negative. Argue about the numerator not being negative. Do not solve the ODE! this is the whole point of using Rayleigh quotient.
Problem We are given an SL problem with boundary conditions and asked to estimate large \(\lambda \) and corresponding eigenfunction.
This is diﬀerent from being asked to estimate the smallest eigenvalue, where we use Rayleigh quotient and trial function. In this one, we instead use WKB. For 1D, just use what is called the physical optics correction method, given by \(\phi \left ( x\right ) \approx \left ( \sigma p\right ) ^{\frac{-1}{4}}\left ( c_{1}\cos \left ( \sqrt{\lambda }\int _{0}^{x}\sqrt{\frac{\sigma \left ( t\right ) }{p\left ( t\right ) }}dt\right ) +c_{2}\sin \left ( \sqrt{\lambda }\int _{0}^{x}\sqrt{\frac{\sigma \left ( t\right ) }{p\left ( t\right ) }}dt\right ) \right ) \). Where \(\sigma ,p\) in this are all known functions (from the problem itself). Notice that \(q\) is not there.
Now use the boundary conditions and solve for one of the constants, should come to be zero. Use the second boundary condition and (remember, the boundary conditions are homogenous), and we get one equation in \(\lambda \) and for non-trivial solution, solve for allowed values of \(\lambda _{n}\). This gives the large eigenvalue estimate. i.e. for \(\lambda _{n}\) when \(n\) is very large. Depending on problem, \(n\) does not have be too large to get good accurate estimate compared with exact solution. See HW7, problem 5.9.2 for example.
Problem We are given 1D heat ODE \(\frac{\partial u}{\partial t}=k\frac{\partial ^{2}u}{\partial x^{2}}\). If B.C. are homogenous, we are done. We know the solution. Separation of variables. If the B.C. is not homogenous, then we have a small problem. We can’t do separation of variables. If you can ﬁnd equilibrium solution, \(u_{r}\left ( x\right ) \), where this solution only needs to satisfy the non-homogenous B.C. then write solution as \(u\left ( x,t\right ) =v\left ( x,t\right ) +u_{r}\left ( x\right ) \), and plug this in the PDE.
This will give \(\frac{\partial v}{\partial t}=k\frac{\partial ^{2}v}{\partial x^{2}}\) but this now satisﬁes the homogenous B.C. This is the case if the original non-homogenous B.C. were Dirichlet. If they were Neumann, then we will get \(\frac{\partial v}{\partial t}=k\frac{\partial ^{2}v}{\partial x^{2}}+extra\) where \(extra\) is extra term that do not vanish because \(u_{r}\left ( x\right ) \) do not vanish in this case after double diﬀerentiating. We solve for \(v\left ( x,t\right ) \) now, since it has homogenous B.C., but it has extra term there, which we treat as new source. We apply initial conditions now also to ﬁnal all the eigenfunction expansion coeﬃcients. We are done. Now we know \(u\left ( x,t\right ) .\)
But if we can ﬁnd the reference function \(u_{r}\) or we do not want to use this method, then we can use another method, called eigenfunction expansion. We assume \(u\left ( x,t\right ) =\sum _{n=1}^{\infty }a_{n}\left ( t\right ) \phi _{n}\left ( x\right ) \) and plug this in the PDE. But now we can’t do term by term diﬀerentiation, since \(\phi _{n}\) are the eigenfunctions that satisﬁes the homogenous B.C., while \(u\left ( x,t\right ) \) has non-homogenous B.C.
So the trick is to use Green formula, to rewrite \(\int _{0}^{L}k\frac{\partial ^{2}u}{\partial x^{2}}\phi _{n}\left ( x\right ) dx\) as \(\int _{0}^{L}ku\frac{\partial ^{2}\phi _{n}\left ( x\right ) }{\partial x^{2}}dx\) plus the contribution from the boundary terms. (this is like doing integration by parts twice, but much easier). Now rewrite \(\frac{\partial ^{2}\phi _{n}\left ( x\right ) }{\partial x^{2}}=-\lambda _{n}\phi _{n}\left ( x\right ) \). See page 355-356, Haberman.
Any linear second order ODE can be converted to S-L form. The S-L form is the following\begin{equation} \left ( p\left ( x\right ) \frac{dy}{dx}\right ) ^{\prime }-q\left ( x\right ) y\left ( x\right ) +\lambda r\left ( x\right ) y\left ( x\right ) =0 \tag{1} \end{equation} Some books write the above with a plus sign \(+q\left ( x\right ) y\left ( x\right ) \) there (for example, Haberman) and some with \(-q\left ( x\right ) y\left ( x\right ) \) for example Boyce/DiPrima. I really never understood why some put a minus sign and some do not. May be I’ll ﬁnd out one day. But for now will use (1) is used as the S-L form. If your book uses the plus sign, the same process works, just ﬂip the sign.
The goal now, is given any general form eigenvalue ODE or second order ODE (with no eigenvalue \(\lambda \)), we want to convert it (rewrite it) in the above form. The second order linear ODE will have this form\begin{equation} a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+\left ( c\left ( x\right ) +\lambda \right ) y=0 \tag{2} \end{equation} The parameter \(\lambda \) in the S-L form, is the eigenvalue. We really only use S-L form for eigenvalue problems, but if the \(\lambda \) is missing, then the form of the input will be \[ a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+c\left ( x\right ) y=0 \] And the process will work the same way.
First, we will show how to convert (2) to (1), and then show few examples. We are given (2), and want to convert it to (1). The ﬁrst step is to convert (2) to standard form\[ y^{\prime \prime }+\frac{b\left ( x\right ) }{a\left ( x\right ) }y^{\prime }+\frac{c\left ( x\right ) +\lambda }{a\left ( x\right ) }y=0 \] Then multiply the above by some unknown \(p\left ( x\right ) \) which is assumed positive. This is called the integrating factor\[ p\left ( x\right ) y^{\prime \prime }+p\left ( x\right ) \frac{b\left ( x\right ) }{a\left ( x\right ) }y^{\prime }+p\left ( x\right ) \frac{c\left ( x\right ) +\lambda }{a\left ( x\right ) }y=0 \] Now, rewrite \(p\left ( x\right ) y^{\prime \prime }\) as \(\left ( p\left ( x\right ) y^{\prime }\right ) ^{\prime }-p^{\prime }\left ( x\right ) y^{\prime }\) in the above\begin{align} \left ( \left ( p\left ( x\right ) y^{\prime }\right ) ^{\prime }-p^{\prime }\left ( x\right ) y^{\prime }\right ) +p\left ( x\right ) \frac{b\left ( x\right ) }{a\left ( x\right ) }y^{\prime }+p\left ( x\right ) \frac{\left ( c\left ( x\right ) +\lambda \right ) }{a\left ( x\right ) }y & =0\nonumber \\ \left ( p\left ( x\right ) y^{\prime }\right ) ^{\prime }+y^{\prime }\left ( p\left ( x\right ) \frac{b\left ( x\right ) }{a\left ( x\right ) }-p^{\prime }\left ( x\right ) \right ) +p\left ( x\right ) \frac{\left ( c\left ( x\right ) +\lambda \right ) }{a\left ( x\right ) }y & =0 \tag{3} \end{align}
Here comes the main trick in all of this. We want to force \(p\left ( x\right ) \frac{b\left ( x\right ) }{a\left ( x\right ) }-p^{\prime }\left ( x\right ) \) to be zero. This implies \[ p\left ( x\right ) =e^{\int \frac{b\left ( x\right ) }{a\left ( x\right ) }dx}\] This comes from just solving the ODE \(p\left ( x\right ) \frac{b\left ( x\right ) }{a\left ( x\right ) }-p^{\prime }\left ( x\right ) =0\). Therefore, if \(p\left ( x\right ) =e^{\int \frac{b\left ( x\right ) }{a\left ( x\right ) }dx}\), then (3) becomes\begin{align} \left ( p\left ( x\right ) y^{\prime }\right ) ^{\prime }+\left ( p\left ( x\right ) \frac{c\left ( x\right ) +\lambda }{a\left ( x\right ) }\right ) y & =0\nonumber \\ \left ( p\left ( x\right ) y^{\prime }\right ) ^{\prime }+\left ( p\left ( x\right ) \frac{c\left ( x\right ) }{a\left ( x\right ) }+\frac{p\left ( x\right ) }{a\left ( x\right ) }\lambda \right ) y & =0\nonumber \\ \left ( p\left ( x\right ) y^{\prime }\right ) ^{\prime }+p\left ( x\right ) \frac{c\left ( x\right ) }{a\left ( x\right ) }y+\frac{p\left ( x\right ) }{a\left ( x\right ) }\lambda y & =0 \tag{4} \end{align}
We are done. Compare (4) to the SL form\[ \left ( p\left ( x\right ) \frac{dy}{dx}\right ) ^{\prime }-q\left ( x\right ) y\left ( x\right ) +\lambda r\left ( x\right ) y\left ( x\right ) =0 \] We see that, given \(a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+\left ( c\left ( x\right ) +\lambda d\right ) y=0\) then
\begin{align*} p\left ( x\right ) & =e^{\int \frac{b}{a}dx}\\ q\left ( x\right ) & =-p\frac{c}{a}\\ r\left ( x\right ) & =\frac{p}{a}d \end{align*}
Let see how this works on some examples
Convert \(y^{\prime \prime }+5\lambda y=0\) to S-L. We see the general form here is \(a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+\left ( c\left ( x\right ) +\lambda d\right ) y=0\), with \(a=1,b=0,c=0,d=5\). Hence \(p\left ( x\right ) =e^{\int \frac{b\left ( x\right ) }{a\left ( x\right ) }dx}=1\), therefore the SL form is\begin{equation} \left ( p\left ( x\right ) \frac{dy}{dx}\right ) ^{\prime }-q\left ( x\right ) y\left ( x\right ) +\lambda r\left ( x\right ) y\left ( x\right ) =0 \tag{1} \end{equation} Where \begin{align*} p\left ( x\right ) & =1\\ q\left ( x\right ) & =-p\left ( x\right ) \frac{c\left ( x\right ) }{a\left ( x\right ) }=0\\ r\left ( x\right ) & =\frac{p\left ( x\right ) }{a\left ( x\right ) }d=5 \end{align*}
Hence (1) becomes\[ \left ( y^{\prime }\right ) ^{\prime }+5\lambda y=0 \] This was easy, since the ODE given was already in SL form.
Convert Bessel ODE to SL \[ x^{2}y^{\prime \prime }+xy^{\prime }+\left ( x^{2}-n^{2}\right ) y=0 \]
Comparing the above to \(a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+\left ( c\left ( x\right ) +\lambda \right ) y=0\), we see that \(a=x^{2},b=x,c=x^{2},d=1,\lambda =-n^{2}\). Hence \[ p\left ( x\right ) =e^{\int \frac{b\left ( x\right ) }{a\left ( x\right ) }dx}=e^{\int \frac{x}{x^{2}}dx}=e^{\ln x}=x \] Therefore the SL form is\begin{equation} \left ( p\left ( x\right ) \frac{dy}{dx}\right ) ^{\prime }-q\left ( x\right ) y\left ( x\right ) +\lambda r\left ( x\right ) y\left ( x\right ) =0 \tag{1} \end{equation} Where \begin{align*} p\left ( x\right ) & =x\\ q\left ( x\right ) & =-p\left ( x\right ) \frac{c\left ( x\right ) }{a\left ( x\right ) }=-x\left ( \frac{x^{2}}{x^{2}}\right ) =-x\\ r\left ( x\right ) & =\frac{p\left ( x\right ) }{a\left ( x\right ) }d=\frac{1}{x} \end{align*}
Therefore (1) becomes\begin{equation} \left ( xy^{\prime }\right ) ^{\prime }+xy-\frac{n^{2}}{x}y\left ( x\right ) =0 \tag{2} \end{equation} Let check if the above gives back \(x^{2}y^{\prime \prime }+xy^{\prime }+\left ( x^{2}-n^{2}\right ) y=0\). Expanding the above gives \[ xy^{\prime \prime }+y^{\prime }+xy-\frac{n^{2}}{x}y=0 \] Multiplying by \(x\)\begin{align*} x^{2}y^{\prime \prime }+xy^{\prime }+x^{2}y-n^{2}y & =0\\ x^{2}y^{\prime \prime }+xy^{\prime }+\left ( x^{2}-n^{2}\right ) y & =0 \end{align*}
This veriﬁes that (2) is the Sturm-Liouville form for Bessel ODE
Convert Legendre ODE to SL\[ \left ( 1-x^{2}\right ) y^{\prime \prime }-2xy^{\prime }+n\left ( n+1\right ) y=0 \] The general form here is \(a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+\left ( c\left ( x\right ) +\lambda \right ) y=0\), with \(a=\left ( 1-x^{2}\right ) ,b=-2x,c=0,d=1\) and \(\lambda =n\left ( n+1\right ) \). Hence \begin{align*} p\left ( x\right ) & =e^{\int \frac{-2x}{\left ( 1-x^{2}\right ) }dx}=e^{\ln \left ( x^{2}-1\right ) }=x^{2}-1\\ q\left ( x\right ) & =-p\frac{c}{a}=0\\ r\left ( x\right ) & =\frac{p}{a}d=\frac{x^{2}-1}{1-x^{2}}=-1 \end{align*}
Hence SL is\begin{align*} \left ( py^{\prime }\right ) ^{\prime }-qy+\lambda ry & =0\\ \left ( \left ( x^{2}-1\right ) y^{\prime }\right ) ^{\prime }-n\left ( n+1\right ) y & =0 \end{align*}
Let check. Expanding the above gives\[ \left ( x^{2}-1\right ) y^{\prime \prime }+2xy^{\prime }-n\left ( n+1\right ) y=0 \] Multiplying both sides by \(-1\) gives \[ \left ( 1-x^{2}\right ) y^{\prime \prime }-2xy^{\prime }+n\left ( n+1\right ) y=0 \] Which is the original form of Legendre ODE
Convert \[ y^{\prime \prime }-\frac{2x}{1-x^{2}}y^{\prime }+\frac{\lambda }{1-x^{2}}y=0 \] Rewrite as\[ \left ( 1-x^{2}\right ) y^{\prime \prime }-2xy^{\prime }+\lambda y=0 \] The general form here is \(a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+\left ( c\left ( x\right ) +\lambda \right ) y=0\), with \(a=\left ( 1-x^{2}\right ) ,b=-2x,c=0,d=1\). Hence \[ \mu \left ( x\right ) =e^{\int \frac{b}{a}dx}=e^{-\int \frac{2x}{1-x^{2}}dx}=e^{-\left ( -\ln \left ( x-1\right ) -\ln \left ( x+1\right ) \right ) }=e^{\ln \left ( x-1\right ) }e^{\ln \left ( x+1\right ) }=\left ( x-1\right ) \left ( x+1\right ) =x^{2}-1 \] Therefore the SL form is\begin{equation} \left ( p\left ( x\right ) \frac{dy}{dx}\right ) ^{\prime }-q\left ( x\right ) y\left ( x\right ) +\lambda r\left ( x\right ) y\left ( x\right ) =0 \tag{1} \end{equation} Where \begin{align*} p\left ( x\right ) & =x^{2}-1\\ q\left ( x\right ) & =-p\left ( x\right ) \frac{c\left ( x\right ) }{a\left ( x\right ) }=0\\ r\left ( x\right ) & =\frac{p\left ( x\right ) }{a\left ( x\right ) }d=\frac{x^{2}-1}{x^{2}-1}=1 \end{align*}
Hence (1) becomes\[ \left ( \left ( x^{2}-1\right ) \frac{dy}{dx}\right ) ^{\prime }+\lambda y\left ( x\right ) =0 \]
Convert \[ \left ( 1-x^{2}\right ) y^{\prime \prime }-xy^{\prime }+\lambda y=0 \] We see the general form here is \(a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+\left ( c\left ( x\right ) +\lambda \right ) y=0\), with \(a=\left ( 1-x^{2}\right ) ,b=-x,c=0,d=1\). Hence \[ p\left ( x\right ) =e^{\int \frac{b}{a}dx}=e^{-\int \frac{x}{1-x^{2}}dx}\] But \(\int \frac{x}{1-x^{2}}dx=-\frac{1}{2}\ln \left \vert 1-x^{2}\right \vert \), therefore \(\mu \left ( x\right ) =e^{\frac{1}{2}\ln \left \vert 1-x^{2}\right \vert }=\sqrt{1-x^{2}}\). Therefore the SL form is\begin{equation} \left ( p\left ( x\right ) \frac{dy}{dx}\right ) ^{\prime }-q\left ( x\right ) y\left ( x\right ) +\lambda r\left ( x\right ) y\left ( x\right ) =0 \tag{1} \end{equation} Where \begin{align*} p\left ( x\right ) & =\sqrt{1-x^{2}}\\ q\left ( x\right ) & =-p\left ( x\right ) \frac{c\left ( x\right ) }{a\left ( x\right ) }=0\\ r\left ( x\right ) & =\frac{p\left ( x\right ) }{a\left ( x\right ) }d=\frac{\sqrt{1-x^{2}}}{\left ( 1-x^{2}\right ) }=\frac{1}{\sqrt{1-x^{2}}} \end{align*}
Hence (1) becomes\[ \left ( \sqrt{1-x^{2}}\frac{dy}{dx}\right ) ^{\prime }+\frac{\lambda }{\sqrt{1-x^{2}}}y\left ( x\right ) =0 \]
Convert\[ a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+c\left ( x\right ) y=0 \] The general form here is \(a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+\left ( c\left ( x\right ) +\lambda \right ) y=0\), with \(a=a\left ( x\right ) ,b=b\left ( x\right ) ,c=c\left ( x\right ) ,d=1\) and \(\lambda =0\). Hence \[ p\left ( x\right ) =e^{\int \frac{b}{a}dx}\] Therefore the SL form is\begin{equation} \left ( p\left ( x\right ) \frac{dy}{dx}\right ) ^{\prime }-q\left ( x\right ) y\left ( x\right ) +\lambda r\left ( x\right ) y\left ( x\right ) =0 \tag{1} \end{equation} Where \begin{align*} p\left ( x\right ) & =e^{\int \frac{b}{a}dx}\\ q\left ( x\right ) & =-p\left ( x\right ) \frac{c\left ( x\right ) }{a\left ( x\right ) }=-e^{\int \frac{b}{a}dx}\frac{c\left ( x\right ) }{a\left ( x\right ) }\\ r\left ( x\right ) & =\frac{p\left ( x\right ) }{a\left ( x\right ) }d=e^{\int \frac{b}{a}dx} \end{align*}
Hence (1) becomes\begin{equation} \left ( e^{\int \frac{b}{a}dx}y^{\prime }\right ) ^{\prime }+e^{\int \frac{b}{a}dx}\frac{c\left ( x\right ) }{a\left ( x\right ) }y\left ( x\right ) =0 \tag{2} \end{equation} We see that the above is the same as \(a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+c\left ( x\right ) y=0\) because\[ \left ( e^{\int \frac{b}{a}dx}y^{\prime }\right ) ^{\prime }=\frac{b}{a}e^{\int \frac{b}{a}dx}y^{\prime }+e^{\int \frac{b}{a}dx}y^{\prime \prime }\] And hence (2) becomes\[ \frac{b}{a}e^{\int \frac{b}{a}dx}y^{\prime }+e^{\int \frac{b}{a}dx}y^{\prime \prime }+e^{\int \frac{b}{a}dx}\frac{c\left ( x\right ) }{a\left ( x\right ) }y\left ( x\right ) =0 \] Canceling\(\ e^{\int \frac{b}{a}dx}\) gives\begin{align*} \frac{b}{a}y^{\prime }+y^{\prime \prime }+\frac{c\left ( x\right ) }{a\left ( x\right ) }y\left ( x\right ) & =0\\ ay^{\prime \prime }+by^{\prime }+cy & =0 \end{align*}
Which is the original ODE.
Convert\[ 3y^{\prime \prime }+2y^{\prime }+5y=0 \] The general form here is \(a\left ( x\right ) y^{\prime \prime }+b\left ( x\right ) y^{\prime }+\left ( c\left ( x\right ) +\lambda \right ) y=0\), with \(a=3,b=2,c=5,d=1\) and \(\lambda =0\). Hence \[ p\left ( x\right ) =e^{\int \frac{b}{a}dx}=e^{\int \frac{2}{3}dx}=e^{\frac{2}{3}x}\] Therefore the SL form is\begin{equation} \left ( p\left ( x\right ) \frac{dy}{dx}\right ) ^{\prime }-q\left ( x\right ) y\left ( x\right ) +\lambda r\left ( x\right ) y\left ( x\right ) =0 \tag{1} \end{equation} Where \begin{align*} p\left ( x\right ) & =e^{\frac{2}{3}x}\\ q\left ( x\right ) & =-p\left ( x\right ) \frac{c\left ( x\right ) }{a\left ( x\right ) }=-e^{\frac{2}{3}x}\left ( \frac{5}{2}\right ) \\ r\left ( x\right ) & =\frac{p\left ( x\right ) }{a\left ( x\right ) }d=\frac{e^{\frac{2}{3}x}}{3} \end{align*}
Hence (1) becomes\begin{equation} \left ( e^{\frac{2}{3}x}y^{\prime }\right ) ^{\prime }+e^{\frac{2}{3}x}\left ( \frac{5}{2}\right ) y\left ( x\right ) =0 \tag{2} \end{equation} The above is the same as \(3y^{\prime \prime }+2y^{\prime }+5y=0\) because\[ \left ( e^{\frac{2}{3}x}y^{\prime }\right ) ^{\prime }=\frac{2}{3}e^{\frac{2}{3}x}y^{\prime }+e^{\frac{2}{3}x}y^{\prime \prime }\] And hence (2) becomes\[ \frac{2}{3}e^{\frac{2}{3}x}y^{\prime }+e^{\frac{2}{3}x}y^{\prime \prime }+e^{\frac{2}{3}x}\left ( \frac{5}{2}\right ) y\left ( x\right ) =0 \] Canceling\(\ e^{\frac{2}{3}x}\) gives\begin{align*} \frac{2}{3}y^{\prime }+y^{\prime \prime }+\frac{5}{2}y\left ( x\right ) & =0\\ 2y^{\prime \prime }+3y^{\prime }+5y & =0 \end{align*}
Which is the original ODE.
\(\blacksquare \) Green function is what is called impulse response in control. But it is more general, and can be used for solving PDE also.
Given a diﬀerential equation with some forcing function on the right side. To solve this, we replace the forcing function with an impulse. The solution of the DE now is called the impulse response, which is the Green’s function of the diﬀerential equation.
Now to ﬁnd the solution to the original problem with the original forcing function, we just convolve the Green function with the original forcing function. Here is an example. Suppose we want to solve \(L\left [ y\left ( t\right ) \right ] =f\left ( t\right ) \) with zero initial conditions. Then we solve \(L\left [ g\left ( t\right ) \right ] =\delta \left ( t\right ) \). The solution is \(g\left ( t\right ) \). Now \(y\left ( t\right ) =g\left ( t\right ) \circledast f\left ( t\right ) \). This is for initial value problem. For example. \(y^{\prime }\left ( t\right ) +kx=e^{at}\), with \(y\left ( 0\right ) =0\). Then we solve \(g^{\prime }\left ( t\right ) +kg=\delta \left ( t\right ) \). The solution is \(g\left ( t\right ) =\left \{ \begin{array} [c]{cc}e^{-kt} & t>0\\ 0 & t<0 \end{array} \right . \), this is for causal system. Hence \(y\left ( t\right ) =g\left ( t\right ) \circledast f\left ( t\right ) \). The nice thing here, is that once we ﬁnd \(g\left ( t\right ) \), we can solve \(y^{\prime }\left ( t\right ) +kx=f\left ( t\right ) \) for any \(f\left ( t\right ) \) by just convolving the Green function (impulse response) with the new \(f\left ( t\right ) \).
\(\blacksquare \) We can think of Green function as an inverse operator. Given \(L\left [ y\left ( t\right ) \right ] =f\left ( t\right ) \), we want to ﬁnd solution \(y\left ( t\right ) =\int _{-\infty }^{\infty }G\left ( t;\tau \right ) f\left ( \tau \right ) d\tau \). So in a sense, \(G\left ( t;\tau \right ) \) is like \(L^{-1}\left [ y\left ( t\right ) \right ] \).
\(\blacksquare \) Need to add notes for Green function for Sturm-Liouville boundary value ODE. Need to be clear on what boundary conditions to use. What is B.C. is not homogeneous?
\(\blacksquare \) Green function properties:
\(\blacksquare \) When solving for \(G\left ( t;\tau \right ) \), in context of 1D, hence two boundary conditions, one at each end, and second order ODE (Sturm-Liouville), we now get two solutions, one for \(t<\tau \) and one for \(t>\tau \).
So we have \(4\) constants of integrations to ﬁnd (this is for second order ODE) not just two constants as normally one would get , since now we have 2 diﬀerent solutions. Two of these constants from the two boundary conditions, and two more come from property of Green function as mentioned above. \(G\left ( t;\tau \right ) =\left \{ \begin{array} [c]{cc}A_{1}y_{1}+A_{2}y_{2} & 0<t<\tau \\ A_{3}y_{1}+A_{4}y_{2} & \tau <t<L \end{array} \right . \)
\(\blacksquare \) Remember that \(u_{c}\left ( t\right ) f\left ( t-c\right ) \Longleftrightarrow e^{-cs}F\left ( s\right ) \) and \(u_{c}\left ( t\right ) f\left ( t\right ) \Longleftrightarrow e^{-cs}\mathcal{L}\left \{ f\left ( t+c\right ) \right \} \). For example, if we are given \(u_{2}\left ( t\right ) t\), then \(\mathcal{L}\left ( u_{2}\left ( t\right ) t\right ) =e^{-2s}\mathcal{L}\left \{ t+2\right \} =e^{-2s}\left ( \frac{1}{s^{2}}+\frac{2}{s}\right ) =e^{-2s}\left ( \frac{1+2s}{s^{2}}\right ) \). Do not do \(u_{c}\left ( t\right ) f\left ( t\right ) \Longleftrightarrow e^{-cs}\mathcal{L}\left \{ f\left ( t\right ) \right \} \) ! That will be a big error. We use this allot when asked to write a piecewise function using Heaviside functions.
\(\blacksquare \) if we have a function \(f\left ( x\right ) \) represented as series (say power series or Fourier series), then we say the series converges to \(f\left ( x\right ) \) uniformly in region \(D\), if given \(\varepsilon >0\), we can number \(N\) which depends only on \(\varepsilon \), such that \(\left \vert f\left ( x\right ) -S_{N}\left ( x\right ) \right \vert <\varepsilon \).
Where here \(S_{N}\left ( x\right ) \) is the partial sum of the series using \(N\) terms. The diﬀerence between uniform convergence and non-uniform convergence, is that with uniform the number \(N\) only depends on \(\varepsilon \) and not on which \(x\) we are trying to approximate \(f\left ( x\right ) \) at. In uniform convergence, the number \(N\) depends on both \(\varepsilon \) and \(x\). So this means at some locations in \(D\) we need much larger \(N\) than in other locations to convergence to \(f\left ( x\right ) \) with same accuracy. Uniform convergence is better. It depends on the basis functions used to approximate \(f\left ( x\right ) \) in the series.
If the function \(f\left ( x\right ) \) is discontinuous at some point, then it is not possible to ﬁnd uniform convergence there. As we get closer and closer to the discontinuity, more and more terms are needed to obtained same approximation away from the discontinuity, hence not uniform convergence. For example, Fourier series approximation of a step function can not be uniformly convergent due to the discontinuity in the step function.
\(\blacksquare \) \(\ \)Geometric series: \begin{align*} \sum _{n=0}^{N}r^{n} & =1+r+r^{2}+r^{3}+\cdots +r^{N}=\frac{1-r^{N-1}}{1-r}\\ \sum _{n=0}^{\infty }r^{n} & =1+r+r^{2}+r^{3}+\cdots =\frac{1}{1-r}\qquad \left \vert r\right \vert <1\\ \sum _{n=0}^{\infty }\left ( -1\right ) ^{n}r^{n} & =1-r+r^{2}-r^{3}+\cdots =\frac{1}{1+r}\qquad \left \vert r\right \vert <1 \end{align*}
\(\blacksquare \) \(\ \)Binomial series:
General binomial is\[ \left ( x+y\right ) ^{n}=x^{n}+nx^{n-1}y+\frac{n\left ( n-1\right ) }{2!}x^{n-2}y^{2}+\frac{n\left ( n-1\right ) \left ( n-2\right ) }{3!}x^{n-3}y^{3}+\cdots \] From the above we can generate all other special cases. For example, \[ \left ( 1+x\right ) ^{n}=1+nx+\frac{n\left ( n-1\right ) x^{2}}{2!}+\frac{n\left ( n-1\right ) \left ( n-2\right ) x^{3}}{3!}+\cdots \] This work for positive and negative \(n\), rational or not. The sum converges when only for \(\left \vert x\right \vert <1\). From this, we can derive the above sums also for the geometric series. For example, for \(n=-1\) the above becomes\begin{align*} \frac{1}{\left ( 1+x\right ) } & =1-x+x^{2}-x^{3}+\cdots \qquad \left \vert x\right \vert <1\\ \frac{1}{\left ( 1-x\right ) } & =1+x+x^{2}+x^{3}+\cdots \qquad \left \vert x\right \vert <1 \end{align*}
For \(\left \vert x\right \vert >1\), we can still ﬁnd series expansion in negative powers of \(x\) as follows\begin{align*} \left ( 1+x\right ) ^{n} & =\left ( x\left ( 1+\frac{1}{x}\right ) \right ) ^{n}\\ & =x^{n}\left ( 1+\frac{1}{x}\right ) ^{n} \end{align*}
And now since \(\left \vert \frac{1}{x}\right \vert <1\), we can use binomial expansion to expand the term \(\left ( 1+\frac{1}{x}\right ) ^{n}\) in the above and obtain a convergent series, since now \(\left \vert \frac{1}{x}\right \vert <1\,.\) This will give the following expansion\begin{align*} \left ( 1+x\right ) ^{n} & =x^{n}\left ( 1+\frac{1}{x}\right ) ^{n}\\ & =x^{n}\left ( 1+n\left ( \frac{1}{x}\right ) +\frac{n\left ( n-1\right ) }{2!}\left ( \frac{1}{x}\right ) ^{2}+\frac{n\left ( n-1\right ) \left ( n-2\right ) }{3!}\left ( \frac{1}{x}\right ) ^{3}+\cdots \right ) \end{align*}
So everything is the same, we just change \(x\) with \(\frac{1}{x}\) and remember to multiply the whole expansion with \(x^{n}\). For example, for \(n=-1\)\begin{align*} \frac{1}{\left ( 1+x\right ) } & =\frac{1}{x\left ( 1+\frac{1}{x}\right ) }=\frac{1}{x}\left ( 1-\frac{1}{x}+\left ( \frac{1}{x}\right ) ^{2}-\left ( \frac{1}{x}\right ) ^{3}+\cdots \right ) \qquad \left \vert x\right \vert >1\\ \frac{1}{\left ( 1-x\right ) } & =\frac{1}{x\left ( 1-\frac{1}{x}\right ) }=\frac{1}{x}\left ( 1+\frac{1}{x}+\left ( \frac{1}{x}\right ) ^{2}+\left ( \frac{1}{x}\right ) ^{3}+\cdots \right ) \qquad \left \vert x\right \vert >1 \end{align*}
These tricks are very useful when working with Laurent series.
\(\blacksquare \) \(\ \)Arithmetic series: \begin{align*} \sum _{n=1}^{N}n & =\frac{1}{2}N\left ( N+1\right ) \\ \sum _{n=1}^{N}a_{n} & =N\left ( \frac{a_{1}+a_{N}}{2}\right ) \end{align*}
i.e. the sum is \(N\) times the arithmetic mean.
\(\blacksquare \) \(\ \)Taylor series: Expanded around \(x=a\) is \[ f\left ( x\right ) =f\left ( a\right ) +\left ( x-a\right ) f^{\prime }\left ( a\right ) +\frac{\left ( x-a\right ) ^{2}f^{\prime \prime }\left ( a\right ) }{2!}+\frac{\left ( x-a\right ) ^{3}f^{\left ( 3\right ) }\left ( a\right ) }{3!}+\cdots +R_{n}\] Where \(R_{n}\) is remainder \(R_{n}=\frac{\left ( x-a\right ) ^{n+1}}{\left ( n+1\right ) !}f^{\left ( n+1\right ) }\left ( x_{0}\right ) \) where \(x_{0}\) is some point between \(x\) and \(a\).
\(\blacksquare \) \(\ \)Maclaurin series: Is just Taylor expanded around zero. i.e. \(a=0\)\[ f\left ( x\right ) =f\left ( 0\right ) +xf^{\prime }\left ( 0\right ) +\frac{x^{2}f^{\prime \prime }\left ( 0\right ) }{2!}+\frac{x^{3}f^{\left ( 3\right ) }\left ( 0\right ) }{3!}+\cdots \] \(\blacksquare \) \(\ \)This diagram shows the diﬀerent convergence of series and the relation between them
The above shows that an absolutely convergent series (\(B\)) is also convergent. Also a uniformly convergent series (\(D\)) is also convergent. But the series \(B\) is absolutely convergent and not uniform convergent. While \(D\) is uniform convergent and not absolutely convergent.
The series \(C\) is both absolutely and uniformly convergent. And ﬁnally the series \(A\) is convergent, but not absolutely (called conditionally convergent). Examples of \(B\) (converges absolutely but not uniformly) is\begin{align*} \sum _{n=0}^{\infty }x^{2}\frac{1}{\left ( 1+x^{2}\right ) ^{n}} & =x^{2}\left ( 1+\frac{1}{1+x^{2}}+\frac{1}{\left ( 1+x^{2}\right ) ^{2}}+\frac{1}{\left ( 1+x^{2}\right ) ^{3}}+\cdots \right ) \\ & =x^{2}+\frac{x^{2}}{1+x^{2}}+\frac{x^{2}}{\left ( 1+x^{2}\right ) ^{2}}+\frac{x^{2}}{\left ( 1+x^{2}\right ) ^{3}}+\cdots \end{align*}
And example of \(D\) (converges uniformly but not absolutely) is\[ \sum _{n=1}^{\infty }\left ( -1\right ) ^{n+1}\frac{1}{x^{2}+n}=\frac{1}{x^{2}+1}-\frac{1}{x^{2}+2}+\frac{1}{x^{3}+3}-\frac{1}{x^{4}+4}+\cdots \] Example of \(A\) (converges but not absolutely) is the alternating harmonic series\[ \sum _{n=1}^{\infty }\left ( -1\right ) ^{n+1}\frac{1}{n}=1-\frac{1}{2}+\frac{1}{3}-\frac{1}{4}+\cdots \] The above converges to \(\ln \left ( 2\right ) \) but absolutely it now becomes the harmonic series and it diverges\[ \sum _{n=1}^{\infty }\frac{1}{n}=1+\frac{1}{2}+\frac{1}{3}+\frac{1}{4}+\cdots \] For uniform convergence, we really need to have an \(x\) in the series and not just numbers, since the idea behind uniform convergence is if the series convergence to within an error tolerance \(\varepsilon \) using the same number of terms independent of the point \(x\) in the region.
\(\blacksquare \) The sequence \(\sum _{n=1}^{\infty }\frac{1}{n^{a}}\) converges for \(a>1\) and diverges for \(a\leq 1\). So \(a=1\) is the ﬂip value. For example\[ 1+\frac{1}{2}+\frac{1}{3}+\frac{1}{4}+\cdots \] Diverges, since \(a=1\), also \(1+\frac{1}{\sqrt{2}}+\frac{1}{\sqrt{3}}+\frac{1}{\sqrt{4}}+\cdots \) diverges, since \(a=\frac{1}{2}\leq 1\). But \(1+\frac{1}{4}+\frac{1}{9}+\frac{1}{16}+\cdots \) converges, where \(a=2\) here and the sum is \(\frac{\pi ^{2}}{6}\).
\(\blacksquare \) Using partial sums. Let \(\sum _{n=0}^{\infty }a_{n}\) be some sequence. The partial sum is \(S_{N}=\sum _{n=0}^{N}a_{n}\). Then\[ \sum _{n=0}^{\infty }a_{n}=\lim _{N\rightarrow \infty }S_{n}\] If \(\lim _{N\rightarrow \infty }S_{n}\) exist and ﬁnite, then we can say that \(\sum _{n=0}^{\infty }a_{n}\) converges. So here we use set up a sequence who terms are partial sum, and them look at what happens in the limit to such a term as \(N\rightarrow \theta \). Need to ﬁnd an example where this method is easier to use to test for convergence than the other method below.
\(\blacksquare \) Given a series, we are allowed to rearrange order of terms only when the series is absolutely convergent. Therefore for the alternating series \(1-\frac{1}{2}+\frac{1}{3}-\frac{1}{4}+\cdots \), do not rearrange terms since this is not absolutely convergent. This means the series sum is independent of the order in which terms are added only when the series is absolutely convergent.
\(\blacksquare \) In an inﬁnite series of complex numbers, the series converges, if the real part of the series and also the complex part of the series, each converges on their own.
\(\blacksquare \) Power series: \(f\left ( z\right ) =\sum _{n=0}^{\infty }a_{n}\left ( z-z_{0}\right ) ^{n}\). This series is centered at \(z_{0}\). Or expanded around \(z_{0}\). This has radius of convergence \(R\) is the series converges for \(\left \vert z-z_{0}\right \vert <R\) and diverges for \(\left \vert z-z_{0}\right \vert >R\).
\(\blacksquare \) Tests for convergence.
\(\blacksquare \) For Laurent series, lets say singularity is at \(z=0\) and \(z=1\). To expand about \(z=0\), get \(f\left ( z\right ) \) to look like \(\frac{1}{1-z}\) and use geometric series for \(\left \vert z\right \vert <1\). To expand about \(z=1\), there are two choices, to the inside and to the outside. For the outside, i.e. \(\left \vert z\right \vert >1\), get \(f\left ( z\right ) \) to have \(\frac{1}{1-\frac{1}{z}}\) form, since this now valid for \(\left \vert z\right \vert >1\).
\(\blacksquare \) Can only use power series \(\sum a_{n}\left ( z-z_{0}\right ) ^{n}\) to expand \(f\left ( z\right ) \) around \(z_{0}\) only if \(f\left ( z\right ) \) is analytic at \(z_{0}\). If \(f\left ( z\right ) \) is not analytic at \(z_{0}\) need to use Laurent series. Think of Laurent series as an extension of power series to handle singularities.
Let us ﬁnd the Laurent series for \(f\left ( z\right ) =\frac{5z-2}{z\left ( z-1\right ) }\). There is a singularity of order \(1\) at \(z=0\) and \(z=1\).
Expansion around \(z=0\). Let\begin{align*} g\left ( z\right ) & =zf\left ( z\right ) \\ & =\frac{5z-2}{\left ( z-1\right ) } \end{align*}
This makes \(g\left ( z\right ) \) analytic around \(z\), since \(g\left ( z\right ) \) do not have a pole at \(z=0\), then it is analytic around \(z=0\) and therefore it has a power series expansion around \(z=0\) given by\begin{equation} g\left ( z\right ) =\sum _{n=0}^{\infty }a_{n}z^{n} \tag{1} \end{equation} Where\[ a_{n}=\frac{1}{n!}\left . g^{\left ( n\right ) }\left ( z\right ) \right \vert _{z=0}\] But \[ g\left ( 0\right ) =2 \] And \begin{align*} g^{\prime }\left ( z\right ) & =\frac{5\left ( z-1\right ) -\left ( 5z-2\right ) }{\left ( z-1\right ) ^{2}}=\frac{-3}{\left ( z-1\right ) ^{2}}\\ g^{\prime }\left ( 0\right ) & =-3 \end{align*}
And \begin{align*} g^{\prime \prime }\left ( z\right ) & =\frac{-3\left ( -2\right ) }{\left ( z-1\right ) ^{3}}=\frac{6}{\left ( z-1\right ) ^{3}}\\ g^{\prime \prime }\left ( 0\right ) & =-6 \end{align*}
And \begin{align*} g^{\prime \prime \prime }\left ( z\right ) & =\frac{6\left ( -3\right ) }{\left ( z-1\right ) ^{4}}=\frac{-18}{\left ( z-1\right ) ^{4}}\\ g^{\prime \prime }\left ( 0\right ) & =-18 \end{align*}
And so on. Therefore, from (1)\begin{align*} g\left ( z\right ) & =g\left ( 0\right ) +g^{\prime }\left ( 0\right ) z+\frac{1}{2!}g^{\prime \prime }\left ( 0\right ) z^{2}+\frac{1}{3!}g^{\prime \prime \prime }\left ( 0\right ) z^{3}+\cdots \\ & =2-3z-\frac{6}{2}z^{2}-\frac{18}{3!}z^{3}-\cdots \\ & =2-3z-3z^{2}-3z^{3}-\cdots \end{align*}
Therefore\begin{align*} f\left ( z\right ) & =\frac{g\left ( z\right ) }{z}\\ & =\frac{2}{z}-3-3z-3z^{2}-\cdots \end{align*}
The residue is \(2\). The above expansion is valid around \(z=0\) up and not including the next singularity, which is at \(z=1\). Now we ﬁnd the expansion of \(f\left ( z\right ) \) around \(z=1\). Let\begin{align*} g\left ( z\right ) & =\left ( z-1\right ) f\left ( z\right ) \\ & =\frac{5z-2}{z} \end{align*}
This makes \(g\left ( z\right ) \) analytic around \(z=1\), since \(g\left ( z\right ) \) do not have a pole at \(z=1\). Therefore it has a power series expansion about \(z=1\) given by\begin{equation} g\left ( z\right ) =\sum _{n=0}^{\infty }a_{n}\left ( z-1\right ) ^{n} \tag{1} \end{equation} Where\[ a_{n}=\frac{1}{n!}\left . g^{\left ( n\right ) }\left ( z\right ) \right \vert _{z=1}\] But \[ g\left ( 1\right ) =3 \] And \begin{align*} g^{\prime }\left ( z\right ) & =\frac{5z-\left ( 5z-2\right ) }{z^{2}}=\frac{2}{z^{2}}\\ g^{\prime }\left ( 1\right ) & =2 \end{align*}
And \begin{align*} g^{\prime \prime }\left ( z\right ) & =\frac{2\left ( -2\right ) }{z^{3}}=\frac{-4}{z^{3}}\\ g^{\prime \prime }\left ( 1\right ) & =-4 \end{align*}
And \begin{align*} g^{\prime \prime \prime }\left ( z\right ) & =\frac{-4\left ( -3\right ) }{z^{4}}=\frac{12}{z^{4}}\\ g^{\prime \prime }\left ( 1\right ) & =12 \end{align*}
And so on. Therefore, from (1)\begin{align*} g\left ( z\right ) & =g\left ( 1\right ) +g^{\prime }\left ( 1\right ) \left ( z-1\right ) +\frac{1}{2!}g^{\prime \prime }\left ( 1\right ) \left ( z-1\right ) ^{2}+\frac{1}{3!}g^{\prime \prime \prime }\left ( 1\right ) \left ( z-1\right ) ^{3}+\cdots \\ & =3+2\left ( z-1\right ) -\frac{4}{2}\left ( z-1\right ) ^{2}+\frac{12}{3!}\left ( z-1\right ) ^{3}-\cdots \\ & =3+2\left ( z-1\right ) -2\left ( z-1\right ) ^{2}+2\left ( z-1\right ) ^{3}-\cdots \end{align*}
Therefore\begin{align*} f\left ( z\right ) & =\frac{g\left ( z\right ) }{z-1}\\ & =\frac{3}{z-1}+2-2\left ( z-1\right ) +2\left ( z-1\right ) ^{2}-2\left ( z-1\right ) ^{3}+\cdots \end{align*}
The residue is \(3\). The above expansion is valid around \(z=1\) up and not including the next singularity, which is at \(z=0\) inside a circle of radius \(1\).
Putting the above two regions together, then we see there is a series expansion of \(f\left ( z\right ) \) that is shared between the two regions, in the shaded region below.
Let check same series in the shared region give same values. Using the series expansion about \(f\left ( 0\right ) \) to ﬁnd \(f\left ( z\right ) \) at point \(z=\frac{1}{2}\), gives \(-2\) when using \(10\) terms in the series. Using series expansion around \(z=1\) to ﬁnd \(f\left ( \frac{1}{2}\right ) \) using \(10\) terms also gives \(-2\). So both series are valid produce same result.
This method is simpler than the above, but it results in diﬀerent regions. It is based on converting the expression in order to use geometric series expansion on it.\[ f\left ( z\right ) =\frac{5z-2}{z\left ( z-1\right ) }\] Since there is a pole at \(z=0\) and at \(z=1\), then we ﬁrst ﬁnd expansion for \(0<\left \vert z\right \vert <1\). To do this, we write the above as\begin{align*} f\left ( z\right ) & =\frac{5z-2}{z}\left ( \frac{1}{z-1}\right ) \\ & =\frac{2-5z}{z}\left ( \frac{1}{1-z}\right ) \end{align*}
And now expand \(\frac{1}{1-z}\) using geometric series, which is valid for \(\left \vert z\right \vert <1\). This gives\begin{align*} f\left ( z\right ) & =\frac{2-5z}{z}\left ( 1+z+z^{2}+z^{3}+\cdots \right ) \\ & =\frac{2}{z}\left ( 1+z+z^{2}+z^{3}+\cdots \right ) -5\left ( 1+z+z^{2}+z^{3}+\cdots \right ) \\ & =\left ( \frac{2}{z}+2+2z+2z^{2}+\cdots \right ) -\left ( 5+5z+5z^{2}+5z^{3}+\cdots \right ) \\ & =\frac{2}{z}-3-3z-3z^{2}-3z^{3}-\cdots \end{align*}
The above is valid for \(0<\left \vert z\right \vert <1\) which agrees with result of method 1.
Now, to ﬁnd expansion for \(\left \vert z\right \vert >1\), we need a term that looks like \(\left ( \frac{1}{1-\frac{1}{z}}\right ) \). Since now it can be expanded for \(\left \vert \frac{1}{z}\right \vert <1\) or \(\left \vert z\right \vert >1\) which is what we want. Therefore, writing \(f\left ( z\right ) \) as\[ f\left ( z\right ) =\frac{5z-2}{z\left ( z-1\right ) }=\frac{5z-2}{z^{2}\left ( 1-\frac{1}{z}\right ) }=\frac{5z-2}{z^{2}}\left ( \frac{1}{1-\frac{1}{z}}\right ) \] But for \(\left \vert \frac{1}{z}\right \vert <1\) the above becomes\begin{align*} f\left ( z\right ) & =\frac{5z-2}{z^{2}}\left ( 1+\frac{1}{z}+\frac{1}{z^{2}}+\frac{1}{z^{3}}+\cdots \right ) \\ & =\frac{5}{z}\left ( 1+\frac{1}{z}+\frac{1}{z^{2}}+\frac{1}{z^{3}}+\cdots \right ) -\frac{2}{z^{2}}\left ( 1+\frac{1}{z}+\frac{1}{z^{2}}+\frac{1}{z^{3}}+\cdots \right ) \\ & =\left ( \frac{5}{z}+\frac{5}{z^{2}}+\frac{5}{z^{3}}+\frac{5}{z^{4}}+\cdots \right ) -\left ( \frac{2}{z^{2}}+\frac{2}{z^{3}}+\frac{2}{z^{4}}+\frac{2}{z^{5}}+\cdots \right ) \\ & =\frac{5}{z}+\frac{3}{z^{3}}+\frac{3}{z^{4}}+\frac{3}{z^{5}}+\cdots \end{align*}
With residue \(5\). The above is valid for \(\left \vert z\right \vert >1\). The following diagram illustrates the result obtained from method 2.
For expansion about \(z=0\), this uses same method as above, giving same series valid for \(\left \vert z\right \vert <1\,.\) This method is a little diﬀerent for those points other than zero. The idea is to replace \(z\) by \(z-z_{0}\) where \(z_{0}\) is the point we want to expand about and do this replacement in \(f\left ( z\right ) \) itself. So for \(z=1\) using this example, we let \(\xi =z-1\) hence \(z=\xi +1\). Then \(f\left ( z\right ) \) becomes
\begin{align*} f\left ( z\right ) & =\frac{5z-2}{z\left ( z-1\right ) }\\ & =\frac{5\left ( \xi +1\right ) -2}{\left ( \xi +1\right ) \left ( \xi \right ) }\\ & =\frac{5\left ( \xi +1\right ) -2}{\xi }\left ( \frac{1}{\xi +1}\right ) \\ & =\frac{5\xi +3}{\xi }\left ( \frac{1}{1+\xi }\right ) \end{align*}
Now we expand \(\frac{1}{1+\xi }\) for \(\left \vert \xi \right \vert <1\) and the above becomes\begin{align*} f\left ( z\right ) & =\frac{5\xi +3}{\xi }\left ( 1-\xi +\xi ^{2}-\xi ^{3}+\xi ^{4}-\cdots \right ) \\ & =\frac{5\xi +3}{\xi }\left ( 1-\xi +\xi ^{2}-\xi ^{3}+\xi ^{4}-\cdots \right ) \\ & =\left ( \frac{5\xi +3}{\xi }-\left ( 5\xi +3\right ) +\left ( 5\xi +3\right ) \xi -\left ( 5\xi +3\right ) \xi ^{2}+\cdots \right ) \\ & =\left ( 5+\frac{3}{\xi }-5\xi -3+5\xi ^{2}+3\xi -5\xi ^{3}-3\xi ^{2}+\cdots \right ) \\ & =\left ( 2+\frac{3}{\xi }-2\xi +2\xi ^{2}-2\xi ^{3}+\cdots \right ) \end{align*}
We now replace \(\xi =z-1\) and the above becomes\[ f\left ( z\right ) =\left ( \frac{3}{\left ( z-1\right ) }+2-2\left ( z-1\right ) +2\left ( z-1\right ) ^{2}-2\left ( z-1\right ) ^{3}+2\left ( z-1\right ) ^{4}-\cdots \right ) \] The above is valid for \(\left \vert \xi \right \vert <1\) or \(\left \vert z-1\,\right \vert <1\) or \(\,-1<\left ( z-1\right ) <1\) or \(0<z<2\). This gives same series and for same region as in method one. But this is little faster as it uses Binomial series short cut to ﬁnd the expansion instead of calculating derivatives as in method one.
Method one and method three give same series and for same regions. Method three uses binomial expansion as short cut and requires one to convert \(f\left ( z\right ) \) to form to allow using Binomial expansion. Method one does not use binomial expansion but requires doing many derivatives to evaluate the terms of the power series. It is more direct method.
Method two also uses binomial expansion, but gives diﬀerent regions that method one and three.
If one is good in diﬀerentiation, method one seems the most direct. Otherwise, the choice is between method two or three as they both use Binomial expansion. Method two seems a little more direct than method three.It also depends what the problem is asking form. If the problem asks to expand around \(z_{0}\) vs. if it is asking to ﬁnd expansion in \(\left \vert z\right \vert >1\) for example, then this decides which method to use.
\(\blacksquare \) Gamma function is deﬁned by \[ \Gamma \left ( x\right ) =\int _{0}^{\infty }t^{x-1}e^{-t}dt\qquad x>0 \] The above is called the Euler representation. Or if we want it deﬁned in complex domain, the above becomes \[ \Gamma \left ( z\right ) =\int _{0}^{\infty }t^{z-1}e^{-t}dt\qquad \operatorname{Re}\left ( z\right ) >0 \] Since the above is deﬁned only for right half plane, there is way to extend this to left half plane, using what is called analytical continuation. More on this below. First, some relations involving \(\Gamma \left ( x\right ) \)\begin{align*} \Gamma \left ( z\right ) & =\left ( z-1\right ) \Gamma \left ( z-1\right ) \qquad \operatorname{Re}\left ( z\right ) >1\\ \Gamma \left ( 1\right ) & =1\\ \Gamma \left ( 2\right ) & =1\\ \Gamma \left ( 3\right ) & =2\\ \Gamma \left ( 4\right ) & =3!\\ \Gamma \left ( n\right ) & =\left ( n-1\right ) !\\ \Gamma \left ( n+1\right ) & =n!\\ \Gamma \left ( \frac{1}{2}\right ) & =\sqrt{\pi }\\ \Gamma \left ( z+1\right ) & =z\Gamma \left ( z\right ) \qquad \text{recursive formula}\\ \Gamma \left ( \bar{z}\right ) & =\overline{\Gamma \left ( z\right ) }\\ \Gamma \left ( n+\frac{1}{2}\right ) & =\frac{1\cdot 3\cdot 5\cdots \left ( 2n-1\right ) }{2^{n}}\sqrt{\pi } \end{align*}
\(\blacksquare \) To extend \(\Gamma \left ( z\right ) \) to the left half plane, i.e. for negative values. Let us deﬁne, using the above recursive formula \[ \bar{\Gamma }\left ( z\right ) =\frac{\Gamma \left ( z+1\right ) }{z}\qquad \operatorname{Re}\left ( z\right ) >-1 \] For example \[ \bar{\Gamma }\left ( -\frac{1}{2}\right ) =\frac{\Gamma \left ( \frac{1}{2}\right ) }{-\frac{1}{2}}=-2\Gamma \left ( \frac{1}{2}\right ) =-2\sqrt{\pi }\] And for \(\operatorname{Re}\left ( z\right ) >-2\) \[ \bar{\Gamma }\left ( -\frac{3}{2}\right ) =\frac{\bar{\Gamma }\left ( -\frac{3}{2}+1\right ) }{-\frac{3}{2}}=\left ( \frac{1}{-\frac{3}{2}}\right ) \bar{\Gamma }\left ( -\frac{1}{2}\right ) =\left ( \frac{1}{-\frac{3}{2}}\right ) \left ( \frac{1}{-\frac{1}{2}}\right ) \Gamma \left ( \frac{1}{2}\right ) =\left ( \frac{1}{-\frac{3}{2}}\right ) \left ( \frac{1}{-\frac{1}{2}}\right ) \sqrt{\pi }=\frac{4}{3}\sqrt{\pi }\] And so on. Notice that for \(x<0\) the functions \(\Gamma \left ( x\right ) \) are not deﬁned for all negative integers \(x=-1,-2,\cdots \) it is also not deﬁned for \(x=0\)
\(\blacksquare \) The above method of extending (or analytical continuation) of the Gamma function to negative values is due to Euler. Another method to extend Gamma is due to Weierstrass. It starts by rewriting from the deﬁnition as follows, where \(a>0\)\begin{align} \Gamma \left ( z\right ) & =\int _{0}^{\infty }t^{z-1}e^{-t}dt\nonumber \\ & =\int _{0}^{a}t^{z-1}e^{-t}dt+\int _{a}^{\infty }t^{z-1}e^{-t}dt \tag{1} \end{align}
Expanding the integrand in the ﬁrst integral using Taylor series gives\begin{align*} \int _{0}^{a}t^{z-1}e^{-t}dt & =\int _{0}^{a}t^{z-1}\left ( 1+\left ( -t\right ) +\frac{\left ( -t\right ) ^{2}}{2!}+\frac{\left ( -t\right ) ^{3}}{3!}+\cdots \right ) dt\\ & =\int _{0}^{a}t^{z-1}\left ( 1+\left ( -t\right ) +\frac{\left ( -t\right ) ^{2}}{2!}+\frac{\left ( -t\right ) ^{3}}{3!}+\cdots \right ) dt\\ & =\int _{0}^{a}t^{z-1}\sum _{n=0}^{\infty }\frac{\left ( -1\right ) ^{n}t^{n}}{n!}dt\\ & =\int _{0}^{a}\sum _{n=0}^{\infty }\frac{\left ( -1\right ) ^{n}t^{n+z-1}}{n!}dt\\ & =\sum _{n=0}^{\infty }\int _{0}^{a}\frac{\left ( -1\right ) ^{n}t^{n+z-1}}{n!}dt\\ & =\sum _{n=0}^{\infty }\frac{\left ( -1\right ) ^{n}}{n!}\int _{0}^{a}t^{n+z-1}dt\\ & =\sum _{n=0}^{\infty }\frac{\left ( -1\right ) ^{n}}{n!}\left [ \frac{t^{n+z}}{n+z}\right ] _{0}^{a}\\ & =\sum _{n=0}^{\infty }\frac{\left ( -1\right ) ^{n}}{n!\left ( n+z\right ) }a^{n+z} \end{align*}
This takes care of the ﬁrst integral in (1). Now, since the lower limits of the second integral in (1) is not zero, then there is no problem integrating it directly. Remember that in the Euler deﬁnition, it had zero in the lower limit, that is why we said there \(\operatorname{Re}\left ( z\right ) >1\). Now can can choose any value for \(a\). Weierstrass choose \(a=1\). Hence (1) becomes\begin{align} \Gamma \left ( z\right ) & =\int _{0}^{a}t^{z-1}e^{-t}dt+\int _{a}^{\infty }t^{z-1}e^{-t}dt\nonumber \\ & =\sum _{n=0}^{\infty }\frac{\left ( -1\right ) ^{n}}{n!\left ( n+z\right ) }+\int _{1}^{\infty }t^{z-1}e^{-t}dt \tag{2} \end{align}
Notice the term \(a^{n+z}\) now is just \(1\) since \(a=1\). The second integral above can now be integrated directly. Let us now verify that Euler continuation \(\bar{\Gamma }\left ( z\right ) \) for say \(z=-\frac{1}{2}\) gives the same result as Weierstrass formula. From above, we found that \(\bar{\Gamma }\left ( z\right ) =-2\sqrt{\pi }\). Equation (2) for \(z=-\frac{1}{2}\) becomes\begin{equation} \bar{\Gamma }\left ( -\frac{1}{2}\right ) =\sum _{n=0}^{\infty }\frac{\left ( -1\right ) ^{n}}{n!\left ( n-\frac{1}{2}\right ) }+\int _{1}^{\infty }t^{-\frac{3}{2}}e^{-t}dt \tag{3} \end{equation} Using the computer \[ \sum _{n=0}^{\infty }\frac{\left ( -1\right ) ^{n}}{n!\left ( n-\frac{1}{2}\right ) }=-2\sqrt{\pi }+2\sqrt{\pi }\left ( 1-\operatorname{erf}\left ( 1\right ) \right ) -2\frac{1}{e}\] And direct integration \[ \int _{1}^{\infty }t^{-\frac{3}{2}}e^{-t}dt=-2\sqrt{\pi }+2\sqrt{\pi }\operatorname{erf}\left ( 1\right ) +\frac{2}{e}\] Hence (3) becomes\begin{align*} \bar{\Gamma }\left ( -\frac{1}{2}\right ) & =\left ( -2\sqrt{\pi }+2\sqrt{\pi }\left ( 1-\operatorname{erf}\left ( 1\right ) \right ) -2\frac{1}{e}\right ) +\left ( -2\sqrt{\pi }+2\sqrt{\pi }\operatorname{erf}\left ( 1\right ) +\frac{2}{e}\right ) \\ & =-2\sqrt{\pi } \end{align*}
Which is the same as using Euler method. Let us check for \(z=-\frac{2}{3}\). We found above that \(\bar{\Gamma }\left ( -\frac{3}{2}\right ) =\frac{4}{3}\sqrt{\pi }\) using Euler method of analytical continuation. Now we will check using Weierstrass method. Equation (2) for \(z=-\frac{3}{2}\) becomes\[ \bar{\Gamma }\left ( -\frac{3}{2}\right ) =\sum _{n=0}^{\infty }\frac{\left ( -1\right ) ^{n}}{n!\left ( n-\frac{3}{2}\right ) }+\int _{1}^{\infty }t^{-\frac{5}{2}}e^{-t}dt \] Using the computer \[ \sum _{n=0}^{\infty }\frac{\left ( -1\right ) ^{n}}{n!\left ( n-\frac{3}{2}\right ) }=\frac{4\sqrt{\pi }}{3}-\frac{4\sqrt{\pi }\left ( 1-\operatorname{erf}\left ( 1\right ) \right ) }{3}+\frac{2}{3e}\] And \[ \int _{1}^{\infty }t^{-\frac{5}{2}}e^{-t}dt=-\frac{4\sqrt{\pi }\operatorname{erf}\left ( 1\right ) }{3}+\frac{4\sqrt{\pi }}{3}-\frac{2}{3e}\] Hence\begin{align*} \bar{\Gamma }\left ( -\frac{3}{2}\right ) & =\left ( \frac{4\sqrt{\pi }}{3}-\frac{4\sqrt{\pi }\left ( 1-\operatorname{erf}\left ( 1\right ) \right ) }{3}+\frac{2}{3e}\right ) +\left ( -\frac{4\sqrt{\pi }\operatorname{erf}\left ( 1\right ) }{3}+\frac{4\sqrt{\pi }}{3}-\frac{2}{3e}\right ) \\ & =\frac{4}{3}\sqrt{\pi } \end{align*}
Which is the same as using the Euler method. Clearly the Euler method for analytical continuation of the Gamma function is simpler to compute.
\(\blacksquare \) Euler reﬂection formula \begin{align*} \Gamma \left ( x\right ) \Gamma \left ( 1-x\right ) & =\int _{0}^{\infty }\frac{t^{x-1}}{1+t}dt\qquad 0<x<1\\ & =\frac{\pi }{\sin \left ( \pi x\right ) } \end{align*}
Where contour integration was used to derive the above. See Mary Boas text book, page 607, second edition, example 5 for full derivation.
\(\blacksquare \) \(\Gamma \left ( z\right ) \) has singularities at \(z=0,-1,-2,\cdots \) and \(\Gamma \left ( 1-z\right ) \) has singularities at \(z=1,2,3,\cdots \) so in the above reﬂection formula, the zeros of \(\sin \left ( \pi x\right ) \) cancel the singularities of \(\Gamma \left ( x\right ) \) when it is written as \[ \Gamma \left ( 1-x\right ) =\frac{\pi }{\Gamma \left ( x\right ) \sin \left ( \pi x\right ) }\]
\(\blacksquare \) \(\frac{1}{\Gamma \left ( z\right ) }\) is entire.
\(\blacksquare \) There are other representations for \(\Gamma \left ( x\right ) \). One that uses products by Euler also is \begin{align*} \Gamma \left ( z\right ) & =\frac{1}{z}\Pi _{n=1}^{\infty }\frac{\left ( 1+\frac{1}{n}\right ) ^{z}}{1+\frac{z}{n}}\\ & =\lim _{n\rightarrow \infty }\frac{n!\left ( n+1\right ) ^{z}}{z\left ( z-1\right ) \cdots \left ( z+n\right ) } \end{align*}
And another due to Weierstrass is \begin{align*} \Gamma \left ( z\right ) & =\frac{e^{-\gamma z}}{z}\Pi _{n=1}^{\infty }\frac{e^{\frac{z}{n}}}{1+\frac{z}{n}}\\ & =e^{-\gamma z}\lim _{n\rightarrow \infty }\frac{n!\exp \left ( z\left ( 1+\frac{1}{2}+\cdots +\frac{1}{n}\right ) \right ) }{z\left ( z+1\right ) \left ( z+2\right ) \cdots \left ( z+n\right ) } \end{align*}
\(\blacksquare \) Given by \(\zeta \left ( s\right ) =\sum _{n=1}^{\infty }\frac{1}{n^{s}}\) for \(\operatorname{Re}\left ( s\right ) >1\). Euler studied this and It was extended to the whole complex plane by Riemann. So the Riemann zeta function refer to the one with the extension to the whole complex plane. Euler only looked at it on the real line. It has pole at \(s=1\). Has trivial zeros at \(s=-2,-4,-6,\cdots \) and all its non trivial zeros are inside the critical strip \(0<s<1\) and they all lie on the critical line \(s=\frac{1}{2}\). \(\zeta \left ( s\right ) \) is also deﬁned by integral formula \[ \zeta \left ( s\right ) =\frac{1}{\Gamma \left ( s\right ) }\int _{0}^{\infty }\frac{1}{e^{t}-1}\frac{t^{s}}{t}dt\qquad \operatorname{Re}\left ( s\right ) >1 \]
\(\blacksquare \) The connection between \(\zeta \left ( s\right ) \) prime numbers is given by the Euler product formula
\begin{align*} \zeta \left ( s\right ) & =\Pi _{p}\frac{1}{1-p^{-s}}\\ & =\left ( \frac{1}{1-2^{-s}}\right ) \left ( \frac{1}{1-3^{-s}}\right ) \left ( \frac{1}{1-5^{-s}}\right ) \left ( \frac{1}{1-7^{-s}}\right ) \cdots \\ & =\left ( \frac{1}{1-\frac{1}{2^{s}}}\right ) \left ( \frac{1}{1-\frac{1}{3^{s}}}\right ) \left ( \frac{1}{1-\frac{1}{5^{s}}}\right ) \left ( \frac{1}{1-\frac{1}{7^{s}}}\right ) \cdots \\ & =\left ( \frac{2^{s}}{2^{s}-1}\right ) \left ( \frac{3^{s}}{3^{s}-1}\right ) \left ( \frac{5^{s}}{5^{s}-1}\right ) \left ( \frac{7^{s}}{7^{s}-1}\right ) \cdots \end{align*}
\(\blacksquare \) \(\zeta \left ( s\right ) \) functional equation is
\[ \zeta \left ( s\right ) =2^{s}\pi ^{s-1}\sin \left ( \frac{\pi s}{2}\right ) \Gamma \left ( 1-s\right ) \zeta \left ( 1-s\right ) \]
\(\blacksquare \) Complex identities \begin{align*} \left \vert z\right \vert ^{2} & =z\bar{z}\\ \overline{\left ( \bar{z}\right ) } & =z\\ \overline{\left ( z_{1}+z_{2}\right ) } & =\bar{z}_{1}+\bar{z}_{2}\\ \left \vert \bar{z}\right \vert & =\left \vert z\right \vert \\ \left \vert z_{1}z_{2}\right \vert & =\left \vert z_{1}\right \vert \left \vert z_{2}\right \vert \\ \operatorname{Re}\left ( z\right ) & =\frac{z+\bar{z}}{2}\\ \operatorname{Im}\left ( z\right ) & =\frac{z+\bar{z}}{2i}\\ \arg \left ( z_{1}z_{2}\right ) & =\arg \left ( z_{1}\right ) +\arg \left ( z_{2}\right ) \end{align*}
\(\blacksquare \) A complex function \(f\left ( z\right ) \) is analytic in a region \(D\) if it is deﬁned and diﬀerentiable at all points in \(D\). One way to check for analyticity is to use the Cauchy Riemann (CR) equations (this is a necessary condition but not suﬃcient). If \(f\left ( z\right ) \) satisﬁes CR everywhere in that region then it is analytic. Let \(f\left ( z\right ) =u\left ( x,y\right ) +iv\left ( x,y\right ) \), then these two equations in Cartesian coordinates are\begin{align*} \frac{\partial u}{\partial x} & =\frac{\partial v}{\partial y}\\ -\frac{\partial u}{\partial y} & =\frac{\partial v}{\partial x} \end{align*}
Sometimes it is easier to use the polar form of these. Let \(f\left ( z\right ) =r\cos \theta +i\sin \theta \), then the equations become\begin{align*} \frac{\partial u}{\partial r} & =\frac{1}{r}\frac{\partial v}{\partial \theta }\\ -\frac{1}{r}\frac{\partial u}{\partial \theta } & =\frac{\partial v}{\partial r} \end{align*}
To remember them, think of the \(r\) as the \(x\) and \(\theta \) as the \(y\).
Let us apply these on \(\sqrt{z}\) to see how it works. Since \(z=re^{i\theta +2n\pi }\) then \(f\left ( z\right ) =\) \(\sqrt{r}e^{i\frac{\theta }{2}+n\pi }\).This is multi-valued function. One value for \(n=0\) and another for \(n=1\). The ﬁrst step is to make it single valued. Choosing \(n=0\) gives the principal value. Then \(f\left ( z\right ) =\sqrt{r}e^{i\frac{\theta }{2}}\). Now we ﬁnd the branch points. \(z=0\) is a branch point. We can pick \(-\pi <\theta <\pi \) and pick the negative real axis as the branch cut (the other branch point being \(-\infty \)). This is one choice.
We could have picked \(0<\theta <2\pi \) and had the positive \(x\) axis as the branch cut, where now the second branch point is \(+\infty \) but in both cases, origin is still part of the branch cut. Let us stick with \(-\pi <\theta <\pi \).
Given all of this, now\(\sqrt{z}=\sqrt{r}e^{i\frac{\theta }{2}}=\sqrt{r}\left ( \cos \left ( \frac{\theta }{2}\right ) +i\sin \left ( \frac{\theta }{2}\right ) \right ) \), hence \(u=\sqrt{r}\cos \left ( \frac{\theta }{2}\right ) \) and \(v=\sqrt{r}\sin \left ( \frac{\theta }{2}\right ) \). Therefore \(\frac{\partial u}{\partial r}=\frac{1}{2}\frac{1}{\sqrt{r}}\cos \left ( \frac{\theta }{2}\right ) ,\) and \(\frac{\partial v}{\partial \theta }=\frac{1}{2}\sqrt{r}\cos \left ( \frac{\theta }{2}\right ) \) and \(\frac{\partial u}{\partial \theta }=-\frac{1}{2}\sqrt{r}\sin \left ( \frac{\theta }{2}\right ) \) and \(\frac{\partial v}{\partial r}=\frac{1}{2}\frac{1}{\sqrt{r}}\sin \left ( \frac{\theta }{2}\right ) \). Applying Cauchy-Riemann above gives \begin{align*} \frac{1}{2}\frac{1}{\sqrt{r}}\cos \left ( \frac{\theta }{2}\right ) & =\frac{1}{r}\frac{1}{2}\sqrt{r}\cos \left ( \frac{\theta }{2}\right ) \\ \frac{1}{2}\frac{1}{\sqrt{r}}\cos \left ( \frac{\theta }{2}\right ) & =\frac{1}{2}\frac{1}{\sqrt{r}}\cos \left ( \frac{\theta }{2}\right ) \end{align*}
Satisﬁed. and for the second equation \begin{align*} -\frac{1}{r}\left ( -\frac{1}{2}\sqrt{r}\sin \left ( \frac{\theta }{2}\right ) \right ) & =\frac{1}{2}\frac{1}{\sqrt{r}}\sin \left ( \frac{\theta }{2}\right ) \\ \frac{1}{2}\frac{1}{\sqrt{r}}\sin \left ( \frac{\theta }{2}\right ) & =\frac{1}{2}\frac{1}{\sqrt{r}}\sin \left ( \frac{\theta }{2}\right ) \end{align*}
so \(\sqrt{z}\) is analytic in the region \(-\pi <\theta <\pi \), and not including branch points and branch cut.
\(\blacksquare \) We can’t just say \(f\left ( z\right ) \) is Analytic and stop. Have to say \(f\left ( z\right ) \) is analytic in a region or at a point. When we say \(f\left ( z\right ) \) analytic at a point, we mean analytic in small region around the point.
If \(f\left ( z\right ) \) is deﬁned only at an isolated point \(z_{0}\) and not deﬁned anywhere around it, then the function can not be analytic at \(z_{0}\) since it is not diﬀerentiable at \(z_{0}\). Also \(f\left ( z\right ) \) is analytic at a point \(z_{0}\) if the power series for \(f\left ( z\right ) \) expanded around \(z_{0}\) converges to \(f\left ( z\right ) \) evaluated at \(z_{0}\). An analytic complex function mean it is inﬁnitely many times diﬀerentiable in the region, which means the limit exist \(\lim _{\Delta z\rightarrow 0}\frac{f\left ( z+\Delta z\right ) -f\left ( z\right ) }{\Delta z}\) and does not depend on direction.
\(\blacksquare \) Before applying the Cauchy Riemann equations, make sure the complex function is ﬁrst made to be single valued.
\(\blacksquare \) Remember that Cauchy Riemann equations as necessary but not suﬃcient condition for function to be analytic. The extra condition needed is that all the partial derivatives are continuous. Need to ﬁnd example where CR is satisﬁed but not the continuity on the partial derivatives. Most of the HW problems just needs the CR but good to keep an eye on this other condition.
\(\blacksquare \) Cauchy-Goursat: If \(f\left ( z\right ) \) is analytic on and inside closed contour \(C\) then \({\displaystyle \oint \limits _{C}} f\left ( z\right ) dz=0\). But remember that if \({\displaystyle \oint \limits _{C}} f\left ( z\right ) dz=0\) then this does not necessarily imply \(f\left ( z\right ) \) is analytic on and inside \(C\). So this is an IF and not an IFF relation. For example \({\displaystyle \oint \limits _{C}} \frac{1}{z^{2}}dz=0\) around unit circle centered at origin, but clearly \(\frac{1}{z^{2}}\) is not analytic everywhere inside \(C\), since it has a singularity at \(z=0\).
proof of Cauchy-Goursat: The proof uses two main ideas. It uses the Cauchy-Riemann equations and also uses Green theorem. Green’s Theorem says\begin{equation} \int _{C}Pdx+Qdy=\int _{D}\left ( \frac{\partial Q}{\partial x}-\frac{\partial P}{\partial y}\right ) dA \tag{1} \end{equation} So Green’s Theorem transforms integration on the boundary \(C\) of region \(D\) by integration over the area inside the boundary \(C\). Let \(f\left ( z\right ) =u+iv\). And since \(z=x+iy\) then \(dz=dx+idy\). Therefore\begin{align}{\displaystyle \oint \limits _{C}} f\left ( z\right ) dz & ={\displaystyle \oint \limits _{C}} \left ( u+iv\right ) \left ( dx+idy\right ) \nonumber \\ & ={\displaystyle \oint \limits _{C}} udx+uidy+ivdx-vdy\nonumber \\ & ={\displaystyle \oint \limits _{C}} \left ( udx-vdy\right ) +i{\displaystyle \oint \limits _{C}} vdx+udy \tag{2} \end{align}
We now apply (1) to each of the two integrals in (3). Hence the ﬁrst integral in (2) becomes\[{\displaystyle \oint \limits _{C}} \left ( udx-vdy\right ) =\int _{D}\left ( -\frac{\partial v}{\partial x}-\frac{\partial u}{\partial y}\right ) dA \] But from CR, we know that \(-\frac{\partial u}{\partial y}=\frac{\partial v}{\partial x}\), hence the above is zero. And the second integral in (2) becomes\[{\displaystyle \oint \limits _{C}} vdx+udy=\int _{D}\left ( \frac{\partial u}{\partial x}-\frac{\partial v}{\partial y}\right ) dA \] But from CR, we know that \(\frac{\partial u}{\partial x}=\frac{\partial v}{\partial y}\), hence the above is zero. Therefore the whole integral in (2) is zero. Therefore \({\displaystyle \oint \limits _{C}} f\left ( z\right ) dz=0\). QED.
\(\blacksquare \) Cauchy residue: If \(f\left ( z\right ) \) is analytic on and inside closed contour \(C\) except at some isolated points \(z_{1},z_{2},\cdots ,z_{N}\) then \({\displaystyle \oint \limits _{C}} f\left ( z\right ) dz=2\pi i\sum _{j=1}^{N}\operatorname{Res}\left ( f\left ( z\right ) \right ) _{z=z_{j}}\). The term \(\operatorname{Res}\left ( f\left ( z\right ) \right ) _{z=z_{j}}\) is the residue of \(f\left ( z\right ) \) at point \(z_{j}\). Use Laurent expansion of \(f\left ( z\right ) \) to ﬁnd residues. See above on methods how to ﬁnd Laurent series.
\(\blacksquare \) Maximum modulus principle: If \(f\left ( z\right ) \) is analytic in some region \(D\) and is not constant inside \(D\), then its maximum value must be on the boundary. Also its minimum on the boundary, as long as \(f\left ( z\right ) \neq 0\) anywhere inside \(D\). In the other hand, if \(f\left ( z\right ) \) happend to have a maximum at some point \(z_{0}\) somewhere inside \(D\), then this implies that \(f\left ( z\right ) \) is constant everywhere and will have the value \(f\left ( z_{0}\right ) \) everywhere. What all this really mean, is that if \(f\left ( z\right ) \) is analytic and not constant in \(D\), then its maximum is on the boundary and not inside.
There is a complicated proof of this. See my notes for Physics 501. Hopefully this will not come up in the exam since I did not study the proof.
\(\blacksquare \) These deﬁnitions from book of Joseph Bak
\(\blacksquare \) Some important formulas.
For example, to solve \({\displaystyle \oint \limits _{C}} \frac{1}{\left ( z+1\right ) \left ( z+2\right ) }dz\) around \(C\) unit circle. Rewriting this as \({\displaystyle \oint \limits _{C}} \frac{\frac{1}{z+2}}{\left ( z-\left ( -1\right ) \right ) }dz\) where now \(f\left ( z\right ) =\frac{1}{z+2}\) and now we can use Cauchy integral formula. So all what we have to do is just evaluate \(\frac{1}{z+2}\) at \(z=-1\), which gives \({\displaystyle \oint \limits _{C}} \frac{1}{\left ( z+1\right ) \left ( z+2\right ) }dz=2\pi i\). This works if \(g\left ( z\right ) \) can be factored into \(\frac{f\left ( z\right ) }{z-z_{0}}\) where \(f\left ( z\right ) \) is analytic on and inside \(C\). This would not work if \(g\left ( z\right ) \) has more than one pole inside \(C\).
where \(\left ( x_{0},y_{0}\right ) \) in the line initial point and \(\left ( x_{1},y_{1}\right ) \) is the line end point. This works for straight lines. Now use the above and rewrite \(z=x+iy\) as \(z\left ( t\right ) =x\left ( t\right ) +iy\left ( t\right ) \) and then plug-in in this \(z\left ( t\right ) \) in \(f\left ( z\right ) \) to obtain \(f\left ( t\right ) \), then the integral becomes \[ \int _{C}f\left ( z\right ) dz=\int _{t=0}^{t=1}f\left ( t\right ) z^{\prime }\left ( t\right ) dt \] And now evaluate this integral using normal integration rules. If the path is a circular arc, then no need to use \(t\), just use \(\theta \). Rewrite \(x=re^{i\theta }\) and use \(\theta \) instead of \(t\) and follow same steps as above.
But if \(h\left ( z_{0}\right ) =0\) then we need to apply L’Hopitals like this. If \(f\left ( z\right ) =\frac{\sin z}{1-z^{4}}\) and we want to ﬁnd residue at \(z=i\). Then do as above, but with extra step, like this \begin{align*} R\left ( i\right ) & =\lim _{z\rightarrow i}\left ( z-i\right ) \frac{\sin z}{1-z^{4}}\\ & =\left ( \lim _{z\rightarrow i}\sin z\right ) \left ( \lim _{z\rightarrow i}\left ( z-i\right ) \frac{1}{1-z^{4}}\right ) \\ & =\sin i\left ( \lim _{z\rightarrow i}\frac{\left ( z-i\right ) }{1-z^{4}}\right ) \qquad \text{Now apply L'Hopitals}\\ & =\sin i\left ( \lim _{z\rightarrow i}\frac{1}{-4z^{3}}\right ) \\ & =\frac{\sin i}{-4i^{3}}\\ & =\frac{1}{4}\sinh \left ( 1\right ) \end{align*}
Now if the pole is not a simple pole or order one,.say of order \(m\), then we ﬁrst multiply \(f\left ( z\right ) \) by \(\left ( z-z_{0}\right ) ^{m}\) then diﬀerentiate the result \(m-1\) times, then divide by \(\left ( m-1\right ) !\), and then evaluate the result at \(z=z_{0}.\) in other words, \[ R\left ( z_{0}\right ) =\lim _{z\rightarrow z_{0}}\frac{1}{\left ( m-1\right ) !}\frac{d^{m-1}}{dz^{m-1}}\left ( \left ( z-z_{0}\right ) ^{m}f\left ( z\right ) \right ) \] For example, if \(f\left ( z\right ) =\frac{z\sin z}{\left ( z-\pi \right ) ^{3}}\) and we want residue at \(z=\pi \). Since order is \(m=3\), then \begin{align*} R\left ( z_{0}\right ) & =\lim _{z\rightarrow \pi }\frac{1}{2!}\frac{d^{2}}{dz^{2}}\left ( \left ( z-\pi \right ) ^{3}\frac{z\sin z}{\left ( z-\pi \right ) ^{3}}\right ) \\ & =\lim _{z\rightarrow \pi }\frac{1}{2}\frac{d^{2}}{dz^{2}}\left ( z\sin z\right ) \\ & =\lim _{z\rightarrow \pi }\frac{1}{2}\left ( -z\sin z+2\cos z\right ) \\ & =-1 \end{align*}
The above methods will work on most of the HW problems I’ve seen so far but If all else fails, try Laurent series, that always works.
Then \[ \frac{df}{f}=\frac{1}{2}\frac{dx}{x}-\frac{3}{2}\frac{dy}{y}\] But \(\frac{dx}{x}\) and \(\frac{dy}{y}\) are the relative errors in \(x\) and \(y\). So if we plug-in \(2\) for \(\frac{dx}{x}\) and \(-2\) for \(\frac{dy}{y}\) we get \(4\%\) is worst relative error in \(f\left ( x,y\right ) \). Notice we used \(-2\%\) relative error for \(y\) and \(+2\%\) relative error for \(x\) since we wanted the worst (largest) relative error. If we wanted the least relative error in \(f\left ( x,y\right ) \), then we will use \(+2\%\) for \(y\) also, which gives \(1-3=-2\) or \(2\%\) relative error in \(f\left ( x,y\right ) \).
\(\blacksquare \) in Mathematica Exp is a symbol. Head[Exp] gives Symbol but in Maple it is not.
In Maple
gives \(\left \{ \Gamma ,\pi ,x,y,z\right \} \) but in Mathematica
gives \(\{e,x,\pi ,z,\text{Gamma},y\}\)
Notice that \(e\) shows up in Mathematica, but not in Maple.
(added December 13, 2018)
The PDE is \begin{equation} \frac{\partial ^{2}\psi }{\partial t^{2}}=c^{2}\frac{\partial ^{2}\psi }{\partial x^{2}} \tag{1} \end{equation} Let \begin{align*} u & =x-ct\\ v & =x+ct \end{align*}
Then\begin{align} \frac{\partial \psi }{\partial t} & =\frac{\partial \psi }{\partial u}\frac{\partial u}{\partial t}+\frac{\partial \psi }{\partial v}\frac{\partial v}{\partial t}\nonumber \\ & =-c\frac{\partial \psi }{\partial u}+c\frac{\partial \psi }{\partial v} \tag{2} \end{align}
And\begin{align} \frac{\partial \psi }{\partial x} & =\frac{\partial \psi }{\partial u}\frac{\partial u}{\partial x}+\frac{\partial \psi }{\partial v}\frac{\partial v}{\partial x}\nonumber \\ & =\frac{\partial \psi }{\partial u}+\frac{\partial \psi }{\partial v} \tag{3} \end{align}
Then, from (2)\begin{align} \frac{\partial ^{2}\psi }{\partial t^{2}} & =-c\left ( \frac{\partial ^{2}\psi }{\partial u^{2}}\frac{\partial u}{\partial t}+\frac{\partial ^{2}\psi }{\partial u\partial v}\frac{\partial v}{\partial t}\right ) +c\left ( \frac{\partial ^{2}\psi }{\partial v^{2}}\frac{\partial v}{\partial t}+\frac{\partial ^{2}\psi }{\partial v\partial u}\frac{\partial u}{\partial t}\right ) \nonumber \\ & =-c\left ( -c\frac{\partial ^{2}\psi }{\partial u^{2}}+c\frac{\partial ^{2}\psi }{\partial u\partial v}\right ) +c\left ( c\frac{\partial ^{2}\psi }{\partial v^{2}}-c\frac{\partial ^{2}\psi }{\partial v\partial u}\right ) \nonumber \\ & =c^{2}\frac{\partial ^{2}\psi }{\partial u^{2}}-c^{2}\frac{\partial ^{2}\psi }{\partial u\partial v}+c^{2}\frac{\partial ^{2}\psi }{\partial v^{2}}-c^{2}\frac{\partial ^{2}\psi }{\partial v\partial u}\nonumber \\ & =c^{2}\frac{\partial ^{2}\psi }{\partial u^{2}}+c^{2}\frac{\partial ^{2}\psi }{\partial v^{2}}-2c^{2}\frac{\partial ^{2}\psi }{\partial v\partial u} \tag{4} \end{align}
And from (3)\begin{align} \frac{\partial ^{2}\psi }{\partial x^{2}} & =\left ( \frac{\partial ^{2}\psi }{\partial u^{2}}\frac{\partial u}{\partial x}+\frac{\partial ^{2}\psi }{\partial u\partial v}\frac{\partial v}{\partial x}\right ) +\left ( \frac{\partial ^{2}\psi }{\partial v^{2}}\frac{\partial v}{\partial x}+\frac{\partial ^{2}\psi }{\partial v\partial u}\frac{\partial u}{\partial x}\right ) \nonumber \\ & =\left ( \frac{\partial ^{2}\psi }{\partial u^{2}}+\frac{\partial ^{2}\psi }{\partial u\partial v}\right ) +\left ( \frac{\partial ^{2}\psi }{\partial v^{2}}+\frac{\partial ^{2}\psi }{\partial v\partial u}\right ) \nonumber \\ & =\frac{\partial ^{2}\psi }{\partial u^{2}}+\frac{\partial ^{2}\psi }{\partial v^{2}}+2\frac{\partial ^{2}\psi }{\partial v\partial u} \tag{5} \end{align}
Substituting (4,5) into (1) gives\begin{align*} -2c^{2}\frac{\partial ^{2}\psi }{\partial v\partial u} & =2c^{2}\frac{\partial ^{2}\psi }{\partial v\partial u}\\ -4c^{2}\frac{\partial ^{2}\psi }{\partial v\partial u} & =0 \end{align*}
Since \(c\neq 0\) then\[ \frac{\partial ^{2}\psi }{\partial v\partial u}=0 \] Integrating w.r.t \(v\) gives\[ \frac{\partial \psi }{\partial u}=f\left ( u\right ) \] Integrating w.r.t \(u\)\[ \psi \left ( x,t\right ) =F\left ( u\right ) +G\left ( v\right ) \] Therefore\begin{equation} \psi \left ( x,t\right ) =F\left ( x-ct\right ) +G\left ( x+ct\right ) \tag{6} \end{equation} The functions \(F\left ( x,t\right ) ,G\left ( x,t\right ) \) are arbitrary functions found from initial and boundary conditions if given. Let initial conditions be\begin{align*} \psi \left ( x,0\right ) & =f_{0}\left ( x\right ) \\ \frac{\partial }{\partial t}\psi \left ( x,0\right ) & =g_{0}\left ( x\right ) \end{align*}
Where the ﬁrst condition above is the shape of the string at time \(t=0\) and the second condition is the initial velocity.
Applying ﬁrst condition to (6) gives\begin{equation} f_{0}\left ( x\right ) =F\left ( x\right ) +G\left ( x\right ) \tag{7} \end{equation} Applying the second condition gives\begin{align} g_{0}\left ( x\right ) & =\left [ \frac{\partial }{\partial t}F\left ( x-ct\right ) \right ] _{t=0}+\left [ \frac{\partial }{\partial t}G\left ( x+ct\right ) \right ] _{t=0}\nonumber \\ & =\left [ \frac{dF\left ( x-ct\right ) }{d\left ( x-ct\right ) }\frac{\partial \left ( x-ct\right ) }{\partial t}\right ] _{t=0}+\left [ \frac{dG\left ( x+ct\right ) }{d\left ( x+ct\right ) }\frac{\partial \left ( x+ct\right ) }{\partial t}\right ] _{t=0}\nonumber \\ & =\left [ -c\frac{dF\left ( x-ct\right ) }{d\left ( x-ct\right ) }\right ] _{t=0}+\left [ c\frac{dG\left ( x+ct\right ) }{d\left ( x+ct\right ) }\right ] _{t=0}\nonumber \\ & =-c\frac{dF\left ( x\right ) }{dx}+c\frac{dG\left ( x\right ) }{dx} \tag{8} \end{align}
Now we have two equations (7,8) and two unknowns \(F,G\) to solve for. But the (8) has derivatives of \(F,G\,\). So to make it easier to solve, we integrate (8) w.r.t. to obtain\begin{equation} \int ^{x}g_{0}\left ( s\right ) ds=-cF\left ( x\right ) +cG\left ( x\right ) \tag{9} \end{equation} So we will use (9) instead of (8) with (7) to solve for \(F,G\). From (7) \begin{equation} F\left ( x\right ) =f_{0}\left ( x\right ) -G\left ( x\right ) \tag{10} \end{equation} Substituting (10) in (9) gives\begin{align} \int ^{x}g_{0}\left ( s\right ) ds & =-c\left ( f_{0}\left ( x\right ) -G\left ( x\right ) \right ) +cG\left ( x\right ) \nonumber \\ & =-cf_{0}\left ( x\right ) +2cG\left ( x\right ) \nonumber \\ G\left ( x\right ) & =\frac{\left ( \int ^{x}g_{0}\left ( s\right ) ds\right ) +cf_{0}\left ( x\right ) }{2c}\nonumber \\ & =\frac{1}{2c}\left ( \int ^{x}g_{0}\left ( s\right ) ds+cf_{0}\left ( x\right ) \right ) \tag{11} \end{align}
Using the above back in (10) gives \(F\left ( x\right ) \) as\begin{equation} F\left ( x\right ) =f_{0}\left ( x\right ) -\frac{1}{2c}\left ( \int ^{x}g_{0}\left ( s\right ) ds+cf_{0}\left ( x\right ) \right ) \tag{12} \end{equation} Using (11,12) in (6) gives the ﬁnal solution\begin{align*} \psi \left ( x,t\right ) & =F\left ( x-ct\right ) +G\left ( x+ct\right ) \\ & =f_{0}\left ( x-ct\right ) -\frac{1}{2c}\left ( \int ^{x-ct}g_{0}\left ( s\right ) ds+cf_{0}\left ( x-ct\right ) \right ) +\frac{1}{2c}\left ( \int ^{x}g_{0}\left ( s\right ) ds+cf_{0}\left ( x\right ) \right ) \\ & =f_{0}\left ( x-ct\right ) -\frac{1}{2c}\int ^{x-ct}g_{0}\left ( s\right ) ds-\frac{1}{2}f_{0}\left ( x-ct\right ) +\frac{1}{2c}\int ^{x+ct}g_{0}\left ( s\right ) ds+\frac{1}{2}f_{0}\left ( x+ct\right ) \\ & =\frac{1}{2}\left ( f_{0}\left ( x-ct\right ) +f_{0}\left ( x-ct\right ) \right ) +\frac{1}{2c}\int _{x-ct}^{x+ct}g_{0}\left ( s\right ) ds \end{align*}
The above is the ﬁnal solution. So if we are given initial position and initial velocity of the string as function of \(x\), we can ﬁnd exact solution to the wave PDE.
Too many references used, but will try to remember to start recording books used from now on. Here is current list