ecterrab

The 2025 Maplesoft Physics Updates and i...

Posted: ecterrab 14540 Product: Maple 2025

April 01 2025

8 8

The Maplesoft Physics Updates, introduced over a decade ago, brought with them an innovative concept: to deliver fixes and new developments continuously, as soon as they enter the development version of the Maple library for the next release. A key aspect of this initiative was prioritizing the resolution of issues reported on MaplePrimes, ensuring that fixes became available to everyone within 24 to 48 hours. Initially focused solely on the Physics package, the scope of the updates quickly expanded to include other parts of the Maple library and the Typesetting system.

This initiative, which I developed outside regular work hours, aimed to enhance the Maple experience—where issues encountered in daily use could be resolved almost immediately, minimizing disruptions and benefiting the entire user community through shared updates.

As of January 1st, I have stepped away from my role at Maplesoft and have been increasingly involved in activities unrelated to Maple. This raises the question of what will happen with the Physics Updates for Maple 2025 and after.

The Physics project remains a unique and personally meaningful endeavor for me. So, for now, I will continue to dedicate some time to these Updates—but only for the Physics package, not for other parts of the library. As before, these fixes and developments will be included in the Physics Updates only after they have been integrated into the development version of Maple’s official library for the next release. In that sense, they will continue to be Maplesoft updates.

On that note, the first release of the Physics Updates for Maple 2025—focused solely on the Physics package—went out today as version 1854. To install it, the first time open Maple 2025 and use the Maplecloud toolbar -> Packages, or else input PackageTools:-Install(5137472255164416). Any next time, just enter Physics:-Version(latest)

As for fixes beyond the Physics package, I understand that Maplesoft is exploring the possibility of offering something similar to what was previously delivered through the Maplesoft Physics Updates.

All the best

PS: to install the last version of the Maplesoft Physics Updates for Maple 2024, open Maple and input Physics:-Version(1852), not 1853.

Edgardo S. Cheb-Terrab
Physics, Differential Equations, and Mathematical Functions
Maplesoft Emmeritus
Research and Education—passionate about all that.

What is new in Physics in Maple 2025...

Posted: ecterrab 14540 Product: Maple

March 27 2025

8 1

What is new in Physics in Maple 2025

This post is, basically, the page of what is new in Physics distributed with Maple 2025. Why post it? Because this time we achieved a result that is a breakthrough/milestone in Computer Algebra: for the first time we can systematically compute Einstein's equations from first principles, using a computer, starting from a Lagrangian for gravity. This is something to celebrate.

Although being able to perform this computation on a computer algebra sheet is relevant mostly for physicists working in general relativity, the developments performed in the tensor manipulation and functional differentiation routines of the Maple system to cover this computation were tremendous, in number of lines of code, efforts and time, resulting in improvements not just fore general relativity but for physics all around, and in an area out of reach of neural network AIs . Other results, like the new routines for linearizing gravity, requested other times here in Mapleprimes, are also part of the relevant novelties.

Lagrange Equations and simplification of tensorial expressions in curved spacetimes

Linearized Gravity

Relative Tensors

New Physics:-Library commands

See Also

Lagrange Equations and simplification of tensorial expressions in curved spacetimes

LagrangeEquations is a Physics command introduced in 2023 taking advantage of the functional differentiation capabilities of the Physics package . This command can handle tensors and vectors of the Physics package as well as derivatives using vectorial differential operators (see d_ and Nabla ), works by performing functional differentiation (see Fundiff ), and handles 1^st, and higher order derivatives of the coordinates in the Lagrangian automatically. LagrangeEquations receives an expression representing a Lagrangian and returns a sequence of Lagrange equations with as many equations as coordinates are indicated. The number of parameters can also be many. For example, in electrodynamics, the "coordinate" is a tensor field , there are then four coordinates, one for each of the values of the index , and there are four parameters .

New in Maple 2025, the "coordinates" can now also be the components of the metric tensor in a curved spacetime, in which case the equations returned are Einstein's equations. Also new, instead of a coordinate or set of them, you can pass the keyword EnergyMomentum, in which case the output is the conserved energy-momentum tensor of the physical model represented by the given Lagrangian L.

Examples

(1)

The model in classical field theory and corresponding field equations, as in previous releases

>

(2)

>

(1/2)*Physics:-d_[mu](Phi(X), [X])*Physics:-d_[`~mu`](Phi(X), [X])-(1/2)*m^2*Phi(X)^2+(1/4)*lambda*Phi(X)^4

(3)

Lagrange's equations

>

(4)

New: The energy-momentum tensor can be computed as the Lagrange equations taking the metric as the coordinate, not equating to 0 the result, but multiplying the variation of the action "(delta S)/(delta g[]^(mu,nu))" by (in flat spacetimes ). For that purpose, you can use the EnergyMomentum keyword. You can optionally indicate the indices to be used in the output as well as their covariant or contravariant character

>

Physics:-EnergyMomentum[mu, nu] = (-(1/4)*lambda*Phi(X)^4+(1/2)*m^2*Phi(X)^2-(1/2)*Physics:-d_[`~beta`](Phi(X), [X])*Physics:-d_[beta](Phi(X), [X]))*Physics:-g_[mu, nu]+Physics:-d_[mu](Phi(X), [X])*Physics:-d_[nu](Phi(X), [X])

(5)

To further compute using the above as the definition for , you can use the Define command

>

${Physics:-Dgamma[mu], Physics:-Psigma[mu], Physics:-d_[mu], Physics:-g_[mu, nu], Physics:-EnergyMomentum[mu, nu], Physics:-LeviCivita[alpha, beta, mu, nu], Physics:-SpaceTimeVector[mu](X)}$

(6)

After which the system knows about the symmetry properties and the components of

>

Physics:-EnergyMomentum[mu, nu] = (-(1/4)*lambda*Phi(X)^4+(1/2)*m^2*Phi(X)^2-(1/2)*Physics:-d_[`~beta`](Phi(X), [X])*Physics:-d_[beta](Phi(X), [X]))*Physics:-g_[mu, nu]+Physics:-d_[mu](Phi(X), [X])*Physics:-d_[nu](Phi(X), [X])

