Physics Courseware Support: Mechanics

Hi
The attached worksheet is the final version that appears in Maple 2023 as Courseware support for Mechanics in the context of Physics courses. Everything below also works in Maple 2022.2 with the last Maplesoft Physics Updates for that release..

What follows is presented as "Topic > Problem > Solution", with typical symbolic problems and how you can solve them on a worksheet. As such, this material does not intend to compete with textbooks nor with teacher's notes but to be a helpful complement, as in "what can computer algebra really do to support the learning activity". Mainly, allow for focusing the logic and thinking while the computer takes care of the intricacies of the algebraic manipulations, that when computing with paper and pencil so frequently take mostly all of our focus.

The material, thus, has 70 solved problems covering all the sort-of-syllabus of hyperlinks below. The presentation uses notation as in textbooks and illustrates different techniques, several not present in help pages. It also shows why it is relevant to have a Vectors package that handles abstract vectors as well as projections using unit vectors, not matrix representations for them. Your feedback about everything you see in the worksheet - suggestions for new topics or problems, or anything else - can be useful and is welcome.

Due to the length of this material (~100 pages), out of the 70 problems, below I left open (visible) the Solution sections of only a few of them, illustrating different things, also new functionality e.g. the first and last ones. That is sufficient to have an idea of what this is about. At the end there is a Maple worksheet with the same contents and a PDF file of the same with all the sections open.



With the best wishes for 2023.

 

  

Explore. While learning, having success is a secondary goal: using your curiosity as a compass is what matters. Things can be done in many different ways, take full permission to make mistakes. Computer algebra can transform the algebraic computation part of physics into interesting discoveries and fun.

  

 

  

The following material assumes knowledge of how to use Maple. If you feel that is not your case, for a compact introduction on reproducing in Maple the computations you do with paper and pencil, see sections 1 to 5 of the Mini-Course: Computer Algebra for Physicists . Also, the presentation assumes an understanding of the subjects and the style is not that of a textbook. Instead, it focuses on conveniently using computer algebra to support the practice and learning process. The selection of topics follows references [1] and [2] at the end. Maple 2023.0 includes Part I. Part II is forthcoming.

NULL

Part I

1. 

Position, velocity and acceleration in Cartesian, cylindrical and spherical coordinates

a. 

The position `#mover(mi("r"),mo("→"))`(t)as a function of time

b. 

The velocity `#mover(mi("v"),mo("→"))`(t)

c. 

The acceleration   `#mover(mi("a"),mo("→"))`(t)

d. 

Deriving these formulas

e. 

Velocity and acceleration in the case of 2-dimensional motion on the x, y plane

1. 

The equations of motion

a. 

A single particle

i. 

The equations of motion - vectorial form

ii. 

The case of constant acceleration

iii. 

Motion under gravitational force close to the Earth's surface

iv. 

Motion under gravitational force not close to the Earth's surface

A. 

Circular motion

B. 

Escape velocity

i. 

Different acceleration in different regions

ii. 

The equations of motion using tensor notation

A. 

Cartesian coordinates

B. 

Curvilinear coordinates

a. 

Many-particle systems

i. 

Center of mass

ii. 

The equations of motion

iii. 

Static: reactions of planes and tensions on cables

a. 

Lagrange equations

i. 

Motion of a pendulum

1. 

Conservation laws

a. 

Work

b. 

Conservation of the total energy of a closed system or a system in a constant external field

c. 

Conservation of the total momentum of a closed system

d. 

Conservation of angular momentum

e. 

Cyclic coordinates

1. 

Integration of the equations of motion

a. 

Motion in one dimension

b. 

Reduced mass

i. 

The two-body problem

ii. 

A many-body problem

a. 

Motion in a central field

b. 

Kepler's problem

1. 

Small Oscillations

a. 

Free oscillations in one dimension

b. 

Forced oscillations

c. 

Oscillations of systems with many degrees of freedom

1. 

Rigid-body motion

a. 

Angular velocity

b. 

Inertia tensor

c. 

Angular momentum of a rigid body

d. 

The equations of motion of a rigid body

1. 

Non-inertial coordinate systems

a. 

Coriolis force and centripetal force

 

Part II (forthcoming)

1. 

The Hamiltonian and equations of motion; Poisson brackets

2. 

Canonical transformations

3. 

The Hamilton-Jacobi equation

 

Position, velocity and acceleration in Cartesian, cylindrical and spherical coordinates

 

 

Load the Physics:-Vectors  package

 

with(Physics:-Vectors)

[`&x`, `+`, `.`, Assume, ChangeBasis, ChangeCoordinates, CompactDisplay, Component, Curl, DirectionalDiff, Divergence, Gradient, Identify, Laplacian, Nabla, Norm, ParametrizeCurve, ParametrizeSurface, ParametrizeVolume, Setup, Simplify, `^`, diff, int]

(1)

 

Depending on the geometry of a problem, it can be convenient to work with either Cartesian or curvilinear coordinates. In an arbitrary reference system, the position in Cartesian coordinates and the basis of unitary vectors`#mover(mi("i"),mo("∧"))`, `#mover(mi("j"),mo("∧"))`, `#mover(mi("k"),mo("∧"))`is given by

r_ = _i*x+_j*y+_k*z

r_ = _i*x+_j*y+_k*z

(2)

NULL

Problem

Rewrite the position vector `#mover(mi("r"),mo("→"))` in cylindrical and spherical coordinates

Solution

   

 

Starting from the position in the Cartesian system, now as functions of the time to allow for differentiation, first note that the Cartesian unit vectors `#mover(mi("i"),mo("∧"))`, `#mover(mi("j"),mo("∧"))`, `#mover(mi("k"),mo("∧"))` do not depend on time, they are constant vectors. So `#mover(mi("r"),mo("→"))`(t) is entered as

 

restart; with(Physics:-Vectors)

 

r_(t) = x(t)*_i+y(t)*_j+z(t)*_k

r_(t) = x(t)*_i+y(t)*_j+z(t)*_k

(20)

Before proceeding further, use a compact display to more clearly visualize the following expressions. When in doubt about the contents behind a given display, input show as shown below.

CompactDisplay((x, y, z, rho, r, theta, phi, _rho, _r, _theta, _phi)(t))

_phi(t)*`will now be displayed as`*_phi

(21)

 

For the velocity and acceleration, note the dot notation for derivatives with respect to t

v_(t) = diff(rhs(r_(t) = x(t)*_i+y(t)*_j+z(t)*_k), t)

v_(t) = (diff(x(t), t))*_i+(diff(y(t), t))*_j+(diff(z(t), t))*_k

(22)

 

show

v_(t) = (diff(x(t), t))*_i+(diff(y(t), t))*_j+(diff(z(t), t))*_k

(23)

a_(t) = diff(rhs(v_(t) = (diff(x(t), t))*_i+(diff(y(t), t))*_j+(diff(z(t), t))*_k), t)

a_(t) = (diff(diff(x(t), t), t))*_i+(diff(diff(y(t), t), t))*_j+(diff(diff(z(t), t), t))*_k

(24)

NULL

The position `#mover(mi("r"),mo("→"))`(t)as a function of time

 

 

Problem

Given the position vector as a function of the time t, rewrite it in cylindrical and spherical coordinates while making the curvilinear unit vectors' time dependency explicit.

Solution

   

 

The velocity `#mover(mi("v"),mo("→"))`(t)

 

 

Problem

Rewrite the velocity `#mover(mi("v"),mo("→"))`(t) = diff(`#mover(mi("r"),mo("→"))`(t), t) in cylindrical and spherical coordinates while making the curvilinear unit vectors' time dependency explicit .

Solution

   

 

The acceleration `#mover(mi("a"),mo("→"))`(t)

 

 

Problem

Rewrite the acceleration `#mover(mi("a"),mo("→"))`(t) = diff(`#mover(mi("r"),mo("→"))`(t), t, t)in cylindrical and spherical components while making the curvilinear unit vectors' time dependency explicit.

Solution

   

 

Deriving these formulas

 

 

All these results for the position `#mover(mi("r"),mo("→"))`, velocity `#mover(mi("v"),mo("→"))` and acceleration `#mover(mi("a"),mo("→"))` are based on the differentiation rules for cylindrical and spherical unit vectors. It is thus instructive to also be able to derive any of these formulas; for that, we need the differentiation rule for the unit vectors. For example, for the spherical ones

 

restart; with(Physics:-Vectors); CompactDisplay((x, y, z, rho, r, theta, phi, _rho, _r, _theta, _phi)(t), quiet)

 

map(%diff = diff, ([_r, _theta, _phi])(t), t)

[%diff(_r(t), t) = (diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t), %diff(_theta(t), t) = -(diff(theta(t), t))*_r(t)+(diff(phi(t), t))*cos(theta(t))*_phi(t), %diff(_phi(t), t) = -(diff(phi(t), t))*_rho(t)]

(38)

The above result contains, in the last equation, the cylindrical radial unit vector `#mover(mi("ρ",fontstyle = "normal"),mo("∧"))`(t); rewrite it in the spherical basis

_rho(t) = ChangeBasis(_rho(t), spherical)

_rho(t) = sin(theta(t))*_r(t)+cos(theta(t))*_theta(t)

(39)

So the differentiation rules for spherical unit vectors, with the result expressed in the spherical system, are

subs(_rho(t) = sin(theta(t))*_r(t)+cos(theta(t))*_theta(t), [%diff(_r(t), t) = (diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t), %diff(_theta(t), t) = -(diff(theta(t), t))*_r(t)+(diff(phi(t), t))*cos(theta(t))*_phi(t), %diff(_phi(t), t) = -(diff(phi(t), t))*_rho(t)])

[%diff(_r(t), t) = (diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t), %diff(_theta(t), t) = -(diff(theta(t), t))*_r(t)+(diff(phi(t), t))*cos(theta(t))*_phi(t), %diff(_phi(t), t) = -(diff(phi(t), t))*(sin(theta(t))*_r(t)+cos(theta(t))*_theta(t))]

(40)

NULL

Problem

With this information at hand, derive, in steps, the expressions for the velocity and acceleration in cylindrical and spherical coordinates

Solution

 

 

We want to compute

%diff(r_(t) = r(t)*_r(t), t)

%diff(r_(t) = r(t)*_r(t), t)

(41)

expand(%diff(r_(t) = r(t)*_r(t), t))

%diff(r_(t), t) = %diff(r(t), t)*_r(t)+r(t)*%diff(_r(t), t)

(42)

Introducing the differentiation rules (40) for the unit vectors

subs([%diff(_r(t), t) = (diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t), %diff(_theta(t), t) = -(diff(theta(t), t))*_r(t)+(diff(phi(t), t))*cos(theta(t))*_phi(t), %diff(_phi(t), t) = -(diff(phi(t), t))*(sin(theta(t))*_r(t)+cos(theta(t))*_theta(t))], %diff(r_(t), t) = %diff(r(t), t)*_r(t)+r(t)*%diff(_r(t), t))

%diff(r_(t), t) = %diff(r(t), t)*_r(t)+r(t)*((diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t))

(43)

Performing the inert (grayed) derivatives

value(%diff(r_(t), t) = %diff(r(t), t)*_r(t)+r(t)*((diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t)))

diff(r_(t), t) = (diff(r(t), t))*_r(t)+r(t)*((diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t))

(44)

In the same way, for the acceleration

%diff(r_(t) = r(t)*_r(t), t, t)

%diff(r_(t) = r(t)*_r(t), t, t)

(45)

expand(%diff(r_(t) = r(t)*_r(t), t, t))

%diff(%diff(r_(t), t), t) = %diff(%diff(r(t), t), t)*_r(t)+2*%diff(r(t), t)*%diff(_r(t), t)+r(t)*%diff(%diff(_r(t), t), t)

(46)

subs([%diff(_r(t), t) = (diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t), %diff(_theta(t), t) = -(diff(theta(t), t))*_r(t)+(diff(phi(t), t))*cos(theta(t))*_phi(t), %diff(_phi(t), t) = -(diff(phi(t), t))*(sin(theta(t))*_r(t)+cos(theta(t))*_theta(t))], %diff(%diff(r_(t), t), t) = %diff(%diff(r(t), t), t)*_r(t)+2*%diff(r(t), t)*%diff(_r(t), t)+r(t)*%diff(%diff(_r(t), t), t))

%diff(%diff(r_(t), t), t) = %diff(%diff(r(t), t), t)*_r(t)+2*%diff(r(t), t)*((diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t))+r(t)*%diff((diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t), t)

(47)

expand(%diff(%diff(r_(t), t), t) = %diff(%diff(r(t), t), t)*_r(t)+2*%diff(r(t), t)*((diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t))+r(t)*%diff((diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t), t))

%diff(%diff(r_(t), t), t) = %diff(%diff(r(t), t), t)*_r(t)+2*%diff(r(t), t)*(diff(theta(t), t))*_theta(t)+2*%diff(r(t), t)*(diff(phi(t), t))*sin(theta(t))*_phi(t)+r(t)*%diff(%diff(theta(t), t), t)*_theta(t)+r(t)*(diff(theta(t), t))*%diff(_theta(t), t)+r(t)*%diff(%diff(phi(t), t), t)*sin(theta(t))*_phi(t)+r(t)*(diff(phi(t), t))*%diff(theta(t), t)*cos(theta(t))*_phi(t)+r(t)*(diff(phi(t), t))*sin(theta(t))*%diff(_phi(t), t)

(48)

subs([%diff(_r(t), t) = (diff(theta(t), t))*_theta(t)+(diff(phi(t), t))*sin(theta(t))*_phi(t), %diff(_theta(t), t) = -(diff(theta(t), t))*_r(t)+(diff(phi(t), t))*cos(theta(t))*_phi(t), %diff(_phi(t), t) = -(diff(phi(t), t))*(sin(theta(t))*_r(t)+cos(theta(t))*_theta(t))], %diff(%diff(r_(t), t), t) = %diff(%diff(r(t), t), t)*_r(t)+2*%diff(r(t), t)*(diff(theta(t), t))*_theta(t)+2*%diff(r(t), t)*(diff(phi(t), t))*sin(theta(t))*_phi(t)+r(t)*%diff(%diff(theta(t), t), t)*_theta(t)+r(t)*(diff(theta(t), t))*%diff(_theta(t), t)+r(t)*%diff(%diff(phi(t), t), t)*sin(theta(t))*_phi(t)+r(t)*(diff(phi(t), t))*%diff(theta(t), t)*cos(theta(t))*_phi(t)+r(t)*(diff(phi(t), t))*sin(theta(t))*%diff(_phi(t), t))

