Physics 229& 100 Homework #7
Name: Nasser Abbasi
1. Consider the problem of a mass m on a spring of spring constant k subject to a force which is zero until t=0 and is equal to Sin[ωt] for t>0. We want to consider cases where ω=ω=
, and cases where ω is slightly different from
. This is a simple mechanical model for an atom illuminated by near-resonant laser light.
a) To exercise your intuition, first sketch what you think the two solutions will look like; use PowerPoint or some other drawing package ( Paint if you have to) , and include the sketch of your expectations here.
For the case when the forcing frequency marches the natural frequency, resonance will occure. The mass will obtain maximum amplitude. When the forcing frequnecy is different from the natural frequency, then the mass will vibrate in sinuosidel manner (harmonic motion). This is sketch using paint (sorry about the bad quality, hard to write/draw using a mouse).
the y-axis is the displacement y(t). We see that when w=x_n, then the displacement will continue to grow.
b) Then use the greens function derived in class to make a plot of the position as a function of time for the case ω=ω=2, m=1. Contrast this with the case in which the force is Sin[0.9*
t] for t>0, i.e. off resonance illumination.
Due to the divide by zero problem when w=2 exactly (resonance), will plot value very near resonance first.
Since =
then k=m
m , hence for w=2, m=1, we get k=4
Notice also there is what is called 'beats' within the exponentional increasing motion. We see this by plotting for large t
Now plot for 0.9 w
Put both plots together
2. Consider a linear molecule as shown above ( like acetylene). Write the Lagrangian for motions along the line in terms of the masses and spring constants. There are two different kinds of bonds, so there are two different spring constants. To make this very clear, copy the graphic into PowerPoint and add text to denote spring constants kHC and kCC , and masses mC and mH, and then re-import the graphic by copying and pasting.
a) Use EulerEquations to generate the equations of motion.
I can also display the above EQM is a more standard format as follows (even though this is not needed to procceed)
b) Try to use DSolve to solve the resulting system of equations.
c) Depending on which version of Mathematica you are running, you will either paralyze the machine or generate quite a lot of gibberish. Now convert the systems of equations to an eigenvalue problem and determine the eigenfrequencies and normal modes of the molecule. This problem is slightly different from the previous ones we discussed because the masses are not all identical, so be careful about the relationship between the eigenvalues and the eigenfrequencies.
Substitute trial solution x(t)= a exp(i w t) into the system of equations
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
Eliminate common factor
![]() |
![]() |
![]() |
![]() |
Now Solve Det[K-w^2 M]==0 to find eigen frequencies. Then for each eigenFrequency, solve for eigenMode V from equation (K-w^2M)*V=0
The matrix I called massStiffnessMatrix above is the same as [ k-w^2 M ]
But first replace w^2 by λ
For each eigenfrequency, sub back into m2 and find eigen vectors. Use NullSpace for this
To see the eigenmodes better, I use some numerical values for mass and stiffness and then do an animation
This below is extra. It is animation of each eigen mode.
Now animate each mode. The first one is not interesting as nothing will move
3. "Shoot the monkey" is a popular physics demonstration in which a projectile is shot from a gun which is aimed directly at a monkey which is in a nearby tree. The monkey tries to avoid being shot by dropping from the tree the moment the projectile is fired. Both the projectile and the monkey are subject to gravity so the poor monkey gets nailed anyway, much to the delight of the freshmen. Make an animation that shows this. Use the projectile example in Plotting&GraphicsExamples.nb as a guide. Model both the projectile and the monkey as a circle or disk.
The equation of motion for the projectile:
Assume the projectile starts from position coordinates (0,h). Assume the initial speed is in the x-direction and
in the y-direction. Hence the position of the project at time t will be given by x=
t and y=h+
t-
For the falling monkey. Let its initial position be given by (a,h), hence at any time t afterwords, its position is given by x=a and y=h-
Now generate the (x,y) coordinates for the projectile
4. Practice in understanding power series solutions: The differential equation y''[x] +y'[x] Sin[x] +y[x]==0 has analytic coefficients and therefor has two independent solutions which can be expressed as a power series. First note that DSolve can not solve the equation. Then construct a solution "by hand" , by explicitly substituting a power series trial solution, expanding the resulting expression as a power series, equating the coefficient of each power of x to zero, and solving the resulting equations. Find the first 15 terms in the series. Remember that the equation must have two constants of integration, which can be taken to be y[0] and y'[0].
Now use a power series as a trial solution and plug in the ode as well
Now apply Series command again on the above to expand the Sin[x] term
Now generate set of algebraic equations to solve for the a[i] coeff.
a[0]+2 a[2]==0 |
2 a[1]+6 a[3]==0 |
3 a[2]+12 a[4]==0 |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
Since we need 2 free parameters (since the ODE needs to free parameters to solve,and these are found from the 2 initial conditions),we generate 2 equations less than the number of coeff.It is better to drop equations from the end than from the beginning.
a[0]+2 a[2]==0 |
2 a[1]+6 a[3]==0 |
3 a[2]+12 a[4]==0 |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
Now solve for coeff a[2],a[3],...,a[15] in terms of a[0] and a[1]
Now plug these coeff. back into the trial series solution
5. Consider the differential equation x y''[x]+2 y'[x] +4 y[x]= 0
a) find the solution using DSolve
b) Construct power series solutions " by hand" by substituting a trial solution of a generalized power series of the form
a[i]
. Find the indicial equation and find the possible values for the exponent b.
We need to choos value for 'b' to make the coeff. of the lowest powers of x vanish.
c) Construct the power series solution corresponding to the largest value of b. Carefully consider the number of equations and the number of unknowns you generate. If there are atypical equations at the end , use Drop to get rid of them. You will probably want a free parameter to start the recursion relation. If you use a free parameter, set it ==1 to get a unique solution.
From part (b) we see that b=0 is the larger of the two solutions of the indicial equation
Now plug the above series into the original DE
Generate the equations to solve for the a's
4 a[0]+2 a[1]==0 |
4 a[1]+6 a[2]==0 |
4 a[2]+12 a[3]==0 |
4 a[3]+20 a[4]==0 |
4 a[4]+30 a[5]==0 |
4 a[5]+42 a[6]==0 |
4 a[6]+56 a[7]==0 |
4 a[7]+72 a[8]==0 |
4 a[8]+90 a[9]==0 |
4 a[9]+110 a[10]==0 |
4 a[10]+132 a[11]==0 |
4 a[11]+156 a[12]==0 |
4 a[12]==0 |
we have 13 equations but 12 unknowns. If a[12]==0 this will cause all other coeff to become zero due to chain reaction. So drop last equation
4 a[0]+2 a[1]==0 |
4 a[1]+6 a[2]==0 |
4 a[2]+12 a[3]==0 |
4 a[3]+20 a[4]==0 |
4 a[4]+30 a[5]==0 |
4 a[5]+42 a[6]==0 |
4 a[6]+56 a[7]==0 |
4 a[7]+72 a[8]==0 |
4 a[8]+90 a[9]==0 |
4 a[9]+110 a[10]==0 |
4 a[10]+132 a[11]==0 |
4 a[11]+156 a[12]==0 |
make a[0] as free parameter. Set it to 1 to solve the equations
4+2 a[1]==0 |
4 a[1]+6 a[2]==0 |
4 a[2]+12 a[3]==0 |
4 a[3]+20 a[4]==0 |
4 a[4]+30 a[5]==0 |
4 a[5]+42 a[6]==0 |
4 a[6]+56 a[7]==0 |
4 a[7]+72 a[8]==0 |
4 a[8]+90 a[9]==0 |
4 a[9]+110 a[10]==0 |
4 a[10]+132 a[11]==0 |
4 a[11]+156 a[12]==0 |
replace a[0] by 1 since that is what we used as free parameter, and now we obtain sol1 finally
d) Now construct another another power series solution using a trial solution of the form
sol1 Log[x] +
a[i]
, where sol1 is the solution from part b), and "b" is the indicial exponent you did not use in part b). As in part b), drop anomalous equations at the end of the list. Be particularly careful about choosing which of the a[i] should be considered as variables and which should be considered as free parameters. When you are done, set the free parameter to 1.
plug into the DE
Now set the coeff. of each power of x to zero
1+4 a[0]==0 |
-6+4 a[1]+2 a[2]==0 |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
Eliminate the last equation
1+4 a[0]==0 |
-6+4 a[1]+2 a[2]==0 |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
![]() |
Now set the free parameter to 1 as requested by the problem. From above we see that the free parameter is a[12]
and now use these coeff's back in the series itself
Now print solutions. append to the first solution the second solution
e) Check your results by computing the series using the MathematicaHandbook function frobenius.
We see they are the same. Verified
f) Now make sure the series solutions you have found are equivalent to some combination of the exact solutions you found in part a). Do this by finding values of C[1] and C[2] in the exact solution which yield a power series that matches your first and second power series solution. Note: EulerGamma is a constant that is approximately 0.57. Collect is a useful tool to collect the terms proportional to Log[x].
Start by just printing again the solutions I found in part (d)
print the exact solution from part (a)
replace EulerGamma by .57
get the terms that are being multiplied by log[x]
CORRECTION:
I solved finding C[1] and C[2] now by inspection.
I Compare the above series solutiuon generated from the exact solution to my series solutions, and looking at the terms that are being multipled by Log[x], I see that =1 hence
C[2]==
And by comparing the coeff. with no 'x' in them, I see that
C[1]=1+=6.12309
Hence these are the C[1] and C[2] that needs to be used in the exact solution generated by Mathematica DSolve result.
6. Make a table to show that Bessel functions of half integer order reduce to trigonometric ( or exponential) functions. These are related to the spherical bessel functions which arise from separating the wave equation in 3 dimensions.
First make some plots to show Bessel functions.
Make table of Bessel functions for half integers side-by-side. We must use odd values of 'n' to get the half part
n | ![]() |
![]() |
1 | ![]() |
![]() |
3 | ![]() |
![]() |
5 | ![]() |
![]() |
7 | ![]() |
![]() |
9 | ![]() |
![]() |
7. Relations between special functions and hypergeometric functions.
a) Make a table for n=1,6 to show that the Legendre functions (1-2x) can be expressed in terms of 2F1[-n,n+1,1,x]. See Mathematical Functions-> Hypergeometric Related in the Help browser for more information on Mathematica's commands for these functions
n | ![]() |
2F1[-n,n+1,1,x] |
1 | 1-2 x | 1-2 x |
2 | ![]() |
![]() |
3 | ![]() |
![]() |
4 | ![]() |
![]() |
5 | ![]() |
![]() |
6 | ![]() |
![]() |
b) Show that the Laguerre polynomials ( important for radial wave functions, often denoted by ) can be expressed in terms of the confluent hypergeometric function U[-n, m+1, x]. Find the constant of proportionality.
Compare both functions. First select m=1
n | ![]() |
HypergeometricU[-n,m+1,x] |
1 | 2-x | -2+x |
2 | ![]() |
![]() |
3 | ![]() |
![]() |
4 | ![]() |
![]() |
5 | ![]() |
![]() |
6 | ![]() |
![]() |
n | ![]() |
HypergeometricU[-n,m+1,x] |
1 | 3-x | -3+x |
2 | ![]() |
![]() |
3 | ![]() |
![]() |
4 | ![]() |
![]() |
5 | ![]() |
![]() |
6 | ![]() |
![]() |
c) Show that the Hermite polynomials can be expressed as U[-n/2, 1/2, ].
First print few Hermite polynomials to see the form
n | HermiteH[n,x] |
0 | 1 |
1 | 2 x |
2 | ![]() |
3 | ![]() |
4 | ![]() |
5 | ![]() |
6 | ![]() |
7 | ![]() |
8 | ![]() |
9 | ![]() |
10 | ![]() |
Now show HermiteH side-by-side with U function to see the mapping
n | HermiteH[n,x] | ![]() |
0 | 1 | 1 |
1 | 2 x | x |
2 | ![]() |
![]() |
3 | ![]() |
![]() |
4 | ![]() |
![]() |
5 | ![]() |
![]() |
6 | ![]() |
![]() |
7 | ![]() |
![]() |
8 | ![]() |
![]() |
9 | ![]() |
![]() |
10 | ![]() |
![]() |
We see that HermiteH[n,x] is just the above HyperU with a constant factor difference. To find the ratio, divide each HermiteH polynomial by the corresponding Hyper polynomial and look at the pattern. Pull out the polynomials from the above table and divide them
We see from above that HermiteH[n,x] is times HypergeometricU[
,
,
]
Lets verify this by regenrate the above table with this factor now added
n | HermiteH[n,x] | ![]() |
0 | 1 | 1 |
1 | 2 x | 2 x |
2 | ![]() |
![]() |
3 | ![]() |
![]() |
4 | ![]() |
![]() |
5 | ![]() |
![]() |
6 | ![]() |
![]() |
7 | ![]() |
![]() |
8 | ![]() |
![]() |
9 | ![]() |
![]() |
10 | ![]() |
![]() |
We see they are the same.
8. In addition to separation of variables problems in cylindrical coordinates, Bessel functions also show up in a number of Fourier-type integrals. Show that dt=π
[x], first by simply using Integrate, and then by expanding in a series of powers of x and integrating term by term and comparing to the series for
Created by Mathematica (January 30, 2007) | ![]() |