4. Lagrangians and Hamiltonians

Legal Notice, Disclaimer, and Book Ordering Information:

Null

Setup

\(Clear["Global`*"]; 
Off[General::spell]; 
Off[General::spell1]; 
    \)

Needs["Calculus`VectorAnalysis`"]

Needs["Graphics`Graphics`"]

<< "D:\Mma\Inits\ch4init.m"

Clear["Global`*"];

4.2 Problems

4.2.2  Orbiting Bodies

Problem 3:  Precessing Ellipse and the Generalized Kepler Problem

Part (a)   Define u[ϕ] = 1/r[ϕ] and assume θ=π/2.  Solve for u[ϕ] when the  potential is given by V = -k/r - b/(r^2).  Express the solution in the form (Subscript[u, 0])(1 + ϵ Cos[ϕ(1 - δ)]) where (Subscript[u, 0])= k m/((ℓ^2) - 2 b m)

eq1 = 
      u''[ϕ] == 
        -u[ϕ] - m f[1/u[ϕ]]/(ℓ u[ϕ])^2

\(forceRule = {V -> \((\(-\((k/#)\)\)\  - \ \((b/#^2)\)&)\)
		, 
        f[1/u[ϕ]] -> \(-\(V'\)[r]\)
			, r -> 1/u[ϕ]}\)

eq2 = eq1 //. forceRule // ExpandAll

initial = {u'[0] == 0, u[0] == u0 (1 + ϵ)}

\(uRule = 
      \((\(\((\((DSolve[Join[{eq2}, initial]
	, {u[ϕ]}, ϕ] // 
                    Flatten)\)
		 /. {Sqrt[x_] -> I\ Sqrt[\(-x\)]})\)

            	 // ExpToTrig\) // Simplify)\)\)

shiftRule = {Sqrt[1 - (2 b m)/ℓ^2] -> 1 - shift};

rhs = u[ϕ] //. uRule //. shiftRule

u1Rule = {
        u[ϕ] -> 
          u0 (1  +  ϵ Cos[(1  -  shift) ϕ])};

lhs = u[ϕ] //. u1Rule

u0Rule = Solve[lhs == rhs, u0] // Flatten

Part (b)  Plot the resuts using the command PolarPlot.

\(Clear[plotOrbit]; 

    plotOrbit[start_, finish_, values_List] := 
	
      PolarPlot[
        \(u[ϕ] //. u1Rule\) //. values
			 // Evaluate, {ϕ, 
          start, finish}]\)

p1 = 
      plotOrbit[0, 6  π, {u0 -> 1, ϵ -> 0.5, shift -> 0.1}];

Part (c)   Use the user-defined function firstOrderPert tofins the first-order perturbative solution in the parameter b.  Show that the result agrees wiht the expansion of the solution derived in Part (a).

\(pert = 
      \((firstOrderPert[eq2, u, {u0 \((1 + ϵ)\), 0}, b, ϕ]

            	 //. u0Rule // Simplify)\)\)

pertSer = Series[pert, {b, 0, 1}] // Normal // ExpandAll

\(exactSer = 
      \(Series[\(u[ϕ] //. uRule\) //. u0Rule, {b, 0, 1}]
			 // 
          Normal\) // Expand\)

pertSer == exactSer