(7)

>

(8)

>

Physics:-EnergyMomentum[mu, nu] = Matrix(%id = 36893488151952536988)

(9)

New: LagrangeEquations takes advantage of the extension of Fundiff to compute functional derivatives in curved spacetimes introduced for Maple 2025, and so it also handles the case of a scalar field in a curved spacetime. Set for instance an arbitrary metric

>

Physics:-g_[mu, nu] = Matrix(%id = 36893488152257720060)

(10)

For the action to be a true scalar in spacetime, the Lagrangian density now needs to be multiplied by the square root of the determinant of the metric

>

(-%g_[determinant])^(1/2)*((1/2)*Physics:-d_[mu](Phi(X), [X])*Physics:-d_[`~mu`](Phi(X), [X])-(1/2)*m^2*Phi(X)^2+(1/4)*lambda*Phi(X)^4)

(11)

New: With the extension of the tensorial simplification algorithms for curved spacetimes, the Lagrange equations can be computed arriving directly to the compact form

>

(12)

Comparing with the result (4) for the same Lagrangian in a flat spacetime, we see the only difference is that the dAlembertian is now expressed in terms of covariant derivatives D_ .

The EnergyMomentum tensor is computed in the same way as when the spacetime is flat

>

Physics:-EnergyMomentum[mu, nu] = (-(1/4)*lambda*Phi(X)^4+(1/2)*m^2*Phi(X)^2-(1/2)*Physics:-d_[`~beta`](Phi(X), [X])*Physics:-d_[beta](Phi(X), [X]))*Physics:-g_[mu, nu]+Physics:-d_[mu](Phi(X), [X])*Physics:-d_[nu](Phi(X), [X])

(13)

General Relativity

New: the most significant development in LagrangeEquations is regarding General Relativity. It can now compute Einstein's equations directly from the Lagrangian, not using tabulated cases, and properly handling several (traditional or not) alternative ways of presenting the Lagrangian.

Einstein's equations concern the case of a curved spacetime with metric as, for instance, the general case of an arbitrary metric set lines above. In the Lagrangian formulation, the coordinates of the problem are the components of the metric , and as in the case of electrodynamics the parameters are the spacetime coordinates . The simplest case is that of Einstein's equation in vacuum, for which the Lagrangian density is expressed in terms of the trace of the Ricci tensor by

>

(14)

Einstein's equations in vacuum:

>

-(1/2)*Physics:-g_[mu, nu]*Physics:-Ricci[alpha, `~alpha`]+Physics:-Ricci[mu, nu] = 0

(15)

where in the above instead of passing as second argument, we passed to get the equations using those free indices. The tensorial equation computed is also the definition of the Einstein tensor

>

Physics:-Einstein[mu, nu] = -(1/2)*Physics:-g_[mu, nu]*Physics:-Ricci[alpha, `~alpha`]+Physics:-Ricci[mu, nu]

(16)

The Lagrangian used to compute Einstein's equations (15) contains first and second derivatives of the metric. To see that, rewrite L in terms of Christoffel symbols

>

(17)

Recalling the definition

>

Physics:-Christoffel[alpha, mu, nu] = (1/2)*Physics:-d_[nu](Physics:-g_[alpha, mu], [X])+(1/2)*Physics:-d_[mu](Physics:-g_[alpha, nu], [X])-(1/2)*Physics:-d_[alpha](Physics:-g_[mu, nu], [X])

(18)

in the two terms containing derivatives of Christoffel symbols contain second order derivatives of . Now, it is always possible to add a total spacetime derivative to without changing Einstein's equations (assuming the variation of the metric in the corresponding boundary integrals vanishes), and in that way, in this particular case of , obtain a Lagrangian involving only 1st order derivatives. The total derivative, expressed using the inert command to see it before the differentiation operation is performed, is

>

(19)

Adding this term to , performing the differentiation operation and simplifying we get

>

(20)

>

(21)

>

(22)

which is a Lagrangian depending only on 1st order derivatives of the metric through Christoffel symbols. As expected, the equations of motion resulting from this Lagrangian are the same Einstein equations computed in (15)

>

-(1/2)*Physics:-Ricci[iota, `~iota`]*Physics:-g_[mu, nu]+Physics:-Ricci[mu, nu] = 0

(23)

To illustrate the new Maple 2025 tensorial simplification capabilities note that is no just ≡ (17) after discarding its two terms involving derivatives of Christoffel symbols. To verify this, split into the terms containing or not derivatives of Christoffel

>

(24)

Comparing, the total derivative TD≡ (19) is not just , but

>

(25)

Things like these, , can now be verified directly with the new tensorial simplification capabilities: take the left-hand side minus the right-hand side, evaluate the inert derivative and simplify to see the equality is true

>

(26)

>

(27)

>

(28)

That said, it is also true that results in the Lagrangian , and since the equations of movement don't depend on the sign of the Lagrangian, for this Lagrangian adding the term happens to be equivalent to just discarding the terms of involving derivatives of Christoffel symbols.

Also new in Maple 2025, due to the extension of Fundiff to compute in curved spacetimes, it is now also possible to compute Einstein's equations from first principles by constructing the action,

>

Int(Int(Int(Int((-%g_[determinant])^(1/2)*Physics:-Ricci[alpha, `~alpha`], x = -infinity .. infinity), y = -infinity .. infinity), z = -infinity .. infinity), t = -infinity .. infinity)

(29)

and equating to zero the functional derivative with respect to the metric. To avoid displaying the resulting large expression, end the input line with ":"

>

Simplifying this result, we get an expression in terms of Christoffel symbols and its derivatives

>

(30)

In this result, we see derivatives of Christoffel symbols, expressed using the D_ command for covariant differentiation. Although, such objects have not the geometrical meaning of a covariant derivative, computationally, they here represent what would be a covariant derivative if the Christoffel symbols were a tensor. For example,

>

(31)

With this computational meaning for the derivatives of Christoffel symbols appearing in (30), rewrite EEC≡(30) in terms of the Ricci and Riemann tensors. For that, consider the definition

>

(32)

Rewrite the noncovariant derivatives in terms of derivatives using the computational representation (31), simplify and isolate one of them