%diff(%diff(r_(t), t), t) = %diff(%diff(r(t), t), t)*_r(t)+2*%diff(r(t), t)*(diff(theta(t), t))*_theta(t)+2*%diff(r(t), t)*(diff(phi(t), t))*sin(theta(t))*_phi(t)+r(t)*%diff(%diff(theta(t), t), t)*_theta(t)+r(t)*(diff(theta(t), t))*(-(diff(theta(t), t))*_r(t)+(diff(phi(t), t))*cos(theta(t))*_phi(t))+r(t)*%diff(%diff(phi(t), t), t)*sin(theta(t))*_phi(t)+r(t)*(diff(phi(t), t))*%diff(theta(t), t)*cos(theta(t))*_phi(t)-r(t)*(diff(phi(t), t))^2*sin(theta(t))*(sin(theta(t))*_r(t)+cos(theta(t))*_theta(t))

(49)

value(%diff(%diff(r_(t), t), t) = %diff(%diff(r(t), t), t)*_r(t)+2*%diff(r(t), t)*(diff(theta(t), t))*_theta(t)+2*%diff(r(t), t)*(diff(phi(t), t))*sin(theta(t))*_phi(t)+r(t)*%diff(%diff(theta(t), t), t)*_theta(t)+r(t)*(diff(theta(t), t))*(-(diff(theta(t), t))*_r(t)+(diff(phi(t), t))*cos(theta(t))*_phi(t))+r(t)*%diff(%diff(phi(t), t), t)*sin(theta(t))*_phi(t)+r(t)*(diff(phi(t), t))*%diff(theta(t), t)*cos(theta(t))*_phi(t)-r(t)*(diff(phi(t), t))^2*sin(theta(t))*(sin(theta(t))*_r(t)+cos(theta(t))*_theta(t)))

diff(diff(r_(t), t), t) = (diff(diff(r(t), t), t))*_r(t)+2*(diff(r(t), t))*(diff(theta(t), t))*_theta(t)+2*(diff(r(t), t))*(diff(phi(t), t))*sin(theta(t))*_phi(t)+r(t)*(diff(diff(theta(t), t), t))*_theta(t)+r(t)*(diff(theta(t), t))*(-(diff(theta(t), t))*_r(t)+(diff(phi(t), t))*cos(theta(t))*_phi(t))+r(t)*(diff(diff(phi(t), t), t))*sin(theta(t))*_phi(t)+r(t)*(diff(phi(t), t))*(diff(theta(t), t))*cos(theta(t))*_phi(t)-r(t)*(diff(phi(t), t))^2*sin(theta(t))*(sin(theta(t))*_r(t)+cos(theta(t))*_theta(t))

(50)

Collect vector components

Physics:-Vectors:-Collect(diff(diff(r_(t), t), t) = (diff(diff(r(t), t), t))*_r(t)+2*(diff(r(t), t))*(diff(theta(t), t))*_theta(t)+2*(diff(r(t), t))*(diff(phi(t), t))*sin(theta(t))*_phi(t)+r(t)*(diff(diff(theta(t), t), t))*_theta(t)+r(t)*(diff(theta(t), t))*(-(diff(theta(t), t))*_r(t)+(diff(phi(t), t))*cos(theta(t))*_phi(t))+r(t)*(diff(diff(phi(t), t), t))*sin(theta(t))*_phi(t)+r(t)*(diff(phi(t), t))*(diff(theta(t), t))*cos(theta(t))*_phi(t)-r(t)*(diff(phi(t), t))^2*sin(theta(t))*(sin(theta(t))*_r(t)+cos(theta(t))*_theta(t)))

diff(diff(r_(t), t), t) = (2*r(t)*(diff(theta(t), t))*(diff(phi(t), t))*cos(theta(t))+2*(diff(r(t), t))*(diff(phi(t), t))*sin(theta(t))+r(t)*(diff(diff(phi(t), t), t))*sin(theta(t)))*_phi(t)+(-r(t)*(diff(phi(t), t))^2*sin(theta(t))^2-r(t)*(diff(theta(t), t))^2+diff(diff(r(t), t), t))*_r(t)+(-(diff(phi(t), t))^2*sin(theta(t))*r(t)*cos(theta(t))+2*(diff(r(t), t))*(diff(theta(t), t))+r(t)*(diff(diff(theta(t), t), t)))*_theta(t)

(51)

NULL

 

Summary

• 

You can express `#mover(mi("r"),mo("→"))`(t), `#mover(mi("v"),mo("→"))`(t) and `#mover(mi("a"),mo("→"))`(t)in any of the Cartesian, cylindrical or spherical systems via three different methods: 1) using the ChangeBasis command 2) differentiating 3) deriving the formulas by differentiating in steps, starting from the differentiation rules for the curvilinear unit vectors.

Velocity and acceleration in the case of 2-dimensional motion on the x, y plane

 


Problem

Derive formulas for velocity and acceleration in the case of 2-dimensional motion on the x, y plane, starting from the general 3-dimensional formulas above, e.g. (44) and (51) in spherical coordinates. Specialize the resulting formulas for the case of circular motion.

Solution

   

 

The equations of motion

 

 

A single particle

 

 

restart; with(Physics:-Vectors); CompactDisplay((r_, p_, F_, L_, N_)(t))

N_(t)*`will now be displayed as`*N_

(62)

The equation of motion of a single particle is Newton's 2^nd law

F_(t) = m*(diff(r_(t), t, t))

F_(t) = m*(diff(diff(r_(t), t), t))

(63)

where diff(`#mover(mi("r"),mo("→"))`(t), t, t) = `#mover(mi("a"),mo("→"))`(t) is the acceleration and m*(diff(`#mover(mi("r"),mo("→"))`(t), t)) = `#mover(mi("p"),mo("→"))`(t) is the linear momentum, so in terms of  `#mover(mi("p"),mo("→"))`

F_(t) = diff(p_(t), t)

F_(t) = diff(p_(t), t)

(64)

We define the angular momentum `#mover(mi("L"),mo("→"))` of a particle, and the torque `#mover(mi("N"),mo("→"))` acting upon it, as

L_(t) = `&x`(r_(t), p_(t))

L_(t) = Physics:-Vectors:-`&x`(r_(t), p_(t))

(65)

N_(t) = `&x`(r_(t), F_(t))

N_(t) = Physics:-Vectors:-`&x`(r_(t), F_(t))

(66)

Differentiating the definition of `#mover(mi("L"),mo("→"))`

diff(L_(t) = Physics:-Vectors:-`&x`(r_(t), p_(t)), t)

diff(L_(t), t) = Physics:-Vectors:-`&x`(diff(r_(t), t), p_(t))+Physics:-Vectors:-`&x`(r_(t), diff(p_(t), t))

(67)

Since diff(`#mover(mi("r"),mo("→"))`(t), t) = `#mover(mi("v"),mo("→"))` is parallel to `#mover(mi("p"),mo("→"))` = m*`#mover(mi("v"),mo("→"))`, the first term in the above cancels, and in the second term, from (64), diff(`#mover(mi("p"),mo("→"))`(t), t) = `#mover(mi("F"),mo("→"))`

eval(diff(L_(t), t) = Physics:-Vectors:-`&x`(diff(r_(t), t), p_(t))+Physics:-Vectors:-`&x`(r_(t), diff(p_(t), t)), [diff(r_(t), t) = 0, diff(p_(t), t) = F_(t)])

diff(L_(t), t) = Physics:-Vectors:-`&x`(r_(t), F_(t))

(68)

from which

subs((rhs = lhs)(diff(L_(t), t) = Physics:-Vectors:-`&x`(r_(t), F_(t))), N_(t) = Physics:-Vectors:-`&x`(r_(t), F_(t)))

N_(t) = diff(L_(t), t)

(69)

NULL

• 

As discussed below , in the case of a closed system, `#mover(mi("F"),mo("→"))` = 0 and these two equations result in

diff(`#mover(mi("p"),mo("→"))`(t), t) = 0, diff(`#mover(mi("L"),mo("→"))`(t), t) = 0

that is, the linear and angular momentum are conserved quantities. Note that diff(`#mover(mi("L"),mo("→"))`(t), t) = 0 does not require that `#mover(mi("F"),mo("→"))` = 0, only that `&x`(`#mover(mi("r"),mo("→"))`, `#mover(mi("F"),mo("→"))`) = 0.

 

The equations of motion - vectorial form

 

 

Problem

Assuming that the acceleration is known as a function of t, compute:

a) The trajectory `#mover(mi("r"),mo("→"))`(t)starting from `#mover(mi("a"),mo("→"))`(t) = diff(`#mover(mi("r"),mo("→"))`(t), t, t)
b) A solution for each of the three Cartesian components

c) A solution for generic initial conditions

Solution

 

 

restart; with(Physics:-Vectors); CompactDisplay((x, y, z, rho, r, theta, phi, _rho, _r, _theta, _phi)(t), quiet)
 

 

a) Let `#mover(mi("r"),mo("→"))`(t) be the position of the particle in a reference system; then, the velocity and acceleration are given by

v_(t) = diff(r_(t), t)

v_(t) = diff(r_(t), t)

(70)

a_(t) = diff(r_(t), t, t)

a_(t) = diff(diff(r_(t), t), t)

(71)

If the acceleration is known as a function of t, the trajectory is computed by integrating (71)

dsolve(a_(t) = diff(diff(r_(t), t), t))

r_(t) = Int(Int(a_(t), t), t)+c__1_*t+c__2_

(72)

where the vectorial integration constants, c__1_ and c__2_, are specified by the initial conditions of the problem (see c) below), typically by the position and velocity at some instant, say %eval(r_(t), `=`(t, t__0)) = r__0_ and %eval(diff(r_(t), t), `=`(t, t__0)) = v__0_.

_______________________________________

 

b) The integration of vectorial equations is also frequently performed after expressing `#mover(mi("r"),mo("→"))`(t), `#mover(mi("v"),mo("→"))`(t) and `#mover(mi("a"),mo("→"))`(t) in a particular system of coordinates. For example, in the Cartesian system (71) has the form

(rhs = lhs)(a_(t) = (diff(diff(x(t), t), t))*_i+(diff(diff(y(t), t), t))*_j+(diff(diff(z(t), t), t))*_k)

(diff(diff(x(t), t), t))*_i+(diff(diff(y(t), t), t))*_j+(diff(diff(z(t), t), t))*_k = a_(t)

(73)

Now suppose that the three components of the acceleration are known as a function of time

subs(a_(t) = a__x(t)*_i+a__y(t)*_j+a__z(t)*_k, (diff(diff(x(t), t), t))*_i+(diff(diff(y(t), t), t))*_j+(diff(diff(z(t), t), t))*_k = a_(t))

(diff(diff(x(t), t), t))*_i+(diff(diff(y(t), t), t))*_j+(diff(diff(z(t), t), t))*_k = a__x(t)*_i+a__y(t)*_j+a__z(t)*_k

(74)

Vectorial equations like this one can be integrated directly, provided that they are expressed in a particular system of coordinates and the unit vectors are constant or known expressions of the time

dsolve((diff(diff(x(t), t), t))*_i+(diff(diff(y(t), t), t))*_j+(diff(diff(z(t), t), t))*_k = a__x(t)*_i+a__y(t)*_j+a__z(t)*_k)

_i*x(t)+_j*y(t)+_k*z(t) = _i*(Int(Int(a__x(t), t), t)+c__3*t+c__4)+_j*(Int(Int(a__y(t), t), t)+c__5*t+c__6)+_k*(Int(Int(a__z(t), t), t)+c__7*t+c__8)

(75)

_______________________________________

 

c) The vectorial initial conditions r__0_ and v__0_, specifying the integration constants {c__3, c__4, c__5, c__6, c__7, c__8}, can also be written in components

x(t__0) = x__0, y(t__0) = y__0, z(t__0) = z__0(t), diff(x(t__0), t__0) = v__x0, diff(y(t__0), t__0) = v__y0, diff(z(t__0), t__0) = v__z0

x(t__0) = x__0, y(t__0) = y__0, z(t__0) = z__0(t), diff(x(t__0), t__0) = v__x0, diff(y(t__0), t__0) = v__y0, diff(z(t__0), t__0) = v__z0

(76)

Passing this information, the system can be solved taking these initial conditions into account -

dsolve([(diff(diff(x(t), t), t))*_i+(diff(diff(y(t), t), t))*_j+(diff(diff(z(t), t), t))*_k = a__x(t)*_i+a__y(t)*_j+a__z(t)*_k, x(t__0) = x__0, y(t__0) = y__0, z(t__0) = z__0(t), diff(x(t__0), t__0) = v__x0, diff(y(t__0), t__0) = v__y0, diff(z(t__0), t__0) = v__z0], ({x, y, z})(t))

{x(t) = Int(Int(a__x(tau), tau = t__0 .. tau), tau = t__0 .. t)+v__x0*t-t__0*v__x0+x__0, y(t) = Int(Int(a__y(tau), tau = t__0 .. tau), tau = t__0 .. t)+v__y0*t-t__0*v__y0+y__0, z(t) = Int(Int(a__z(tau), tau = t__0 .. tau), tau = t__0 .. t)+v__z0*t+c__4}

(77)

Note that a vectorial equation is also always equivalent to a system of equations, one for each of the components, with or without initial conditions:

convert((diff(diff(x(t), t), t))*_i+(diff(diff(y(t), t), t))*_j+(diff(diff(z(t), t), t))*_k = a__x(t)*_i+a__y(t)*_j+a__z(t)*_k, setofequations)

{diff(diff(x(t), t), t) = a__x(t), diff(diff(y(t), t), t) = a__y(t), diff(diff(z(t), t), t) = a__z(t)}

(78)

dsolve({diff(diff(x(t), t), t) = a__x(t), diff(diff(y(t), t), t) = a__y(t), diff(diff(z(t), t), t) = a__z(t)}, {x, y, z})

{x(t) = Int(Int(a__x(t), t), t)+c__7*t+c__8, y(t) = Int(Int(a__y(t), t), t)+c__5*t+c__6, z(t) = Int(Int(a__z(t), t), t)+c__3*t+c__4}

(79)

 

NULL

 

The case of constant acceleration

 

 

Problem

Starting from the vectorial equation (72) for `#mover(mi("r"),mo("→"))`(t), derive the formula for constant acceleration

