User:Egm6322.s14.team2/report2

Report 2 PEA 2 Spring 2014 Team 2

Problems are given in lecture notes Sec66a

R2.1

Problem Statement

The atomic-resolution TEM figure on [p.66b-16 ] is not obvious (not easy to see and to understand); search the literature (web, pdf papers, pdf reports, etc.) for better (good quality) figure(s) that show the dislocation at atomic level with a burger vector. for example, start with an image web search for “dislocation atomic TEM burger vector”.

Find

(a): A clear image that shows the atomic dislocation Burger vector

Solution

On our honor, we did this assignment on our own, without looking at the solutions in previous semesters or other online solutions.

From looking through the experimental results of atomic scale imaging of dislocations it seems that a clear illustration of the burgers vector is difficult to find.

An imageof atomic scale edge dislocation is available from a paper published the journal Nature .

There is also a diagram of an atomic scale edge dislocation with a burgers vector notated on the image found on a website about solid mechanics. The image is small however.

R2.2 Roots of Associated Laguerre Polynomials

Problem Statement

Plot the associated Laguerre polynomials for $n=3$ and $\alpha =0,-1,-2,-3,-4$. Use MATLAB to find the roots of $L_{n}^{(-1)}(x)$ for $n=7,8$, and verify the numbers in the table on p.66a-7

Note: If you already have the plots made in report R1, you simply copy and paste these plots in report R2; the task then is only to find the roots of these polynomials.

Find

(a): Associated Laguerre polynomial plots

(b): Roots of Laguerre polynomials

(c): Verify numbers in the table on p.66-7

Solution

On our honor, we did this assignment on our own, without looking at the solutions in previous semesters or other online solutions.

Part (a): Associated Laguerre Polynomial Plots

The plot was created using the Python scripting language.

Figure 1: Associated Laguerre Polynomial $L_{3}^{(\alpha )}$ ([1])

Part (b): Roots of Associated Laguerre Polynomials

The roots of the associated Laguerre polynomial were computed using Wolfram Alpha. The syntax for creating associated Laguerre polynomials in Wolfram Alpha: LaguerreL[n,k,x] where n is the degree, k is the family of Laguerre polynomials, and x is the independent variable.

The Laguerre polynomial given by,

$L_{3}^{(0)}={\frac {1}{6}}(6-18x+9x^{2}-x^{3})$

(Eq 2-1)

Has three roots. The approximate roots of the polynomial are:
$x\approx 0.41577$
$x\approx 2.2943$
$x\approx 6.2899$

The Laguerre polynomial given by,

$L_{3}^{(-1)}={\frac {1}{6}}(-6x+6x^{2}-x^{3})$

(Eq 2-2)

Has three roots. The roots of the polynomial are:
$x=0$
$x=3-{\sqrt {(}}3)$
$x=3+{\sqrt {(}}3)$

The Laguerre polynomial given by,

$L_{3}^{(-2)}={\frac {1}{6}}(3x^{2}-x^{3})$

(Eq 2-3)

Has three roots. The roots of the polynomial are:
$x=0$
$x=0$
$x=3$

The Laguerre polynomial given by,

$L_{3}^{(-3)}=-{\frac {x^{3}}{6}}$

(Eq 2-4)

Has three roots. The roots of the polynomial are:
$x=0$
$x=0$
$x=0$

The Laguerre polynomial given by,

$L_{3}^{(-4)}={\frac {1}{6}}(-6-6x-3x^{2}-x^{3})$

(Eq 2-5)

Has three roots. The roots of the polynomial are:
$x\approx -1.5961$
$x\approx -0.7020-1.8073i$
$x\approx -0.7020+1.8073i$

Part (c): Two Sets of Roots of Associated Laguerre Polynomials For $n=7,8$

In this section we are asked to compute the roots of the associated Laguerre polynomials with $\alpha =-1{\text{ and }}n=7,8$ . We are then asked to compare these results to previously tabulated results.

The Laguerre polynomial given by,

$L_{7}^{(-1)}={\frac {-5040x+15120x^{2}-12600x^{3}+4200x^{4}-630x^{5}+42x^{6}-x^{7}}{5040}}$

(Eq 2-6)

The following roots:

Associated Laguerre Polynomial Roots of $L_{7}^{(-1)}$
Root #	Root
$x_{1}$	0
$x_{2}$	0.5277
$x_{3}$	1.796
$x_{4}$	3.877
$x_{5}$	6.919
$x_{6}$	11.235
$x_{7}$	17.65

The Laguerre polynomial given by,

$L_{8}^{(-1)}={\frac {-40320x+141120x^{2}-141120x^{3}+58800x^{4}-11760x^{5}+1176x^{6}-56x^{7}+x^{8}}{40320}}$

(Eq 2-7)

The following roots:

Associated Laguerre Polynomial Roots of $L_{8}^{(-1)}$
Root #	Root
$x_{1}$	0
$x_{2}$	0.4610
$x_{3}$	1.564
$x_{4}$	3.352
$x_{5}$	5.916
$x_{6}$	9.421
$x_{7}$	14.19
$x_{8}$	21.09

R2.3: Comparison of Laguerre Polynomial Power Series Representations

Given Information

(Eq 3-1) was given for the value of P_N(u) and (Eq 3-2) through (Eq 3-4) were given for the Laguerre polynomials written with the same variables as (Eq 3-1) for comparisons.

$P_{N}(u)=u+\sum _{n=2}^{N}{\dfrac {(-1)^{n-1}(N-1)(N-2)...(N-n+1)}{(n-1)!n!}}u^{n}$