>

(33)

>

(34)

>

(35)

Analogously, derive an expression to rewrite derivatives of Christoffel symbols using the Riemann tensor

>

(36)

>

(37)

>

(38)

>

(39)

Substitute these two equations, in sequence, into Einstein's equations EEC≡(30)

>

(40)

Simplify to arrive at the traditional compact form of Einstein's equations

>

(1/2)*Physics:-Ricci[chi, `~chi`]*Physics:-g_[`~alpha`, `~beta`]-Physics:-Ricci[`~alpha`, `~beta`] = 0

(41)

Linearized Gravity

Generally speaking, linearizing gravity is about discarding in Einstein's field equations the terms that are quadratic in the metric and its derivatives, an approximation valid when the gravitational field is weak (the deviation from a flat Minkowski spacetime is small). Linearizing gravity is used, e.g. in the study of gravitational waves. In the context of Maple's Physics , the formulation of linearized gravity can be done using the general relativity tensors that come predefined in Physics plus a new in Maple 2025 Physics:-Library:-Linearize command.

In what follows it is shown how to linearize the Ricci tensor and through it Einstein 's equations. To compare results, see for instance the Wikipedia page for Linearized gravity. Start setting coordinates, you could use Cartesian, spherical, cylindrical, or define your own.

>

restart; with(Physics); Setup(coordinates = cartesian)

(42)

The default metric when Physics is loaded is the Minkowski metric, representing a flat (no curvature) spacetime

>

Physics:-g_[mu, nu] = Matrix(%id = 36893488151987392132)

(43)

The weakly perturbed metric

Suppose you want to define a small perturbation around this metric. For that purpose, define a perturbation tensor , that in the general case depends on the coordinates and is not diagonal, the only requirement is that it is symmetric (to have it diagonal, change symmetric by diagonal; to have it constant, change by )

>

h[mu, nu] = Matrix(%id = 36893488152568729716)

(44)

In the above it is understood that , so that quadratic or higher powers of it or its derivatives can be approximated to 0 (discarded). Define the components of accordingly

>

${Physics:-Dgamma[mu], Physics:-Psigma[mu], Physics:-d_[mu], Physics:-g_[mu, nu], h[mu, nu], Physics:-LeviCivita[alpha, beta, mu, nu], Physics:-SpaceTimeVector[mu](X)}$

(45)

Define also a tensor representing the unperturbed Minkowski metric

>

eta[mu, nu] = Matrix(%id = 36893488151987392132)

(46)

>

${Physics:-Dgamma[mu], Physics:-Psigma[mu], Physics:-d_[mu], eta[mu, nu], Physics:-g_[mu, nu], h[mu, nu], Physics:-LeviCivita[alpha, beta, mu, nu], Physics:-SpaceTimeVector[mu](X)}$

(47)

The weakly perturbed metric is given by

>

(48)

Make this be the definition of the metric

>

Physics:-g_[mu, nu] = Matrix(%id = 36893488152066920564)

${Physics:-D_[mu], Physics:-Dgamma[mu], Physics:-Psigma[mu], Physics:-Ricci[mu, nu], Physics:-Riemann[mu, nu, alpha, beta], Physics:-Weyl[mu, nu, alpha, beta], Physics:-d_[mu], eta[mu, nu], Physics:-g_[mu, nu], Physics:-gamma_[i, j], h[mu, nu], Physics:-Christoffel[mu, nu, alpha], Physics:-Einstein[mu, nu], Physics:-LeviCivita[alpha, beta, mu, nu], Physics:-SpaceTimeVector[mu](X)}$

(49)

Linearizing the Ricci tensor

The linearized form of the Ricci tensor is computed by introducing this weakly perturbed metric (48) in the expression of the Ricci tensor as a function of the metric. This can be accomplished in different ways, the simpler being to use the conversion network between tensors , but for illustration purposes, showing steps one at time, a substitution of definitions one into the other one is used

>

(50)

>

Physics:-Christoffel[`~alpha`, mu, nu] = (1/2)*Physics:-g_[`~alpha`, `~beta`]*(Physics:-d_[nu](Physics:-g_[beta, mu], [X])+Physics:-d_[mu](Physics:-g_[beta, nu], [X])-Physics:-d_[beta](Physics:-g_[mu, nu], [X]))

(51)

>

(52)

Introducing (48)≡, and also the inert form of the Ricci tensor to facilitate simplification some steps below,

>

(53)

This expression contains several terms quadratic in the small perturbation and its derivatives. The new in Maple 2025 routine to filter out those terms is Physics:-Library:-Linearize , which requires specifying the symbol representing the small quantities

>

(54)

Important: in this result, is the flat Minkowski metric, not the perturbed metric . However, in the context of a linearized formulation, raises and lowers tensor indices the same way as . Hence, to further simplify contracted products of in (54) , it is practical to reintroduce representing that Minkowski metric and simplify using the internal algorithms for a flat metric

>

Physics:-g_[mu, nu] = Matrix(%id = 36893488152066895868)

(55)

To proceed simplifying, replace in the expression (54) for the Ricci tensor the intermediate Minkowski by

>

(56)

Simplifying, results in the linearized form of the Ricci tensor shown in the Wikipedia page for Linearized gravity.

>

%Ricci[mu, nu] = -(1/2)*Physics:-d_[mu](Physics:-d_[nu](h[tau, `~tau`], [X]), [X])-(1/2)*Physics:-dAlembertian(h[mu, nu], [X])+(1/2)*Physics:-d_[nu](Physics:-d_[tau](h[mu, `~tau`], [X]), [X])+(1/2)*Physics:-d_[mu](Physics:-d_[tau](h[nu, `~tau`], [X]), [X])

(57)

Linearizing Einstein's equations

Einstein's equations are the components of Einstein's tensor , whose definition in terms of the Ricci tensor is

>

Physics:-Einstein[mu, nu] = Physics:-Ricci[mu, nu]-(1/2)*Physics:-g_[mu, nu]*Physics:-Ricci[alpha, `~alpha`]

(58)

Compute the trace directly from the linearized form (57) of the Ricci tensor,

>