Solution

   

 

Motion under gravitational force close to the Earth's surface

 

 

Problem

Derive a formula for motion under gravitational force close to the Earth's surface

Solution

   

 

Motion under gravitational force not close to the Earth's surface

 

 

The problem of two particles of masses m__1 and m__2 gravitationally attracted to each other, discarding relativistic effects, is formulated by Newton's law of gravity: the particles attract each other - so both move - with a force along the line that joins the particles and whose magnitude is proportional to 1/r^2, where r represents the distance between the particles (this problem is treated in general form  in the more advanced sections).

 

Problem

As a specific case, consider the problem of a particle of mass `≪`(m, M), where M is earth's mass, moving not close to the surface (if compared with the radius of earth).

Solution

   

 

Circular motion

 

 

Problem

Determine the angular velocity diff(phi(t), t) in the case of circular motion and show it is constant.

Solution

   

 

Escape velocity

 

 

Problem

Determine the velocity that a particle of mass m should have at Earth's surface in order to escape the planet's gravitational attraction.

Solution

   

 

Different acceleration in different regions

 

 

Problem

Suppose a particle is moving along the x axis according to

x(t) = t^3-8*t^2+18*t+3

a) Determine the regions where the motion has positive and negative acceleration. Compute the position at proc (t) options operator, arrow; infinity end proc.

b) Compute the velocity v__x(t)corresponding to x(t) = t^3-8*t^2+18*t+3, starting - not from this expression for x(t) but from the acceleration a__x(t) = diff(x(t), t, t)

Solution

   

 

The equations of motion using tensor notation

 

 

Using vector notation to formulate the equations of motion of a particle in Cartesian coordinates is relatively simple. However, for certain problems it may be advantageous to use curvilinear coordinates and / or tensor notation.

 

restart; with(Physics); with(Vectors)

NULL

Cartesian coordinates

   

Curvilinear coordinates

   

 

Many-particle systems

 

 

Center of Mass

 

Given a system of n particles of masses m__i with positions `#mover(mi("r"),mo("→"))`[i] in some frame of reference K, the center of mass of the system is defined as

 `#mover(mi("R"),mo("→"))` = (Sum(m[i]*`#mover(mi("r"),mo("→"))`[i], i = 1 .. n))/(Sum(m[i], i = 1 .. n))

The velocity of the center of mass is thus

  `#mover(mi("V"),mo("→"))` = diff(`#mover(mi("R"),mo("→"))`(t), t) and diff(`#mover(mi("R"),mo("→"))`(t), t) = (Sum(m[i]*(diff(`#mover(mi("r"),mo("→"))`[i](t), t)), i = 1 .. n))/(Sum(m[i], i = 1 .. n)) 

Problem

Consider a system of particles viewed from two systems of reference, K and K', that move with respect to each other at a constant velocity `#mover(mi("V"),mo("→"))` measured in K. Show that:

a) When `#mover(mi("V"),mo("→"))` is the velocity of the center of mass, the total momentum "(P )'" measured in K'  is equal to 0.

b) The relation between `#mover(mi("P"),mo("→"))` and the velocity `#mover(mi("V"),mo("→"))` of the center of mass, both measured in K, is the same as the relation `#mover(mi("p"),mo("→"))` = m*`#mover(mi("v"),mo("→"))` between the momentum, velocity and mass of a single particle of mass mu = Sum(m[i], i = 1 .. n).

Solution

   

 

The equations of motion

 

 

Problem
Show that, for a system of particles with total mass mu = Sum(m[i], i = 1 .. n), Newton's 2nd law for each particle `#mover(mi("F"),mo("→"))`[i] = m[i]*(diff(`#mover(mi("r"),mo("→"))`[i](t), t, t)) implies that `#mscripts(mi("F"),mi("ext"),none(),none(),mo("→"),none(),none())` = mu*(diff(`#mover(mi("R"),mo("→"))`(t), t, t)), where `#mover(mi("R"),mo("→"))` is the center of mass and `#mscripts(mi("F"),mi("ext"),none(),none(),mo("→"),none(),none())` is the external force applied to the system (it excludes the force that the particles exercise on each other).

Solution

   

 

Problem

Show that :

a) The total linear momentum `#mover(mi("P"),mo("→"))` satisfies diff(`#mover(mi("P"),mo("→"))`(t), t) = `#mscripts(mi("F"),mi("ext"),none(),none(),mo("→"),none(),none())`

b) The total torque `#mover(mi("N"),mo("→"))` = diff(`#mover(mi("L"),mo("→"))`(t), t) satisfies `#mover(mi("N"),mo("→"))` = Sum(`&x`(`#mover(mi("r"),mo("→"))`[i], `#mover(mi("f"),mo("→"))`[i, ext]), i = 1 .. n)

Solution

   

 

Static: reactions of planes and tensions on cables

 

 

Problem

A bar AB of weight w and length L has one extreme on a horizontal plane and the other on a vertical place, and is kept in that position by two cables AD and BC. The bar forms an angle alpha with the horizontal plane and its projection BC over this plane forms an angle beta with the vertical plane. The cable BC is on the same vertical plane as the bar.

 

Determine the reactions of the planes at A and B as well as the tensions on the cables.

Solution

   

 

Lagrange equations

 

 

restart; with(Physics:-Vectors); CompactDisplay((r_, v_)(t))

v_(t)*`will now be displayed as`*v_

(208)

 

In the case of a closed system, or a system in a constant external field, the equations of motions can also be derived from the knowledge of the kinetic energy T and the potential energy U . For this purpose, construct the Lagrange function L = T-U and derive the equations of motion as the Lagrange equations for L.

 

For closed systems, the potential energy U(`#mover(mi("r"),mo("→"))`[i])is related to the force acting on each particle by the equation `#mover(mi("F"),mo("→"))`[i] = -`∇__i`(U(`$`(r_[i], `=`(i, 1 .. n)))). Formally, "`∇__i`≡diff(L,r_[i]" is the gradient operator in the basis onto which `#mover(mi("r"),mo("→"))`[i] is projected, and with respect to its coordinates - say in Cartesian basis and coordinates "(∇)[i]=(i)*(`%diff`(,)+j*(`%diff`(,y))+k*(`%diff`(,z))".


The kinetic energy - say T - of a single particle is given by

T := (1/2)*m*v_(t)^2

(1/2)*m*Physics:-`^`(v_(t), 2)

(209)

Since the kinetic energy T is additive, a system of n particles has

T := %sum((1/2)*m*v_[i](t)^2, i = 1 .. n)

%sum((1/2)*m*Physics:-`^`(v_[i](t), 2), i = 1 .. n)

(210)

where `#mover(mi("v"),mo("→"))`[i] is the velocity of the ith particle. Generally speaking, the potential energy U(`$`(r_[i], `=`(i, 1 .. n))) of the system is a function of the coordinates `#mover(mi("r"),mo("→"))`[i] of the n particles, and the Lagrangian is defined asNULL

L = T-U(`$`(r_[i], `=`(i, 1 .. n)))

L = %sum((1/2)*m*Physics:-`^`(v_[i](t), 2), i = 1 .. n)-U(`$`(r_[i], i = 1 .. n))

(211)

The potential energy U(`#mover(mi("r"),mo("→"))`[i])is related to the force acting on each particle by the equation `#mover(mi("F"),mo("→"))`[i] = -`∇__i`(U(`$`(r_[i], `=`(i, 1 .. n)))). Formally, "`∇__i`≡diff(L,r_[i]" is the gradient operator in the basis onto which `#mover(mi("r"),mo("→"))`[i] is projected. Knowing the Lagrangian, we can derive the (Lagrange) equations of motion  as

%diff(%diff(L, v_[i]), t) = %diff(L, r_[i])

%diff(%diff(L, v_[i]), t) = %diff(L, r_[i])

(212)

NULL

Motion of a pendulum

 

 

Problem

Determine the Lagrangian and equation of motion of a plane pendulum with a mass m at its extremity and a suspension point which:

a) Moves uniformly over a vertical circumference with a constant frequency "omega."

 

b) Oscillates horizontally on the plane of the pendulum according to x = cos(omega*t).

c) Is fixed (omega = 0). Integrate the equation of motion for small oscillations.

Solution

   

 

Conservation laws

 

Work

 

By definition, the work performed by a force `#mover(mi("F"),mo("→"))` to move a particle an infinitesimal amount d*`#mover(mi("r"),mo("→"))` is Typesetting[delayDotProduct](`#mover(mi("F"),mo("→"))`.d, `#mover(mi("r"),mo("→"))`, true). Thus, the work to move it from `#mover(mi("A"),mo("→"))` to `#mover(mi("B"),mo("→"))` along some path C is

 

"((∫)[A]^(B) (F ) . ⅆr)[path=C]"

Problem

A particle is submitted to a force whose Cartesian components are given by" F[x]=a x^3+b x y^2+c z,  ""F[y]=a y^3+b x y^2,  "F[z] = c*z. Calculate the work done by this force when moving the particle along a straight line from the origin to a point (x[0], y[0], z[0]).

Solution

   

 

Conservation of the total energy of a closed system or system in a constant external field

 

 

Problem

Consider a closed system, or a system in a constant external field, for which the total force acting on the ith particle of the system can be derived from a potential, `#mscripts(mi("F"),mi("i"),none(),none(),mo("→"),none(),none())` = -`∇i`*U. Show that the total energy of the system is conserved.

Solution

   

 

Problem

 Consider a system of n particles in two reference systems K and K'  that move relative to each other with constant velocity `#mover(mi("V"),mo("→"))`. Show that the relation between the energies of the system, E and E', is given by

E(x) = diff(E(x), x)+(1/2)*LinearAlgebra[Norm](`#mover(mi("V"),mo("→"))`)^2*(Sum(m[a], a = 1 .. n))+`#mover(mi("V"),mo("→"))`.`#mover(mi("P'"),mo("→"))`

Solution

   

 

Conservation of the total momentum of a closed system

 


The conservation of the total momentum of a closed system of one particle is clear: if the particle does not interact with anything external, the force acting on it is zero, and from Newton's 2nd law `#mover(mi("F"),mo("→"))` = diff(`#mover(mi("p"),mo("→"))`(t), t) follows diff(`#mover(mi("p"),mo("→"))`(t), t) = 0.

 

For a closed system of many particles, while the total force acting on the system is equal to 0, there can be internal forces different from zero acting on each particle due to its interaction with the other particles. These internal forces, however, produce no acceleration of the system; in general, they cancel each other out due to Newton's 3rd law.

 

Problem

Consider a system of n particles measured in two frames of reference K and K'  that move relative to each other with velocity `#mover(mi("V"),mo("→"))`. Show that the system's momenta `#mover(mi("P"),mo("→"))` and "(P)^( ')" are related by

"P=(P)^( ')+V *(∑)m[a]".

Solution

   

 

Problem

A particle of mass m moving with velocity `#mover(mi("v"),mo("→"))`[1] leaves a half-space in which its potential energy is a constant U[1] and enters another in which its potential energy is a different constant U[2]. Determine the change in direction of motion of the particle; that is, sin(`θ__1`)/sin(`θ__2`) where `θ__1` and `θ__2` are the angles between an axis perpendicular to the separating plane and the momentum `#mover(mi("p"),mo("→"))` in the regions 1 and 2.

Solution

   

 

Conservation of angular momentum

 


Like the conservation of linear momentum, the conservation of the total angular momentum of a closed system of one particle is natural: if the particle does not interact with anything external, the force acting on it is zero and therefore its torque `#mover(mi("N"),mo("→"))` = `&x`(`#mover(mi("r"),mo("→"))`, `#mover(mi("F"),mo("→"))`) and `&x`(`#mover(mi("r"),mo("→"))`, `#mover(mi("F"),mo("→"))`) = 0. From diff(`#mover(mi("L"),mo("→"))`(t), t) = `#mover(mi("N"),mo("→"))` it follows that diff(`#mover(mi("L"),mo("→"))`(t), t) = 0, that is, the angular momentum `#mover(mi("L"),mo("→"))` is conserved.

 

Problem

a) Express the Cartesian components of the angular momentum `#mover(mi("L"),mo("→"))`, as well as its norm, in cylindrical and spherical coordinates.

b) Rewrite `#mover(mi("L"),mo("→"))` in cylindrical coordinates and using the cylindrical orthonormal basis of unit vectors, then do the same using spherical coordinates and the spherical basis.

Solution

   

 

Problem

Consider a system of n particles measured in two frames of reference K and K' whose origins have distance `#mover(mi("A"),mo("→"))` from each other. Show that the momenta `#mover(mi("L"),mo("→"))` and "((L))^( ')" of the system are related by

 "(L)=((L))^( ')+(A) *P"

where `#mover(mi("P"),mo("→"))` = Sum(`#mover(mi("p"),mo("→"))`[a], a = 1 .. n) is the total momentum of the system as seen from K.

Solution

   

 

Problem

a) Consider a closed system of n particles, and two frames of reference K and K'  that move relative to each other with a constant velocity `#mover(mi("V"),mo("→"))`. Show that the momenta `#mover(mi("L"),mo("→"))` and "((L))^( ')"respectively measured in K and K' are related by

`#mover(mi("L"),mo("→"))` = (Sum(m[a], a = 1 .. n))*`&x`(`#mover(mi("R"),mo("→"))`, `#mover(mi("V"),mo("→"))`)+`&x`(`#mover(mi("A"),mo("→"))`, `#mover(mi("P'"),mo("→"))`)+`#mover(mi("L'"),mo("→"))`

where `#mover(mi("A"),mo("→"))` is the distance from the origin of K to the origin of K',  "R=(∑)m[a] (r)[a]/(∑)m[a]" is the position of the center of mass as seen from K, and `#mover(mi("P'"),mo("→"))` and `#mover(mi("L'"),mo("→"))` are the total linear and angular momenta measured in K'.
b) Show that when the origin of K' is the center of mass `#mover(mi("R"),mo("→"))`, this formula reduces to

`#mover(mi("L"),mo("→"))` = `#mover(mi("L'"),mo("→"))`+`&x`(`#mover(mi("R"),mo("→"))`, `#mover(mi("P"),mo("→"))`),

where `#mover(mi("P"),mo("→"))` = Sum(m[a]*`#mover(mi("v"),mo("→"))`[a], a = 1 .. n) is the total linear momentum in K.

Solution

 

 

restart; with(Physics:-Vectors)

 