(Eq 3-1)

$L_{N}^{(\alpha )}(u)=\sum _{n=0}^{N}{\dfrac {(\alpha +n+1)_{N-n}}{(N-n)!n!}}(-u)^{n}$

(Eq 3-2)

$L_{N}^{(\alpha )}(u)=\sum _{n=0}^{N}(-1)^{n}{N+\alpha \choose N-n}{\dfrac {1}{n!}}u^{n}$

(Eq 3-3)

$L_{N}^{(\alpha )}(u)=\sum _{n=0}^{N}(-u)^{n}{N \choose n}\prod _{k=n+1}^{N}(\alpha +k)$

(Eq 3-4)

Problem Statement

Compare (Eq 3-1) to the various finite power series expressions for the associated Laguerre polynomials in (Eq 3-2) through (Eq 3-4). For (Eq 3-2), read about the Pochhammer symbol in King 2003, p.60 ( there is a misprint in the name Pochhammer in King 2003; "hh" not "kh"). Read also other refs on the Pochhammer symbol such as Wikipedia, MathWorld,etc. The problem statement is also given in R2.3: sec.66a, Pb-66-6.

Solution

On our honor, we did this assignment on our own, without looking at the solutions in previous semesters or other online solutions.

Part A

The falling factorial term in (Eq 3-1) can be rewritten with a Pochhammer symbol.

$P_{N}(u)=u+\sum _{n=2}^{N}{(-1)^{n-1}{\frac {(N)_{n}}{N\,n!\,(n-1)!}}u^{n}},N\geq 2$

(Eq 3-5)

Where,

$(N)_{n}=N(N-1)(N-2)...(N-n+1)$

(Eq 3-6)

The Pochhammer symbol can be related to a binomial coefficient by the equation:

${\frac {(N)_{n}}{n!}}={\binom {N}{n}}$

(Eq 3-7)

Plugging (Eq 3-7) into (Eq 3-5) results in:

$P_{N}(u)=u+\sum _{n=2}^{N}{\binom {N}{n}}{\frac {(-1)^{n-1}}{N(n-1)!}}u^{n}$

(Eq 3-8)

The binomial coefficient has the following closed form for integer arguments:

${\binom {N}{n}}=\left\{{\begin{matrix}{\frac {N!}{n!\,(N-n)!}}&0\leq n\leq N\\0&{\text{otherwise}}\end{matrix}}\right.$

(Eq 3-9)

Substituting the closed form in (Eq 3-9) into (Eq 3-8) gives:

$P_{N}(u)=u+\sum _{n=2}^{N}{\frac {(N-1)!\,(-1)^{n-1}}{(N-n)!\,(n-1)!\,n!}}u^{n}$

(Eq 3-10)

Factoring out a $-1$ and changing the reference for the $(-1)^{n-1}$ to $(-1)^{n}$ and making $(u)^{n}$ to $(-u)^{n}$ gives:
I do not believe this is correct because you are factoring out a single negative 1, not (-1)^n. This equation is always positive where n is an integer.

$P_{N}(u)=u+\sum _{n=2}^{N}{\frac {(N-1)!\,(-1)^{n}}{(N-n)!\,(n-1)!\,n!}}(-u)^{n}$

(Eq 3-11)

This is what I get...

$P_{N}(u)=u-\sum _{n=2}^{N}{\frac {(N-1)!\,(-u)^{n}}{(N-n)!\,(n-1)!\,n!}}$

(Eq 3-11)

Taking out the first term two terms of (Eq 3-2) (the first term is 0) with $\alpha =-1$ gives:

$L_{N}^{(-1)}(u)=u+\sum _{n=2}^{N}{\frac {(N-1)!\,(-1)^{n}}{(N-n)!\,(n-1)!\,n!}}(-u)^{n}$

(Eq 3-12)

And when I change the range of LN, I get -1 as the coefficient of u.

$L_{N}^{(-1)}(u)=-u+\sum _{n=2}^{N}{\frac {(N-1)!}{(N-n)!(n-1)!n!}}(-u)^{n}$

(Eq 3-12)

It should be noted that the summation index in (Eq 3-12) changes from $n=0$ to $n=1$. This is because when $n=0$ the binomial coefficient is zero and cannot be represented by the factorial expression. This term must be brought outside of the sum, but since the term is zero, all that results is an index change.

Part B

For the the case where $\alpha =-1$, (Eq 3-3) has the form:

$L_{N}^{(-1)}(u)=\sum _{n=0}^{N}{\binom {N-1}{N-n}}{\frac {(-1)^{n}}{n!}}(-u)^{n}$

(Eq 3-13)

Using the expression (Eq 3-7) for the binomial term in (Eq 3-13) yields:

${\binom {N-1}{N-n}}=\left\{{\begin{matrix}{\frac {(N-1)!}{(N-n)!\,(n-1)!}}&0\leq N-n\leq N-1\\0&{\text{otherwise}}\end{matrix}}\right.$

(Eq 3-14)

Plugging (Eq 3-14) into (Eq 3-13) gives:

$L_{N}^{(-1)}(u)=\sum _{n=1}^{N}{\frac {(N-1)!\,(-1)^{n}}{(N-n)!\,(n-1)!\,n!}}(-u)^{n}$

(Eq 3-15)

The first term of (Eq 3-13) is brought outside the summation. To reproduce (Eq 3-12)

$L_{N}^{(-1)}(u)=u+\sum _{n=2}^{N}{\frac {(N-1)!\,(-1)^{n}}{(N-n)!\,(n-1)!\,n!}}(-u)^{n}$

(Eq 3-16)

Part C

Making $\alpha =-1$ for (Eq 3-4) gives

$L_{N}^{(\alpha )}(u)=\sum _{n=0}^{N}(-u)^{n}{N \choose n}\prod _{k=n+1}^{N}(-1+k)$

(Eq 3-17)

$\prod _{k=n+1}^{N}(k-1)={\frac {(N-1)!}{n!}}$

(Eq 3-18)

Substituting the closed form of (Eq 3-9) and (Eq 3-18) into (Eq 3-17) gives:

$L_{N}^{(-1)}(u)=\sum _{n=0}^{N}{\dfrac {N!(N-1)!}{(N-n)!(n)!n!}}(-u)^{n}$

(Eq 3-19)

R2.4:Differential Element in Spherical Coordinates

Given Information

The following equation is given for ds²:

$ds^{2}=\sum _{i=1}^{3}(dx_{i})^{2}=h_{1}^{2}dr^{2}+h_{2}^{2}d{\theta }^{2}+h_{3}^{2}d{\varphi }^{2}$

(Eq 4-1)

x_i is written in spherical coordinates in (Eq 4-2) through (Eq 4-4).

$x_{1}=rcos{\theta }cos{\varphi }$

(Eq 4-2)

$x_{2}=rcos{\theta }sin{\varphi }$

(Eq 4-3)

$x_{3}=rsin{\theta }$

(Eq 4-4)

Problem Statement

Show that the infinitesimal length ds² in (Eq 4-1) can be written in spherical coordinates, (Eq 4-5), as given in Pb.R*7.3 on pp.39-[1,3] of sec.39 of the notes.

$ds^{2}=dr^{2}+r^{2}d{\theta }^{2}+r^{2}cos^{2}{\theta }d{\varphi }^{2}$

(Eq 4-5)

Solution

On our honor, we did this assignment on our own, without looking at the solutions in previous semesters or other online solutions.

The derivatives of x_i are given in (Eq 4-6) through (Eq 4-8).

$dx_{1}=cos{\theta }cos{\varphi }dr-rsin{\theta }cos{\varphi }d{\theta }-rcos{\theta }sin{\varphi }d{\varphi }$

(Eq 4-6)

$dx_{2}=cos{\theta }sin{\varphi }dr-rsin{\theta }sin{\varphi }d{\theta }+rcos{\theta }cos{\varphi }d{\varphi }$

(Eq 4-7)

$dx_{3}=sin{\theta }dr+rcos{\theta }d{\theta }$

(Eq 4-8)

(Eq 4-9) was developed by substituting the squared derivatives in the summation of (Eq 4-1) and combining like terms.

$ds^{2}=(sin^{2}{\theta }(cos^{2}{\varphi }+sin^{2}{\varphi })+cos^{2}{\theta })dr^{2}+(2rsin{\theta }cos{\theta }-2rsin{\theta }cos{\theta }(sin^{2}{\varphi }+cos^{2}{\varphi }))drd{\theta }+(r^{2}(sin^{2}{\theta }(sin^{2}{\varphi }+cos^{2}{\varphi })+cos^{2}{\theta })d{\theta }^{2}+(r^{2}cos^{2}{\theta }(sin^{2}{\varphi }+cos^{2}{\varphi }))d{\varphi }^{2}$

(Eq 4-9)

Using the Pythagorean Trigonometric Identity, (Eq 4-9) reduces to (Eq 4-10).

$ds^{2}=dr^{2}+r^{2}d{\theta }^{2}+r^{2}cos^{2}d{\varphi }^{2}$

(Eq 4-10)

The magnitude of the tangent vector h are then identified as follows:

$h_{1}=1$

(Eq 4-11)

$h_{2}=r$

(Eq 4-12)

$h_{3}=rcos{\theta }$

(Eq 4-13)

The definition of the Laplace operator, ${\triangle }u$, in cartesian coordinates:

${\triangle }u=\sum _{i=1}^{3}{\frac {{\partial }u}{{\partial }x_{i}}}$

(Eq 4-14)

The Laplace operator can be expressed in spherical coordinates from equation (2) from the class notes 39-2,

${\triangle }u={\frac {1}{h_{1}h_{2}h_{3}}}({\frac {\partial }{{\partial }r}}({\frac {h_{1}h_{2}h_{3}}{h_{1}^{2}}}{\frac {{\partial }u}{{\partial }r}})+{\frac {\partial }{{\partial }{\theta }}}({\frac {h_{1}h_{2}h_{3}}{h_{2}^{2}}}{\frac {{\partial }u}{{\partial }{\theta }}})+{\frac {\partial }{{\partial }{\varphi }}}({\frac {h_{1}h_{2}h_{3}}{h_{3}^{2}}}{\frac {{\partial }u}{{\partial }{\varphi }}}))$

(Eq 4-15)

After substituting the values of h from (Eq 4-11) through (Eq 4-13), the Laplacian takes the form:

${\triangle }u={\frac {1}{r^{2}cos{\theta }}}({\frac {\partial }{{\partial }r}}(r^{2}cos{\theta }{\frac {{\partial }u}{{\partial }r}})+{\frac {\partial }{{\partial }{\theta }}}(cos{\theta }{\frac {{\partial }u}{{\partial }{\theta }}})+{\frac {\partial }{{\partial }{\varphi }}}({\frac {1}{cos{\theta }}}{\frac {{\partial }u}{{\partial }{\varphi }}}))$

(Eq 4-16)

Factoring the $cos{\theta }$ from the first term, the partial derivative with respect to r, and the $sec{\theta }$ from the third term, the partial derivative with respect to ${\varphi }$, since both values are fixed in the partial derivatives, the Laplacian takes the form:

${\triangle }u={\frac {1}{r^{2}}}{\frac {\partial }{{\partial }r}}(r^{2}{\frac {{\partial }u}{{\partial }r}})+{\frac {1}{r^{2}cos{\theta }}}{\frac {\partial }{{\partial }{\theta }}}(cos{\theta }{\frac {{\partial }u}{{\partial }{\theta }}})+{\frac {1}{r^{2}cos^{2}{\theta }}}{\frac {{\partial }^{2}u}{{\partial }{{\varphi }^{2}}}}$

(Eq 4-17)

R2.5: Using Variation of Parameters To Find Second Solution To The Legendre Differential Equation

Given Information

The following equation is given for P₂(x), the Legendre Polynomial of degree n=2:

$P_{2}(x)={\dfrac {1}{2}}(3x^{2}-1)$

(Eq 5-1)

Problem Statement

Using variation of parameters, show that Equation 5-2 is a second solution to the Legendre differential equation in Equation 5-3 as given in Pb.R*6.11 on pp.37-4 in sec.37 of the notes.

$Q_{2}(x)={\dfrac {1}{4}}(3x^{2}-1)log({\dfrac {1+x}{1-x}})-{\dfrac {3}{2}}x$

(Eq 5-2)

Solution

On our honor, we did this assignment on our own, without looking at the solutions in previous semesters or other online solutions.

The Legendre differential equation of degree n=2 is given in equation 5-3.

$(1-x^{2})y''-2xy'+6y=0$

(Eq 5-3)

Using reduction of order, some function, v, times P₂(x) and the first two derivatives of the product are substituted into the Legendre differential equation in equation 5-7. Furthermore, since P₂(x) is a solution to the Legendre equation the coefficient of the first order term, v, reduces to zero. Since the lowest order term cancels out, the second order differential equation is reduced to a first order differential equation where w=v' and w'=v'' in Equation 5-8 which can be solved through separation of variables.

$R_{2}(x)=v({\dfrac {3}{2}}x^{2}-{\dfrac {1}{2}})$

(Eq 5-4)

$R'_{2}(x)=v'({\dfrac {3}{2}}x^{2}-{\dfrac {1}{2}})+v(3x)$

(Eq 5-5)

$R''_{2}(x)=v''({\dfrac {3}{2}}x^{2}-{\dfrac {1}{2}})+6v'x+3v$

(Eq 5-6)

$(1-x^{2})({\dfrac {3}{2}}x^{2}-{\dfrac {1}{2}})v''+(6x(1-x^{2})-2x({\dfrac {3}{2}}x^{2}-{\dfrac {1}{2}}))v'+(3(1-x^{2})-2x(3x)+6({\dfrac {3}{2}}x^{2}-{\dfrac {1}{2}}))v=0$

(Eq 5-7)

${\dfrac {dw}{w}}={\dfrac {(6x(1-x^{2})-2x({\dfrac {3}{2}}x^{2}-{\dfrac {1}{2}}))dx}{(1-x^{2})({\dfrac {3}{2}}x^{2}-{\dfrac {1}{2}})}}$

(Eq 5-7)

Through partial fraction decomposition, the integral is

$\int {\dfrac {dw}{w}}=\int {\dfrac {-6x}{{\dfrac {3}{2}}x^{2}dx-{\dfrac {1}{2}}}}+\int {\dfrac {2x}{1-x^{2}}}dx$

(Eq 5-8)

Through integration

$log(w)=-2log({\dfrac {3}{2}}x^{2}-{\dfrac {1}{2}})-log(1-x^{2})$

(Eq 5-9)

Using properties of the natural logarithm and then inverting it using the exponential function, equation 5-9 becomes

$w={\dfrac {1}{({\dfrac {3}{2}}x^{2}-{\dfrac {1}{2}})^{2}(1-x^{2})}}$

(Eq 5-10)

Substituting v' back into Equation 5-10 and again using partial fraction decomposition

$\int {dv}={\frac {1}{2}}\int {\dfrac {dx}{(1-x)}}+{\frac {1}{2}}\int {\dfrac {dx}{(1+x)}}+3\int {\dfrac {(3x^{2}+1)}{(3x^{2}-1)^{2}}}dx$

(Eq 5-11)

Using integration by parts

$\int {udv}=uv-\int {vdu}$

(Eq 5-12)

where u and v are defined as follows:

$dv=x(3x^{2}-1)^{-2}$

(Eq 5-13)

$u={\frac {(3x^{2}+1)}{x}}$

(Eq 5-14)

$v={\frac {-1}{2}}log(1-x)+{\frac {1}{2}}log(1+x)-3({\dfrac {(3x^{2}+1)}{6x(3x^{2}-1)}}-\int {\frac {dx}{6x^{2}}})$

(Eq 5-15)

Using the properties of the natural logarithm to combine the first two terms and performing the simple integration on the last term:

$v={\dfrac {1}{2}}log({\dfrac {1+x}{1-x}})-{\dfrac {3x}{3x^{2}-1}}$

(Eq 5-16)

Substituting v into Equation 5-4 then identifying the new solution as Q₂

$Q_{2}(x)={\dfrac {1}{4}}(3x^{2}-1)log({\dfrac {1+x}{1-x}})-{\dfrac {3x}{2}}$

(Eq 5-17)

R2.6: Derive Lapace's Equation in Cylindrical Coordinates

Given Information

The following equations are given for cylindrical coordinates:

$x=r\cos \theta =\xi _{1}\,\cos {\xi _{2}}$

(Eq 6-1)

$y=r\sin \theta =\xi _{1}\,\sin {\xi _{2}}$

(Eq 6-2)

$z=\xi _{3}$

(Eq 6-3)

Problem Statement

(a) Find

$\left\{dx_{i}\right\}=\left\{dx_{1},dx_{2},dx_{3}\right\}$

(Eq 6-4)

in terms of

$\left\{\xi _{j}\right\}=\left\{\xi _{1},\xi _{2},\xi _{3}\right\}$

(Eq 6-5)

and

$\left\{d\xi _{k}\right\}=\left\{d\xi _{1},d\xi _{2},d\xi _{3}\right\}$

(Eq 6-6)

(b) Find

$ds^{2}=\sum _{i}(dx_{i})^{2}=\sum _{k}(h_{k})^{2}(d\xi _{k})^{2}$

(Eq 6-7)

Identify $\{h_{i}\}$ in terms of $\{\xi _{i}\}$

(c) Find $\Delta u$ in terms of cylindrical coordinates

(d) Use separation of variable to find the separated equations of Laplace's equation in cylindrical coordinates and compare to the Bessel equation:

$x^{2}{\frac {d^{2}y}{dx^{2}}}+x{\frac {dy}{dx}}+(x^{2}-\alpha ^{2})y=0$

(Eq 6-8)

Solution

On our honor, we consulted the solution of Dr. Gaely of Concordia College (http://wwwp.cord.edu/faculty/gealy/physics315/SepVarsCyl.pdf) to motivate the change of variables in the solution of part d. We improved upon the referenced solution by explicitly showing the functional arguments as a result of the change of variables.

Part a

The total differential of $x$ can be written as:

$dx={\frac {\partial x}{\partial \xi _{1}}}d\xi _{1}+{\frac {\partial x}{\partial \xi _{2}}}d\xi _{2}+{\frac {\partial x}{\partial \xi _{3}}}d\xi _{3}$

(Eq 6-9)

Similar expressions can be written for $y$ and $z$.

$dy={\frac {\partial y}{\partial \xi _{1}}}d\xi _{1}+{\frac {\partial y}{\partial \xi _{2}}}d\xi _{2}+{\frac {\partial y}{\partial \xi _{3}}}d\xi _{3}$

(Eq 6-10)

$dz={\frac {\partial z}{\partial \xi _{1}}}d\xi _{1}+{\frac {\partial z}{\partial \xi _{2}}}d\xi _{2}+{\frac {\partial z}{\partial \xi _{3}}}d\xi _{3}$

(Eq 6-11)

Plugging in the expressions from (Eq 6-1), (Eq 6-2) and (Eq 6-3):

$dx=\cos {\xi _{2}}\,d\xi _{1}-\xi _{1}\,\sin {\xi _{2}}\,d\xi _{2}$

(Eq 6-12)

$dy=\sin {\xi _{2}}\,d\xi _{1}+\xi _{1}\,\cos {\xi _{2}}\,d\xi _{2}$

(Eq 6-13)

$dz=d\xi _{3}$

(Eq 6-14)

Part b

Plugging in (Eq 6-12), (Eq 6-13), and (Eq 6-14) into (Eq 6-7).

$ds^{2}=(h_{1})^{2}(\cos {\xi _{2}}\,d\xi _{1}-\xi _{1}\sin {\xi _{2}}\,d\xi _{2})^{2}+(h_{2})^{2}(\sin {\xi _{2}}\,d\xi _{1}+\xi _{1}\cos {\xi _{2}}\,d\xi _{2})^{2}+(h_{3})^{2}(d\xi _{3})^{2}$

(Eq 6-15)

Expanding (Eq 6-14) out, and applying the following trig identity:

$\sin ^{2}{\xi _{2}}+\cos ^{2}{\xi _{2}}=1$

(Eq 6-16)

Results in the following expression:

$ds^{2}=(d\xi _{1})^{2}+(\xi _{1})^{2}(d\xi _{2})^{2}+(d\xi _{3})^{2}$

(Eq 6-17)

Comparing (Eq 6-17) to the left hand side of (Eq 6-7):

$h_{1}=1$

(Eq 6-18)

$h_{2}=\xi _{1}$

(Eq 6-19)

$h_{3}=1$

(Eq 6-20)

Part c

The Laplacian operator for a scalar function $u$ in general curvilinear coordinates is defined as:

$\Delta u={\frac {1}{h_{1}h_{2}h_{3}}}\left[{\frac {\partial }{\partial \xi _{1}}}\left({\frac {h_{2}h_{3}}{h_{1}}}{\frac {\partial }{\partial \xi _{1}}}\right)+{\frac {\partial }{\partial \xi _{2}}}\left({\frac {h_{1}h_{3}}{h_{2}}}{\frac {\partial }{\partial \xi _{2}}}\right)+{\frac {\partial }{\partial \xi _{3}}}\left({\frac {h_{1}h_{2}}{h_{3}}}{\frac {\partial }{\partial \xi _{3}}}\right)\right]u$

(Eq 6-21)

Plugging in (Eq 6-18), (Eq 6-19) and (Eq 6-20).

$\Delta u={\frac {1}{\xi _{1}}}\left[{\frac {\partial }{\partial \xi _{1}}}\left(\xi _{1}{\frac {\partial }{\partial \xi _{1}}}\right)+{\frac {\partial }{\partial \xi _{2}}}\left({\frac {1}{\xi _{1}}}{\frac {\partial }{\partial \xi _{2}}}\right)+{\frac {\partial }{\partial \xi _{3}}}\left(\xi _{1}{\frac {h_{1}h_{2}}{h_{3}}}{\frac {\partial }{\partial \xi _{3}}}\right)\right]u$

(Eq 6-22)

$\Delta u={\frac {1}{\xi _{1}}}{\frac {\partial u}{\partial \xi _{1}}}+{\frac {\partial ^{2}u}{\partial \xi _{1}^{2}}}+{\frac {1}{\xi _{1}^{2}}}{\frac {\partial ^{2}u}{\partial \xi _{2}^{2}}}+{\frac {\partial ^{2}u}{\partial \xi _{3}^{2}}}$

(Eq 6-23)

$\Delta u={\frac {1}{r}}{\frac {\partial u}{\partial r}}+{\frac {\partial ^{2}u}{\partial r^{2}}}+{\frac {1}{r^{2}}}{\frac {\partial ^{2}u}{\partial \theta ^{2}}}+{\frac {\partial ^{2}u}{\partial z^{2}}}$

(Eq 6-24)

Part d

Laplace's equation in cylindrical coordinates is:

$\Delta u=0$

(Eq 6-25)

${\frac {1}{r}}{\frac {\partial u}{\partial r}}+{\frac {\partial ^{2}u}{\partial r^{2}}}+{\frac {1}{r^{2}}}{\frac {\partial ^{2}u}{\partial \theta ^{2}}}+{\frac {\partial ^{2}u}{\partial z^{2}}}=0$

(Eq 6-26)

Assume a solution to (Eq 6-26) has the following form:

$u(r,\theta ,z)=R(r)\Theta (\theta )Z(z)$

(Eq 6-27)

Plugging (Eq 6-27) into (Eq 6-26):

$\Theta (\theta )Z(z){\frac {1}{r}}{\frac {dR(r)}{dr}}+\Theta (\theta )Z(z){\frac {d^{2}R(r)}{dr^{2}}}+R(r)Z(z){\frac {1}{r^{2}}}{\frac {d^{2}\Theta (\theta )}{d\theta ^{2}}}+R(r)\Theta (\theta ){\frac {d^{2}Z(z)}{dz^{2}}}=0$

(Eq 6-28)

Dividing by $u$

${\frac {1}{rR(r)}}{\frac {dR(r)}{dr}}+{\frac {1}{R(r)}}{\frac {d^{2}R(r)}{dr^{2}}}+{\frac {1}{\Theta (\theta )}}{\frac {1}{r^{2}}}{\frac {d^{2}\Theta (\theta )}{d\theta ^{2}}}+{\frac {1}{Z(z)}}(\theta ){\frac {d^{2}Z(z)}{dz^{2}}}=0$

(Eq 6-29)

The terms in z are then "seperated" and set equal to a constant: $\lambda ^{2}$.

${\frac {1}{rR(r)}}{\frac {dR(r)}{dr}}+{\frac {1}{R(r)}}{\frac {d^{2}R(r)}{dr^{2}}}+{\frac {1}{\Theta (\theta )}}{\frac {1}{r^{2}}}{\frac {d^{2}\Theta (\theta )}{d\theta ^{2}}}+\lambda ^{2}=0$

(Eq 6-30)

Multiplying by $r^{2}$:

${\frac {r}{R(r)}}{\frac {dR(r)}{dr}}+{\frac {r^{2}}{R(r)}}{\frac {d^{2}R(r)}{dr^{2}}}+{\frac {1}{\Theta (\theta )}}{\frac {d^{2}\Theta (\theta )}{d\theta ^{2}}}+(\lambda r)^{2}=0$

(Eq 6-31)

Now the $\theta$ term can be separated and set to $-\alpha ^{2}$.

${\frac {r}{R(r)}}{\frac {dR(r)}{dr}}+{\frac {r^{2}}{R(r)}}{\frac {d^{2}R(r)}{dr^{2}}}+-\alpha ^{2}+(\lambda r)^{2}=0$

(Eq 6-32)

Making a change of variable :

$\rho (r)=\lambda r$

(Eq 6-33)

${\frac {dR(r)}{dr}}=\lambda {\frac {dR(\rho )}{d\rho }}$

(Eq 6-34)

${\frac {d^{2}R(r)}{dr^{2}}}=\lambda ^{2}{\frac {d^{2}R(\rho )}{d\rho ^{2}}}$

(Eq 6-35)

Plugging the above change of variable into (Eq 6-32) and multiplying by $R(\rho )$ gives:

$\rho ^{2}{\frac {d^{2}R(\rho )}{dr^{2}}}+\rho {\frac {dR(\rho )}{dr}}+(\rho ^{2}-\alpha ^{2})R(\rho )=0$

(Eq 6-36)

Which is Bessel's equation given in (Eq 6-8).

R2.7:Testing Gauss-Laguerre Quadrature

Given Information

The following equations are given:

$e^{-2}\int _{2}^{\infty }{\frac {1}{x\log ^{2}(x)}}\,dx$

(Eq 7-1)

$J_{0}(z)={\frac {1}{\pi }}\int _{0}^{\pi }\cos(-z\sin t)\,dt={\frac {2}{\pi }}\int _{0}^{\infty }\sin(z\cosh t)\,dt$

(Eq 7-2)

(Eq 7-1) is motivated by quadrature tests written by Dr. John Burkardt found here:

http://people.sc.fsu.edu/~%20jburkardt/c_src/test_int_laguerre/test_int_laguerre.html

Problem Statement

(a) Find the values of the integral of the function (7-1) by other means (exact or approximate), then evaluate these integrals with the Gauss-Laguerre quadrature for comparison.

(b) Find and appropriate transformation to show the equality of (7-2). Find the value of $J_{0}(1/2)$ in some table and use the Gauss-Laguerre quadrature to find an approximate value of the integral to compare.

Solution

On our honor, we did this assignment on our own, without looking at the solutions in previous semesters or other online solutions.

Part a

$e^{-2}\int _{2}^{\infty }{\frac {1}{x\log ^{2}(x)}}\,dx=e^{-2}{\frac {1}{\log {2}}}$

Gauss-Laguerre quadrature is a Gaussian quadrature that integrates a function $g(x)$ over the interval $[0,\infty ]$ with the weighting factor $e^{-x}$.

$\int _{0}^{+\infty }e^{-x}g(x)\,dx\approx \sum _{i=1}^{n}w_{i}\,g(x_{i})$

(Eq 7-3)

Since the integral in (Eq 7-3) has the weighting term $e^{-x}$ in it, a transformation can be applied to integrate an arbitrary function of interest.

$f(x)=e^{-x}\left(e^{x}f(x)\right)=e^{-x}g(x)$

(Eq 7-4)

The integration points $x_{i}$ are found to be the roots of the Laguerre polynomial $P_{n}$, and the weights $w_{i}$ are found by the following expression :

$w_{i}={\frac {x_{i}}{(n+1)^{2}[L_{n+1}(x_{i})]^{2}}}.$

(Eq 7-5)

Since the integral in (Eq 7-1) has a lower limit of 2 instead of zero, a change of variable is needed.

$x=y-a$

(Eq 7-6)

(Eq 7-3) then becomes :

$\int _{a}^{+\infty }e^{-y}g(y)\,dy=\int _{0}^{+\infty }e^{-(x+a)}g(x+a)\,dx\approx e^{-a}\sum _{i=1}^{n}w_{i}\,g(x_{i}+a)$

(Eq 7-7)

A numerical python (NumPy) script was written to test the quadrature for the integral in (Eq 7-1). The source code is given below:

Source:

#################################################################################
# Name: gausslagtest.py
# Author: Cameron Stewart
# Date: 2/5/2014
# Purpose: This script uses Gauss-Laguerre quad  to integrate the integral:
#     exp(-2) Integral (2 to infinity) (1/(x*log(x)^2)
#################################################################################
import numpy as np
import matplotlib.pyplot as plt
#################################################################################
# Function name: TestFunction
# Input: a value, or array of values of x
# Output: g(x) at given points
# f(x)=exp(-2)/(x*log(x)*log(x); g(x)=f(x)*(e^x)
# Function f(x) that is desired to be integrated this function is multiplied by
# e^x to account for the weighting function e^-x in G-L Quad
#################################################################################
def TestFunction(x):
    y=np.exp(x)*1/(x*np.log(x)*np.log(x))
    return y

#################################################################################
# Uses the numpy Laguerre Gauss quadrature routine:
# numpy.polynomial.laguerre.laggauss(deg)  to find roots and weight points.
# Output is tuple of 1d arrays for the weights and integration points.  
# Note: Not tested for deg>100
#################################################################################

weight=(np.polynomial.laguerre.laggauss(4))
# weight[0] contains the roots/integration points
# weight[1] contains the weights at said points

# COV: y=x+a to account for lower limit of 
# integration
a=2

func=TestFunction(weight[0]+a)
integ=weight[1]*func
integ=np.exp(-a)*np.sum(integ)
print('n=4')
print(np.exp(-2)*integ)

weight=(np.polynomial.laguerre.laggauss(16))
func=TestFunction(weight[0]+a)
integ=weight[1]*func
integ=np.exp(-a)*np.sum(integ)
print('n=16')
print(np.exp(-2)*integ)

weight=(np.polynomial.laguerre.laggauss(100))

func=TestFunction(weight[0]+a)
integ=weight[1]*func
integ=np.exp(-a)*np.sum(integ)
print('n=100')
print(np.exp(-2)*integ)

# Exact value computed from analytic solution
exact=np.exp(-2)*1./(np.log(2))
print('exact')
print(exact)

#################################################################################
# Plot to show the convergence with increasing integration points
#################################################################################
t=np.linspace(2,100,99)
error=[]
for i in range(2, 101):
    weight=(np.polynomial.laguerre.laggauss(i))
    func=TestFunction(weight[0]+a)
    integ=weight[1]*func
    integ=np.exp(-a)*np.sum(integ)
    error.append((np.abs(np.exp(-2)*integ-exact))\
                     /exact*100.)
plt.plot(t,error)
plt.title('Convergence test for Gauss-Laguerre quadrature: logarithmic integral')
plt.xlabel(r'Integration points: $n$')
plt.ylabel(r'Error %')
plt.show()

Script terminal output:

n=4
0.145107525093
n=16
0.162236347114
n=100
0.172558662198
exact
0.195247541983

Plot generated by script:

Figure 7.1: Convergence of ( 7-1) with increasing integration points

The above script uses NumPy's Laguerre module quadrature routine: numpy.polynomial.laguerre.laggauss(deg). This routine solves the roots of the Laguerre polynomials by first finding the eigenvalues of the Companion Matrix, and then applying an iteration of Newton's Method to improve accuracy.

Documentation

source

The above script uses ( 7-5) to calculate ( 7-1).

$e^{-2}\int _{2}^{+\infty }e^{-y}\left(e^{y}{\frac {1}{y\log ^{2}(y)}}\right)\,dy\approx e^{-2}\left(e^{-2}\sum _{i=1}^{n}w_{i}\,g(x_{i}+2)\right)$

(Eq 7-8)

Where,

$g(x_{i}+2)=e^{x_{i}+2}{\frac {1}{(x_{i}+2)\log ^{2}(x_{i}+2)}}$

(Eq 7-9)

The script above compares favourably with values calculated by Dr. Burkardt .

Table 7.1 Comparison of quadrature calculations

n	Python Script	Dr. Burkardt
4	0.145107525093	0.145108
16	0.162236347114	0.162236

As can be seen in Figure 1, even at 100 integration points there is still significant error of around 12 percent.

Part b

The equivalence of ( 7-2) is not shown here. A suitable transformation has been elusive. It is the Author's belief that this relation was originally developed via contour integration, and is not trivially shown

( 7-2) can be evaluated at Wolfram Alpha.

$J_{0}(1/2)=0.9384698072408129...$

(Eq 7-10)

To integrate ( 7-2) by Gauss-Laguerre quadrature, ( 7-3) has the following form:

${\frac {2}{\pi }}\int _{0}^{\infty }\sin(z\cosh t)\,dt\approx {\frac {2}{\pi }}\sum _{i=1}^{n}w_{i}\,g(t_{i})$

(Eq 7-11)

Where,

$f(t)=e^{t}\sin {\left({\frac {\cosh {t}}{2}}\right)}=e^{t}\sin {\left({\frac {e^{t}+e^{-t}}{4}}\right)}$

(Eq 7-12)

A NumPy script similar to the one used in part a was used to integrate ( 7-2). A figure similar to Figure 7-1 was made testing the percent error as the number of integration points are increased.

Source:

#################################################################################                                    
# Name: gausslagbessel.py                                                                                            
# Author: Cameron Stewart                                                                                            
# Date: 2/5/2014                                                                                                     
# Purpose: This script uses Gauss-Laguerre quad  to integrate the integral:                                          
#     2/pi Integral (0 to infinity) sin( 0.5 cosh(x) )                                                               
#################################################################################                                    
import numpy as np
import matplotlib.pyplot as plt
from scipy.special import jn
#################################################################################                                    
# Function name: TestFunction                                                                                        
# Input: a value, or array of values of x                                                                            
# Output: g(x) at given points                                                                                       
# f(x)=sin( 0.5 * cosh(x) )                                                                                          
# Function f(x) that is desired to be integrated this function is multiplied by                                      
# e^x to account for the weighting function e^-x in G-L Quad                                                         
#################################################################################                                    
def TestFunction(x):
    y=np.exp(x)*np.sin(0.5*np.cosh(x))
    return y

#################################################################################                                    
# Uses the numpy Laguerre Gauss quadrature routine:                                                                  
# numpy.polynomial.laguerre.laggauss(deg)  to find roots and weight points.                                          
# Output is tuple of 1d arrays for the weights and integration points.                                               
# Note: Not tested for deg>100                                                                                       
#################################################################################                                    
# Exact value computed from analytic solution                                                                        
exact=jn(0,0.5)
print('exact')
print(exact)
#################################################################################                                    
# Plot to show the convergence with increasing integration points                                                    
#################################################################################                                    
t=np.linspace(2,100,99)
a=0
error=[]
for i in range(2, 101):
    weight=(np.polynomial.laguerre.laggauss(i))
    func=TestFunction(weight[0]+a)
    integ=weight[1]*func
    integ=np.exp(-a)*np.sum(integ)
    error.append((np.abs((2./np.pi)*integ-exact))\
                     /exact*100.)
plt.plot(t,error)
plt.title('Convergence test for Gauss-Laguerre quadrature: Bessel integral')
plt.xlabel(r'Integration points: $n$')
plt.ylabel(r'Error %')
plt.ylim(0,5999)
plt.show()

Script terminal output:

exact
0.938469807241

Plot generated by script:

Figure 7.2: Percent error of ( 7-2) with increasing integration points

Figure 7-2 does not show any convergence to the correct solution. Observing (Eq 7-12) there are two apparent red flags. The first is that $f(t)$ is unbounded as $t\to \infty$. The second red flag is that $f(t)$ oscillates with an unbounded frequency as $t\to \infty$. Both of these characteristics make a polynomial extrapolation as used in Gaussian quadrature schemes problematic. With a limited number of sampling integration points, it is not expected that Gauss-Laguerre quadrature will yield accurate results for equations with these characteristics.

Contributing Team Members

• Cameron Stewart Solved problems 6 and 7, Reviewed problem 4
• Elizabeth Bartlett Solved problems 3, 4 and 5, Reviewed 2, 6, and 7.
• Kevin Frost Solved problem 3
• Christopher Neal Solved Problems 1 and 2

Problem Assignments
Problem #	Solved by	Reviewed by
1	Neal, Christopher	all
2	Neal, Christopher	all
3	Bartlett, Elizabeth	all
4	Bartlett, Elizabeth; Frost, Kevin	Stewart. Cameron
5	Bartlett, Elizabeth	Stewart, Cameron
6	Stewart, Cameron	Bartlett, Elizabeth
7	Stewart, Cameron	Bartlett, Elizabeth