%Ricci[mu, nu]*Physics:-g_[`~mu`, `~nu`] = (-(1/2)*Physics:-d_[mu](Physics:-d_[nu](h[tau, `~tau`], [X]), [X])-(1/2)*Physics:-dAlembertian(h[mu, nu], [X])+(1/2)*Physics:-d_[nu](Physics:-d_[tau](h[mu, `~tau`], [X]), [X])+(1/2)*Physics:-d_[mu](Physics:-d_[tau](h[nu, `~tau`], [X]), [X]))*Physics:-g_[`~mu`, `~nu`]

(59)

>

(60)

The linearized Einstein equations are constructed reproducing the definition (58) using (57) and (60)

>

(61)

which is the same formula shown in the Wikipedia page for Linearized gravity.

You can now redefine the general introduced in (44) in different ways (see discussion in the Wikipedia page), or, depending on the case, just substitute your preferred gauge in this formula (61) for the general case. For example, the condition for the Harmonic gauge also known as Lorentz gauge reduces the linearized field equations to their simplest form

>

Physics:-d_[mu](h[`~mu`, nu], [X]) = (1/2)*Physics:-d_[nu](h[alpha, `~alpha`], [X])

(62)

>

(63)

>

%Ricci[mu, nu]-(1/2)*Physics:-g_[mu, nu]*%Ricci[alpha, `~alpha`] = -(1/2)*Physics:-dAlembertian(h[mu, nu], [X])+(1/4)*Physics:-dAlembertian(h[alpha, `~alpha`], [X])*Physics:-g_[mu, nu]

(64)

Relative Tensors

In General Relativity, the context of a curved spacetime, it is sometimes necessary to work with relative tensors, for which the transformation rule under a transformation of coordinates involves powers of the determinant of the transformation - see Chapter 4 of "Lovelock, D., and Rund, H. Tensors, Differential Forms and Variational Principles, Dover, 1989." Physics in Maple 2025 includes a complete, new implementation of relative tensors.

To indicate that a tensor being defined is relative pass its relative weight. For example, set a curved spacetime,

>

restart; with(Physics); g_[sc]

$`Default differentiation variables for d_, D_ and dAlembertian are:`*{X = (r, theta, phi, t)}$

Physics:-g_[mu, nu] = Matrix(%id = 36893488152573138332)

(65)

Define now two tensors of one index, one of them being relative

>

${Physics:-D_[mu], Physics:-Dgamma[mu], Physics:-Psigma[mu], Physics:-Ricci[mu, nu], Physics:-Riemann[mu, nu, alpha, beta], T[mu], Physics:-Weyl[mu, nu, alpha, beta], Physics:-d_[mu], Physics:-g_[mu, nu], Physics:-gamma_[i, j], Physics:-Christoffel[mu, nu, alpha], Physics:-Einstein[mu, nu], Physics:-LeviCivita[alpha, beta, mu, nu], Physics:-SpaceTimeVector[mu](X)}$

(66)

>

${Physics:-D_[mu], Physics:-Dgamma[mu], Physics:-Psigma[mu], R[mu], Physics:-Ricci[mu, nu], Physics:-Riemann[mu, nu, alpha, beta], T[mu], Physics:-Weyl[mu, nu, alpha, beta], Physics:-d_[mu], Physics:-g_[mu, nu], Physics:-gamma_[i, j], Physics:-Christoffel[mu, nu, alpha], Physics:-Einstein[mu, nu], Physics:-LeviCivita[alpha, beta, mu, nu], Physics:-SpaceTimeVector[mu](X)}$

(67)

Transformation of Coordinates

Consider a transformation of coordinates, from spherical to where

>

TR := r = (1+m/(2*rho))^2*rho

r = (1+(1/2)*m/rho)^2*rho

(68)

The transformed components of and are, respectively,

>

Vector[column](%id = 36893488151820263660)

(69)

>

Vector[column](%id = 36893488151823155676)

(70)

where, when comparing both results, we see that the transformed components for are all multiplied by with and is the determinant of the transformation:

>

Matrix(%id = 36893488151964220700)

(71)

>

J = (1/4)*(-m^2+4*rho^2)/rho^2

(72)

Relative weight

The relative weight of a scalar, tensor or tensorial expression can be computed using the Physics:-Library:-GetRelativeWeight command. For the two tensors and used above,

>

(73)

>

(74)

The relative weight of a tensor does not depend on the covariant or contravariant character of its indices

>

(75)

The LeviCivita tensor is a special case, has its relative weight defined when Physics is loaded, and because in a curved spacetime it is not a tensor its relative weight depends on the covariant or contravariant character of its indices

>

(76)

>

(77)

The relative weight w of a product is equal to the sum of relative weights of each factor

>

(78)

>

(79)

The relative weight w of a power is equal to the relative weight of the base multiplied by the power

>

1/(R[mu]*R[`~mu`])

(80)

(81)

The relative weight w of a sum is equal to the relative weight of one of its terms and exists if all the terms have the same w.

>

(82)

>

(83)

The relative weight of any determinant is always equal to 2

>

(84)

>

(85)

Relative Term in covariant derivatives

When computing the covariant derivative of a relative scalar, tensor or tensorial expression that has non-zero relative weight w, a relative term is added, that can be computed using the Physics:-Library:-GetRelativeWeight command.

>

(86)

>

(87)

Consequently,

>

(88)

To understand this zero value on the right-hand side, express the left-hand side in terms of d_

>

(89)

evaluate the inert %d_

>

(90)

The factor in parentheses is equal to , where the covariant derivative of the metric is equal to zero, so

>

(91)

Consider the covariant derivative of and defined in (66) and (67)

>

(92)

>

(93)

The corresponding covariant derivatives

>

(94)

>

%D_[mu](T[`~mu`](X)) = 2*T[`~mu`](X)*Physics:-d_[mu](r, [X])/r+Physics:-d_[mu](theta, [X])*cos(theta)*T[`~mu`](X)/sin(theta)+Physics:-d_[mu](T[`~mu`](X), [X])

(95)

>

(96)

>