a) The momentum `#mover(mi("L"),mo("→"))` of the system measured in K is given by

L_ = Sum(`&x`(r_[a], p_[a]), a = 1 .. n)

L_ = Sum(Physics:-Vectors:-`&x`(r_[a], p_[a]), a = 1 .. n)

(322)

The right-hand side of the above can be expressed in terms of `#mover(mi("P'"),mo("→"))` and `#mover(mi("L'"),mo("→"))` using  `#mover(mi("p"),mo("→"))`[a] = m[a]*`#mover(mi("v"),mo("→"))`[a], where `#mover(mi("v"),mo("→"))`[a] = `#mover(mi("v'"),mo("→"))`[a]+`#mover(mi("V"),mo("→"))` :

p_[a] = (diff(p(x), x))*_[a]+m[a]*V_

p_[a] = V_*m[a]+`p'_`[a]

(323)

subs(p_[a] = V_*m[a]+`p'_`[a], L_ = Sum(Physics[Vectors][`&x`](r_[a], p_[a]), a = 1 .. n))

L_ = Sum(Physics:-Vectors:-`&x`(r_[a], V_*m[a]+`p'_`[a]), a = 1 .. n)

(324)

ee := expand(L_ = Sum(Physics[Vectors][`&x`](r_[a], V_*m[a]+`p'_`[a]), a = 1 .. n))

L_ = -Physics:-Vectors:-`&x`(V_, Sum(m[a]*r_[a], a = 1 .. n))+Sum(Physics:-Vectors:-`&x`(r_[a], `p'_`[a]), a = 1 .. n)

(325)

The term with Sum(m[a]*`#mover(mi("r"),mo("→"))`[a], a = 1 .. n) can be expressed in terms of the position vector of the center of mass `#mover(mi("R"),mo("→"))` (copy the subexpression from (325) , paste, then edit)

subs(Sum(m[a]*`#mover(mi("r"),mo("→"))`[a], a = 1 .. n) = (Sum(m[a], a = 1 .. n))*R_, L_ = -Physics[Vectors][`&x`](V_, Sum(m[a]*r_[a], a = 1 .. n))+Sum(Physics[Vectors][`&x`](r_[a], `p'_`[a]), a = 1 .. n))

L_ = -Physics:-Vectors:-`&x`(V_, (Sum(m[a], a = 1 .. n))*R_)+Sum(Physics:-Vectors:-`&x`(r_[a], `p'_`[a]), a = 1 .. n)

(326)

To express Sum(`&x`(`#mover(mi("r"),mo("→"))`[a], `#mover(mi("p'"),mo("→"))`[a]), a = 1 .. n) in terms of `#mover(mi("L'"),mo("→"))` and `#mover(mi("A"),mo("→"))`, introduce the relation between `#mover(mi("r"),mo("→"))`[a] and `#mover(mi("r'"),mo("→"))`[a]

r_[a] = A_+(diff(r(x), x))*_[a]

r_[a] = A_+`r'_`[a]

(327)

subs(r_[a] = A_+`r'_`[a], L_ = -Physics[Vectors][`&x`](V_, (Sum(m[a], a = 1 .. n))*R_)+Sum(Physics[Vectors][`&x`](r_[a], `p'_`[a]), a = 1 .. n))

L_ = -(Sum(m[a], a = 1 .. n))*Physics:-Vectors:-`&x`(V_, R_)+Sum(Physics:-Vectors:-`&x`(A_+`r'_`[a], `p'_`[a]), a = 1 .. n)

(328)

expand(L_ = -(Sum(m[a], a = 1 .. n))*Physics[Vectors][`&x`](V_, R_)+Sum(Physics[Vectors][`&x`](A_+`r'_`[a], `p'_`[a]), a = 1 .. n))

L_ = -(Sum(m[a], a = 1 .. n))*Physics:-Vectors:-`&x`(V_, R_)+Physics:-Vectors:-`&x`(A_, Sum(`p'_`[a], a = 1 .. n))+Sum(Physics:-Vectors:-`&x`(`r'_`[a], `p'_`[a]), a = 1 .. n)

(329)

On the right-hand side, two of the sums represent `#mover(mi("P'"),mo("→"))` and `#mover(mi("L'"),mo("→"))` (copy the sum subexpressions, paste into the next line, then edit)

subs(Sum(`&x`(`#mover(mi("r'"),mo("→"))`[a], `#mover(mi("p'"),mo("→"))`[a]), a = 1 .. n) = (diff(L(x), x))*_, Sum(`#mover(mi("p'"),mo("→"))`[a], a = 1 .. n) = (diff(P(x), x))*_, L_ = -(Sum(m[a], a = 1 .. n))*Physics[Vectors][`&x`](V_, R_)+Physics[Vectors][`&x`](A_, Sum(`p'_`[a], a = 1 .. n))+Sum(Physics[Vectors][`&x`](`r'_`[a], `p'_`[a]), a = 1 .. n))

L_ = -(Sum(m[a], a = 1 .. n))*Physics:-Vectors:-`&x`(V_, R_)+Physics:-Vectors:-`&x`(A_, `P'_`)+`L'_`

(330)

This is already the desired result.

_______________________________________

 

b) If the origin of K' is the center of mass `#mover(mi("R"),mo("→"))`, then `#mover(mi("A"),mo("→"))` = `#mover(mi("R"),mo("→"))`

subs(A_ = R_, L_ = -(Sum(m[a], a = 1 .. n))*Physics[Vectors][`&x`](V_, R_)+Physics[Vectors][`&x`](A_, `P'_`)+`L'_`)

L_ = -(Sum(m[a], a = 1 .. n))*Physics:-Vectors:-`&x`(V_, R_)+Physics:-Vectors:-`&x`(R_, `P'_`)+`L'_`

(331)

`#mover(mi("P"),mo("→"))` and `#mover(mi("P'"),mo("→"))` were related in "P_ = V_*Sum(m[a],a = 1 .. n)+`P'_`≡ "`#mover(mi("P"),mo("→"))` = `#mover(mi("V"),mo("→"))`*(Sum(m[a], a = 1 .. n))+`#mover(mi("P'"),mo("→"))`.  This relation can be used to express `#mover(mi("L"),mo("→"))` in terms of  `#mover(mi("P"),mo("→"))` instead of `#mover(mi("P'"),mo("→"))`

simplify(L_ = -(Sum(m[a], a = 1 .. n))*Physics[Vectors][`&x`](V_, R_)+Physics[Vectors][`&x`](R_, `P'_`)+`L'_`, {P_ = V_*(Sum(m[a], a = 1 .. n))+`P'_`}, {(diff(P(x), x))*_})

L_ = -(Sum(m[a], a = 1 .. n))*Physics:-Vectors:-`&x`(V_, R_)+Physics:-Vectors:-`&x`(R_, -V_*(Sum(m[a], a = 1 .. n))+P_)+`L'_`

(332)

expand(L_ = -(Sum(m[a], a = 1 .. n))*Physics[Vectors][`&x`](V_, R_)+Physics[Vectors][`&x`](R_, -V_*(Sum(m[a], a = 1 .. n))+P_)+`L'_`)

L_ = Physics:-Vectors:-`&x`(R_, P_)+`L'_`

(333)

NULL

which is the result we were looking for.

 

 

Cyclic coordinates

 

 

Any generalized coordinate q__i which does not appear explicitly in the Lagrangian is called cyclic. To any cyclic coordinate corresponds a conserved quantity. From

%diff(%diff(L, diff(q[i](t), t)), t) = %diff(L, q[i])

when q__i is cyclic, the right-hand side is 0 and so the quantity %diff(L, diff(q[i](t), t)) is conserved.

 

Problem

The Lagrangian describing the movement of a particle in a central field has phi as a cyclic coordinate. Using cylindrical coordinates, show that the corresponding conserved quantity is the z component of the angular momentum `#mover(mi("L"),mo("→"))`

Solution

   

 

Integration of the equations of motion

 

 

Motion in one dimension

 

Problem
For a closed system, or any system where the total energy E = T+U is conserved, show the following:

a) The trajectory in implicit form can always be computed directly from E.

b) The turning points, if any, can be computed directly from U.

Solution

   

 

Problem

Determine the period of oscillations of a pendulum of mass m and length l in a gravitational field as a function of the amplitude of the oscillations

Solution

   

 

Problem

Integrate the equations of motion for a particle of mass m moving in a field whose potential energy is U = A*abs(x)^n.

Solution

   

Reduced mass

 

The two-body problem

 

 

Problem
Show that by placing the origin of the reference system at the center of mass  `#mover(mi("R"),mo("→"))` = (Sum(m[i]*`#mover(mi("r"),mo("→"))`[i], i = 1 .. n))/(Sum(m[i], i = 1 .. n)), the problem of two particles that interact with each other can be reduced to the problem of a single particle of mass mu = m__1*m__2/(m__1+m__2), herein called reduced mass, in an external field U(LinearAlgebra[Norm](r__2_-r__1_)).

Solution

   

 

A many-body problem

 


Problem

A system consists of one particle of mass M and n particles of equal masses m.

a) Show, in steps, that eliminating the motion of the center of mass reduces the problem to one involving only n particles.

b) Show that when n = 1the result of a) becomes the result obtained for the previous two-body problem, equation `≡`(L = (1/2)*Physics[Vectors][Norm](diff(`ℝ_`(t), t))^2*mu-U(Physics[Vectors][Norm](`ℝ_`(t))), L) = (1/2)*mu*LinearAlgebra[Norm](diff(`#mover(mi("ℝ"),mo("→"))`(t), t))^2-U.

Solution

 

 

restart; with(Physics:-Vectors)

 

a) Let R__M_ represent the position vector of the particle of mass M and `#mover(mi("r"),mo("→"))`[a] those of the particles of mass m. The Lagrangian L = T-U is given by

CompactDisplay((r__M_, r_, `ℝ_`)(t))

`ℝ_`(t)*`will now be displayed as`*`ℝ_`

(386)

Physics:-Assume(m > 0, M > 0)

{m::(RealRange(Open(0), infinity))}, {M::(RealRange(Open(0), infinity))}

(387)

L = (1/2)*M*(diff(r__M_(t), t))^2+(1/2)*(Sum(m*(diff(r_[a](t), t))^2, a = 1 .. n))-U

L = (1/2)*M*Physics:-`^`(diff(r__M_(t), t), 2)+(1/2)*(Sum(m*Physics:-`^`(diff(r_[a](t), t), 2), a = 1 .. n))-U

(388)

Expand the formal powers of vectors to express them in terms of their Norm

expand(L = (1/2)*M*Physics[`^`](diff(r__M_(t), t), 2)+(1/2)*(Sum(m*Physics[`^`](diff(r_[a](t), t), 2), a = 1 .. n))-U)

L = (1/2)*M*Physics:-Vectors:-Norm(diff(r__M_(t), t))^2+(1/2)*m*(Sum(Physics:-Vectors:-Norm(diff(r_[a](t), t))^2, a = 1 .. n))-U

(389)

As done in the two-body problem, introduce the relative vector positions of the n particles with respect to the particle of mass M

`ℝ_`[a](t) = r_[a](t)-r__M_(t)

`ℝ_`[a](t) = r_[a](t)-r__M_(t)

(390)

and place the origin of the reference system at the center of mass,

(M*r__M_(t)+Sum(m*r_[a](t), a = 1 .. n))/(M+Sum(m, a = 1 .. n)) = 0

(M*r__M_(t)+Sum(m*r_[a](t), a = 1 .. n))/(M+Sum(m, a = 1 .. n)) = 0

(391)

expand((M*r__M_(t)+Sum(m*r_[a](t), a = 1 .. n))/(M+Sum(m, a = 1 .. n)) = 0)

M*r__M_(t)/(m*n+M)+m*(Sum(r_[a](t), a = 1 .. n))/(m*n+M) = 0

(392)

Using this equation and (390)`#mover(mi("r"),mo("→"))`[a] = `#mover(mi("ℝ"),mo("→"))`[a]-R__M_, fontstyle = "normal", Typesetting:-msemantics = "function" is sufficient to eliminate R__M_ and `#mover(mi("r"),mo("→"))`[a] from the Lagrangian (389). In steps, eliminating `#mover(mi("r"),mo("→"))`[a] from the equation for the center of mass

simplify(M*r__M_(t)/(m*n+M)+m*(Sum(r_[a](t), a = 1 .. n))/(m*n+M) = 0, {`ℝ_`[a](t) = r_[a](t)-r__M_(t)}, {r_[a](t)})

(M*r__M_(t)+m*(Sum(`ℝ_`[a](t)+r__M_(t), a = 1 .. n)))/(m*n+M) = 0

(393)

expand((M*r__M_(t)+m*(Sum(`ℝ_`[a](t)+r__M_(t), a = 1 .. n)))/(m*n+M) = 0)

M*r__M_(t)/(m*n+M)+m*(Sum(`ℝ_`[a](t), a = 1 .. n))/(m*n+M)+m*r__M_(t)*n/(m*n+M) = 0

(394)

To eliminate `#msub(mi("R"),mi("M_"))` , use

simplify(isolate(M*r__M_(t)/(m*n+M)+m*(Sum(`ℝ_`[a](t), a = 1 .. n))/(m*n+M)+m*r__M_(t)*n/(m*n+M) = 0, r__M_(t)))

r__M_(t) = -m*(Sum(`ℝ_`[a](t), a = 1 .. n))/(m*n+M)

(395)

Likewise, to eliminate R__M_ use

isolate(`ℝ_`[a](t) = r_[a](t)-r__M_(t), r_[a](t))

r_[a](t) = `ℝ_`[a](t)+r__M_(t)

(396)

Substitute, sequentially, these two equations into the Lagrangian (389)

subs(r_[a](t) = `ℝ_`[a](t)+r__M_(t), r__M_(t) = -m*(Sum(`ℝ_`[a](t), a = 1 .. n))/(m*n+M), L = (1/2)*M*Physics[Vectors][Norm](diff(r__M_(t), t))^2+(1/2)*m*(Sum(Physics[Vectors][Norm](diff(r_[a](t), t))^2, a = 1 .. n))-U)

L = (1/2)*M*Physics:-Vectors:-Norm(diff(-m*(Sum(`ℝ_`[a](t), a = 1 .. n))/(m*n+M), t))^2+(1/2)*m*(Sum(Physics:-Vectors:-Norm(diff(`ℝ_`[a](t)-m*(Sum(`ℝ_`[a](t), a = 1 .. n))/(m*n+M), t))^2, a = 1 .. n))-U

(397)

To verify by eye each step for correctness, one can manipulate this expression surgically using subsindets . First expand only the Norms 

subsindets(L = (1/2)*M*Physics[Vectors][Norm](diff(-m*(Sum(`ℝ_`[a](t), a = 1 .. n))/(m*n+M), t))^2+(1/2)*m*(Sum(Physics[Vectors][Norm](diff(`ℝ_`[a](t)-m*(Sum(`ℝ_`[a](t), a = 1 .. n))/(m*n+M), t))^2, a = 1 .. n))-U, specfunc(Norm), expand)

L = (1/2)*M*m^2*Physics:-Vectors:-Norm(Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2/(m*n+M)^2+(1/2)*m*(Sum(Physics:-Vectors:-Norm(diff(`ℝ_`[a](t), t))^2-2*m*Physics:-Vectors:-`.`(diff(`ℝ_`[a](t), t), Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))/(m*n+M)+m^2*Physics:-Vectors:-Norm(Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2/(m*n+M)^2, a = 1 .. n))-U

(398)

This result is correct. Next expand only the Sums 

subsindets(L = (1/2)*M*m^2*Physics[Vectors][Norm](Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2/(m*n+M)^2+(1/2)*m*(Sum(Physics[Vectors][Norm](diff(`ℝ_`[a](t), t))^2-2*m*Physics[Vectors][`.`](diff(`ℝ_`[a](t), t), Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))/(m*n+M)+m^2*Physics[Vectors][Norm](Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2/(m*n+M)^2, a = 1 .. n))-U, specfunc(Sum), expand)

L = (1/2)*M*m^2*Physics:-Vectors:-Norm(Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2/(m*n+M)^2+(1/2)*m*(m^2*n^2*(Sum(Physics:-Vectors:-Norm(diff(`ℝ_`[a](t), t))^2, a = 1 .. n))/(m*n+M)^2+2*M*m*n*(Sum(Physics:-Vectors:-Norm(diff(`ℝ_`[a](t), t))^2, a = 1 .. n))/(m*n+M)^2+M^2*(Sum(Physics:-Vectors:-Norm(diff(`ℝ_`[a](t), t))^2, a = 1 .. n))/(m*n+M)^2-m^2*Physics:-Vectors:-Norm(Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2*n/(m*n+M)^2-2*M*m*Physics:-Vectors:-Norm(Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2/(m*n+M)^2)-U

(399)

This result also verifies visually. Collect the terms polynomial in Norm then Sum and normalize the coefficients

collect(L = (1/2)*M*m^2*Physics[Vectors][Norm](Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2/(m*n+M)^2+(1/2)*m*(m^2*n^2*(Sum(Physics[Vectors][Norm](diff(`ℝ_`[a](t), t))^2, a = 1 .. n))/(m*n+M)^2+2*M*m*n*(Sum(Physics[Vectors][Norm](diff(`ℝ_`[a](t), t))^2, a = 1 .. n))/(m*n+M)^2+M^2*(Sum(Physics[Vectors][Norm](diff(`ℝ_`[a](t), t))^2, a = 1 .. n))/(m*n+M)^2-m^2*Physics[Vectors][Norm](Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2*n/(m*n+M)^2-2*M*m*Physics[Vectors][Norm](Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2/(m*n+M)^2)-U, [Norm, Sum], normal)

L = -(1/2)*m^2*Physics:-Vectors:-Norm(Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2/(m*n+M)+(1/2)*m*(Sum(Physics:-Vectors:-Norm(diff(`ℝ_`[a](t), t))^2, a = 1 .. n))-U

(400)

This Lagrangian involves only the n relative position vectors `#mover(mi("ℝ"),mo("→"))`[a], achieving the reduction of the original problem of n+1 particles to a problem of only n particles.

_______________________________________

 

b) To obtain the result `≡`(L = (1/2)*Physics[Vectors][Norm](diff(`ℝ_`(t), t))^2*mu-U(Physics[Vectors][Norm](`ℝ_`(t))), L) = (1/2)*LinearAlgebra[Norm](diff(`#mover(mi("ℝ"),mo("→"))`(t), t))^2*mu-U for the two-body problem, substitute n = 1 into `≡`(L = -(1/2)*m^2*Physics[Vectors][Norm](Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2/(m*n+M)+(1/2)*m*(Sum(Physics[Vectors][Norm](diff(`ℝ_`[a](t), t))^2, a = 1 .. n))-U, L) = -(1/2)*m^2*LinearAlgebra[Norm](Sum(diff(`#mover(mi("ℝ"),mo("→"))`[a](t), t), a = 1 .. n))^2/(m*n+M)+(1/2)*m*(Sum(LinearAlgebra[Norm](diff(`#mover(mi("ℝ",fontstyle = "normal"),mo("→"))`[a](t), t))^2, a = 1 .. n))-U

subs(n = 1, L = -(1/2)*m^2*Physics[Vectors][Norm](Sum(diff(`ℝ_`[a](t), t), a = 1 .. n))^2/(m*n+M)+(1/2)*m*(Sum(Physics[Vectors][Norm](diff(`ℝ_`[a](t), t))^2, a = 1 .. n))-U)

L = -(1/2)*m^2*Physics:-Vectors:-Norm(Sum(diff(`ℝ_`[a](t), t), a = 1 .. 1))^2/(M+m)+(1/2)*m*(Sum(Physics:-Vectors:-Norm(diff(`ℝ_`[a](t), t))^2, a = 1 .. 1))-U

(401)

Expand only the Sums and collect the Norm

subsindets(L = -(1/2)*m^2*Physics[Vectors][Norm](Sum(diff(`ℝ_`[a](t), t), a = 1 .. 1))^2/(M+m)+(1/2)*m*(Sum(Physics[Vectors][Norm](diff(`ℝ_`[a](t), t))^2, a = 1 .. 1))-U, specfunc(Sum), expand)

L = -(1/2)*m^2*Physics:-Vectors:-Norm(diff(`ℝ_`[1](t), t))^2/(M+m)+(1/2)*m*Physics:-Vectors:-Norm(diff(`ℝ_`[1](t), t))^2-U

(402)

collect(L = -(1/2)*m^2*Physics[Vectors][Norm](diff(`ℝ_`[1](t), t))^2/(M+m)+(1/2)*m*Physics[Vectors][Norm](diff(`ℝ_`[1](t), t))^2-U, Norm, normal)

L = (1/2)*M*m*Physics:-Vectors:-Norm(diff(`ℝ_`[1](t), t))^2/(M+m)-U

(403)

Comparing with the definition of reduced mass (384)

m[1]*m[2]/(m[1]+m[2]) = mu

m[1]*m[2]/(m[1]+m[2]) = mu

(404)

we see the substitutions to transform (403) into (384) are

M = m[1], m = m[2], r_[1] = R_

M = m[1], m = m[2], r_[1] = R_

(405)

Substitute them, simultaneously (enclose the sequence of equations into a list), then simplify taking m[1]*m[2]/(m[1]+m[2]) = mu into account

simplify(subs([M = m[1], m = m[2], r_[1] = R_], L = (1/2)*M*m*Physics[Vectors][Norm](diff(`ℝ_`[1](t), t))^2/(M+m)-U), {m[1]*m[2]/(m[1]+m[2]) = mu})

L = (1/2)*Physics:-Vectors:-Norm(diff(`ℝ_`[1](t), t))^2*mu-U

(406)

NULL

This is the same as (385), the reduced Lagrangian for the two-body problem.

 

Motion in a central field

 

A one-body problem in a central field, is about the motion of a single particle of mass m in a field where the potential energy, and so the magnitude of the force, depends only on the distance between the particle and a fixed point, frequently chosen as the origin of the reference system. As seen above, the two-body problem , is reducible to a one-body problem in a central field.

 

Problem

The angular momentum `#mover(mi("L"),mo("→"))` = `&x`(`#mover(mi("r"),mo("→"))`, `#mover(mi("p"),mo("→"))`) of a particle that moves in a central field is conserved, so `#mover(mi("r"),mo("→"))` evolves in time on a fixed plane perpendicular to the constant `#mover(mi("L"),mo("→"))`. Show that the surface swept per second by the position vector `#mover(mi("r"),mo("→"))` is constant (Kepler's second law).

Solution

   

 

Problem
Starting from the constancy of the energy E and the angular momentum `#mover(mi("L"),mo("→"))`, compute the equations of movement and integrate them according to:

a) using differential elimination techniques to uncouple the system of equations of movement that involve both of rho(t) and phi(t) 

b) interactively, one step at a time, uncouple the equations of movement eliminating phi from the problem, resulting in an implicit solution `≡`(t, t(rho)).

c) eliminate t from the problem to obtain an equation for diff(phi(rho), rho), whose solution is the trajectory as phi(rho)

Solution

   

 

Kepler's problem

 

 

Problem

Show that, when U(rho) = -alpha/rho, with "alpha>0,"the solution (442) for the motion in a central field, part c) of the previous problem,

phi(rho) = Int(L__z/(rho^2*(2*m*(E-U(rho))-L__z^2/rho^2)^(1/2)), rho)+c__1

becomes the equation of a conic section

p/rho = 1+xi*cos(phi) 

where p = L__z^2/(m*alpha), xi = sqrt(1+2*E*L__z^2/(m*alpha^2)).

Solution

   

 

Small Oscillations

 

Free oscillations in one dimension

 

 

Problem

Consider the case of 1-dimensional motion where the acting force opposes the motion as a function of the position, `#mover(mi("F"),mo("→"))` = -k*x*`#mover(mi("i"),mo("∧"))`. This is the case, for example, of a spring, the more you stretch it (the bigger x), the more it opposes the stretching in the opposite direction -`#mover(mi("i"),mo("∧"))`. Write the equation of motion as Newton's 2nd law, then the Lagrangian and Lagrange equation for it, and integrate the equation for generic initial conditions

Solution

   

 

Forced oscillations

 

 

Problem

Consider oscillations in 1 dimension of a system on which an external force `#mover(mi("F"),mo("→"))`[ext](t) acts. For the oscillations to be small, `#mover(mi("F"),mo("→"))`[ext](t) must produce only small displacements. The total force is  `#mover(mi("F"),mo("→"))` = -k*x*`#mover(mi("i"),mo("∧"))`+`#mover(mi("F"),mo("→"))`[ext](t).

a) Write the equation of motion as Newton's 2nd law, then write the Lagrangian and Lagrange equations.

b) Integrate the equation of motion for generic initial conditions.

c) Specialize the solution computed in b) for LinearAlgebra[Norm](`#mover(mi("F"),mo("→"))`[ext](t)) = f__0*cos(lambda*t+beta) to obtain

x(t) = (a*cos(omega*t+alpha)-f__0*omega*cos(lambda*t+beta)/(m*(lambda^2-omega^2))-f__0*cos(-omega*t+beta)/(2*m*(lambda+omega))+f__0*cos(omega*t+beta)/((2*(lambda-omega))*m))/omega

for some constants a, b, alpha, beta.

d) Show that a solution for the case considered in c), that is, LinearAlgebra[Norm](`#mover(mi("F"),mo("→"))`[ext](t)) = f__0*cos(lambda*t+beta), can be computed manually to get

x(t) = a*cos(omega*t+alpha)+f__0*cos(lambda*t+beta)/(m*(-lambda^2+omega^2))

for some other constants a, b, alpha, beta.

e) Specialize the solution of item d) in the case of resonance, when lambda = omega, by taking limits, thus obtaining

x(t) = a*cos(omega*t+alpha)+f__0*t*sin(omega*t+beta)/(2*omega*m)

f) Show that the solution computed taking limits in e) can be compute directly by using dsolve and specializing the integration constants c__1 and c__2 that appear when solving the underlying differential equation.

Solution

   

 

Oscillations of systems with many degrees of freedom

 

 

Problem
Formulate the equations of motion for the free oscillations of a system with n degrees of freedom as

m[a, c]*(diff(x[a](t), t, t))+`k__a,c`*x[a] = 0