%D_[mu](R[`~mu`](X)) = 2*R[`~mu`](X)*Physics:-d_[mu](r, [X])/r+Physics:-d_[mu](theta, [X])*cos(theta)*R[`~mu`](X)/sin(theta)+Physics:-d_[mu](R[`~mu`](X), [X])-Physics:-Christoffel[`~nu`, mu, nu]*R[`~mu`](X)

(97)

where in the above we see the additional (relative) term

>

(98)

New Physics:-Library commands

ConvertToF , Linearize , GetRelativeTerm , GetRelativeWeight.

Examples

•

ConvertToF receives an algebraic expression involving tensors and/or tensor functions and rewrites them in terms of the tensor of name F when that is possible. This routine is similar, however more general than the standard convert which only handles the existing conversion network for the tensors of General Relativity in that ConvertToF also uses any tensor definition you introduce using Define , expressing a tensor in terms of others.

Load any curved spacetime metric automatically setting the coordinates

restart; with(Physics); g_[sc]

Physics:-g_[mu, nu] = Matrix(%id = 36893488152573203868)

(99)

For example, rewrite the Christoffel symbols in terms of the metric g_ ; this works as in previous releases

>

Physics:-Christoffel[mu, alpha, beta] = (1/2)*Physics:-d_[beta](Physics:-g_[alpha, mu], [X])+(1/2)*Physics:-d_[alpha](Physics:-g_[beta, mu], [X])-(1/2)*Physics:-d_[mu](Physics:-g_[alpha, beta], [X])

(100)

Define a representing the 4D electromagnetic potential as a function of the coordinates and representing the electromagnetic field tensors

>

${A[mu], Physics:-D_[mu], Physics:-Dgamma[mu], Physics:-Psigma[mu], Physics:-Ricci[mu, nu], Physics:-Riemann[mu, nu, alpha, beta], Physics:-Weyl[mu, nu, alpha, beta], Physics:-d_[mu], Physics:-g_[mu, nu], Physics:-gamma_[i, j], Physics:-Christoffel[mu, nu, alpha], Physics:-Einstein[mu, nu], Physics:-LeviCivita[alpha, beta, mu, nu], Physics:-SpaceTimeVector[mu](X)}$

(101)

>

${A[mu], Physics:-D_[mu], Physics:-Dgamma[mu], F[mu, nu], Physics:-Psigma[mu], Physics:-Ricci[mu, nu], Physics:-Riemann[mu, nu, alpha, beta], Physics:-Weyl[mu, nu, alpha, beta], Physics:-d_[mu], Physics:-g_[mu, nu], Physics:-gamma_[i, j], Physics:-Christoffel[mu, nu, alpha], Physics:-Einstein[mu, nu], Physics:-LeviCivita[alpha, beta, mu, nu], Physics:-SpaceTimeVector[mu](X)}$

(102)

Rewrite the following expression in terms of the electromagnetic potential

>

(103)

In the example above, the output is similar to this other one

>

(104)

The rewriting, however, works also with tensorial expressions

>

(105)

>

(106)

•	Linearize receives a tensorial expression T and an indication of the small quantities h in T , and discards terms quadratic or of higher order in h. For an example of this new routine in action, see the section Linearized Gravity above.

•	GetRelativeTerm and GetRelativeWeight are illustrated in the section Relative Tensors above.

	See Also
	Index of New Maple 2025 Features , Physics , Computer Algebra for Theoretical Physics, The Physics project, The Physics Updates

Download New_in_Physics_2025.mw

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

Computing Einstein's equations using var...

Posted: ecterrab 14540 Product: Maple

March 27 2025

8 1

Computing Einstein's equations using variational principles

Edgardo S. Cheb-Terrab

Freddy Baudine

One of the biggest challenges for a computer algebra system is to compute Einstein's equations from first principles, by equating to zero the functional derivative of the Action for gravity, without using human-shortcuts or tricks. Developments during 2024, appearing as new in Maple 2025, broke through that barrier and now the LagrangeEquations Physics command, introduced in 2023, can perform that elusive computation of Einstein's equations, not using tabulated cases, properly handling several (traditional or not) alternative ways of presenting the Lagrangian (the integrand in the Action), taking advantage of the functional differentiation capabilities of the Physics package .

The computation can be performed in one call to Physics:-LagrangeEquations, or in steps interactively using the Physics:-Fundiff and Physics:-Simplify commands that include newly implemented capabilities for simplifying tensorial expressions in curved spacetimes. This is an exciting breakthrough/milestone in Computer Algebra, also in an area out of reach of neural network AIs.

Einstein's equations are a system of second order nonlinear coupled partial differential equations for the 10 components of the spacetime metric

>

Physics:-g_[mu, nu] = Matrix(%id = 36893488152225885228)

(1)

In the Lagrangian formulation, the coordinates of the problem are the 10 components of the metric,

>

(2)

and the parameters of the variational problem are the spacetime coordinates .

The traditional Lagrangian

The simplest case is that of Einstein's equation in vacuum, for which the Lagrangian density is expressed in terms of the trace of the Ricci tensor and the determinant of the metric by

>

(3)

What was almost a dream in previous years, Einstein's equations can now be computed in one simple instruction as the Lagrange equations for this Lagrangian, in the traditional compact form, taking the components of the metric tensor as the coordinates (for an interactive, step by step computation, see further below)

>

-(1/2)*Physics:-g_[mu, nu]*Physics:-Ricci[alpha, `~alpha`]+Physics:-Ricci[mu, nu] = 0

(4)

The tensorial equation computed is also the definition of the Einstein tensor

>

Physics:-Einstein[mu, nu] = -(1/2)*Physics:-g_[mu, nu]*Physics:-Ricci[alpha, `~alpha`]+Physics:-Ricci[mu, nu]

(5)

A Lagrangian depending only on first derivatives of the metric

The Lagrangian used to compute Einstein's equations (4) contains first and second derivatives of the metric; the latter make the computation significantly more complicated. The second order derivatives, however, can be removed from the formulation. To see that, rewrite L in terms of Christoffel symbols

>

(6)

Recalling the definition

>

Physics:-Christoffel[alpha, mu, nu] = (1/2)*Physics:-d_[nu](Physics:-g_[alpha, mu], [X])+(1/2)*Physics:-d_[mu](Physics:-g_[alpha, nu], [X])-(1/2)*Physics:-d_[alpha](Physics:-g_[mu, nu], [X])