where x[a] represents the displacement of the a^th generalized coordinate q__a, the index a runs from 1 to n and there is an implicit sum over repeated indices (Einstein's convention).

Solution

 

 

restart; with(Physics); with(Vectors); CompactDisplay(x(t))

x(t)*`will now be displayed as`*x

(524)

 

Denoting the generalized coordinates by q__a, the potential energy U(q__a) can be expanded in series around the minimums q__a0. Since the movement consist only in small displacements x__a = q__a-q__a0 around q__a0, it is sufficient to keep terms in the expansion up to order 2, resulting in an expression analogous to U = (1/2)*k*x^2 of the 1-dimensional case:

U = (1/2)*(Sum(Sum(k[a, b]*x[a](t)*x[b](t), i = 1 .. n), j = 1 .. n))

U = (1/2)*(Sum(Sum(k[a, b]*x[a](t)*x[b](t), i = 1 .. n), j = 1 .. n))

(525)

Likewise, for the kinetic energy,

T = (1/2)*(Sum(Sum(A[a, b](q__0)*(diff(x[a](t), t))*(diff(x[b](t), t)), a = 1 .. n), b = 1 .. n))

T = (1/2)*(Sum(Sum(A[a, b](q__0)*(diff(x[a](t), t))*(diff(x[b](t), t)), a = 1 .. n), b = 1 .. n))

(526)

where we take the A[a, b] at the minimums `#msub(mi("q"),mn("0"))` and denote them as `m__i,j`

subs(A[a, b](q__0) = m[a, b], T = (1/2)*(Sum(Sum(A[a, b](q__0)*(diff(x[a](t), t))*(diff(x[b](t), t)), a = 1 .. n), b = 1 .. n)))

T = (1/2)*(Sum(Sum(m[a, b]*(diff(x[a](t), t))*(diff(x[b](t), t)), a = 1 .. n), b = 1 .. n))

(527)

Both `k__a,b` and `m__a,b` can be split into symmetric and antisymmetric parts, with the antisymmetric parts canceling out in view of the symmetric character of x[a]*x[b] in U and (diff(x[a](t), t))*(diff(x[b](t), t)) in T. Therefore, we can take `k__a,b` and `m__a,b` as symmetric without any loss of generality.

 

Before proceeding, note the similarity in notation between the three formulas (525) to (527) for T and U and tensor notation. In T and U the x__a describe independent displacements, so one can think of `#msub(mi("x"),mi("a"))` as a tensor in an Euclidean space of displacements of generic abstract dimension n, with KroneckerDelta  as the metric. It is then simpler to write the Lagrangian using tensor notation with a generic type of index (that admits and abstract n-dimension). For this purpose, introduce lowercaselatin indices from a to h to represent generic indices, and when necessary use KroneckerDelta  as the metric.

Setup(genericindices = lowercaselatin_ah)

[genericindices = lowercaselatin_ah]

(528)

Now introduce the tensors while making sure to indicate `m__a,b` and `k__a,b` are symmetric (passing x__a together is not a problem, it has only one index)

Define(x[a], m[a, b], k[a, b], symmetric)

{Physics:-Dgamma[mu], Physics:-Psigma[mu], Physics:-d_[mu], Physics:-g_[mu, nu], k[a, b], m[a, b], x[a], Physics:-LeviCivita[alpha, beta, mu, nu]}

(529)

The Lagrangian L = T-U can now be written as

L := 1/2*(m[a, b]*(diff(x[a](t), t))*(diff(x[b](t), t))-k[a, b]*x[a](t)*x[b](t))

(1/2)*m[a, b]*(diff(x[a](t), t))*(diff(x[b](t), t))-(1/2)*k[a, b]*x[a](t)*x[b](t)

(530)

where Einstein's summation rule for repeated indices is used. Einstein's rule is taken into account by the system when differentiating, computing products and simplifying tensor indices. The simplest way to compute the Lagrange equations is to use the LagrangeEquations  command

LagrangeEquations(L, x)

k[a, c]*x[a](t)+m[a, c]*(diff(diff(x[a](t), t), t)) = 0

(531)

The same result can be computed via %diff(%diff(L, diff(x[c](t), t)), t) = %diff(L, x[c](t)). Although not necessary, enclose the operation with forward quotes to delay its evaluation in order to see what is being computed

'%diff(%diff(L, diff(x[c](t), t)), t) = %diff(L, x[c](t))'

%diff(%diff(L, Physics:-Vectors:-diff(x[c](t), t)), t) = %diff(L, x[c](t))

(532)

value(%diff(%diff(L, Physics[Vectors][diff](x[c](t), t)), t) = %diff(L, x[c](t)))

(1/2)*m[a, b]*(diff(diff(x[a](t), t), t))*Physics:-KroneckerDelta[b, c]+(1/2)*m[a, b]*(diff(diff(x[b](t), t), t))*Physics:-KroneckerDelta[a, c] = -(1/2)*k[a, b]*x[a](t)*Physics:-KroneckerDelta[b, c]-(1/2)*k[a, b]*x[b](t)*Physics:-KroneckerDelta[a, c]

(533)

Simplifying tensor indices,

Simplify((1/2)*m[a, b]*(diff(diff(x[a](t), t), t))*Physics[KroneckerDelta][b, c]+(1/2)*m[a, b]*(diff(diff(x[b](t), t), t))*Physics[KroneckerDelta][a, c] = -(1/2)*k[a, b]*x[a](t)*Physics[KroneckerDelta][b, c]-(1/2)*k[a, b]*x[b](t)*Physics[KroneckerDelta][a, c])

m[a, c]*(diff(diff(x[a](t), t), t)) = -k[a, c]*x[a](t)

(534)

NULL

which is the same as (531).

``

 

Rigid-body motion

 

 

A rigid body is one where (in approximation) the distances between the body's parts remain unchanged. In what follows, for simplicity, the body is considered as discrete set of particles; the formulas for a continuous body can be obtained from those by replacing the masses m__i of each particle by rho(`#mover(mi("r"),mo("→"))`)*dV, where rho(`#mover(mi("r"),mo("→"))`) is the mass density as a function of the position and dV is the volume element, whose integration represent the body's volume.

 

This problem is systematically treated by using two reference systems: an inertial one K, where the observer is, and another one K', rigidly attached to the body, that moves with it and thus it is typically non-inertial. It is customary (not necessary) to place the origin `#mover(mi("R"),mo("→"))`(t) of K'  at the body's center of mass .

 

A rigid body is thus a system with six degrees of freedom: three indicating the position `#mover(mi("R"),mo("→"))`(t) of the center of mass plus three angles specifying the orientation of the axes of K' with respect to those of K.

 

Angular velocity

 

 

Problem
a) Show, using graphs, that the velocity `#mover(mi("v"),mo("→"))` of a point P of a body, measured in an inertial reference system K, can be written as

`#mover(mi("v"),mo("→"))` = `#mover(mi("V"),mo("→"))`+`&x`(`#mover(mi("Ω",fontstyle = "normal"),mo("→"))`, `#mover(mi("r'"),mo("→"))`)

where `#mover(mi("V"),mo("→"))` is the velocity of the origin of K', a frame of reference attached to the body's center of mass, `#mover(mi("Ω",fontstyle = "normal"),mo("→"))` = (diff(Phi(t), t))*`#mover(mi("Φ",fontstyle = "normal"),mo("∧"))` is the body's angular velocity (its instantaneous counter-clockwise rotation speed around some axis in the direction of a unit vector`#mover(mi("Φ",fontstyle = "normal"),mo("∧"))`) and `#mover(mi("r'"),mo("→"))`  is the distance from the center of mass (origin of K') to the point P.

b) Derive algebraically the same result of a) , using the fact that vectors are defined up to parallel translation and so `#mover(mi("r"),mo("→"))` and `#mover(mi("r'"),mo("→"))` are related by a rotation matrix omega[a, b] which, as all rotation matrices, is orthogonal.

Solution

   

 

Inertia tensor

 

 

Problem

a) Show that using `#mover(mi("v"),mo("→"))` = `#mover(mi("V"),mo("→"))`+`&x`(`#mover(mi("Ω",fontstyle = "normal"),mo("→"))`, `#mover(mi("r'"),mo("→"))`), derived in Angular velocity  for the velocity `#mover(mi("v"),mo("→"))` of a point P of a rigid body in terms of `#mover(mi("Ω",fontstyle = "normal"),mo("→"))` and the position `#mover(mi("r'"),mo("→"))` of P viewed from the center of mass `#mover(mi("R"),mo("→"))`, the kinetic energy of the rigid body can be written in terms of the positions `#mover(mi("r'"),mo("→"))` of the n particles (not their velocities), the velocity of the center of mass `#mover(mi("V"),mo("→"))` and angular velocity `#mover(mi("Ω",fontstyle = "normal"),mo("→"))` as

T = (1/2)*LinearAlgebra[Norm](`#mover(mi("V"),mo("→"))`)^2*mu+Sum((1/2)*m[i]*LinearAlgebra[Norm](`#mover(mi("Ω",fontstyle = "normal"),mo("→"))`)^2*LinearAlgebra[Norm](`#mover(mi("r'"),mo("→"))`[i])^2-(1/2)*m[i]*(`#mover(mi("Ω",fontstyle = "normal"),mo("→"))`.`#mover(mi("r'"),mo("→"))`[i])^2, i = 1 .. n)

b) Use tensor notation to show that this result can be rewritten as

T = (1/2)*LinearAlgebra[Norm](`#mover(mi("V"),mo("→"))`(t))^2*mu+(1/2)*`𝕀`[a, b]*Omega[b]*Omega[a]

where

"`𝕀`[a,b]=(∑)m[i] (r[i,c]^2 delta[a,b]-`r'`[i,a] `r'`[i,b])"

is the Inertia tensor, `Ω__a` represents the components of the vector `#mover(mi("Ω",fontstyle = "normal"),mo("→"))` and `r'`[i, a] represents the components of the position vector `#mover(mi("r'"),mo("→"))`[i] of the i^th particle.

Solution

   

 

NULL

Problem

Determine the Inertia tensor corresponding to a triatomic molecule that has the form of an isosceles triangle with two masses m[1] in the extremes of the base and a mass m[2] in the third vertex. The distance between the two masses m[1] is equal to a, and the height of the triangle is equal to h.

Solution

   

 

Angular momentum of a rigid body

 

 

In the section related to the conservation of angular momentum , the solution to the second Problem  shows that the value of the angular momentum `#mover(mi("L"),mo("→"))` of a system of particles depends on the origin of the frame of reference. In this section, it is assumed that the origin is at the center of mass, so Sum(m[i]*`#mover(mi("r"),mo("→"))`[i], i = 1 .. n) = 0.

 

Problem
Show, using tensor notation, that in the K' system whose origin is at the center of mass, the components (diff(L(x), x))[a] of the angular momentum of a rigid body can be expressed in terms of the inertia tensor `𝕀`[a, b] and the components of the angular velocity Omega[b] as

(diff(L(x), x))[a] = `𝕀`[a, b]*Omega[b]

Solution

   

 

The equations of motion of a rigid body

 

 

A rigid body is a system with six degrees of freedom: three indicating the position `#mover(mi("R"),mo("→"))`(t) plus three angles specifying the orientation of the axes of K' with respect to those of K. As discussed in the equations of motion  for many-particle systems , the two vectorial equations of motion are

value(%diff(P_(t), t) = F_[ext])

diff(P_(t), t) = F_[ext]

(615)

diff(L_(t), t) = N_(t)

diff(L_(t), t) = N_(t)

(616)

N_ = Sum(Physics[Vectors][`&x`](r_[i](t), f_[i, ext]), i = 1 .. n)

N_ = Sum(Physics:-Vectors:-`&x`(r_[i](t), f_[i, ext]), i = 1 .. n)

(617)

NULL

where `#mover(mi("P"),mo("→"))` is the total momentum, `#mover(mi("F"),mo("→"))`[ext] is the total external force acting upon the body, `#mover(mi("f"),mo("→"))`[i, ext] is the external force acting upon the i^th particle, `#mover(mi("L"),mo("→"))`the total angular momentum and `#mover(mi("N"),mo("→"))` is the total torque.

 

Problem
Show that the equations of movement of a rigid body can be computed as the Lagrange equations for `#mover(mi("R"),mo("→"))` and `#mover(mi("Φ",fontstyle = "normal"),mo("→"))` from the Lagrangian

L = (1/2)*mu*LinearAlgebra[Norm](`#mover(mi("V"),mo("→"))`)^2+(1/2)*`𝕀`[a, b]*Omega[b]*Omega[a]+U(`#mover(mi("R"),mo("→"))`, `#mover(mi("Φ",fontstyle = "normal"),mo("→"))`)

where mu = Sum(m[i], i = 1 .. n) is the total mass, `𝕀`[a, b] is the inertia tensor, U(`#mover(mi("R"),mo("→"))`, `#mover(mi("Φ",fontstyle = "normal"),mo("→"))`) is the potential energy for the external force `#mover(mi("F"),mo("→"))`[ext], and  `#mover(mi("Ω",fontstyle = "normal"),mo("→"))` = diff(`#mover(mi("Φ",fontstyle = "normal"),mo("→"))`(t), t).

Solution

 

 

restart; with(Physics:-Vectors); CompactDisplay((R_, V_, Phi_, Omega)(t))

Omega(t)*`will now be displayed as`*Omega

(618)

The required derivation is easy by expressing the Lagrangian in tensor notation from the beginning. In what follows, however, with the purpose of illustrating different techniques, Lagrange's equations are computed using vectorial notation, only switching to tensor notation at the time of expressing the time derivative of the angular momentum.

 

The kinetic energy T in vectorial form is derived in this problem  for the inertia tensor

subs((diff(r(x), x))*_[i](t) = (diff(r(x), x))*_[i], T = (1/2)*(Sum(m[i]*Physics[Vectors][Norm](Omega_(t))^2*Physics[Vectors][Norm](`r'_`[i](t))^2-m[i]*Physics[Vectors][`.`](Omega_(t), `r'_`[i](t))^2, i = 1 .. n))+(1/2)*Physics[Vectors][Norm](V_(t))^2*mu)

T = (1/2)*(Sum(m[i]*Physics:-Vectors:-Norm(Omega_(t))^2*Physics:-Vectors:-Norm(`r'_`[i])^2-m[i]*Physics:-Vectors:-`.`(Omega_(t), `r'_`[i])^2, i = 1 .. n))+(1/2)*Physics:-Vectors:-Norm(V_(t))^2*mu

(619)

Adding the potential energy as a function of the coordinates `#mover(mi("R"),mo("→"))` (location of the center of mass) and `#mover(mi("Φ",fontstyle = "normal"),mo("→"))` (the rotation axis, so that `#mover(mi("Ω",fontstyle = "normal"),mo("→"))` = diff(`#mover(mi("Φ",fontstyle = "normal"),mo("→"))`(t), t)), the Lagrangian in vectorial form is given by

L = rhs(T = (1/2)*(Sum(m[i]*Physics[Vectors][Norm](Omega_(t))^2*Physics[Vectors][Norm](`r'_`[i])^2-m[i]*Physics[Vectors][`.`](Omega_(t), `r'_`[i])^2, i = 1 .. n))+(1/2)*Physics[Vectors][Norm](V_(t))^2*mu)+U(R_(t), Phi_(t))

L = (1/2)*(Sum(m[i]*Physics:-Vectors:-Norm(Omega_(t))^2*Physics:-Vectors:-Norm(`r'_`[i])^2-m[i]*Physics:-Vectors:-`.`(Omega_(t), `r'_`[i])^2, i = 1 .. n))+(1/2)*Physics:-Vectors:-Norm(V_(t))^2*mu+U(R_(t), Phi_(t))

(620)

The first equation of movement, for the total momentum `#mover(mi("P"),mo("→"))`(t) is derived as the Lagrange equation for this Lagrangian

%diff(%diff(L, V_(t)), t) = %diff(L, R_(t))

%diff(%diff(L, V_(t)), t) = %diff(L, R_(t))

(621)

subs(L = (1/2)*(Sum(m[i]*Physics[Vectors][Norm](Omega_(t))^2*Physics[Vectors][Norm](`r'_`[i])^2-m[i]*Physics[Vectors][`.`](Omega_(t), `r'_`[i])^2, i = 1 .. n))+(1/2)*Physics[Vectors][Norm](V_(t))^2*mu+U(R_(t), Phi_(t)), %diff(%diff(L, V_(t)), t) = %diff(L, R_(t)))

%diff(%diff((1/2)*(Sum(m[i]*Physics:-Vectors:-Norm(Omega_(t))^2*Physics:-Vectors:-Norm(`r'_`[i])^2-m[i]*Physics:-Vectors:-`.`(Omega_(t), `r'_`[i])^2, i = 1 .. n))+(1/2)*Physics:-Vectors:-Norm(V_(t))^2*mu+U(R_(t), Phi_(t)), V_(t)), t) = %diff((1/2)*(Sum(m[i]*Physics:-Vectors:-Norm(Omega_(t))^2*Physics:-Vectors:-Norm(`r'_`[i])^2-m[i]*Physics:-Vectors:-`.`(Omega_(t), `r'_`[i])^2, i = 1 .. n))+(1/2)*Physics:-Vectors:-Norm(V_(t))^2*mu+U(R_(t), Phi_(t)), R_(t))

(622)

value(%)

mu*(diff(V_(t), t)) = (D[1](U))(R_(t), Phi_(t))

(623)

On the right-hand side, (D[1](U))(`#mover(mi("R"),mo("→"))`, `#mover(mi("Φ",fontstyle = "normal"),mo("→"))`)is*the*derivative*ofU(`#mover(mi("R"),mo("→"))`, `#mover(mi("Φ",fontstyle = "normal"),mo("→"))`) with respect to its first argument `#mover(mi("R"),mo("→"))`, equivalent to the gradient taking `#mover(mi("R"),mo("→"))` as the coordinates, and mu*(diff(`#mover(mi("V"),mo("→"))`(t), t)) = diff(`#mover(mi("P"),mo("→"))`(t), t) is the total momentum

subs(rhs(mu*(diff(V_(t), t)) = (D[1](U))(R_(t), Phi_(t))) = F_[ext](t), diff(V_(t), t) = (diff(P_(t), t))/mu, mu*(diff(V_(t), t)) = (D[1](U))(R_(t), Phi_(t)))

diff(P_(t), t) = F_[ext](t)

(624)

Note however that in Maple the Vectors:-Gradient  command computes the gradient with respect to Cartesian, cylindrical or spherical coordinates, not an arbitrary vector `#mover(mi("R"),mo("→"))`. If more precision is required, the dependency of the potential energy could be expressed in terms of the norm of `#mover(mi("R"),mo("→"))`, as in

U(Norm(R_))

U(Physics:-Vectors:-Norm(R_))

(625)

in which case differentiating with respect to `#mover(mi("R"),mo("→"))` is equivalent to the definition of directional derivative

(%diff = diff)(U(Physics[Vectors][Norm](R_)), `#mover(mi("R"),mo("→"))`)

%diff(U(Physics:-Vectors:-Norm(R_)), R_) = (D(U))(Physics:-Vectors:-Norm(R_))*R_/Physics:-Vectors:-Norm(R_)

(626)

The same computation, this time with respect to the coordinates `#mover(mi("Φ",fontstyle = "normal"),mo("→"))` where `#mover(mi("Ω",fontstyle = "normal"),mo("→"))` = diff(`#mover(mi("Φ",fontstyle = "normal"),mo("→"))`(t), t) is the corresponding velocity,

%diff(%diff(L, Omega_(t)), t) = %diff(L, Phi_(t))

%diff(%diff(L, Omega_(t)), t) = %diff(L, Phi_(t))

(627)

subs(L = (1/2)*(Sum(m[i]*Physics[Vectors][Norm](Omega_(t))^2*Physics[Vectors][Norm](`r'_`[i])^2-m[i]*Physics[Vectors][`.`](Omega_(t), `r'_`[i])^2, i = 1 .. n))+(1/2)*Physics[Vectors][Norm](V_(t))^2*mu+U(R_(t), Phi_(t)), %diff(%diff(L, Omega_(t)), t) = %diff(L, Phi_(t)))

%diff(%diff((1/2)*(Sum(m[i]*Physics:-Vectors:-Norm(Omega_(t))^2*Physics:-Vectors:-Norm(`r'_`[i])^2-m[i]*Physics:-Vectors:-`.`(Omega_(t), `r'_`[i])^2, i = 1 .. n))+(1/2)*Physics:-Vectors:-Norm(V_(t))^2*mu+U(R_(t), Phi_(t)), Omega_(t)), t) = %diff((1/2)*(Sum(m[i]*Physics:-Vectors:-Norm(Omega_(t))^2*Physics:-Vectors:-Norm(`r'_`[i])^2-m[i]*Physics:-Vectors:-`.`(Omega_(t), `r'_`[i])^2, i = 1 .. n))+(1/2)*Physics:-Vectors:-Norm(V_(t))^2*mu+U(R_(t), Phi_(t)), Phi_(t))

(628)

subs(sum = Sum, value(%diff(%diff((1/2)*(Sum(m[i]*Physics[Vectors][Norm](Omega_(t))^2*Physics[Vectors][Norm](`r'_`[i])^2-m[i]*Physics[Vectors][`.`](Omega_(t), `r'_`[i])^2, i = 1 .. n))+(1/2)*Physics[Vectors][Norm](V_(t))^2*mu+U(R_(t), Phi_(t)), Omega_(t)), t) = %diff((1/2)*(Sum(m[i]*Physics[Vectors][Norm](Omega_(t))^2*Physics[Vectors][Norm](`r'_`[i])^2-m[i]*Physics[Vectors][`.`](Omega_(t), `r'_`[i])^2, i = 1 .. n))+(1/2)*Physics[Vectors][Norm](V_(t))^2*mu+U(R_(t), Phi_(t)), Phi_(t))))

(1/2)*(Sum(2*m[i]*Physics:-Vectors:-Norm(`r'_`[i])^2*(diff(Omega_(t), t))-2*m[i]*Physics:-Vectors:-`.`(diff(Omega_(t), t), `r'_`[i])*`r'_`[i], i = 1 .. n)) = (D[2](U))(R_(t), Phi_(t))

(629)

Switching to tensor notation as done in the other problems, the left-hand side can be rewritten as the time derivative of the components (diff(L(x), x))[a] of the angular momentum .

with(Physics)

Setup(spaceindices = lowercase_ah, generic = lowercase_is, tensors = {Omega[a], (diff(r(x), x))[i, a]})

[genericindices = lowercaselatin_is, spaceindices = lowercaselatin_ah, tensors = {Physics:-Dgamma[mu], Omega[a], Physics:-Psigma[mu], Physics:-d_[mu], Physics:-g_[mu, nu], Physics:-gamma_[a, b], `r'`[i, a], Physics:-LeviCivita[alpha, beta, mu, nu]}]

(630)

Introducing tensor components for the vectors of (629)

Typesetting[delayDotProduct]((diff(Omega_(t), t)).(diff(r(x), x)), _[i], true) = (diff(Omega[b](t), t))*(diff(r(x), x))[i, b], diff(Omega_(t), t) = KroneckerDelta[a, b]*(diff(Omega[b](t), t)), Norm((diff(r(x), x))*_[i])^2 = (diff(r(x), x))[i, c]^2, (diff(r(x), x))*_[i] = (diff(r(x), x))[i, a], (D[2](U))(R_(t), Phi_(t)) = Component((D[2](U))(R_(t), Phi_(t)), a)

Physics:-Vectors:-`.`(diff(Omega_(t), t), `r'_`[i]) = (diff(Omega[b](t), t))*`r'`[i, b], diff(Omega_(t), t) = Physics:-KroneckerDelta[a, b]*(diff(Omega[b](t), t)), Physics:-Vectors:-Norm(`r'_`[i])^2 = `r'`[i, c]^2, `r'_`[i] = `r'`[i, a], (D[2](U))(R_(t), Phi_(t)) = Physics:-Vectors:-Component((D[2](U))(R_(t), Phi_(t)), a)

(631)

eval((1/2)*(Sum(2*m[i]*Physics[Vectors][Norm](`r'_`[i])^2*(diff(Omega_(t), t))-2*m[i]*Physics[Vectors][`.`](diff(Omega_(t), t), `r'_`[i])*`r'_`[i], i = 1 .. n)) = (D[2](U))(R_(t), Phi_(t)), [Physics[Vectors][`.`](diff(Omega_(t), t), `r'_`[i]) = (diff(Omega[b](t), t))*`r'`[i, b], diff(Omega_(t), t) = Physics[KroneckerDelta][a, b]*(diff(Omega[b](t), t)), Physics[Vectors][Norm](`r'_`[i])^2 = `r'`[i, c]^2, `r'_`[i] = `r'`[i, a], (D[2](U))(R_(t), Phi_(t)) = Physics[Vectors][Component]((D[2](U))(R_(t), Phi_(t)), a)])

(1/2)*(Sum(2*m[i]*`r'`[i, c]^2*Physics:-KroneckerDelta[a, b]*(diff(Omega[b](t), t))-2*m[i]*(diff(Omega[b](t), t))*`r'`[i, b]*`r'`[i, a], i = 1 .. n)) = Physics:-Vectors:-Component((D[2](U))(R_(t), Phi_(t)), a)

(632)

expand((1/2)*(Sum(2*m[i]*`r'`[i, c]^2*Physics[KroneckerDelta][a, b]*(diff(Omega[b](t), t))-2*m[i]*(diff(Omega[b](t), t))*`r'`[i, b]*`r'`[i, a], i = 1 .. n)) = Physics[Vectors][Component]((D[2](U))(R_(t), Phi_(t)), a))

Physics:-KroneckerDelta[a, b]*(diff(Omega[b](t), t))*(Sum(m[i]*`r'`[i, c]^2, i = 1 .. n))-(diff(Omega[b](t), t))*(Sum(m[i]*`r'`[i, a]*`r'`[i, b], i = 1 .. n)) = Physics:-Vectors:-Component((D[2](U))(R_(t), Phi_(t)), a)

(633)

collect(Physics[KroneckerDelta][a, b]*(diff(Omega[b](t), t))*(Sum(m[i]*`r'`[i, c]^2, i = 1 .. n))-(diff(Omega[b](t), t))*(Sum(m[i]*`r'`[i, a]*`r'`[i, b], i = 1 .. n)) = Physics[Vectors][Component]((D[2](U))(R_(t), Phi_(t)), a), diff(Omega[b](t), t), `@`(factor, combine))

(Sum(-m[i]*(-`r'`[i, c]^2*Physics:-KroneckerDelta[a, b]+`r'`[i, a]*`r'`[i, b]), i = 1 .. n))*(diff(Omega[b](t), t)) = Physics:-Vectors:-Component((D[2](U))(R_(t), Phi_(t)), a)

(634)

Introducing the expression for the inertia tensor,

`𝕀`[a, b] = Sum(-m[i]*(-`r'`[i, c]^2*Physics[KroneckerDelta][a, b]+`r'`[i, a]*`r'`[i, b]), i = 1 .. n)

`𝕀`[a, b] = Sum(-m[i]*(-`r'`[i, c]^2*Physics:-KroneckerDelta[a, b]+`r'`[i, a]*`r'`[i, b]), i = 1 .. n)

(635)

subs((rhs = lhs)(`𝕀`[a, b] = Sum(-m[i]*(-`r'`[i, c]^2*Physics[KroneckerDelta][a, b]+`r'`[i, a]*`r'`[i, b]), i = 1 .. n)), (Sum(-m[i]*(-`r'`[i, c]^2*Physics[KroneckerDelta][a, b]+`r'`[i, a]*`r'`[i, b]), i = 1 .. n))*(diff(Omega[b](t), t)) = Physics[Vectors][Component]((D[2](U))(R_(t), Phi_(t)), a))

`𝕀`[a, b]*(diff(Omega[b](t), t)) = Physics:-Vectors:-Component((D[2](U))(R_(t), Phi_(t)), a)

(636)

From the result (614) for the problem of representing the angular momentum of a rigid body in terms of the inertia tensor  the left-hand side is the derivative of the components of the angular momentum

`L'`[a] = `𝕀`[a, b]*Omega[b]

`L'`[a] = `𝕀`[a, b]*Omega[b]

(637)

subs(Omega[b] = Omega[b](t), (diff(L(x), x))[a] = (diff(L(x), x))[a](t), `L'`[a] = `𝕀`[a, b]*Omega[b])

`L'`[a](t) = `𝕀`[a, b]*Omega[b](t)

(638)

diff(`L'`[a](t) = `𝕀`[a, b]*Omega[b](t), t)

diff(`L'`[a](t), t) = `𝕀`[a, b]*(diff(Omega[b](t), t))

(639)

subs((rhs = lhs)(diff(`L'`[a](t), t) = `𝕀`[a, b]*(diff(Omega[b](t), t))), `𝕀`[a, b]*(diff(Omega[b](t), t)) = Physics[Vectors][Component]((D[2](U))(R_(t), Phi_(t)), a))

diff(`L'`[a](t), t) = Physics:-Vectors:-Component((D[2](U))(R_(t), Phi_(t)), a)

(640)

The right-hand side is the a^th component of the variation of the potential energy U(`#mover(mi("R"),mo("→"))`, `#mover(mi("Φ",fontstyle = "normal"),mo("→"))`) with respect to `#mover(mi("Φ",fontstyle = "normal"),mo("→"))`. A graphics analysis of (D[2](U))(`#mover(mi("R"),mo("→"))`, `#mover(mi("Φ",fontstyle = "normal"),mo("→"))`) and algebraic derivation as done in the problem of the section Angular velocity  results in

(D[2](U))(R_(t), Phi_(t)) = `&x`(R_(t), F_[ext])

(D[2](U))(R_(t), Phi_(t)) = Physics:-Vectors:-`&x`(R_(t), F_[ext])

(641)

subs((D[2](U))(R_(t), Phi_(t)) = Physics[Vectors][`&x`](R_(t), F_[ext]), diff(`L'`[a](t), t) = Physics[Vectors][Component]((D[2](U))(R_(t), Phi_(t)), a))

diff(`L'`[a](t), t) = Physics:-Vectors:-Component(Physics:-Vectors:-`&x`(R_(t), F_[ext]), a)

(642)

As shown in this problem  of the section on the equations of motion for many-particle systems , the right-hand side of this result of that is the total torque `#mover(mi("N"),mo("→"))`(t), resulting in the second equation of movement of a rigid body, for the time derivative of the angular momentum

subs(rhs(diff(`L'`[a](t), t) = Physics[Vectors][Component](Physics[Vectors][`&x`](R_(t), F_[ext]), a)) = N_[a](t), diff(`L'`[a](t), t) = Physics[Vectors][Component](Physics[Vectors][`&x`](R_(t), F_[ext]), a))

diff(`L'`[a](t), t) = N_[a](t)

(643)

NULL

Non-inertial coordinate systems

 

 