(7)

in the two terms containing derivatives of Christoffel symbols contain second order derivatives of .

Now, it is always possible to add a total spacetime derivative to without changing Einstein's equations (assuming the variation of the metric in the corresponding boundary integrals vanishes), and in that way, in this particular case of , obtain a Lagrangian involving only 1st order derivatives. The total derivative, expressed using the inert command to see it before the differentiation operation is performed, is

>

(8)

Adding this term to , performing the differentiation operation and simplifying we get

>

(9)

>

(10)

>

(11)

which is a Lagrangian depending only on 1st order derivatives of the metric through Christoffel symbols. As expected, the equations of motion resulting from this Lagrangian are the same Einstein equations computed in (4); this computation can also now be performed in a single call

>

-(1/2)*Physics:-Ricci[iota, `~iota`]*Physics:-g_[mu, nu]+Physics:-Ricci[mu, nu] = 0

(12)

Simplification capabilities developed to cover these computations

To illustrate the new Maple 2025 tensorial simplification capabilities note that is no just ≡ (6) after discarding its two terms involving derivatives of Christoffel symbols. To verify this, split into the terms containing or not derivatives of Christoffel

>

(13)

Comparing, the total derivative TD≡ (8) is not just , but

>

(14)

Things like these, , can now be verified directly with the new tensorial simplification capabilities: take the left-hand side minus the right-hand side, evaluate the inert derivative and simplify to see the equality is true

>

(15)

>

(16)

>

(17)

That said, it is also true that results in the Lagrangian , and since the equations of movement don't depend on the sign of the Lagrangian, for this Lagrangian adding the term happens to be equivalent to just discarding the terms of involving derivatives of Christoffel symbols.

Derivation step by step using functional differentiation (variational principle)

Also new in Maple 2025, due to the extension of Fundiff to compute in curved spacetimes, it is now possible to compute Einstein's equations from first principles by constructing the action,

>

Int(Int(Int(Int((-%g_[determinant])^(1/2)*Physics:-Ricci[alpha, `~alpha`], x = -infinity .. infinity), y = -infinity .. infinity), z = -infinity .. infinity), t = -infinity .. infinity)

(18)

and equating to zero the functional derivative with respect to the metric. To avoid displaying the resulting large expression, end the input line with ":"

>

Simplifying this result, we get an expression in terms of Christoffel symbols and its derivatives

>

(19)

In this result, we see derivatives of Christoffel symbols, expressed using the D_ command for covariant differentiation. Although, such objects have not the geometrical meaning of a covariant derivative, computationally, they here represent what would be a covariant derivative if the Christoffel symbols were a tensor. For example,

>

(20)

With this computational meaning for the derivatives of Christoffel symbols appearing in (19), rewrite EEC≡(19) in terms of the Ricci and Riemann tensors. For that, consider the definition

>

(21)

Rewrite the noncovariant derivatives in terms of derivatives using the computational representation (20), simplify and isolate one of them

>

(22)

>

(23)

>

(24)

Analogously, derive an expression to rewrite derivatives of Christoffel symbols using the Riemann tensor

>

(25)

>

(26)

>

(27)

>

(28)

Substitute these two equations, in sequence, into Einstein's equations EEC≡(19)

>

(29)

Simplify to arrive at the traditional compact form of Einstein's equations

>

(1/2)*Physics:-Ricci[chi, `~chi`]*Physics:-g_[`~alpha`, `~beta`]-Physics:-Ricci[`~alpha`, `~beta`] = 0

(30)

Download Einstein_Equations_From_First_Principles.mw

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

Inert arithmetic operators %*, %+, ... i...

Posted: ecterrab 14540 Product: Maple

July 23 2023

7 0

Hi
It's been years since expressions like A %* B %+ C involving inert arithmetic operators used in infix form are correctly understood (parsed) when written on a 1D-Math input line. The idea is simple: have the operators %. %*, %+, %-, %^, %/ work on input as infix operators the same way their active forms: ., *, +, -, ^, / do. This useful functionality, however, remained elusive when using 2D-Math input notation, so one would have to resort to using the functional form of the operators. E.g., input the above as `%+`(`%*`(A, B) ,C), which for me is really ugly. Besides being a bit demoralizing: we do all this fuzz about how great computer algebra and 2D-Math input notation is, and then input things in that way …

So this is to mention that this elusive functionality of inert arithmetic operators used in infix form within a 2D-Math input line now works. The novelty is present in the latest Maplesoft Physics Updates for Maple 2023, which is version 1490. As usual, to install the Updates open a Maple worksheet and input Physics:-Version(latest).

Here is an image (worksheet at the end) showing the new thing

The implementation is pretty new; reports of anything related to these inert operators not displayed/working as you'd expect are much appreciated.

Download Inert_arithmetic_operators_in_2D_Math.mw

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

Physics Courseware Support: Mechanics...

Posted: ecterrab 14540 Product: Maple

February 09 2023

12 1

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.

Part I

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

a.	The position as a function of time

b.	The velocity

c.	The acceleration

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