When describing the motion of a particle as seen from a non-inertial reference system (e.g. a rotating planet, like the Earth), we also see "acceleration" that is not due to any force but instead to the fact that the reference system itself is accelerated.

 

Problem
Consider a non-inertial reference system K' which moves with non-constant translational velocity `#mover(mi("V"),mo("→"))`(t) with regards to an inertial reference system K__0. 

a) Show that the Lagrangian L' in K' is given by

diff(L(x), x) = (1/2)*m*`#mover(mi("v'"),mo("→"))`-m*`#mover(mi("W"),mo("→"))`.`#mover(mi("r'"),mo("→"))`-U

where `#mover(mi("W"),mo("→"))` = diff(`#mover(mi("V"),mo("→"))`(t), t) is the translational acceleration of the frame K'  as seen from K__0.

 

b) Show that the Lagrange equation derived from this Lagrangian in the frame K'  is

m*(diff(`#mover(mi("v'"),mo("→"))`(t), t)) = -Nabla(U)-m*`#mover(mi("W"),mo("→"))`

Solution

   

 

Coriolis force and centripetal force

 

 

Problem

Consider a second non-inertial frame of reference J whose origin coincides with that of K', but which rotates relative to K' with variable angular velocity `#mover(mi("Ω",fontstyle = "normal"),mo("→"))`(t). Denote the position vector and velocity in the non-inertial frame J as `#mover(mi("r"),mo("→"))` and `#mover(mi("v"),mo("→"))`,

a) Show that the Lagrangian L in the non-inertial frame J is given by

L = (1/2)*m*`#mover(mi("v"),mo("→"))`^2+m*`#mover(mi("v"),mo("→"))`.`&x`(`#mover(mi("Ω",fontstyle = "normal"),mo("→"))`, `#mover(mi("r"),mo("→"))`)+(1/2)*m*`&x`(`#mover(mi("Ω",fontstyle = "normal"),mo("→"))`, `#mover(mi("r"),mo("→"))`)^2-m*(`#mover(mi("W"),mo("→"))`.`#mover(mi("r"),mo("→"))`)-U

b) Show that the Lagrange equation derived from this Lagrangian in the frame J  is

m*(diff(`#mover(mi("v"),mo("→"))`(t), t)) = -Nabla(U)-m*`#mover(mi("W"),mo("→"))`+m*`&x`(`#mover(mi("r"),mo("→"))`, diff(`#mover(mi("Ω",fontstyle = "normal"),mo("→"))`(t), t))+2*m*`&x`(`#mover(mi("v"),mo("→"))`, `#mover(mi("Ω",fontstyle = "normal"),mo("→"))`)+m*`&x`(`#mover(mi("Ω",fontstyle = "normal"),mo("→"))`, `&x`(`#mover(mi("r"),mo("→"))`, `#mover(mi("Ω",fontstyle = "normal"),mo("→"))`))

where 2*m*`&x`(`#mover(mi("v"),mo("→"))`, `#mover(mi("Ω",fontstyle = "normal"),mo("→"))`) is the Coriolis force and m*`&x`(`#mover(mi("Ω",fontstyle = "normal"),mo("→"))`, `&x`(`#mover(mi("r"),mo("→"))`, `#mover(mi("Ω",fontstyle = "normal"),mo("→"))`)) is the centrifugal force.

 

Solution

 

 

restart; with(Physics:-Vectors)

 

a) The starting point is the Lagrangian in the frame K__0. Denoting vectors and the Lagrangian in K__0 with the suffix 0, L__0 is given by

CompactDisplay((r__0_, v__0_, (diff(r(x), x))*_, (diff(v(x), x))*_, r_, v_, V_, W_, Omega_)(t))

Omega_(t)*`will now be displayed as`*Omega_

(658)

L__0 = (1/2)*m*v__0_(t)^2-U(r__0_(t))

L__0 = (1/2)*m*Physics:-`^`(v__0_(t), 2)-U(r__0_(t))

(659)

The velocities of the particle in the frames K__0 and K' are related by

v__0_(t) = (diff(v(x), x))*_(t)+V_(t)

v__0_(t) = `v'_`(t)+V_(t)

(660)

where `#mover(mi("V"),mo("→"))` is the translational velocity of K' viewed from K__0. Inserting this relation (660) into L__0 gives diff(L(x), x), the result of the previous problem

`L'` = (1/2)*m*Physics[Vectors][Norm](`v'_`(t))^2-m*Physics[Vectors][`.`](`r'_`(t), W_(t))-U(`r'_`(t))

`L'` = (1/2)*m*Physics:-Vectors:-Norm(`v'_`(t))^2-m*Physics:-Vectors:-`.`(`r'_`(t), W_(t))-U(`r'_`(t))

(661)

In turn, the velocities `#mover(mi("v'"),mo("→"))` and `#mover(mi("v"),mo("→"))` in the frames K' and J are related by

(diff(v(x), x))*_(t) = v_(t)+`&x`(Omega_(t), r_(t))

`v'_`(t) = v_(t)+Physics:-Vectors:-`&x`(Omega_(t), r_(t))

(662)

where `#mover(mi("Ω",fontstyle = "normal"),mo("→"))` is the angular velocity of the frame J viewed from K'. Also, since K' and J have the same origin, `#mover(mi("r'"),mo("→"))` = `#mover(mi("r"),mo("→"))` and the Lagrangian L in J  is

subs(`v'_`(t) = v_(t)+Physics[Vectors][`&x`](Omega_(t), r_(t)), (diff(r(x), x))*_ = r_, diff(L(x), x) = L(x), `L'` = (1/2)*m*Physics[Vectors][Norm](`v'_`(t))^2-m*Physics[Vectors][`.`](`r'_`(t), W_(t))-U(`r'_`(t)))

L = (1/2)*m*Physics:-Vectors:-Norm(v_(t)+Physics:-Vectors:-`&x`(Omega_(t), r_(t)))^2-m*Physics:-Vectors:-`.`(r_(t), W_(t))-U(r_(t))

(663)

expand(L = (1/2)*m*Physics[Vectors][Norm](v_(t)+Physics[Vectors][`&x`](Omega_(t), r_(t)))^2-m*Physics[Vectors][`.`](r_(t), W_(t))-U(r_(t)))

L = (1/2)*m*Physics:-Vectors:-Norm(v_(t))^2+m*Physics:-Vectors:-`.`(v_(t), Physics:-Vectors:-`&x`(Omega_(t), r_(t)))+(1/2)*m*Physics:-Vectors:-Norm(Omega_(t))^2*Physics:-Vectors:-Norm(r_(t))^2-(1/2)*m*Physics:-Vectors:-`.`(Omega_(t), r_(t))^2-m*Physics:-Vectors:-`.`(r_(t), W_(t))-U(r_(t))

(664)

Two of these terms can be regrouped as `&x`(`#mover(mi("Ω",fontstyle = "normal"),mo("→"))`, `#mover(mi("r"),mo("→"))`)^2

`&x`(Omega_(t), r_(t))^2

Physics:-`^`(Physics:-Vectors:-`&x`(Omega_(t), r_(t)), 2)

(665)

expand(Physics[`^`](Physics[Vectors][`&x`](Omega_(t), r_(t)), 2)) = Physics[`^`](Physics[Vectors][`&x`](Omega_(t), r_(t)), 2)

Physics:-Vectors:-Norm(Omega_(t))^2*Physics:-Vectors:-Norm(r_(t))^2-Physics:-Vectors:-`.`(Omega_(t), r_(t))^2 = Physics:-`^`(Physics:-Vectors:-`&x`(Omega_(t), r_(t)), 2)

(666)

simplify(L = (1/2)*m*Physics[Vectors][Norm](v_(t))^2+m*Physics[Vectors][`.`](v_(t), Physics[Vectors][`&x`](Omega_(t), r_(t)))+(1/2)*m*Physics[Vectors][Norm](Omega_(t))^2*Physics[Vectors][Norm](r_(t))^2-(1/2)*m*Physics[Vectors][`.`](Omega_(t), r_(t))^2-m*Physics[Vectors][`.`](r_(t), W_(t))-U(r_(t)), {Physics[Vectors][Norm](Omega_(t))^2*Physics[Vectors][Norm](r_(t))^2-Physics[Vectors][`.`](Omega_(t), r_(t))^2 = Physics[`^`](Physics[Vectors][`&x`](Omega_(t), r_(t)), 2)})

L = (1/2)*m*Physics:-Vectors:-Norm(v_(t))^2+m*Physics:-Vectors:-`.`(v_(t), Physics:-Vectors:-`&x`(Omega_(t), r_(t)))-m*Physics:-Vectors:-`.`(r_(t), W_(t))+(1/2)*m*Physics:-`^`(Physics:-Vectors:-`&x`(Omega_(t), r_(t)), 2)-U(r_(t))

(667)

which is the expected result.

_______________________________________

 

b) To compute Lagrange's equation in vectorial form, one can use the standard formula as in the previous problem

 

%diff(%diff(L = (1/2)*m*Physics[Vectors][Norm](v_(t))^2+m*Physics[Vectors][`.`](v_(t), Physics[Vectors][`&x`](Omega_(t), r_(t)))-m*Physics[Vectors][`.`](r_(t), W_(t))+(1/2)*m*Physics[`^`](Physics[Vectors][`&x`](Omega_(t), r_(t)), 2)-U(r_(t)), v_(t)), t)-%diff(L = (1/2)*m*Physics[Vectors][Norm](v_(t))^2+m*Physics[Vectors][`.`](v_(t), Physics[Vectors][`&x`](Omega_(t), r_(t)))-m*Physics[Vectors][`.`](r_(t), W_(t))+(1/2)*m*Physics[`^`](Physics[Vectors][`&x`](Omega_(t), r_(t)), 2)-U(r_(t)), r_(t))

%diff(%diff(L = (1/2)*m*Physics:-Vectors:-Norm(v_(t))^2+m*Physics:-Vectors:-`.`(v_(t), Physics:-Vectors:-`&x`(Omega_(t), r_(t)))-m*Physics:-Vectors:-`.`(r_(t), W_(t))+(1/2)*m*Physics:-`^`(Physics:-Vectors:-`&x`(Omega_(t), r_(t)), 2)-U(r_(t)), v_(t)), t)-%diff(L = (1/2)*m*Physics:-Vectors:-Norm(v_(t))^2+m*Physics:-Vectors:-`.`(v_(t), Physics:-Vectors:-`&x`(Omega_(t), r_(t)))-m*Physics:-Vectors:-`.`(r_(t), W_(t))+(1/2)*m*Physics:-`^`(Physics:-Vectors:-`&x`(Omega_(t), r_(t)), 2)-U(r_(t)), r_(t))

(668)

 

Evaluate these derivatives and replace diff(`#mover(mi("r"),mo("→"))`(t), t) = `#mover(mi("v"),mo("→"))` 

subs(diff(r_(t), t) = v_(t), value(%diff(%diff(L = (1/2)*m*Physics[Vectors][Norm](v_(t))^2+m*Physics[Vectors][`.`](v_(t), Physics[Vectors][`&x`](Omega_(t), r_(t)))-m*Physics[Vectors][`.`](r_(t), W_(t))+(1/2)*m*Physics[`^`](Physics[Vectors][`&x`](Omega_(t), r_(t)), 2)-U(r_(t)), v_(t)), t)-%diff(L = (1/2)*m*Physics[Vectors][Norm](v_(t))^2+m*Physics[Vectors][`.`](v_(t), Physics[Vectors][`&x`](Omega_(t), r_(t)))-m*Physics[Vectors][`.`](r_(t), W_(t))+(1/2)*m*Physics[`^`](Physics[Vectors][`&x`](Omega_(t), r_(t)), 2)-U(r_(t)), r_(t))))

0 = m*(diff(v_(t), t))+m*(Physics:-Vectors:-`&x`(diff(Omega_(t), t), r_(t))+Physics:-Vectors:-`&x`(Omega_(t), v_(t)))-m*Physics:-Vectors:-`&x`(v_(t), Omega_(t))+m*W_(t)-m*Physics:-Vectors:-`&x`(Physics:-Vectors:-`&x`(Omega_(t), r_(t)), Omega_(t))+(D(U))(r_(t))

(669)

Isolating m*(diff(`#mover(mi("v"),mo("→"))`(t), t))  and collecting vector products we get the expected result

isolate(0 = m*(diff(v_(t), t))+m*(Physics[Vectors][`&x`](diff(Omega_(t), t), r_(t))+Physics[Vectors][`&x`](Omega_(t), v_(t)))-m*Physics[Vectors][`&x`](v_(t), Omega_(t))+m*W_(t)-m*Physics[Vectors][`&x`](Physics[Vectors][`&x`](Omega_(t), r_(t)), Omega_(t))+(D(U))(r_(t)), m*(diff(v_(t), t)))

m*(diff(v_(t), t)) = -m*(Physics:-Vectors:-`&x`(diff(Omega_(t), t), r_(t))-Physics:-Vectors:-`&x`(v_(t), Omega_(t)))+m*Physics:-Vectors:-`&x`(v_(t), Omega_(t))-m*W_(t)+m*Physics:-Vectors:-`&x`(Physics:-Vectors:-`&x`(Omega_(t), r_(t)), Omega_(t))-(D(U))(r_(t))

(670)

collect(m*(diff(v_(t), t)) = -m*(Physics[Vectors][`&x`](diff(Omega_(t), t), r_(t))-Physics[Vectors][`&x`](v_(t), Omega_(t)))+m*Physics[Vectors][`&x`](v_(t), Omega_(t))-m*W_(t)+m*Physics[Vectors][`&x`](Physics[Vectors][`&x`](Omega_(t), r_(t)), Omega_(t))-(D(U))(r_(t)), `&x`)

m*(diff(v_(t), t)) = m*Physics:-Vectors:-`&x`(Physics:-Vectors:-`&x`(Omega_(t), r_(t)), Omega_(t))-m*Physics:-Vectors:-`&x`(diff(Omega_(t), t), r_(t))+2*m*Physics:-Vectors:-`&x`(v_(t), Omega_(t))-m*W_(t)-(D(U))(r_(t))

(671)

NULL

 

 

Part II (forthcoming)

The Hamiltonian and equations of motion; Poisson brackets

   

Canonical transformations

   

The Hamilton-Jacobi equation

   

References

   

 

NULL

Download Mechanics_(2022).mw

Download Mechanics_(2022).pdf
 

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft


Please Wait...