[`&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 vectorsis given by

(2)

Problem

Rewrite the position vector 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 do not depend on time, they are constant vectors. So is entered as

(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.

(21)

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

(22)

(23)

(24)

The position 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

Problem

Rewrite the velocity in cylindrical and spherical coordinates while making the curvilinear unit vectors' time dependency explicit .

	Solution

The acceleration

Problem

Rewrite the acceleration in cylindrical and spherical components while making the curvilinear unit vectors' time dependency explicit.

	Solution

Deriving these formulas

All these results for the position , velocity and acceleration 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)

(38)

The above result contains, in the last equation, the cylindrical radial unit vector ; rewrite it in the spherical basis

(39)

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

(40)

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

(41)

(42)

Introducing the differentiation rules (40) for the unit vectors

(43)

Performing the inert (grayed) derivatives

(44)

In the same way, for the acceleration

(45)

%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)

%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)

(48)

(49)

(50)

Collect vector components

(51)

Summary

•	You can express and 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 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))

(62)

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

(63)

where is the acceleration and is the linear momentum, so in terms of

(64)

We define the angular momentum of a particle, and the torque acting upon it, as

(65)

(66)

Differentiating the definition of

(67)

Since is parallel to , the first term in the above cancels, and in the second term, from (64),

(68)

from which

(69)

•	As discussed below , in the case of a closed system, and these two equations result in

that is, the linear and angular momentum are conserved quantities. Note that does not require that , only that .

The equations of motion - vectorial form

Problem

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

a) The trajectory starting from
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 be the position of the particle in a reference system; then, the velocity and acceleration are given by

(70)

(71)

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

(72)

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

_______________________________________

b) The integration of vectorial equations is also frequently performed after expressing , and in a particular system of coordinates. For example, in the Cartesian system (71) has the form

(73)

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

(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

_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 and , specifying the integration constants ${c__3, c__4, c__5, c__6, c__7, c__8}$ , can also be written in components

(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:

${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)

The case of constant acceleration

Problem

Starting from the vectorial equation (72) for , 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 and 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 , where 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 , 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 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

a) Determine the regions where the motion has positive and negative acceleration. Compute the position at .

b) Compute the velocity corresponding to , starting - not from this expression for but from the acceleration

	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.

	Cartesian coordinates

	Curvilinear coordinates

Many-particle systems

Center of Mass

Given a system of n particles of masses with positions 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 measured in K. Show that:

a) When is the velocity of the center of mass, the total momentum measured in K' is equal to 0.

b) The relation between and the velocity of the center of mass, both measured in K, is the same as the relation 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 2^nd law for each particle implies that , where is the center of mass and 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 satisfies

b) The total torque 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 with the horizontal plane and its projection BC over this plane forms an angle 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))

(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 and derive the equations of motion as the Lagrange equations for L.

For closed systems, the potential energy is related to the force acting on each particle by the equation . Formally, is the gradient operator in the basis onto which is projected, and with respect to its coordinates - say in Cartesian basis and coordinates .

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

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

(209)

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

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

(210)

where is the velocity of the i^th particle. Generally speaking, the potential energy of the system is a function of the coordinates of the n particles, and the Lagrangian is defined as

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

(211)

The potential energy is related to the force acting on each particle by the equation . Formally, is the gradient operator in the basis onto which is projected. Knowing the Lagrangian, we can derive the (Lagrange) equations of motion as

(212)

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

b) Oscillates horizontally on the plane of the pendulum according to .

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

	Solution

Conservation laws

Work

By definition, the work performed by a force to move a particle an infinitesimal amount is . Thus, the work to move it from to along some path C is

"((∫)[A]^(B) (F ) . &DifferentialD;r)[path=C]"

Problem

A particle is submitted to a force whose Cartesian components are given by. Calculate the work done by this force when moving the particle along a straight line from the origin to a point ().

	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 i^th particle of the system can be derived from a potential, . 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 . 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 2^nd law follows .

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 3^rdlaw.

Problem

Consider a system of n particles measured in two frames of reference K and K' that move relative to each other with velocity . Show that the system's momenta and are related by

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

	Solution

Problem

A particle of mass m moving with velocity leaves a half-space in which its potential energy is a constant and enters another in which its potential energy is a different constant . Determine the change in direction of motion of the particle; that is, sin(`θ__1`)/sin(`θ__2`) where and are the angles between an axis perpendicular to the separating plane and the momentum 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 . From it follows that , that is, the angular momentum is conserved.

Problem

a) Express the Cartesian components of the angular momentum , as well as its norm, in cylindrical and spherical coordinates.

b) Rewrite 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 from each other. Show that the momenta and of the system are related by

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 . Show that the momenta and 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 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 and are the total linear and angular momenta measured in K'.
b) Show that when the origin of K' is the center of mass , this formula reduces to

,

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

Solution

a) The momentum of the system measured in K is given by

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 and using , where :

(323)

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

(324)

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 (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 and , introduce the relation between and

(327)

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)

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 and (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 , then

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

(331)

and were related in `#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 in terms of instead of

$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)

(333)

which is the result we were looking for.

Cyclic coordinates

Any generalized coordinate 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 is cyclic, the right-hand side is 0 and so the quantity is conserved.

Problem

The Lagrangian describing the movement of a particle in a central field has as a cyclic coordinate. Using cylindrical coordinates, show that the corresponding conserved quantity is the z component of the angular momentum

	Solution

Integration of the equations of motion

Motion in one dimension

Problem
For a closed system, or any system where the total energy 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 .

	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 .

	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 the result of a) becomes the result obtained for the previous two-body problem, equation .

Solution

a) Let represent the position vector of the particle of mass M and those of the particles of mass m. The Lagrangian is given by

(386)

(387)

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

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

(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

(391)

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)≡ is sufficient to eliminate and from the Lagrangian (389). In steps, eliminating 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, {`&Ropf;_`[a](t) = r_[a](t)-r__M_(t)}, {r_[a](t)})$

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

(393)

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

(394)

To eliminate , use

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

(395)

Likewise, to eliminate use

(396)

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

L = (1/2)*M*Physics:-Vectors:-Norm(diff(-m*(Sum(`&Ropf;_`[a](t), a = 1 .. n))/(m*n+M), t))^2+(1/2)*m*(Sum(Physics:-Vectors:-Norm(diff(`&Ropf;_`[a](t)-m*(Sum(`&Ropf;_`[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

L = (1/2)*M*m^2*Physics:-Vectors:-Norm(Sum(diff(`&Ropf;_`[a](t), t), a = 1 .. n))^2/(m*n+M)^2+(1/2)*m*(Sum(Physics:-Vectors:-Norm(diff(`&Ropf;_`[a](t), t))^2-2*m*Physics:-Vectors:-`.`(diff(`&Ropf;_`[a](t), t), Sum(diff(`&Ropf;_`[a](t), t), a = 1 .. n))/(m*n+M)+m^2*Physics:-Vectors:-Norm(Sum(diff(`&Ropf;_`[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

(399)

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

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

(400)

This Lagrangian involves only the n relative position vectors , achieving the reduction of the original problem of particles to a problem of only n particles.

_______________________________________

b) To obtain the result for the two-body problem, substitute into

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

(401)

Expand only the Sums and collect the Norm

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

(402)

L = (1/2)*M*m*Physics:-Vectors:-Norm(diff(`&Ropf;_`[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

(404)

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

(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(`&Ropf;_`[1](t), t))^2/(M+m)-U), {m[1]*m[2]/(m[1]+m[2]) = mu})$

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

(406)

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 of a particle that moves in a central field is conserved, so evolves in time on a fixed plane perpendicular to the constant . Show that the surface swept per second by the position vector is constant (Kepler's second law).

	Solution

Problem
Starting from the constancy of the energy E and the angular momentum , 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 and

b) interactively, one step at a time, uncouple the equations of movement eliminating from the problem, resulting in an implicit solution .

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

	Solution

Kepler's problem

Problem

Show that, when , with 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

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, . 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 . Write the equation of motion as Newton's 2^nd 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 acts. For the oscillations to be small, must produce only small displacements. The total force is .

a) Write the equation of motion as Newton's 2^nd 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 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 .

d) Show that a solution for the case considered in c), that is, , 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 .

e) Specialize the solution of item d) in the case of resonance, when , 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 and 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

where represents the displacement of the generalized coordinate , 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))

(524)

Denoting the generalized coordinates by , the potential energy can be expanded in series around the minimums . Since the movement consist only in small displacements around , it is sufficient to keep terms in the expansion up to order 2, resulting in an expression analogous to 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))

(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))

(526)

where we take the at the minimums and denote them as

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 and can be split into symmetric and antisymmetric parts, with the antisymmetric parts canceling out in view of the symmetric character of in U and in T. Therefore, we can take and 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 describe independent displacements, so one can think of 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.

(528)

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

${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 can now be written as

(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

(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, Physics:-Vectors:-diff(x[c](t), t)), t) = %diff(L, x[c](t))

(532)

(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,

(534)

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 of each particle by , where is the mass density as a function of the position and 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 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 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 of a point P of a body, measured in an inertial reference system K, can be written as

where is the velocity of the origin of K', a frame of reference attached to the body's center of mass, is the body's angular velocity (its instantaneous counter-clockwise rotation speed around some axis in the direction of a unit vector) and 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 and are related by a rotation matrix which, as all rotation matrices, is orthogonal.

	Solution

Inertia tensor

Problem

a) Show that using , derived in Angular velocity for the velocity of a point P of a rigid body in terms of and the position of P viewed from the center of mass , the kinetic energy of the rigid body can be written in terms of the positions of the n particles (not their velocities), the velocity of the center of mass and angular velocity 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

where

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

is the Inertia tensor, represents the components of the vector and represents the components of the position vector of the particle.

	Solution

Problem

Determine the Inertia tensor corresponding to a triatomic molecule that has the form of an isosceles triangle with two masses in the extremes of the base and a mass in the third vertex. The distance between the two masses 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 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 of the angular momentum of a rigid body can be expressed in terms of the inertia tensor and the components of the angular velocity as

	Solution

The equations of motion of a rigid body

A rigid body is a system with six degrees of freedom: three indicating the position 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

(615)

(616)

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

(617)

where is the total momentum, is the total external force acting upon the body, is the external force acting upon the particle, the total angular momentum and is the total torque.

Problem
Show that the equations of movement of a rigid body can be computed as the Lagrange equations for and from the Lagrangian

where mu = Sum(m[i], i = 1 .. n) is the total mass, is the inertia tensor, is the potential energy for the external force , and .

Solution

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

(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

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 (location of the center of mass) and (the rotation axis, so that ), the Lagrangian in vectorial form is given by

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 is derived as the Lagrange equation for this Lagrangian

(621)

(622)

(623)

On the right-hand side, with respect to its first argument , equivalent to the gradient taking as the coordinates, and is the total momentum

(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 . If more precision is required, the dependency of the potential energy could be expressed in terms of the norm of , as in

(625)

in which case differentiating with respect to is equivalent to the definition of directional derivative

(626)

The same computation, this time with respect to the coordinates where is the corresponding velocity,

(627)

(628)

(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 of the angular momentum .

$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)

(631)

(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)

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)

(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,

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

(635)

(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

(637)

(638)

(639)

(640)

The right-hand side is the component of the variation of the potential energy with respect to . A graphics analysis of and algebraic derivation as done in the problem of the section Angular velocity results in

(641)

(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 , resulting in the second equation of movement of a rigid body, for the time derivative of the angular momentum

(643)

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 with regards to an inertial reference system .

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 is the translational acceleration of the frame K' as seen from .

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

	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 . Denote the position vector and velocity in the non-inertial frame J as and ,

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

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

where is the Coriolis force and is the centrifugal force.

Solution

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

(658)

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

(659)

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

(660)

where is the translational velocity of K' viewed from . Inserting this relation (660) into gives , 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))

(661)

In turn, the velocities and in the frames K' and J are related by

(662)

where is the angular velocity of the frame J viewed from K'. Also, since K' and J have the same origin, and the Lagrangian in J is

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)

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

(665)

(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

(668)

Evaluate these derivatives and replace

(669)

Isolating and collecting vector products we get the expected result

(670)

(671)

Part II (forthcoming)

	The Hamiltonian and equations of motion; Poisson brackets

	Canonical transformations

	The Hamilton-Jacobi equation

	References

Download Mechanics_(2022).mw

Download Mechanics_(2022).pdf

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

E-Mail Address:
Password:
Remember Me:	Automatically sign in on future visits

E-Mail Address:
Password:
Remember Me:	Automatically sign in on future visits

Ask a Question

Create a Post

14540 Reputation

24 Badges

MaplePrimes Activity

The 2025 Maplesoft Physics Updates and i...

What is new in Physics in Maple 2025...

Computing Einstein's equations using var...

Inert arithmetic operators %*, %+, ... i...

Physics Courseware Support: Mechanics...

Save this setting as your default sorting preference?

Ask a Question

Create a Post

Generating PDF…

Save this setting as your default sorting preference?
Note: You can change your preference any time in your account settings
Don't show this again