MaplePrimes Posts

MaplePrimes Posts are for sharing your experiences, techniques and opinions about Maple, MapleSim and related products, as well as general interests in math and computing.

Latest Post
  • Latest Posts Feed
  • I don't have the latest Maple, and I'm sure this isn't in the latest version. 

    One thing that has been an annoyance for all time, and it gets me time and time again, is not having a global degrees or radians setting. 

    Of course it needs to be a setting option, otherwise it would break many older worksheets. 

    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.

       So, someone just posted spam here. Cleverly designed, the title was an HTTP link that sent the reader to a .in website. Accessing the actual post was impossible, and marking it as spam was impossible as well. Debugging tools had to be used.

       The design of this website should be improved to make this impossible. An unalterable icon should be present to connect the reader to the post without relying on the user-defined title.

    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 1st, 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 A[mu](x, y, z, t), there are then four coordinates, one for each of the values of the index mu, and there are four parameters x, y, z, t.

     

    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

     

    with(Physics)

    Setup(mathematicalnotation = true, coordinates = cartesian)

    [coordinatesystems = {X}, mathematicalnotation = true]

    (1)

    The lambda*Phi^4 model in classical field theory and corresponding field equations, as in previous releases

    CompactDisplay(Phi(X))

    Phi(x, y, z, t)*`will now be displayed as`*Phi

    (2)

    L := (1/2)*d_[mu](Phi(X))*d_[mu](Phi(X))-(1/2)*m^2*Phi(X)^2+(1/4)*lambda*Phi(X)^4

    (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

    LagrangeEquations(L, Phi)

    Phi(X)^3*lambda-Phi(X)*m^2-Physics:-dAlembertian(Phi(X), [X]) = 0

    (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 2/sqrt(-%g_[determinant]) (in flat spacetimes sqrt(-%g_[determinant]) = 1). 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

    LagrangeEquations(L, EnergyMomentum[mu, nu])

    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 T[mu, nu], you can use the Define  command

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

    `Defined objects with tensor properties`

     

    {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 T[mu, nu]

    EnergyMomentum[definition]

    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)

    Library:-IsTensorialSymmetric(EnergyMomentum[mu, nu])

    true

    (8)

    EnergyMomentum[]

    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

    g_[arb]

    _______________________________________________________

     

    `Setting `*lowercaselatin_is*` letters to represent `*space*` indices`

     

    `The arbitrary metric in coordinates `*[x, y, z, t]

     

    `Signature: `(`- - - +`)

     

    _______________________________________________________

     

    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

    L := sqrt(-%g_[determinant])*L

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

    LagrangeEquations(L, Phi)

    Phi(X)^3*lambda-Phi(X)*m^2-Physics:-D_[kappa](Physics:-d_[`~kappa`](Phi(X), [X]), [X]) = 0

    (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

    LagrangeEquations(L, EnergyMomentum[mu, nu])

    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 g[mu, nu] 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 g[mu, nu], and as in the case of electrodynamics the parameters are the spacetime coordinates X^alpha. 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

    L := sqrt(-%g_[determinant])*Ricci[alpha, `~alpha`]

    (-%g_[determinant])^(1/2)*Physics:-Ricci[alpha, `~alpha`]

    (14)

    Einstein's equations in vacuum:

    LagrangeEquations(L, g_[mu, nu])

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

    (15)

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

    Einstein[definition]

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

    (16)

    The Lagrangian L 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

    L__C := convert(L, Christoffel)

    (-%g_[determinant])^(1/2)*Physics:-g_[`~alpha`, `~lambda`]*(Physics:-d_[nu](Physics:-Christoffel[`~nu`, alpha, lambda], [X])-Physics:-d_[lambda](Physics:-Christoffel[`~nu`, alpha, nu], [X])+Physics:-Christoffel[`~beta`, alpha, lambda]*Physics:-Christoffel[`~nu`, beta, nu]-Physics:-Christoffel[`~beta`, alpha, nu]*Physics:-Christoffel[`~nu`, beta, lambda])

    (17)

    Recalling the definition

    Christoffel[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 L[C] the two terms containing derivatives of Christoffel symbols contain second order derivatives of g[mu, nu]. Now, it is always possible to add a total spacetime derivative to L[C] 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 L[C], 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

    TD := %d_[alpha](g_[`~mu`, `~nu`]*sqrt(-%g_[determinant])*(g_[`~alpha`, mu]*Christoffel[`~beta`, nu, beta]-Christoffel[`~alpha`, mu, nu]))

    %d_[alpha](Physics:-g_[`~mu`, `~nu`]*(-%g_[determinant])^(1/2)*(Physics:-g_[mu, `~alpha`]*Physics:-Christoffel[`~beta`, beta, nu]-Physics:-Christoffel[`~alpha`, mu, nu]))

    (19)

    Adding this term to L[C], performing the `∂` differentiation operation and simplifying we get

    L__1 := L__C+TD

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

    (20)

    L__1 := eval(L__1, %d_ = d_)

    (-%g_[determinant])^(1/2)*Physics:-g_[`~alpha`, `~lambda`]*(Physics:-d_[nu](Physics:-Christoffel[`~nu`, alpha, lambda], [X])-Physics:-d_[lambda](Physics:-Christoffel[`~nu`, alpha, nu], [X])+Physics:-Christoffel[`~beta`, alpha, lambda]*Physics:-Christoffel[`~nu`, beta, nu]-Physics:-Christoffel[`~beta`, alpha, nu]*Physics:-Christoffel[`~nu`, beta, lambda])+Physics:-d_[alpha](Physics:-g_[`~mu`, `~nu`], [X])*(-%g_[determinant])^(1/2)*(Physics:-g_[mu, `~alpha`]*Physics:-Christoffel[`~beta`, beta, nu]-Physics:-Christoffel[`~alpha`, mu, nu])-(1/2)*Physics:-g_[`~mu`, `~nu`]*(Physics:-g_[mu, `~alpha`]*Physics:-Christoffel[`~beta`, beta, nu]-Physics:-Christoffel[`~alpha`, mu, nu])*%g_[determinant]*Physics:-g_[`~kappa`, `~lambda`]*Physics:-d_[alpha](Physics:-g_[kappa, lambda], [X])/(-%g_[determinant])^(1/2)+Physics:-g_[`~mu`, `~nu`]*(-%g_[determinant])^(1/2)*(Physics:-g_[mu, `~alpha`]*Physics:-d_[alpha](Physics:-Christoffel[`~beta`, beta, nu], [X])-Physics:-d_[alpha](Physics:-Christoffel[`~alpha`, mu, nu], [X]))

    (21)

    L__1 := Simplify(L__1)

    (Physics:-Christoffel[alpha, beta, kappa]*Physics:-Christoffel[`~beta`, `~alpha`, `~kappa`]-Physics:-Christoffel[alpha, beta, `~alpha`]*Physics:-Christoffel[`~beta`, kappa, `~kappa`])*(-%g_[determinant])^(1/2)

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

    LagrangeEquations(L__1, g_[mu, nu])

    -(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 `≡`(L[1], (Physics[Christoffel][alpha, beta, kappa]*Physics[Christoffel][`~beta`, `~alpha`, `~kappa`]-Physics[Christoffel][alpha, beta, `~alpha`]*Physics[Christoffel][`~beta`, kappa, `~kappa`])*(-%g_[determinant])^(1/2)) is no just L[C] ≡ (17) after discarding its two terms involving derivatives of Christoffel symbols. To verify this, split L[C] into the terms containing or not derivatives of Christoffel

    L__22, L__11 := selectremove(has, expand(L__C), d_)

    (-%g_[determinant])^(1/2)*Physics:-g_[`~alpha`, `~lambda`]*Physics:-d_[nu](Physics:-Christoffel[`~nu`, alpha, lambda], [X])-(-%g_[determinant])^(1/2)*Physics:-g_[`~alpha`, `~lambda`]*Physics:-d_[lambda](Physics:-Christoffel[`~nu`, alpha, nu], [X]), (-%g_[determinant])^(1/2)*Physics:-g_[`~alpha`, `~lambda`]*Physics:-Christoffel[`~beta`, alpha, lambda]*Physics:-Christoffel[`~nu`, beta, nu]-(-%g_[determinant])^(1/2)*Physics:-g_[`~alpha`, `~lambda`]*Physics:-Christoffel[`~beta`, alpha, nu]*Physics:-Christoffel[`~nu`, beta, lambda]

    (24)

    Comparing, the total derivative TD≡ (19) is not just -L[22], but

    TD = -L__22-2*L__11

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

    (25)

    Things like these, TD = -L__22-2*L__11, 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

    (lhs-rhs)(%d_[alpha](Physics[g_][`~mu`, `~nu`]*(-%g_[determinant])^(1/2)*(Physics[g_][mu, `~alpha`]*Physics[Christoffel][`~beta`, beta, nu]-Physics[Christoffel][`~alpha`, mu, nu])) = -(-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*Physics[d_][nu](Physics[Christoffel][`~nu`, alpha, lambda], [X])+(-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*Physics[d_][lambda](Physics[Christoffel][`~nu`, alpha, nu], [X])-2*(-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*Physics[Christoffel][`~beta`, alpha, lambda]*Physics[Christoffel][`~nu`, beta, nu]+2*(-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*Physics[Christoffel][`~beta`, alpha, nu]*Physics[Christoffel][`~nu`, beta, lambda])

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

    (26)

    eval(%d_[alpha](Physics[g_][`~mu`, `~nu`]*(-%g_[determinant])^(1/2)*(Physics[g_][mu, `~alpha`]*Physics[Christoffel][`~beta`, beta, nu]-Physics[Christoffel][`~alpha`, mu, nu]))+(-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*Physics[d_][nu](Physics[Christoffel][`~nu`, alpha, lambda], [X])-(-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*Physics[d_][lambda](Physics[Christoffel][`~nu`, alpha, nu], [X])+2*(-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*Physics[Christoffel][`~beta`, alpha, lambda]*Physics[Christoffel][`~nu`, beta, nu]-2*(-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*Physics[Christoffel][`~beta`, alpha, nu]*Physics[Christoffel][`~nu`, beta, lambda], %d_ = d_)

    (-%g_[determinant])^(1/2)*(Physics:-g_[mu, `~alpha`]*Physics:-Christoffel[`~beta`, beta, nu]-Physics:-Christoffel[`~alpha`, mu, nu])*Physics:-d_[alpha](Physics:-g_[`~mu`, `~nu`], [X])-(1/2)*(Physics:-g_[mu, `~alpha`]*Physics:-Christoffel[`~beta`, beta, nu]-Physics:-Christoffel[`~alpha`, mu, nu])*Physics:-g_[`~mu`, `~nu`]*%g_[determinant]*Physics:-g_[`~kappa`, `~lambda`]*Physics:-d_[alpha](Physics:-g_[kappa, lambda], [X])/(-%g_[determinant])^(1/2)+(-%g_[determinant])^(1/2)*(Physics:-g_[mu, `~alpha`]*Physics:-d_[alpha](Physics:-Christoffel[`~beta`, beta, nu], [X])-Physics:-d_[alpha](Physics:-Christoffel[`~alpha`, mu, nu], [X]))*Physics:-g_[`~mu`, `~nu`]+(-%g_[determinant])^(1/2)*Physics:-g_[`~alpha`, `~lambda`]*Physics:-d_[nu](Physics:-Christoffel[`~nu`, alpha, lambda], [X])-(-%g_[determinant])^(1/2)*Physics:-g_[`~alpha`, `~lambda`]*Physics:-d_[lambda](Physics:-Christoffel[`~nu`, alpha, nu], [X])+2*(-%g_[determinant])^(1/2)*Physics:-g_[`~alpha`, `~lambda`]*Physics:-Christoffel[`~beta`, alpha, lambda]*Physics:-Christoffel[`~nu`, beta, nu]-2*(-%g_[determinant])^(1/2)*Physics:-g_[`~alpha`, `~lambda`]*Physics:-Christoffel[`~beta`, alpha, nu]*Physics:-Christoffel[`~nu`, beta, lambda]

    (27)

    Simplify((-%g_[determinant])^(1/2)*(Physics[g_][mu, `~alpha`]*Physics[Christoffel][`~beta`, beta, nu]-Physics[Christoffel][`~alpha`, mu, nu])*Physics[d_][alpha](Physics[g_][`~mu`, `~nu`], [X])-(1/2)*(Physics[g_][mu, `~alpha`]*Physics[Christoffel][`~beta`, beta, nu]-Physics[Christoffel][`~alpha`, mu, nu])*Physics[g_][`~mu`, `~nu`]*%g_[determinant]*Physics[g_][`~kappa`, `~lambda`]*Physics[d_][alpha](Physics[g_][kappa, lambda], [X])/(-%g_[determinant])^(1/2)+(-%g_[determinant])^(1/2)*(Physics[g_][mu, `~alpha`]*Physics[d_][alpha](Physics[Christoffel][`~beta`, beta, nu], [X])-Physics[d_][alpha](Physics[Christoffel][`~alpha`, mu, nu], [X]))*Physics[g_][`~mu`, `~nu`]+(-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*Physics[d_][nu](Physics[Christoffel][`~nu`, alpha, lambda], [X])-(-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*Physics[d_][lambda](Physics[Christoffel][`~nu`, alpha, nu], [X])+2*(-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*Physics[Christoffel][`~beta`, alpha, lambda]*Physics[Christoffel][`~nu`, beta, nu]-2*(-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*Physics[Christoffel][`~beta`, alpha, nu]*Physics[Christoffel][`~nu`, beta, lambda])

    0

    (28)

    That said, it is also true that TD = -L[22]-2*L[11] results in the Lagrangian L[1] = -L[11], and since the equations of movement don't depend on the sign of the Lagrangian, for this Lagrangian `≡`(L[C], (-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*(Physics[d_][nu](Physics[Christoffel][`~nu`, alpha, lambda], [X])-Physics[d_][lambda](Physics[Christoffel][`~nu`, alpha, nu], [X])+Physics[Christoffel][`~beta`, alpha, lambda]*Physics[Christoffel][`~nu`, beta, nu]-Physics[Christoffel][`~beta`, alpha, nu]*Physics[Christoffel][`~nu`, beta, lambda])) adding the term TD happens to be equivalent to just discarding the terms of L__C 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,

    S := Intc(L, X)

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

    EE__unsimplified := Fundiff(S, g_[alpha, beta]) = 0

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

    EEC := Simplify(EE__unsimplified)

    (1/4)*(2*Physics:-Christoffel[chi, iota, kappa]*Physics:-Christoffel[`~iota`, `~chi`, `~kappa`]-2*Physics:-Christoffel[chi, iota, `~chi`]*Physics:-Christoffel[`~iota`, kappa, `~kappa`]-2*Physics:-D_[iota](Physics:-Christoffel[chi, `~chi`, `~iota`], [X])+2*Physics:-D_[chi](Physics:-Christoffel[`~chi`, iota, `~iota`], [X]))*Physics:-g_[`~alpha`, `~beta`]+(1/4)*(Physics:-Christoffel[`~alpha`, chi, `~beta`]+Physics:-Christoffel[`~beta`, chi, `~alpha`])*Physics:-Christoffel[`~chi`, iota, `~iota`]-(1/4)*Physics:-Christoffel[`~beta`, chi, iota]*Physics:-Christoffel[`~chi`, `~alpha`, `~iota`]-(1/2)*Physics:-Christoffel[chi, iota, `~alpha`]*Physics:-Christoffel[`~iota`, `~beta`, `~chi`]+(1/2)*Physics:-Christoffel[chi, `~alpha`, `~beta`]*Physics:-Christoffel[iota, `~chi`, `~iota`]-(1/4)*Physics:-Christoffel[`~alpha`, chi, iota]*Physics:-Christoffel[`~chi`, `~beta`, `~iota`]+(1/4)*Physics:-D_[chi](Physics:-Christoffel[`~alpha`, `~beta`, `~chi`], [X])+(1/4)*Physics:-D_[chi](Physics:-Christoffel[`~beta`, `~alpha`, `~chi`], [X])-(1/2)*Physics:-D_[chi](Physics:-Christoffel[`~chi`, `~alpha`, `~beta`], [X])-(1/4)*Physics:-D_[`~alpha`](Physics:-Christoffel[`~beta`, chi, `~chi`], [X])+(1/2)*Physics:-D_[`~beta`](Physics:-Christoffel[chi, `~alpha`, `~chi`], [X])-(1/4)*Physics:-D_[`~beta`](Physics:-Christoffel[`~alpha`, chi, `~chi`], [X]) = 0

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

    "`▿`[chi](GAMMA[]^(alpha,beta,chi)) :"

    % = expand(%)

    Physics:-D_[chi](Physics:-Christoffel[`~alpha`, `~beta`, `~chi`], [X]) = Physics:-d_[chi](Physics:-Christoffel[`~alpha`, `~beta`, `~chi`], [X])+Physics:-Christoffel[`~alpha`, chi, mu]*Physics:-Christoffel[`~mu`, `~beta`, `~chi`]+Physics:-Christoffel[`~beta`, chi, mu]*Physics:-Christoffel[`~alpha`, `~chi`, `~mu`]+Physics:-Christoffel[`~chi`, chi, mu]*Physics:-Christoffel[`~alpha`, `~beta`, `~mu`]

    (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

    Ricci[definition]

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

    (32)

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

    convert(Physics[Ricci][mu, nu] = Physics[d_][alpha](Physics[Christoffel][`~alpha`, mu, nu], [X])-Physics[d_][nu](Physics[Christoffel][`~alpha`, mu, alpha], [X])+Physics[Christoffel][`~beta`, mu, nu]*Physics[Christoffel][`~alpha`, beta, alpha]-Physics[Christoffel][`~beta`, mu, alpha]*Physics[Christoffel][`~alpha`, nu, beta], D_)

    Physics:-Ricci[mu, nu] = Physics:-D_[alpha](Physics:-Christoffel[`~alpha`, mu, nu], [X])-Physics:-Christoffel[`~alpha`, alpha, kappa]*Physics:-Christoffel[`~kappa`, mu, nu]+Physics:-Christoffel[`~kappa`, alpha, mu]*Physics:-Christoffel[`~alpha`, kappa, nu]+Physics:-Christoffel[`~kappa`, alpha, nu]*Physics:-Christoffel[`~alpha`, mu, kappa]-Physics:-D_[nu](Physics:-Christoffel[`~alpha`, alpha, mu], [X])-Physics:-Christoffel[`~lambda`, mu, nu]*Physics:-Christoffel[`~alpha`, alpha, lambda]+Physics:-Christoffel[`~beta`, mu, nu]*Physics:-Christoffel[`~alpha`, alpha, beta]-Physics:-Christoffel[`~beta`, alpha, mu]*Physics:-Christoffel[`~alpha`, beta, nu]

    (33)

    Simplify(Physics[Ricci][mu, nu] = D_[alpha](Physics[Christoffel][`~alpha`, mu, nu], [X])-Physics[Christoffel][`~alpha`, alpha, kappa]*Physics[Christoffel][`~kappa`, mu, nu]+Physics[Christoffel][`~kappa`, alpha, mu]*Physics[Christoffel][`~alpha`, kappa, nu]+Physics[Christoffel][`~kappa`, alpha, nu]*Physics[Christoffel][`~alpha`, mu, kappa]-D_[nu](Physics[Christoffel][`~alpha`, alpha, mu], [X])-Physics[Christoffel][`~lambda`, mu, nu]*Physics[Christoffel][`~alpha`, alpha, lambda]+Physics[Christoffel][`~beta`, mu, nu]*Physics[Christoffel][`~alpha`, alpha, beta]-Physics[Christoffel][`~beta`, alpha, mu]*Physics[Christoffel][`~alpha`, beta, nu])

    Physics:-Ricci[mu, nu] = Physics:-Christoffel[alpha, beta, mu]*Physics:-Christoffel[`~beta`, nu, `~alpha`]-Physics:-Christoffel[beta, mu, nu]*Physics:-Christoffel[alpha, `~alpha`, `~beta`]+Physics:-D_[alpha](Physics:-Christoffel[`~alpha`, mu, nu], [X])-Physics:-D_[nu](Physics:-Christoffel[alpha, mu, `~alpha`], [X])

    (34)

    C_to_Ricci := isolate(Physics[Ricci][mu, nu] = Physics[Christoffel][alpha, beta, mu]*Physics[Christoffel][`~beta`, nu, `~alpha`]-Physics[Christoffel][beta, mu, nu]*Physics[Christoffel][alpha, `~alpha`, `~beta`]+D_[alpha](Physics[Christoffel][`~alpha`, mu, nu], [X])-D_[nu](Physics[Christoffel][alpha, mu, `~alpha`], [X]), D_[alpha](Christoffel[`~alpha`, mu, nu]))

    Physics:-D_[alpha](Physics:-Christoffel[`~alpha`, mu, nu], [X]) = -Physics:-Christoffel[alpha, beta, mu]*Physics:-Christoffel[`~beta`, nu, `~alpha`]+Physics:-Christoffel[beta, mu, nu]*Physics:-Christoffel[alpha, `~alpha`, `~beta`]+Physics:-Ricci[mu, nu]+Physics:-D_[nu](Physics:-Christoffel[alpha, mu, `~alpha`], [X])

    (35)

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

    Riemann[`~alpha`, beta, mu, nu, definition]

    Physics:-Riemann[`~alpha`, beta, mu, nu] = Physics:-d_[mu](Physics:-Christoffel[`~alpha`, beta, nu], [X])-Physics:-d_[nu](Physics:-Christoffel[`~alpha`, beta, mu], [X])+Physics:-Christoffel[`~alpha`, upsilon, mu]*Physics:-Christoffel[`~upsilon`, beta, nu]-Physics:-Christoffel[`~alpha`, upsilon, nu]*Physics:-Christoffel[`~upsilon`, beta, mu]

    (36)

    convert(Physics[Riemann][`~alpha`, beta, mu, nu] = Physics[d_][mu](Physics[Christoffel][`~alpha`, beta, nu], [X])-Physics[d_][nu](Physics[Christoffel][`~alpha`, beta, mu], [X])+Physics[Christoffel][`~alpha`, upsilon, mu]*Physics[Christoffel][`~upsilon`, beta, nu]-Physics[Christoffel][`~alpha`, upsilon, nu]*Physics[Christoffel][`~upsilon`, beta, mu], D_)

    Physics:-Riemann[`~alpha`, beta, mu, nu] = Physics:-D_[mu](Physics:-Christoffel[`~alpha`, beta, nu], [X])+Physics:-Christoffel[`~kappa`, mu, nu]*Physics:-Christoffel[`~alpha`, beta, kappa]-Physics:-Christoffel[`~alpha`, kappa, mu]*Physics:-Christoffel[`~kappa`, beta, nu]+Physics:-Christoffel[`~kappa`, beta, mu]*Physics:-Christoffel[`~alpha`, kappa, nu]-Physics:-D_[nu](Physics:-Christoffel[`~alpha`, beta, mu], [X])-Physics:-Christoffel[`~lambda`, mu, nu]*Physics:-Christoffel[`~alpha`, beta, lambda]-Physics:-Christoffel[`~lambda`, beta, nu]*Physics:-Christoffel[`~alpha`, lambda, mu]+Physics:-Christoffel[`~alpha`, lambda, nu]*Physics:-Christoffel[`~lambda`, beta, mu]+Physics:-Christoffel[`~alpha`, mu, upsilon]*Physics:-Christoffel[`~upsilon`, beta, nu]-Physics:-Christoffel[`~alpha`, nu, upsilon]*Physics:-Christoffel[`~upsilon`, beta, mu]

    (37)

    Simplify(Physics[Riemann][`~alpha`, beta, mu, nu] = D_[mu](Physics[Christoffel][`~alpha`, beta, nu], [X])+Physics[Christoffel][`~kappa`, mu, nu]*Physics[Christoffel][`~alpha`, beta, kappa]-Physics[Christoffel][`~alpha`, kappa, mu]*Physics[Christoffel][`~kappa`, beta, nu]+Physics[Christoffel][`~kappa`, beta, mu]*Physics[Christoffel][`~alpha`, kappa, nu]-D_[nu](Physics[Christoffel][`~alpha`, beta, mu], [X])-Physics[Christoffel][`~lambda`, mu, nu]*Physics[Christoffel][`~alpha`, beta, lambda]-Physics[Christoffel][`~lambda`, beta, nu]*Physics[Christoffel][`~alpha`, lambda, mu]+Physics[Christoffel][`~alpha`, lambda, nu]*Physics[Christoffel][`~lambda`, beta, mu]+Physics[Christoffel][`~alpha`, mu, upsilon]*Physics[Christoffel][`~upsilon`, beta, nu]-Physics[Christoffel][`~alpha`, nu, upsilon]*Physics[Christoffel][`~upsilon`, beta, mu])

    Physics:-Riemann[`~alpha`, beta, mu, nu] = -Physics:-Christoffel[`~alpha`, kappa, mu]*Physics:-Christoffel[`~kappa`, beta, nu]+Physics:-Christoffel[`~kappa`, beta, mu]*Physics:-Christoffel[`~alpha`, kappa, nu]+Physics:-D_[mu](Physics:-Christoffel[`~alpha`, beta, nu], [X])-Physics:-D_[nu](Physics:-Christoffel[`~alpha`, beta, mu], [X])

    (38)

    C_to_Riemann := isolate(Physics[Riemann][`~alpha`, beta, mu, nu] = -Physics[Christoffel][`~alpha`, kappa, mu]*Physics[Christoffel][`~kappa`, beta, nu]+Physics[Christoffel][`~kappa`, beta, mu]*Physics[Christoffel][`~alpha`, kappa, nu]+D_[mu](Physics[Christoffel][`~alpha`, beta, nu], [X])-D_[nu](Physics[Christoffel][`~alpha`, beta, mu], [X]), D_[mu](Christoffel[`~alpha`, beta, nu]))

    Physics:-D_[mu](Physics:-Christoffel[`~alpha`, beta, nu], [X]) = Physics:-Christoffel[`~alpha`, kappa, mu]*Physics:-Christoffel[`~kappa`, beta, nu]-Physics:-Christoffel[`~kappa`, beta, mu]*Physics:-Christoffel[`~alpha`, kappa, nu]+Physics:-Riemann[`~alpha`, beta, mu, nu]+Physics:-D_[nu](Physics:-Christoffel[`~alpha`, beta, mu], [X])

    (39)

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

    Substitute(C_to_Riemann, C_to_Ricci, EEC)

    -(1/4)*Physics:-Christoffel[`~beta`, psi, `~alpha`]*Physics:-Christoffel[`~psi`, omega, `~omega`]+(1/4)*Physics:-Christoffel[`~beta`, psi, `~omega`]*Physics:-Christoffel[`~psi`, omega, `~alpha`]+(1/4)*Physics:-Christoffel[`~alpha`, lambda, mu]*Physics:-Christoffel[`~lambda`, `~beta`, `~mu`]-(1/4)*Physics:-Christoffel[`~alpha`, lambda, `~mu`]*Physics:-Christoffel[`~lambda`, mu, `~beta`]+(1/4)*Physics:-Christoffel[`~beta`, nu, sigma]*Physics:-Christoffel[`~sigma`, `~alpha`, `~nu`]-(1/4)*Physics:-Christoffel[`~beta`, sigma, `~nu`]*Physics:-Christoffel[`~sigma`, nu, `~alpha`]+(1/2)*Physics:-Christoffel[omicron, zeta, `~beta`]*Physics:-Christoffel[`~zeta`, `~alpha`, `~omicron`]-(1/2)*Physics:-Christoffel[omicron, zeta, `~omicron`]*Physics:-Christoffel[`~zeta`, `~alpha`, `~beta`]-Physics:-Ricci[`~alpha`, `~beta`]-(1/2)*Physics:-Christoffel[`~tau`, `~alpha`, `~beta`]*Physics:-Christoffel[`~upsilon`, tau, upsilon]-(1/2)*Physics:-Christoffel[chi, iota, `~alpha`]*Physics:-Christoffel[`~iota`, `~beta`, `~chi`]+(1/4)*(Physics:-Christoffel[`~alpha`, chi, `~beta`]+Physics:-Christoffel[`~beta`, chi, `~alpha`])*Physics:-Christoffel[`~chi`, iota, `~iota`]+(1/2)*Physics:-Christoffel[`~upsilon`, tau, `~beta`]*Physics:-Christoffel[`~tau`, upsilon, `~alpha`]-(1/4)*Physics:-Christoffel[`~beta`, chi, iota]*Physics:-Christoffel[`~chi`, `~alpha`, `~iota`]+(1/2)*Physics:-Christoffel[chi, `~alpha`, `~beta`]*Physics:-Christoffel[iota, `~chi`, `~iota`]-(1/4)*Physics:-Christoffel[`~alpha`, chi, iota]*Physics:-Christoffel[`~chi`, `~beta`, `~iota`]-(1/4)*Physics:-D_[`~rho1`](Physics:-Christoffel[`~alpha`, rho1, `~beta`], [X])+(1/4)*Physics:-D_[`~mu`](Physics:-Christoffel[`~alpha`, mu, `~beta`], [X])+(1/4)*Physics:-D_[`~nu`](Physics:-Christoffel[`~beta`, nu, `~alpha`], [X])-(1/2)*Physics:-D_[`~beta`](Physics:-Christoffel[`~upsilon`, upsilon, `~alpha`], [X])-(1/4)*Physics:-D_[`~omega`](Physics:-Christoffel[`~beta`, omega, `~alpha`], [X])+(1/2)*Physics:-D_[`~beta`](Physics:-Christoffel[alpha9, `~alpha`, `~alpha9`], [X])-(1/4)*Physics:-Christoffel[`~alpha`, rho, `~beta`]*Physics:-Christoffel[`~rho`, rho1, `~rho1`]+(1/4)*Physics:-Christoffel[`~alpha`, rho, `~rho1`]*Physics:-Christoffel[`~rho`, rho1, `~beta`]+(1/2)*Physics:-Christoffel[alpha10, `~alpha`, `~beta`]*Physics:-Christoffel[alpha9, `~alpha10`, `~alpha9`]-(1/2)*Physics:-Christoffel[alpha9, alpha10, `~alpha`]*Physics:-Christoffel[`~alpha10`, `~alpha9`, `~beta`]+(1/4)*(2*Physics:-Christoffel[chi, iota, kappa]*Physics:-Christoffel[`~iota`, `~chi`, `~kappa`]-2*Physics:-Christoffel[chi, iota, `~chi`]*Physics:-Christoffel[`~iota`, kappa, `~kappa`]-2*Physics:-Christoffel[alpha4, alpha5, alpha6]*Physics:-Christoffel[`~alpha6`, `~alpha4`, `~alpha5`]+2*Physics:-Christoffel[alpha4, alpha6, `~alpha5`]*Physics:-Christoffel[`~alpha6`, alpha5, `~alpha4`]-2*Physics:-D_[`~alpha5`](Physics:-Christoffel[alpha4, alpha5, `~alpha4`], [X])+2*Physics:-Christoffel[`~alpha1`, alpha1, alpha3]*Physics:-Christoffel[`~alpha3`, alpha2, `~alpha2`]-2*Physics:-Christoffel[`~alpha1`, alpha3, `~alpha2`]*Physics:-Christoffel[`~alpha3`, alpha1, alpha2]+2*Physics:-Ricci[alpha2, `~alpha2`]+2*Physics:-D_[`~alpha2`](Physics:-Christoffel[`~alpha1`, alpha1, alpha2], [X]))*Physics:-g_[`~alpha`, `~beta`] = 0

    (40)

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

    Simplify(-(1/4)*Physics[Christoffel][`~beta`, psi, `~alpha`]*Physics[Christoffel][`~psi`, omega, `~omega`]+(1/4)*Physics[Christoffel][`~beta`, psi, `~omega`]*Physics[Christoffel][`~psi`, omega, `~alpha`]+(1/4)*Physics[Christoffel][`~alpha`, lambda, mu]*Physics[Christoffel][`~lambda`, `~beta`, `~mu`]-(1/4)*Physics[Christoffel][`~alpha`, lambda, `~mu`]*Physics[Christoffel][`~lambda`, mu, `~beta`]+(1/4)*Physics[Christoffel][`~beta`, nu, sigma]*Physics[Christoffel][`~sigma`, `~alpha`, `~nu`]-(1/4)*Physics[Christoffel][`~beta`, sigma, `~nu`]*Physics[Christoffel][`~sigma`, nu, `~alpha`]+(1/2)*Physics[Christoffel][omicron, zeta, `~beta`]*Physics[Christoffel][`~zeta`, `~alpha`, `~omicron`]-(1/2)*Physics[Christoffel][omicron, zeta, `~omicron`]*Physics[Christoffel][`~zeta`, `~alpha`, `~beta`]-(1/2)*Physics[Christoffel][`~tau`, `~alpha`, `~beta`]*Physics[Christoffel][`~upsilon`, tau, upsilon]-(1/2)*Physics[Christoffel][chi, iota, `~alpha`]*Physics[Christoffel][`~iota`, `~beta`, `~chi`]+(1/4)*(Physics[Christoffel][`~alpha`, chi, `~beta`]+Physics[Christoffel][`~beta`, chi, `~alpha`])*Physics[Christoffel][`~chi`, iota, `~iota`]+(1/2)*Physics[Christoffel][`~upsilon`, tau, `~beta`]*Physics[Christoffel][`~tau`, upsilon, `~alpha`]-(1/4)*Physics[Christoffel][`~beta`, chi, iota]*Physics[Christoffel][`~chi`, `~alpha`, `~iota`]+(1/2)*Physics[Christoffel][chi, `~alpha`, `~beta`]*Physics[Christoffel][iota, `~chi`, `~iota`]-(1/4)*Physics[Christoffel][`~alpha`, chi, iota]*Physics[Christoffel][`~chi`, `~beta`, `~iota`]-(1/4)*Physics[Christoffel][`~alpha`, rho, `~beta`]*Physics[Christoffel][`~rho`, rho1, `~rho1`]+(1/4)*Physics[Christoffel][`~alpha`, rho, `~rho1`]*Physics[Christoffel][`~rho`, rho1, `~beta`]+(1/2)*Physics[Christoffel][alpha10, `~alpha`, `~beta`]*Physics[Christoffel][alpha9, `~alpha10`, `~alpha9`]-(1/2)*Physics[Christoffel][alpha9, alpha10, `~alpha`]*Physics[Christoffel][`~alpha10`, `~alpha9`, `~beta`]+(1/4)*(2*Physics[Christoffel][chi, iota, kappa]*Physics[Christoffel][`~iota`, `~chi`, `~kappa`]-2*Physics[Christoffel][chi, iota, `~chi`]*Physics[Christoffel][`~iota`, kappa, `~kappa`]-2*Physics[Christoffel][alpha4, alpha5, alpha6]*Physics[Christoffel][`~alpha6`, `~alpha4`, `~alpha5`]+2*Physics[Christoffel][alpha4, alpha6, `~alpha5`]*Physics[Christoffel][`~alpha6`, alpha5, `~alpha4`]-2*D_[`~alpha5`](Physics[Christoffel][alpha4, alpha5, `~alpha4`], [X])+2*Physics[Christoffel][`~alpha1`, alpha1, alpha3]*Physics[Christoffel][`~alpha3`, alpha2, `~alpha2`]-2*Physics[Christoffel][`~alpha1`, alpha3, `~alpha2`]*Physics[Christoffel][`~alpha3`, alpha1, alpha2]+2*Physics[Ricci][alpha2, `~alpha2`]+2*D_[`~alpha2`](Physics[Christoffel][`~alpha1`, alpha1, alpha2], [X]))*Physics[g_][`~alpha`, `~beta`]-(1/4)*D_[`~rho1`](Physics[Christoffel][`~alpha`, rho1, `~beta`], [X])+(1/4)*D_[`~mu`](Physics[Christoffel][`~alpha`, mu, `~beta`], [X])+(1/4)*D_[`~nu`](Physics[Christoffel][`~beta`, nu, `~alpha`], [X])-(1/2)*D_[`~beta`](Physics[Christoffel][`~upsilon`, upsilon, `~alpha`], [X])-(1/4)*D_[`~omega`](Physics[Christoffel][`~beta`, omega, `~alpha`], [X])+(1/2)*D_[`~beta`](Physics[Christoffel][alpha9, `~alpha`, `~alpha9`], [X])-Physics[Ricci][`~alpha`, `~beta`] = 0)

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

    `Systems of spacetime coordinates are:`*{X = (x, y, z, t)}

     

    _______________________________________________________

     

    [coordinatesystems = {X}]

    (42)

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

    g_[]

    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 h[mu, nu], 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 delta[i, j](X) by delta[i, j])

    h[mu, nu] = Matrix(4, proc (i, j) options operator, arrow; delta[i, j](X) end proc, shape = symmetric)

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

    (44)

    In the above it is understood that `≪`(`δ__i,j`, 1), so that quadratic or higher powers of it or its derivatives can be approximated to 0 (discarded). Define the components of `h__μ,ν` accordingly

    "Define(?)"

    `Defined objects with tensor properties`

     

    {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] = rhs(?)"

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

    (46)

    "Define(?)"

    `Defined objects with tensor properties`

     

    {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

    g_[mu, nu] = eta[mu, nu]+h[mu, nu]

    Physics:-g_[mu, nu] = eta[mu, nu]+h[mu, nu]

    (48)

    Make this be the definition of the metric

    Define(Physics[g_][mu, nu] = eta[mu, nu]+h[mu, nu])

    _______________________________________________________

     

    `Coordinates: `[x, y, z, t]*`. Signature: `(`- - - +`)

     

    _______________________________________________________

     

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

     

    _______________________________________________________

     

    `Setting `*lowercaselatin_is*` letters to represent `*space*` indices`

     

    `Defined objects with tensor properties`

     

    {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

    Ricci[definition]

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

    (50)

    Christoffel[`~alpha`, mu, nu, definition]

    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)

    Substitute(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])), Physics[Ricci][mu, nu] = Physics[d_][alpha](Physics[Christoffel][`~alpha`, mu, nu], [X])-Physics[d_][nu](Physics[Christoffel][`~alpha`, mu, alpha], [X])+Physics[Christoffel][`~beta`, mu, nu]*Physics[Christoffel][`~alpha`, beta, alpha]-Physics[Christoffel][`~beta`, mu, alpha]*Physics[Christoffel][`~alpha`, nu, beta])

    Physics:-Ricci[mu, nu] = Physics:-d_[alpha]((1/2)*Physics:-g_[`~alpha`, `~kappa`]*(Physics:-d_[nu](Physics:-g_[kappa, mu], [X])+Physics:-d_[mu](Physics:-g_[kappa, nu], [X])-Physics:-d_[kappa](Physics:-g_[mu, nu], [X])), [X])-Physics:-d_[nu]((1/2)*Physics:-g_[`~alpha`, `~tau`]*(Physics:-d_[mu](Physics:-g_[tau, alpha], [X])+Physics:-d_[alpha](Physics:-g_[tau, mu], [X])-Physics:-d_[tau](Physics:-g_[alpha, mu], [X])), [X])+(1/4)*Physics:-g_[`~beta`, `~iota`]*(Physics:-d_[nu](Physics:-g_[iota, mu], [X])+Physics:-d_[mu](Physics:-g_[iota, nu], [X])-Physics:-d_[iota](Physics:-g_[mu, nu], [X]))*Physics:-g_[`~alpha`, `~lambda`]*(Physics:-d_[beta](Physics:-g_[lambda, alpha], [X])+Physics:-d_[alpha](Physics:-g_[lambda, beta], [X])-Physics:-d_[lambda](Physics:-g_[alpha, beta], [X]))-(1/4)*Physics:-g_[`~beta`, `~omega`]*(Physics:-d_[mu](Physics:-g_[omega, alpha], [X])+Physics:-d_[alpha](Physics:-g_[omega, mu], [X])-Physics:-d_[omega](Physics:-g_[alpha, mu], [X]))*Physics:-g_[`~alpha`, `~chi`]*(Physics:-d_[nu](Physics:-g_[chi, beta], [X])+Physics:-d_[beta](Physics:-g_[chi, nu], [X])-Physics:-d_[chi](Physics:-g_[beta, nu], [X]))

    (52)

    Introducing (48)g[mu, nu] = eta[mu, nu]+h[mu, nu], and also the inert form of the Ricci tensor to facilitate simplification some steps below,

    Substitute(Physics[g_][mu, nu] = eta[mu, nu]+h[mu, nu], Ricci = %Ricci, Physics[Ricci][mu, nu] = Physics[d_][alpha]((1/2)*Physics[g_][`~alpha`, `~kappa`]*(Physics[d_][nu](Physics[g_][kappa, mu], [X])+Physics[d_][mu](Physics[g_][kappa, nu], [X])-Physics[d_][kappa](Physics[g_][mu, nu], [X])), [X])-Physics[d_][nu]((1/2)*Physics[g_][`~alpha`, `~tau`]*(Physics[d_][mu](Physics[g_][tau, alpha], [X])+Physics[d_][alpha](Physics[g_][tau, mu], [X])-Physics[d_][tau](Physics[g_][alpha, mu], [X])), [X])+(1/4)*Physics[g_][`~beta`, `~iota`]*(Physics[d_][nu](Physics[g_][iota, mu], [X])+Physics[d_][mu](Physics[g_][iota, nu], [X])-Physics[d_][iota](Physics[g_][mu, nu], [X]))*Physics[g_][`~alpha`, `~lambda`]*(Physics[d_][beta](Physics[g_][lambda, alpha], [X])+Physics[d_][alpha](Physics[g_][lambda, beta], [X])-Physics[d_][lambda](Physics[g_][alpha, beta], [X]))-(1/4)*Physics[g_][`~beta`, `~omega`]*(Physics[d_][mu](Physics[g_][omega, alpha], [X])+Physics[d_][alpha](Physics[g_][omega, mu], [X])-Physics[d_][omega](Physics[g_][alpha, mu], [X]))*Physics[g_][`~alpha`, `~chi`]*(Physics[d_][nu](Physics[g_][chi, beta], [X])+Physics[d_][beta](Physics[g_][chi, nu], [X])-Physics[d_][chi](Physics[g_][beta, nu], [X])))

    %Ricci[mu, nu] = (1/2)*Physics:-d_[alpha](eta[`~alpha`, `~kappa`]+h[`~alpha`, `~kappa`], [X])*(Physics:-d_[nu](eta[kappa, mu]+h[kappa, mu], [X])+Physics:-d_[mu](eta[kappa, nu]+h[kappa, nu], [X])-Physics:-d_[kappa](eta[mu, nu]+h[mu, nu], [X]))+(1/2)*(eta[`~alpha`, `~kappa`]+h[`~alpha`, `~kappa`])*(Physics:-d_[alpha](Physics:-d_[nu](eta[kappa, mu]+h[kappa, mu], [X]), [X])+Physics:-d_[alpha](Physics:-d_[mu](eta[kappa, nu]+h[kappa, nu], [X]), [X])-Physics:-d_[alpha](Physics:-d_[kappa](eta[mu, nu]+h[mu, nu], [X]), [X]))-(1/2)*Physics:-d_[nu](eta[`~alpha`, `~tau`]+h[`~alpha`, `~tau`], [X])*(Physics:-d_[mu](eta[alpha, tau]+h[alpha, tau], [X])+Physics:-d_[alpha](eta[mu, tau]+h[mu, tau], [X])-Physics:-d_[tau](eta[alpha, mu]+h[alpha, mu], [X]))-(1/2)*(eta[`~alpha`, `~tau`]+h[`~alpha`, `~tau`])*(Physics:-d_[mu](Physics:-d_[nu](eta[alpha, tau]+h[alpha, tau], [X]), [X])+Physics:-d_[alpha](Physics:-d_[nu](eta[mu, tau]+h[mu, tau], [X]), [X])-Physics:-d_[nu](Physics:-d_[tau](eta[alpha, mu]+h[alpha, mu], [X]), [X]))+(1/4)*(eta[`~beta`, `~iota`]+h[`~beta`, `~iota`])*(Physics:-d_[nu](eta[iota, mu]+h[iota, mu], [X])+Physics:-d_[mu](eta[iota, nu]+h[iota, nu], [X])-Physics:-d_[iota](eta[mu, nu]+h[mu, nu], [X]))*(eta[`~alpha`, `~lambda`]+h[`~alpha`, `~lambda`])*(Physics:-d_[beta](eta[alpha, lambda]+h[alpha, lambda], [X])+Physics:-d_[alpha](eta[beta, lambda]+h[beta, lambda], [X])-Physics:-d_[lambda](eta[alpha, beta]+h[alpha, beta], [X]))-(1/4)*(eta[`~beta`, `~omega`]+h[`~beta`, `~omega`])*(Physics:-d_[mu](eta[alpha, omega]+h[alpha, omega], [X])+Physics:-d_[alpha](eta[mu, omega]+h[mu, omega], [X])-Physics:-d_[omega](eta[alpha, mu]+h[alpha, mu], [X]))*(eta[`~alpha`, `~chi`]+h[`~alpha`, `~chi`])*(Physics:-d_[nu](eta[beta, chi]+h[beta, chi], [X])+Physics:-d_[beta](eta[chi, nu]+h[chi, nu], [X])-Physics:-d_[chi](eta[beta, nu]+h[beta, nu], [X]))

    (53)

    This expression contains several terms quadratic in the small perturbation `h__μ,ν` 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

    Library:-Linearize(%Ricci[mu, nu] = (1/2)*Physics[d_][alpha](eta[`~alpha`, `~kappa`]+h[`~alpha`, `~kappa`], [X])*(Physics[d_][nu](eta[kappa, mu]+h[kappa, mu], [X])+Physics[d_][mu](eta[kappa, nu]+h[kappa, nu], [X])-Physics[d_][kappa](eta[mu, nu]+h[mu, nu], [X]))+(1/2)*(eta[`~alpha`, `~kappa`]+h[`~alpha`, `~kappa`])*(Physics[d_][alpha](Physics[d_][nu](eta[kappa, mu]+h[kappa, mu], [X]), [X])+Physics[d_][alpha](Physics[d_][mu](eta[kappa, nu]+h[kappa, nu], [X]), [X])-Physics[d_][alpha](Physics[d_][kappa](eta[mu, nu]+h[mu, nu], [X]), [X]))-(1/2)*Physics[d_][nu](eta[`~alpha`, `~tau`]+h[`~alpha`, `~tau`], [X])*(Physics[d_][mu](eta[alpha, tau]+h[alpha, tau], [X])+Physics[d_][alpha](eta[mu, tau]+h[mu, tau], [X])-Physics[d_][tau](eta[alpha, mu]+h[alpha, mu], [X]))-(1/2)*(eta[`~alpha`, `~tau`]+h[`~alpha`, `~tau`])*(Physics[d_][mu](Physics[d_][nu](eta[alpha, tau]+h[alpha, tau], [X]), [X])+Physics[d_][alpha](Physics[d_][nu](eta[mu, tau]+h[mu, tau], [X]), [X])-Physics[d_][nu](Physics[d_][tau](eta[alpha, mu]+h[alpha, mu], [X]), [X]))+(1/4)*(eta[`~beta`, `~iota`]+h[`~beta`, `~iota`])*(Physics[d_][nu](eta[iota, mu]+h[iota, mu], [X])+Physics[d_][mu](eta[iota, nu]+h[iota, nu], [X])-Physics[d_][iota](eta[mu, nu]+h[mu, nu], [X]))*(eta[`~alpha`, `~lambda`]+h[`~alpha`, `~lambda`])*(Physics[d_][beta](eta[alpha, lambda]+h[alpha, lambda], [X])+Physics[d_][alpha](eta[beta, lambda]+h[beta, lambda], [X])-Physics[d_][lambda](eta[alpha, beta]+h[alpha, beta], [X]))-(1/4)*(eta[`~beta`, `~omega`]+h[`~beta`, `~omega`])*(Physics[d_][mu](eta[alpha, omega]+h[alpha, omega], [X])+Physics[d_][alpha](eta[mu, omega]+h[mu, omega], [X])-Physics[d_][omega](eta[alpha, mu]+h[alpha, mu], [X]))*(eta[`~alpha`, `~chi`]+h[`~alpha`, `~chi`])*(Physics[d_][nu](eta[beta, chi]+h[beta, chi], [X])+Physics[d_][beta](eta[chi, nu]+h[chi, nu], [X])-Physics[d_][chi](eta[beta, nu]+h[beta, nu], [X])), h)

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

    (54)

    Important: in this result, `η__μ,ν` is the flat Minkowski metric, not the perturbed metric g[mu, nu]. However, in the context of a linearized formulation, `η__μ,ν` raises and lowers tensor indices the same way as g[mu, nu]. Hence, to further simplify contracted products of eta[mu, nu] in (54) , it is practical to reintroduce `g__μ,ν`representing that Minkowski metric and simplify using the internal algorithms for a flat metric

    g_[min]

    _______________________________________________________

     

    `The Minkowski metric in coordinates `*[x, y, z, t]

     

    `Signature: `(`- - - +`)

     

    _______________________________________________________

     

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

    (55)

    To proceed simplifying, replace in the expression (54) for the Ricci tensor the intermediate Minkowski `η__μ,ν`by `#msub(mi("g"),mrow(mi("μ",fontstyle = "normal"),mo(","),mi("ν",fontstyle = "normal")))`

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

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

    (56)

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

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

    %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

    Einstein[definition]

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

    (58)

    Compute the trace R[alpha, `~alpha`] directly from the linearized form (57) of the Ricci tensor,

    g_[mu, nu]*(%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]))

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

    Simplify(%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`])

    %Ricci[nu, `~nu`] = -Physics:-dAlembertian(h[alpha, `~alpha`], [X])+Physics:-d_[alpha](Physics:-d_[tau](h[`~alpha`, `~tau`], [X]), [X])

    (60)

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

    (%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]))-(1/2)*g_[mu, nu]*(%Ricci[nu, `~nu`] = -Physics[dAlembertian](h[alpha, `~alpha`], [X])+Physics[d_][alpha](Physics[d_][tau](h[`~alpha`, `~tau`], [X]), [X]))

    %Ricci[mu, nu]-(1/2)*Physics:-g_[mu, nu]*%Ricci[alpha, `~alpha`] = -(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])-(1/2)*Physics:-g_[mu, nu]*(-Physics:-dAlembertian(h[alpha, `~alpha`], [X])+Physics:-d_[alpha](Physics:-d_[tau](h[`~alpha`, `~tau`], [X]), [X]))

    (61)

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


    You can now redefine the general h[mu, nu] 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

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

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

    (62)

    Substitute(Physics[d_][mu](h[`~mu`, nu], [X]) = (1/2)*Physics[d_][nu](h[alpha, `~alpha`], [X]), %Ricci[mu, nu]-(1/2)*Physics[g_][mu, nu]*%Ricci[alpha, `~alpha`] = -(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])-(1/2)*Physics[g_][mu, nu]*(-Physics[dAlembertian](h[alpha, `~alpha`], [X])+Physics[d_][alpha](Physics[d_][tau](h[`~alpha`, `~tau`], [X]), [X])))

    %Ricci[mu, nu]-(1/2)*Physics:-g_[mu, nu]*%Ricci[alpha, `~alpha`] = -(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]((1/2)*Physics:-d_[mu](h[lambda, `~lambda`], [X]), [X])+(1/2)*Physics:-d_[mu]((1/2)*Physics:-d_[nu](h[kappa, `~kappa`], [X]), [X])-(1/2)*Physics:-g_[mu, nu]*(-Physics:-dAlembertian(h[alpha, `~alpha`], [X])+Physics:-d_[alpha]((1/2)*Physics:-d_[`~alpha`](h[beta, `~beta`], [X]), [X]))

    (63)

    Simplify(%Ricci[mu, nu]-(1/2)*Physics[g_][mu, nu]*%Ricci[alpha, `~alpha`] = -(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]((1/2)*Physics[d_][mu](h[lambda, `~lambda`], [X]), [X])+(1/2)*Physics[d_][mu]((1/2)*Physics[d_][nu](h[kappa, `~kappa`], [X]), [X])-(1/2)*Physics[g_][mu, nu]*(-Physics[dAlembertian](h[alpha, `~alpha`], [X])+Physics[d_][alpha]((1/2)*Physics[d_][`~alpha`](h[beta, `~beta`], [X]), [X])))

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

    _______________________________________________________

     

    `Systems of spacetime coordinates are:`*{X = (r, theta, phi, t)}

     

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

     

    `Setting `*lowercaselatin_is*` letters to represent `*space*` indices`

     

    `The Schwarzschild metric in coordinates `*[r, theta, phi, t]

     

    `Parameters: `[m]

     

    `Signature: `(`- - - +`)

     

    _______________________________________________________

     

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

    (65)

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

    Define(T[mu])

    `Defined objects with tensor properties`

     

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

    Define(R[mu], relativeweight = 1)

    `Defined objects with tensor properties`

     

    {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 r, theta, phi, t to  rho, theta, phi, t where

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

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

    (68)

    The transformed components of T[mu] and R[mu] are, respectively,

    TransformCoordinates(TR, T[mu], [rho, theta, phi, t])

    Vector[column](%id = 36893488151820263660)

    (69)

    TransformCoordinates(TR, R[mu], [rho, theta, phi, t])

    Vector[column](%id = 36893488151823155676)

    (70)

    where, when comparing both results, we see that the transformed components for R[mu] are all multiplied by J^n with n = 1 and J is the determinant of the transformation:

    J__matrix := simplify(VectorCalculus:-Jacobian([rhs(TR), theta, phi, t], [rho, theta, phi, t]))

    Matrix(%id = 36893488151964220700)

    (71)

    J = LinearAlgebra:-Determinant(J__matrix)

    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 T[mu] and R[mu] used above,

    Library:-GetRelativeWeight(T[mu])

    0

    (73)

    Library:-GetRelativeWeight(R[mu])

    1

    (74)

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

    Library:-GetRelativeWeight(R[`~mu`])

    1

    (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

    Library:-GetRelativeWeight(LeviCivita[alpha, beta, mu, nu])

    -1

    (76)

    "Library:-GetRelativeWeight(LeviCivita[~alpha, ~beta, ~mu,~nu])   "

    1

    (77)

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

    R[mu]^2

    R[mu]*R[`~mu`]

    (78)

    Library:-GetRelativeWeight(R[mu]*R[`~mu`])

    2

    (79)

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

    1/R[mu]^2

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

    (80)

    Library:-GetRelativeWeight(1/(R[mu]*R[`~mu`]))

    -2

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

    "R[~mu] + LeviCivita[~alpha, ~beta, ~mu,~nu]*T[alpha] T[beta] T[nu]"

    T[alpha]*T[beta]*T[nu]*Physics:-LeviCivita[`~alpha`, `~beta`, `~mu`, `~nu`]+R[`~mu`]

    (82)

    Library:-GetRelativeWeight(T[alpha]*T[beta]*T[nu]*Physics[LeviCivita][`~alpha`, `~beta`, `~mu`, `~nu`]+R[`~mu`])

    1

    (83)

    The relative weight of any determinant is always equal to 2

    %g_[determinant]

    %g_[determinant]

    (84)

    Library:-GetRelativeWeight(%g_[determinant])

    2

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

    g__det := %g_[:-determinant]

    %g_[determinant]

    (86)

    Library:-GetRelativeTerm(g__det, mu)

    -2*Physics:-Christoffel[`~nu`, mu, nu]*%g_[determinant]

    (87)
      

    Consequently,

    (%D_[mu] = D_[mu])(g__det)

    %D_[mu](%g_[determinant]) = 0

    (88)
      

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

    convert(%D_[mu](%g_[determinant]) = 0, d_)

    %d_[mu](%g_[determinant], [X])-2*Physics:-Christoffel[`~alpha`, alpha, mu]*%g_[determinant] = 0

    (89)
      

    evaluate the inert %d_

    factor(eval(%d_[mu](%g_[determinant], [X])-2*Physics[Christoffel][`~alpha`, alpha, mu]*%g_[determinant] = 0, %d_ = d_))

    %g_[determinant]*(Physics:-g_[`~alpha`, `~nu`]*Physics:-d_[mu](Physics:-g_[alpha, nu], [X])-2*Physics:-Christoffel[`~alpha`, alpha, mu]) = 0

    (90)
      

    The factor in parentheses is equal to `#mrow(msubsup(mi("g"),none(),mrow(mi("α",fontstyle = "normal"),mo(","),mi("ν",fontstyle = "normal"))),msub(mi("▿"),mi("μ",fontstyle = "normal")),mo("⁡"),mfenced(msub(mi("g"),mrow(mi("α",fontstyle = "normal"),mo(","),mi("ν",fontstyle = "normal")))))`, where the covariant derivative of the metric is equal to zero, so

    Simplify(%g_[determinant]*(Physics[g_][`~alpha`, `~nu`]*Physics[d_][mu](Physics[g_][alpha, nu], [X])-2*Physics[Christoffel][`~alpha`, alpha, mu]) = 0)

    0 = 0

    (91)
      

    Consider the covariant derivative of T[mu] and R[mu] defined in (66) and (67)

    Library:-GetRelativeWeight(T[mu])

    0

    (92)

    Library:-GetRelativeWeight(R[mu])

    1

    (93)

    The corresponding covariant derivatives

    (%D_[mu] = D_[mu])(T[mu](X))

    %D_[mu](T[`~mu`](X)) = Physics:-D_[mu](T[`~mu`](X), [X])

    (94)

    expand(%D_[mu](T[`~mu`](X)) = D_[mu](T[`~mu`](X), [X]))

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

    (%D_[mu] = D_[mu])(R[mu](X))

    %D_[mu](R[`~mu`](X)) = Physics:-D_[mu](R[`~mu`](X), [X])

    (96)

    expand(%D_[mu](R[`~mu`](X)) = D_[mu](R[`~mu`](X), [X]))

    %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

    Library:-GetRelativeTerm(R[`~mu`](X), mu)

    -Physics:-Christoffel[`~nu`, mu, nu]*R[`~mu`](X)

    (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

    Christoffel[mu, alpha, beta] = Library:-ConvertToF(Christoffel[mu, alpha, beta], g_)

    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 A[mu] representing the 4D electromagnetic potential as a function of the coordinates X and F[mu, nu] representing the electromagnetic field tensors

    Define(A[mu] = A[mu](X), quiet)

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

    Define(F[mu, nu] = d_[mu](A[nu])-d_[nu](A[mu]))

    `Defined objects with tensor properties`

     

    {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 A[mu]

    F[mu, nu] = Library:-ConvertToF(F[mu, nu], A)

    F[mu, nu] = Physics:-d_[mu](A[nu], [X])-Physics:-d_[nu](A[mu], [X])

    (103)

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

    F[definition]

    F[mu, nu] = Physics:-d_[mu](A[nu], [X])-Physics:-d_[nu](A[mu], [X])

    (104)

    The rewriting, however, works also with tensorial expressions

    F[mu, nu]*A[mu]*A[nu]

    F[mu, nu]*A[`~mu`]*A[`~nu`]

    (105)

    Library:-ConvertToF(F[mu, nu]*A[`~mu`]*A[`~nu`], A)

    (Physics:-d_[mu](A[nu], [X])-Physics:-d_[nu](A[mu], [X]))*A[`~mu`]*A[`~nu`]

    (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

    NULL


     

    Download New_in_Physics_2025.mw

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

    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 g[mu, nu]

    with(Physics); Setup(coordinates = cartesian, metric = arbitrary)

    `Systems of spacetime coordinates are:`*{X = (x, y, z, t)}

     

    _______________________________________________________

     

    `Setting `*lowercaselatin_is*` letters to represent `*space*` indices`

     

    `The arbitrary metric in coordinates `*[x, y, z, t]

     

    `Signature: `(`- - - +`)

     

    _______________________________________________________

     

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

     

    _______________________________________________________; "_noterminate"

    (1)

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

    CompactDisplay(g_[])

    f__1(x, y, z, t)*`will now be displayed as`*f__1

     

    f__10(x, y, z, t)*`will now be displayed as`*f__10

     

    f__2(x, y, z, t)*`will now be displayed as`*f__2

     

    f__3(x, y, z, t)*`will now be displayed as`*f__3

     

    f__4(x, y, z, t)*`will now be displayed as`*f__4

     

    f__5(x, y, z, t)*`will now be displayed as`*f__5

     

    f__6(x, y, z, t)*`will now be displayed as`*f__6

     

    f__7(x, y, z, t)*`will now be displayed as`*f__7

     

    f__8(x, y, z, t)*`will now be displayed as`*f__8

     

    f__9(x, y, z, t)*`will now be displayed as`*f__9

    (2)

    and the parameters of the variational problem are the spacetime coordinates X^alpha.

     

    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

    L := sqrt(-%g_[determinant])*Ricci[alpha, `~alpha`]

    (-%g_[determinant])^(1/2)*Physics:-Ricci[alpha, `~alpha`]

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

    LagrangeEquations(L, g_[mu, nu])

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

    Einstein[definition]

    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 L 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

    L__C := convert(L, Christoffel)

    (-%g_[determinant])^(1/2)*Physics:-g_[`~alpha`, `~lambda`]*(Physics:-d_[nu](Physics:-Christoffel[`~nu`, alpha, lambda], [X])-Physics:-d_[lambda](Physics:-Christoffel[`~nu`, alpha, nu], [X])+Physics:-Christoffel[`~beta`, alpha, lambda]*Physics:-Christoffel[`~nu`, beta, nu]-Physics:-Christoffel[`~beta`, alpha, nu]*Physics:-Christoffel[`~nu`, beta, lambda])

    (6)

    Recalling the definition

    Christoffel[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 L[C] the two terms containing derivatives of Christoffel symbols contain second order derivatives of g[mu, nu].

     

    Now, it is always possible to add a total spacetime derivative to L[C] 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 L[C], 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

    TD := %d_[alpha](g_[`~mu`, `~nu`]*sqrt(-%g_[determinant])*(g_[`~alpha`, mu]*Christoffel[`~beta`, nu, beta]-Christoffel[`~alpha`, mu, nu]))

    %d_[alpha](Physics:-g_[`~mu`, `~nu`]*(-%g_[determinant])^(1/2)*(Physics:-g_[mu, `~alpha`]*Physics:-Christoffel[`~beta`, beta, nu]-Physics:-Christoffel[`~alpha`, mu, nu]))

    (8)

    Adding this term to L[C], performing the `∂` differentiation operation and simplifying we get

    L__1 := L__C+TD

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

    (9)

    L__1 := eval(L__1, %d_ = d_)

    (-%g_[determinant])^(1/2)*Physics:-g_[`~alpha`, `~lambda`]*(Physics:-d_[nu](Physics:-Christoffel[`~nu`, alpha, lambda], [X])-Physics:-d_[lambda](Physics:-Christoffel[`~nu`, alpha, nu], [X])+Physics:-Christoffel[`~beta`, alpha, lambda]*Physics:-Christoffel[`~nu`, beta, nu]-Physics:-Christoffel[`~beta`, alpha, nu]*Physics:-Christoffel[`~nu`, beta, lambda])+Physics:-d_[alpha](Physics:-g_[`~mu`, `~nu`], [X])*(-%g_[determinant])^(1/2)*(Physics:-g_[mu, `~alpha`]*Physics:-Christoffel[`~beta`, beta, nu]-Physics:-Christoffel[`~alpha`, mu, nu])-(1/2)*Physics:-g_[`~mu`, `~nu`]*(Physics:-g_[mu, `~alpha`]*Physics:-Christoffel[`~beta`, beta, nu]-Physics:-Christoffel[`~alpha`, mu, nu])*%g_[determinant]*Physics:-g_[`~kappa`, `~lambda`]*Physics:-d_[alpha](Physics:-g_[kappa, lambda], [X])/(-%g_[determinant])^(1/2)+Physics:-g_[`~mu`, `~nu`]*(-%g_[determinant])^(1/2)*(Physics:-g_[mu, `~alpha`]*Physics:-d_[alpha](Physics:-Christoffel[`~beta`, beta, nu], [X])-Physics:-d_[alpha](Physics:-Christoffel[`~alpha`, mu, nu], [X]))

    (10)

    L__1 := Simplify(L__1)

    (Physics:-Christoffel[alpha, beta, kappa]*Physics:-Christoffel[`~beta`, `~alpha`, `~kappa`]-Physics:-Christoffel[alpha, beta, `~alpha`]*Physics:-Christoffel[`~beta`, kappa, `~kappa`])*(-%g_[determinant])^(1/2)

    (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

    LagrangeEquations(L__1, g_[mu, nu])

    -(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 `≡`(L[1], (Physics[Christoffel][alpha, beta, kappa]*Physics[Christoffel][`~beta`, `~alpha`, `~kappa`]-Physics[Christoffel][alpha, beta, `~alpha`]*Physics[Christoffel][`~beta`, kappa, `~kappa`])*(-%g_[determinant])^(1/2)) is no just L[C] ≡ (6) after discarding its two terms involving derivatives of Christoffel symbols. To verify this, split L[C] into the terms containing or not derivatives of Christoffel

    L__22, L__11 := selectremove(has, expand(L__C), d_)

    (-%g_[determinant])^(1/2)*Physics:-d_[nu](Physics:-Christoffel[`~nu`, alpha, lambda], [X])*Physics:-g_[`~alpha`, `~lambda`]-(-%g_[determinant])^(1/2)*Physics:-d_[lambda](Physics:-Christoffel[`~nu`, alpha, nu], [X])*Physics:-g_[`~alpha`, `~lambda`], (-%g_[determinant])^(1/2)*Physics:-Christoffel[`~beta`, alpha, lambda]*Physics:-Christoffel[`~nu`, beta, nu]*Physics:-g_[`~alpha`, `~lambda`]-(-%g_[determinant])^(1/2)*Physics:-Christoffel[`~beta`, alpha, nu]*Physics:-Christoffel[`~nu`, beta, lambda]*Physics:-g_[`~alpha`, `~lambda`]

    (13)

    Comparing, the total derivative TD≡ (8) is not just -L[22], but

    TD = -L__22-2*L__11

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

    (14)

    Things like these, TD = -L__22-2*L__11, 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

    (lhs-rhs)(%d_[alpha](Physics[g_][`~mu`, `~nu`]*(-%g_[determinant])^(1/2)*(Physics[g_][mu, `~alpha`]*Physics[Christoffel][`~beta`, beta, nu]-Physics[Christoffel][`~alpha`, mu, nu])) = -(-%g_[determinant])^(1/2)*Physics[d_][nu](Physics[Christoffel][`~nu`, alpha, lambda], [X])*Physics[g_][`~alpha`, `~lambda`]+(-%g_[determinant])^(1/2)*Physics[d_][lambda](Physics[Christoffel][`~nu`, alpha, nu], [X])*Physics[g_][`~alpha`, `~lambda`]-2*(-%g_[determinant])^(1/2)*Physics[Christoffel][`~beta`, alpha, lambda]*Physics[Christoffel][`~nu`, beta, nu]*Physics[g_][`~alpha`, `~lambda`]+2*(-%g_[determinant])^(1/2)*Physics[Christoffel][`~beta`, alpha, nu]*Physics[Christoffel][`~nu`, beta, lambda]*Physics[g_][`~alpha`, `~lambda`])

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

    (15)

    eval(%d_[alpha](Physics[g_][`~mu`, `~nu`]*(-%g_[determinant])^(1/2)*(Physics[g_][mu, `~alpha`]*Physics[Christoffel][`~beta`, beta, nu]-Physics[Christoffel][`~alpha`, mu, nu]))+(-%g_[determinant])^(1/2)*Physics[d_][nu](Physics[Christoffel][`~nu`, alpha, lambda], [X])*Physics[g_][`~alpha`, `~lambda`]-(-%g_[determinant])^(1/2)*Physics[d_][lambda](Physics[Christoffel][`~nu`, alpha, nu], [X])*Physics[g_][`~alpha`, `~lambda`]+2*(-%g_[determinant])^(1/2)*Physics[Christoffel][`~beta`, alpha, lambda]*Physics[Christoffel][`~nu`, beta, nu]*Physics[g_][`~alpha`, `~lambda`]-2*(-%g_[determinant])^(1/2)*Physics[Christoffel][`~beta`, alpha, nu]*Physics[Christoffel][`~nu`, beta, lambda]*Physics[g_][`~alpha`, `~lambda`], %d_ = d_)

    (-%g_[determinant])^(1/2)*(Physics:-g_[mu, `~alpha`]*Physics:-Christoffel[`~beta`, beta, nu]-Physics:-Christoffel[`~alpha`, mu, nu])*Physics:-d_[alpha](Physics:-g_[`~mu`, `~nu`], [X])-(1/2)*(Physics:-g_[mu, `~alpha`]*Physics:-Christoffel[`~beta`, beta, nu]-Physics:-Christoffel[`~alpha`, mu, nu])*Physics:-g_[`~mu`, `~nu`]*%g_[determinant]*Physics:-g_[`~kappa`, `~lambda`]*Physics:-d_[alpha](Physics:-g_[kappa, lambda], [X])/(-%g_[determinant])^(1/2)+(-%g_[determinant])^(1/2)*(Physics:-g_[mu, `~alpha`]*Physics:-d_[alpha](Physics:-Christoffel[`~beta`, beta, nu], [X])-Physics:-d_[alpha](Physics:-Christoffel[`~alpha`, mu, nu], [X]))*Physics:-g_[`~mu`, `~nu`]+(-%g_[determinant])^(1/2)*Physics:-d_[nu](Physics:-Christoffel[`~nu`, alpha, lambda], [X])*Physics:-g_[`~alpha`, `~lambda`]-(-%g_[determinant])^(1/2)*Physics:-d_[lambda](Physics:-Christoffel[`~nu`, alpha, nu], [X])*Physics:-g_[`~alpha`, `~lambda`]+2*(-%g_[determinant])^(1/2)*Physics:-Christoffel[`~beta`, alpha, lambda]*Physics:-Christoffel[`~nu`, beta, nu]*Physics:-g_[`~alpha`, `~lambda`]-2*(-%g_[determinant])^(1/2)*Physics:-Christoffel[`~beta`, alpha, nu]*Physics:-Christoffel[`~nu`, beta, lambda]*Physics:-g_[`~alpha`, `~lambda`]

    (16)

    Simplify((-%g_[determinant])^(1/2)*(Physics[g_][mu, `~alpha`]*Physics[Christoffel][`~beta`, beta, nu]-Physics[Christoffel][`~alpha`, mu, nu])*Physics[d_][alpha](Physics[g_][`~mu`, `~nu`], [X])-(1/2)*(Physics[g_][mu, `~alpha`]*Physics[Christoffel][`~beta`, beta, nu]-Physics[Christoffel][`~alpha`, mu, nu])*Physics[g_][`~mu`, `~nu`]*%g_[determinant]*Physics[g_][`~kappa`, `~lambda`]*Physics[d_][alpha](Physics[g_][kappa, lambda], [X])/(-%g_[determinant])^(1/2)+(-%g_[determinant])^(1/2)*(Physics[g_][mu, `~alpha`]*Physics[d_][alpha](Physics[Christoffel][`~beta`, beta, nu], [X])-Physics[d_][alpha](Physics[Christoffel][`~alpha`, mu, nu], [X]))*Physics[g_][`~mu`, `~nu`]+(-%g_[determinant])^(1/2)*Physics[d_][nu](Physics[Christoffel][`~nu`, alpha, lambda], [X])*Physics[g_][`~alpha`, `~lambda`]-(-%g_[determinant])^(1/2)*Physics[d_][lambda](Physics[Christoffel][`~nu`, alpha, nu], [X])*Physics[g_][`~alpha`, `~lambda`]+2*(-%g_[determinant])^(1/2)*Physics[Christoffel][`~beta`, alpha, lambda]*Physics[Christoffel][`~nu`, beta, nu]*Physics[g_][`~alpha`, `~lambda`]-2*(-%g_[determinant])^(1/2)*Physics[Christoffel][`~beta`, alpha, nu]*Physics[Christoffel][`~nu`, beta, lambda]*Physics[g_][`~alpha`, `~lambda`])

    0

    (17)

    That said, it is also true that TD = -L[22]-2*L[11] results in the Lagrangian L[1] = -L[11], and since the equations of movement don't depend on the sign of the Lagrangian, for this Lagrangian `≡`(L[C], (-%g_[determinant])^(1/2)*Physics[g_][`~alpha`, `~lambda`]*(Physics[d_][nu](Physics[Christoffel][`~nu`, alpha, lambda], [X])-Physics[d_][lambda](Physics[Christoffel][`~nu`, alpha, nu], [X])+Physics[Christoffel][`~beta`, alpha, lambda]*Physics[Christoffel][`~nu`, beta, nu]-Physics[Christoffel][`~beta`, alpha, nu]*Physics[Christoffel][`~nu`, beta, lambda])) adding the term TD happens to be equivalent to just discarding the terms of L__C 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,

    S := Intc(L, X)

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

    EE__unsimplified := Fundiff(S, g_[alpha, beta]) = 0

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

    EEC := Simplify(EE__unsimplified)

    (1/4)*(2*Physics:-Christoffel[chi, iota, kappa]*Physics:-Christoffel[`~iota`, `~chi`, `~kappa`]-2*Physics:-Christoffel[chi, iota, `~chi`]*Physics:-Christoffel[`~iota`, kappa, `~kappa`]-2*Physics:-D_[iota](Physics:-Christoffel[chi, `~chi`, `~iota`], [X])+2*Physics:-D_[chi](Physics:-Christoffel[`~chi`, iota, `~iota`], [X]))*Physics:-g_[`~alpha`, `~beta`]+(1/4)*(Physics:-Christoffel[`~alpha`, chi, `~beta`]+Physics:-Christoffel[`~beta`, chi, `~alpha`])*Physics:-Christoffel[`~chi`, iota, `~iota`]-(1/4)*Physics:-Christoffel[`~beta`, chi, iota]*Physics:-Christoffel[`~chi`, `~alpha`, `~iota`]-(1/2)*Physics:-Christoffel[chi, iota, `~alpha`]*Physics:-Christoffel[`~iota`, `~beta`, `~chi`]+(1/2)*Physics:-Christoffel[chi, `~alpha`, `~beta`]*Physics:-Christoffel[iota, `~chi`, `~iota`]-(1/4)*Physics:-Christoffel[`~alpha`, chi, iota]*Physics:-Christoffel[`~chi`, `~beta`, `~iota`]+(1/4)*Physics:-D_[chi](Physics:-Christoffel[`~alpha`, `~beta`, `~chi`], [X])+(1/4)*Physics:-D_[chi](Physics:-Christoffel[`~beta`, `~alpha`, `~chi`], [X])-(1/2)*Physics:-D_[chi](Physics:-Christoffel[`~chi`, `~alpha`, `~beta`], [X])-(1/4)*Physics:-D_[`~alpha`](Physics:-Christoffel[`~beta`, chi, `~chi`], [X])+(1/2)*Physics:-D_[`~beta`](Physics:-Christoffel[chi, `~alpha`, `~chi`], [X])-(1/4)*Physics:-D_[`~beta`](Physics:-Christoffel[`~alpha`, chi, `~chi`], [X]) = 0

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

    "`▿`[chi](GAMMA[]^(alpha,beta,chi)) :"

    % = expand(%)

    Physics:-D_[chi](Physics:-Christoffel[`~alpha`, `~beta`, `~chi`], [X]) = Physics:-d_[chi](Physics:-Christoffel[`~alpha`, `~beta`, `~chi`], [X])+Physics:-Christoffel[`~alpha`, chi, mu]*Physics:-Christoffel[`~mu`, `~beta`, `~chi`]+Physics:-Christoffel[`~beta`, chi, mu]*Physics:-Christoffel[`~alpha`, `~chi`, `~mu`]+Physics:-Christoffel[`~chi`, chi, mu]*Physics:-Christoffel[`~alpha`, `~beta`, `~mu`]

    (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

    Ricci[definition]

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

    (21)

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

    convert(Physics[Ricci][mu, nu] = Physics[d_][alpha](Physics[Christoffel][`~alpha`, mu, nu], [X])-Physics[d_][nu](Physics[Christoffel][`~alpha`, mu, alpha], [X])+Physics[Christoffel][`~beta`, mu, nu]*Physics[Christoffel][`~alpha`, beta, alpha]-Physics[Christoffel][`~beta`, mu, alpha]*Physics[Christoffel][`~alpha`, nu, beta], D_)

    Physics:-Ricci[mu, nu] = Physics:-D_[alpha](Physics:-Christoffel[`~alpha`, mu, nu], [X])-Physics:-Christoffel[`~alpha`, alpha, kappa]*Physics:-Christoffel[`~kappa`, mu, nu]+Physics:-Christoffel[`~kappa`, alpha, mu]*Physics:-Christoffel[`~alpha`, kappa, nu]+Physics:-Christoffel[`~kappa`, alpha, nu]*Physics:-Christoffel[`~alpha`, mu, kappa]-Physics:-D_[nu](Physics:-Christoffel[`~alpha`, alpha, mu], [X])-Physics:-Christoffel[`~lambda`, mu, nu]*Physics:-Christoffel[`~alpha`, alpha, lambda]+Physics:-Christoffel[`~beta`, mu, nu]*Physics:-Christoffel[`~alpha`, alpha, beta]-Physics:-Christoffel[`~beta`, alpha, mu]*Physics:-Christoffel[`~alpha`, beta, nu]

    (22)

    Simplify(Physics[Ricci][mu, nu] = D_[alpha](Physics[Christoffel][`~alpha`, mu, nu], [X])-Physics[Christoffel][`~alpha`, alpha, kappa]*Physics[Christoffel][`~kappa`, mu, nu]+Physics[Christoffel][`~kappa`, alpha, mu]*Physics[Christoffel][`~alpha`, kappa, nu]+Physics[Christoffel][`~kappa`, alpha, nu]*Physics[Christoffel][`~alpha`, mu, kappa]-D_[nu](Physics[Christoffel][`~alpha`, alpha, mu], [X])-Physics[Christoffel][`~lambda`, mu, nu]*Physics[Christoffel][`~alpha`, alpha, lambda]+Physics[Christoffel][`~beta`, mu, nu]*Physics[Christoffel][`~alpha`, alpha, beta]-Physics[Christoffel][`~beta`, alpha, mu]*Physics[Christoffel][`~alpha`, beta, nu])

    Physics:-Ricci[mu, nu] = Physics:-Christoffel[alpha, beta, mu]*Physics:-Christoffel[`~beta`, nu, `~alpha`]-Physics:-Christoffel[beta, mu, nu]*Physics:-Christoffel[alpha, `~alpha`, `~beta`]+Physics:-D_[alpha](Physics:-Christoffel[`~alpha`, mu, nu], [X])-Physics:-D_[nu](Physics:-Christoffel[alpha, mu, `~alpha`], [X])

    (23)

    C_to_Ricci := isolate(Physics[Ricci][mu, nu] = Physics[Christoffel][alpha, beta, mu]*Physics[Christoffel][`~beta`, nu, `~alpha`]-Physics[Christoffel][beta, mu, nu]*Physics[Christoffel][alpha, `~alpha`, `~beta`]+D_[alpha](Physics[Christoffel][`~alpha`, mu, nu], [X])-D_[nu](Physics[Christoffel][alpha, mu, `~alpha`], [X]), D_[alpha](Christoffel[`~alpha`, mu, nu]))

    Physics:-D_[alpha](Physics:-Christoffel[`~alpha`, mu, nu], [X]) = -Physics:-Christoffel[alpha, beta, mu]*Physics:-Christoffel[`~beta`, nu, `~alpha`]+Physics:-Christoffel[beta, mu, nu]*Physics:-Christoffel[alpha, `~alpha`, `~beta`]+Physics:-Ricci[mu, nu]+Physics:-D_[nu](Physics:-Christoffel[alpha, mu, `~alpha`], [X])

    (24)

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

    Riemann[`~alpha`, beta, mu, nu, definition]

    Physics:-Riemann[`~alpha`, beta, mu, nu] = Physics:-d_[mu](Physics:-Christoffel[`~alpha`, beta, nu], [X])-Physics:-d_[nu](Physics:-Christoffel[`~alpha`, beta, mu], [X])+Physics:-Christoffel[`~alpha`, upsilon, mu]*Physics:-Christoffel[`~upsilon`, beta, nu]-Physics:-Christoffel[`~alpha`, upsilon, nu]*Physics:-Christoffel[`~upsilon`, beta, mu]

    (25)

    convert(Physics[Riemann][`~alpha`, beta, mu, nu] = Physics[d_][mu](Physics[Christoffel][`~alpha`, beta, nu], [X])-Physics[d_][nu](Physics[Christoffel][`~alpha`, beta, mu], [X])+Physics[Christoffel][`~alpha`, upsilon, mu]*Physics[Christoffel][`~upsilon`, beta, nu]-Physics[Christoffel][`~alpha`, upsilon, nu]*Physics[Christoffel][`~upsilon`, beta, mu], D_)

    Physics:-Riemann[`~alpha`, beta, mu, nu] = Physics:-D_[mu](Physics:-Christoffel[`~alpha`, beta, nu], [X])+Physics:-Christoffel[`~kappa`, mu, nu]*Physics:-Christoffel[`~alpha`, beta, kappa]-Physics:-Christoffel[`~alpha`, kappa, mu]*Physics:-Christoffel[`~kappa`, beta, nu]+Physics:-Christoffel[`~kappa`, beta, mu]*Physics:-Christoffel[`~alpha`, kappa, nu]-Physics:-D_[nu](Physics:-Christoffel[`~alpha`, beta, mu], [X])-Physics:-Christoffel[`~lambda`, mu, nu]*Physics:-Christoffel[`~alpha`, beta, lambda]-Physics:-Christoffel[`~lambda`, beta, nu]*Physics:-Christoffel[`~alpha`, lambda, mu]+Physics:-Christoffel[`~alpha`, lambda, nu]*Physics:-Christoffel[`~lambda`, beta, mu]+Physics:-Christoffel[`~alpha`, mu, upsilon]*Physics:-Christoffel[`~upsilon`, beta, nu]-Physics:-Christoffel[`~alpha`, nu, upsilon]*Physics:-Christoffel[`~upsilon`, beta, mu]

    (26)

    Simplify(Physics[Riemann][`~alpha`, beta, mu, nu] = D_[mu](Physics[Christoffel][`~alpha`, beta, nu], [X])+Physics[Christoffel][`~kappa`, mu, nu]*Physics[Christoffel][`~alpha`, beta, kappa]-Physics[Christoffel][`~alpha`, kappa, mu]*Physics[Christoffel][`~kappa`, beta, nu]+Physics[Christoffel][`~kappa`, beta, mu]*Physics[Christoffel][`~alpha`, kappa, nu]-D_[nu](Physics[Christoffel][`~alpha`, beta, mu], [X])-Physics[Christoffel][`~lambda`, mu, nu]*Physics[Christoffel][`~alpha`, beta, lambda]-Physics[Christoffel][`~lambda`, beta, nu]*Physics[Christoffel][`~alpha`, lambda, mu]+Physics[Christoffel][`~alpha`, lambda, nu]*Physics[Christoffel][`~lambda`, beta, mu]+Physics[Christoffel][`~alpha`, mu, upsilon]*Physics[Christoffel][`~upsilon`, beta, nu]-Physics[Christoffel][`~alpha`, nu, upsilon]*Physics[Christoffel][`~upsilon`, beta, mu])

    Physics:-Riemann[`~alpha`, beta, mu, nu] = -Physics:-Christoffel[`~alpha`, kappa, mu]*Physics:-Christoffel[`~kappa`, beta, nu]+Physics:-Christoffel[`~kappa`, beta, mu]*Physics:-Christoffel[`~alpha`, kappa, nu]+Physics:-D_[mu](Physics:-Christoffel[`~alpha`, beta, nu], [X])-Physics:-D_[nu](Physics:-Christoffel[`~alpha`, beta, mu], [X])

    (27)

    C_to_Riemann := isolate(Physics[Riemann][`~alpha`, beta, mu, nu] = -Physics[Christoffel][`~alpha`, kappa, mu]*Physics[Christoffel][`~kappa`, beta, nu]+Physics[Christoffel][`~kappa`, beta, mu]*Physics[Christoffel][`~alpha`, kappa, nu]+D_[mu](Physics[Christoffel][`~alpha`, beta, nu], [X])-D_[nu](Physics[Christoffel][`~alpha`, beta, mu], [X]), D_[mu](Christoffel[`~alpha`, beta, nu]))

    Physics:-D_[mu](Physics:-Christoffel[`~alpha`, beta, nu], [X]) = Physics:-Christoffel[`~alpha`, kappa, mu]*Physics:-Christoffel[`~kappa`, beta, nu]-Physics:-Christoffel[`~kappa`, beta, mu]*Physics:-Christoffel[`~alpha`, kappa, nu]+Physics:-Riemann[`~alpha`, beta, mu, nu]+Physics:-D_[nu](Physics:-Christoffel[`~alpha`, beta, mu], [X])

    (28)

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

    Substitute(C_to_Riemann, C_to_Ricci, EEC)

    (1/4)*(2*Physics:-Christoffel[chi, iota, kappa]*Physics:-Christoffel[`~iota`, `~chi`, `~kappa`]-2*Physics:-Christoffel[chi, iota, `~chi`]*Physics:-Christoffel[`~iota`, kappa, `~kappa`]-2*Physics:-Christoffel[alpha4, alpha5, alpha6]*Physics:-Christoffel[`~alpha6`, `~alpha4`, `~alpha5`]+2*Physics:-Christoffel[alpha4, alpha6, `~alpha5`]*Physics:-Christoffel[`~alpha6`, alpha5, `~alpha4`]-2*Physics:-D_[`~alpha5`](Physics:-Christoffel[alpha4, alpha5, `~alpha4`], [X])+2*Physics:-Christoffel[`~alpha1`, alpha1, alpha3]*Physics:-Christoffel[`~alpha3`, alpha2, `~alpha2`]-2*Physics:-Christoffel[`~alpha1`, alpha3, `~alpha2`]*Physics:-Christoffel[`~alpha3`, alpha1, alpha2]+2*Physics:-Ricci[alpha2, `~alpha2`]+2*Physics:-D_[`~alpha2`](Physics:-Christoffel[`~alpha1`, alpha1, alpha2], [X]))*Physics:-g_[`~alpha`, `~beta`]-Physics:-Ricci[`~alpha`, `~beta`]+(1/2)*Physics:-Christoffel[alpha10, `~alpha`, `~beta`]*Physics:-Christoffel[alpha9, `~alpha10`, `~alpha9`]-(1/2)*Physics:-Christoffel[alpha9, alpha10, `~alpha`]*Physics:-Christoffel[`~alpha10`, `~alpha9`, `~beta`]-(1/4)*Physics:-Christoffel[`~alpha`, chi, iota]*Physics:-Christoffel[`~chi`, `~beta`, `~iota`]-(1/4)*Physics:-Christoffel[`~beta`, chi, iota]*Physics:-Christoffel[`~chi`, `~alpha`, `~iota`]+(1/2)*Physics:-Christoffel[chi, `~alpha`, `~beta`]*Physics:-Christoffel[iota, `~chi`, `~iota`]+(1/4)*(Physics:-Christoffel[`~alpha`, chi, `~beta`]+Physics:-Christoffel[`~beta`, chi, `~alpha`])*Physics:-Christoffel[`~chi`, iota, `~iota`]-(1/2)*Physics:-Christoffel[chi, iota, `~alpha`]*Physics:-Christoffel[`~iota`, `~beta`, `~chi`]-(1/4)*Physics:-Christoffel[`~alpha`, rho, `~beta`]*Physics:-Christoffel[`~rho`, rho1, `~rho1`]+(1/4)*Physics:-Christoffel[`~alpha`, rho, `~rho1`]*Physics:-Christoffel[`~rho`, rho1, `~beta`]-(1/2)*Physics:-Christoffel[omicron, zeta, `~omicron`]*Physics:-Christoffel[`~zeta`, `~alpha`, `~beta`]+(1/2)*Physics:-Christoffel[omicron, zeta, `~beta`]*Physics:-Christoffel[`~zeta`, `~alpha`, `~omicron`]-(1/4)*Physics:-Christoffel[`~beta`, psi, `~alpha`]*Physics:-Christoffel[`~psi`, omega, `~omega`]+(1/4)*Physics:-Christoffel[`~beta`, psi, `~omega`]*Physics:-Christoffel[`~psi`, omega, `~alpha`]+(1/2)*Physics:-Christoffel[`~tau`, upsilon, `~alpha`]*Physics:-Christoffel[`~upsilon`, tau, `~beta`]-(1/2)*Physics:-Christoffel[`~tau`, `~alpha`, `~beta`]*Physics:-Christoffel[`~upsilon`, tau, upsilon]+(1/4)*Physics:-Christoffel[`~beta`, nu, sigma]*Physics:-Christoffel[`~sigma`, `~alpha`, `~nu`]-(1/4)*Physics:-Christoffel[`~beta`, sigma, `~nu`]*Physics:-Christoffel[`~sigma`, nu, `~alpha`]-(1/4)*Physics:-Christoffel[`~alpha`, lambda, `~mu`]*Physics:-Christoffel[`~lambda`, mu, `~beta`]+(1/4)*Physics:-Christoffel[`~alpha`, lambda, mu]*Physics:-Christoffel[`~lambda`, `~beta`, `~mu`]-(1/2)*Physics:-D_[`~beta`](Physics:-Christoffel[`~upsilon`, upsilon, `~alpha`], [X])+(1/4)*Physics:-D_[`~nu`](Physics:-Christoffel[`~beta`, nu, `~alpha`], [X])+(1/4)*Physics:-D_[`~mu`](Physics:-Christoffel[`~alpha`, mu, `~beta`], [X])-(1/4)*Physics:-D_[`~omega`](Physics:-Christoffel[`~beta`, omega, `~alpha`], [X])+(1/2)*Physics:-D_[`~beta`](Physics:-Christoffel[alpha9, `~alpha`, `~alpha9`], [X])-(1/4)*Physics:-D_[`~rho1`](Physics:-Christoffel[`~alpha`, rho1, `~beta`], [X]) = 0

    (29)

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

    Simplify(-(1/2)*D_[`~beta`](Physics[Christoffel][`~upsilon`, upsilon, `~alpha`], [X])+(1/4)*D_[`~nu`](Physics[Christoffel][`~beta`, nu, `~alpha`], [X])+(1/4)*D_[`~mu`](Physics[Christoffel][`~alpha`, mu, `~beta`], [X])-(1/4)*D_[`~omega`](Physics[Christoffel][`~beta`, omega, `~alpha`], [X])+(1/2)*D_[`~beta`](Physics[Christoffel][alpha9, `~alpha`, `~alpha9`], [X])-(1/4)*D_[`~rho1`](Physics[Christoffel][`~alpha`, rho1, `~beta`], [X])-Physics[Ricci][`~alpha`, `~beta`]-(1/2)*Physics[Christoffel][alpha9, alpha10, `~alpha`]*Physics[Christoffel][`~alpha10`, `~alpha9`, `~beta`]-(1/4)*Physics[Christoffel][`~alpha`, chi, iota]*Physics[Christoffel][`~chi`, `~beta`, `~iota`]-(1/4)*Physics[Christoffel][`~beta`, chi, iota]*Physics[Christoffel][`~chi`, `~alpha`, `~iota`]+(1/2)*Physics[Christoffel][chi, `~alpha`, `~beta`]*Physics[Christoffel][iota, `~chi`, `~iota`]+(1/4)*(Physics[Christoffel][`~alpha`, chi, `~beta`]+Physics[Christoffel][`~beta`, chi, `~alpha`])*Physics[Christoffel][`~chi`, iota, `~iota`]-(1/2)*Physics[Christoffel][chi, iota, `~alpha`]*Physics[Christoffel][`~iota`, `~beta`, `~chi`]-(1/4)*Physics[Christoffel][`~alpha`, rho, `~beta`]*Physics[Christoffel][`~rho`, rho1, `~rho1`]+(1/4)*Physics[Christoffel][`~alpha`, rho, `~rho1`]*Physics[Christoffel][`~rho`, rho1, `~beta`]-(1/2)*Physics[Christoffel][omicron, zeta, `~omicron`]*Physics[Christoffel][`~zeta`, `~alpha`, `~beta`]+(1/2)*Physics[Christoffel][omicron, zeta, `~beta`]*Physics[Christoffel][`~zeta`, `~alpha`, `~omicron`]-(1/4)*Physics[Christoffel][`~beta`, psi, `~alpha`]*Physics[Christoffel][`~psi`, omega, `~omega`]+(1/4)*Physics[Christoffel][`~beta`, psi, `~omega`]*Physics[Christoffel][`~psi`, omega, `~alpha`]+(1/2)*Physics[Christoffel][`~tau`, upsilon, `~alpha`]*Physics[Christoffel][`~upsilon`, tau, `~beta`]-(1/2)*Physics[Christoffel][`~tau`, `~alpha`, `~beta`]*Physics[Christoffel][`~upsilon`, tau, upsilon]+(1/4)*Physics[Christoffel][`~beta`, nu, sigma]*Physics[Christoffel][`~sigma`, `~alpha`, `~nu`]-(1/4)*Physics[Christoffel][`~beta`, sigma, `~nu`]*Physics[Christoffel][`~sigma`, nu, `~alpha`]-(1/4)*Physics[Christoffel][`~alpha`, lambda, `~mu`]*Physics[Christoffel][`~lambda`, mu, `~beta`]+(1/4)*Physics[Christoffel][`~alpha`, lambda, mu]*Physics[Christoffel][`~lambda`, `~beta`, `~mu`]+(1/4)*(2*Physics[Christoffel][chi, iota, kappa]*Physics[Christoffel][`~iota`, `~chi`, `~kappa`]-2*Physics[Christoffel][chi, iota, `~chi`]*Physics[Christoffel][`~iota`, kappa, `~kappa`]-2*Physics[Christoffel][alpha4, alpha5, alpha6]*Physics[Christoffel][`~alpha6`, `~alpha4`, `~alpha5`]+2*Physics[Christoffel][alpha4, alpha6, `~alpha5`]*Physics[Christoffel][`~alpha6`, alpha5, `~alpha4`]-2*D_[`~alpha5`](Physics[Christoffel][alpha4, alpha5, `~alpha4`], [X])+2*Physics[Christoffel][`~alpha1`, alpha1, alpha3]*Physics[Christoffel][`~alpha3`, alpha2, `~alpha2`]-2*Physics[Christoffel][`~alpha1`, alpha3, `~alpha2`]*Physics[Christoffel][`~alpha3`, alpha1, alpha2]+2*Physics[Ricci][alpha2, `~alpha2`]+2*D_[`~alpha2`](Physics[Christoffel][`~alpha1`, alpha1, alpha2], [X]))*Physics[g_][`~alpha`, `~beta`]+(1/2)*Physics[Christoffel][alpha10, `~alpha`, `~beta`]*Physics[Christoffel][alpha9, `~alpha10`, `~alpha9`] = 0)

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

    (30)

    NULL


     

    Download Einstein_Equations_From_First_Principles.mw

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

    I’m absolutely delighted to announce the launch of Maple 2025!

    Although you see a new release every year, new features take anything from a few fast-paced weeks to develop, to months of careful cultivation.

    Working on so many features in parallel, each with varying time scales, isn't easy! We have to fastidiously manage and track our work.

    So it's easy to lose ourselves in the daily minutiae of software development. To help us maintain perspective, we constantly ask ourselves questions like:

    • What user problem are we solving and how often does this problem occur?
    • Can we validate our proposed solution with preliminary user feedback?
    • Is this a solution to a problem that doesn't exist and will never exist, or are we pre-empting a future need?
    • Are we offering value to our users?

    Given the answers, we course-correct to make sure we stay on track for our central mission - to make you happy, and to keep you coming back year-after-year.

    With Maple 2025, I think we've smashed that goal. We have many new features that'll appeal to many different types of users - from students, educators and mathematicians, to engineers, scientists and technical professionals

    Let me walk you through some of my personal highlights.

    It’ll be difficult for anyone to miss this - Maple 2025 has a new interface! It’s a ribbon-based UI that look clean and contemporary, and helps you find and discover tools more quickly than before.

    You have large, meaningful icons.

    Items are logically grouped.

    The ribbons is contextual. If you click on a plot, you get new tabs for interacting with and drawing on the plot.

    A new Education tab collects pedagogical resources that were scattered around the interface in prior releases.

    This is the biggest visual overhaul to Maple in many years. We hope you like it! 

    We also appreciate that changes in look and feel can be divisive. Please rest assured that we will refine and finesse the interface with each successive release; your comments and suggestions are most welcome.

    The new interface is available on Windows and Linux, and as a technology preview on Mac.

    The right arrow key on my keyboard is wearing out…and it’s all because of Maple. I’m knee deep in Maple nearly every day entering equations, and I’m always using right-arrow to move the cursor. It gets kind of tedious!

    This anecdote reflects some investigative work we did. We comprehensively examined our internal library of thousands of Maple worksheets and discovered that these three input patterns are extremely common.

    Previously, you’d use the right-arrow key to move the cursor out of the exponential, division or subscript.

    Now, in Maple 2025, when you

    • type ^, /, or enter a literal subscript with a double-underscore,
    • followed by a number or symbol
    • and then input another operator (such as +)

    the operator is automatically inserted on the baseline (except when y = 1).

    Of course, you can also make the cursor return to or stay in the exponent or denominator with a simple keystroke, when that is what is needed.

    This is one of those little quality of life refinements that I’m very fond of - it’s a little visual and usability dopamine hit.

    The sum command (and its typeset form) now indexes into vectors without you needing to spam unevaluation quotes all over your expression.

    Maple 2024

    Maple 2025

    We’ve been integrating units deeper into the Maple system, release after release. Much of this is driven by our engineering users.

    A few releases ago, we made int(numeric) compatible with units. With Maple 2025, you can now numerically differentiate  expressions and procedures that have units.

    I’m a grizzled thermodynamics hack, so here’s an example in which I calculate the specific heat capacity of water by differentiating enthalpy with respect to temperature (and then confirm the result with the built-in value):

    This is in addition to many other improvements to the units experience.

    Although this is a part of Maple that I don’t touch often (my colleague Karishma takes point on the education side), I REALLY wish I’d had this when I was struggling with math.

    You can now automatically generate unlimited variants of the same problem for students to solve with the Try Another feature, which has been added to Maple’s Check My Work tools (another feature I really could have used!). This is available for many common math principles, including factorization, simplification, integration and more.

    This is just one of the improvements in Maple 2025 for teaching and learning.

     If you’ve ever found yourself going back and forth (and back and forth) between two large, almost identical-looking Maple expressions, trying to figure out how they are different, you’re going to love this one.  ExpressionTools is a new package that lets you compare the differences between two expressions.

    I really like the use of color to highlight differences. Less squinting at the screen!

    You can now run Maple Flow worksheets from Maple (you don’t need Maple Flow installed to do this). You can send parameters into the Flow worksheet and extract your desired results.

    This means you can use the entire flexibility of Maple to analysis and manipulate your Flow worksheet. You could, for example:

    • Attach a Flow worksheet to a Maple workbook and create an interactive application
    • Carry out parameter studies of a Flow worksheet by evaluating it over many parameter sets in Maple
    • Create an Excel interface for a Flow worksheet using the Maple add-in for Excel

    Simplify is one of those functions that literally tens of thousands of people use each day. Every time we make an incremental improvement, the cumulative benefits across our entire user base are significant.

    We’ve refined simplify in a number of critical ways. For example, simplify now recognizes when exponentials can be profitably converted to hyperbolic trig functions:

    The analysis of many scientific phenomena result in Laplace transforms that do not have a symbolic inverse which can be expressed in terms of elementary functions. This includes applications in heat transfer, fluid mechanics, fractional diffusion processes, control systems and electrical transmission.

    For example, this monster Laplace transform results from an analysis of voltage on a transmission line:

    You can now numerically invert this transform courtesy of an enhancement to inttrans:-invlaplace - a fast quadrature method.

    I’ve saved what I think has the most future potential for last.

    I’m sure nearly all of you have experimented with the various AI tools. They’re an inevitable part of our present and future, whether we're comfortable with it or not.

    This is something we've been mulling over for some time.

    • In Maple 2019, the DeepLearning package made its debut. This package provides tools for machine learning, supporting operations such as classification and regression using neural networks.
    • In Maple 2024, we introduced an AI-powered formula lookup feature.

    In Maple 2025, we’re giving you an early-stage technology preview of AI-powered document generation.

    You can automatically generate worksheet content by prompting an AI, and then gradually refine the content

    If you’re an educator, you might want some content that describes applications of calculus. So you might ask the AI “How do I derive the formula for the area of a circle” by entering your prompt into this text box:

    This is the worksheet content that may be returned:

    If you’re structural engineer who wants to know how to calculate the hardness of concrete, you might ask the AI: “How do I calculate the compressive strength of slow hardening concrete as a function of time? Use the CEB-FIP Model Code 90. Include a worked example with Maple code”.

    This worksheet content that could be generated (note the live Maple code):

    We’re labelling AI-generated worksheet content as a technology preview. You might see

    • text that might be misleading (but sounds plausible)
    • code that doesn’t work (but looks plausible)
    • or different results each time you click “Generate Document”

    For the moment, I would not rely on AI-generated worksheet content without realistic expectations, a healthy dose of scepticism and a modicum of detached analysis. But AI models are rapidly growing in robustness, and we want to position ourselves to best exploit their future potential. The next few years will be VERY exciting.

    We can never cover everything in a short blog post like this. So if you want to know more, head on over to the What’s New pages for Maple 2025!

    Hi MaplePrimes,
    I made a quick code to calculate Aliquot Sequences.

    The code works.

    I was inspired by a Numberphile YouTube video.

    see

    aliquot_sequence_cut_1.mw

    aliquot_sequence_cut_1.pdf

    This should be helpful.

    Regards,

    Matt

    LIMITS

    Limits in maths are defined as the values that a function approaches the output for the given input values. Limits play a vital role in calculus and mathematical analysis and are used to define integrals, derivatives, and continuity. It is used in the analysis process, and it always concerns the behavior of the function at a particular point. The limit of a sequence is further generalized in the concept of the limit of a topological net and related to the limit and direct limit in the theory category. Generally, the integrals are classified into two types namely, definite and indefinite integrals. For definite integrals, the upper limit and lower limits are defined properly. Whereas indefinite integrals are expressed without limits, and it will have an arbitrary constant while integrating the function.

    Sometimes we can't work something out directly ... but we can see what it should be as we get closer and closer!

    Example 1

    "restart;  f(x):=(|x|-3)/(x-3);"

    proc (x) options operator, arrow, function_assign; (abs(x)-3)/(x-3) end proc

    (1)

    plot(f(x), x = -10 .. 10, discont = true, color = "Green")

     

    f(3)

    Error, (in f) numeric exception: division by zero

     

    Now 0/0 is a difficulty! We don't really know the value of 0/0 (it is "indeterminate"), so we need another way of answering this.

    So instead of trying to work it out for x=3 let's try approaching it closer and closer:

    f(3.01)

    1.000000000

    (2)

    f(3.0000001)

    1.000000000

    (3)

    f(2.9999999)

    1.000000000

    (4)

    Limit(f(x), x = 3)

    Limit((abs(x)-3)/(x-3), x = 3)

    (5)

    limit(f(x), x = 3)

    1

    (6)

    limit(f(x), x = 3, left)

    1

    (7)

    limit(f(x), x = 3, right)

    1

    (8)

    Example 2

    Sometimes some functions are not continuous. That is, they appear to be approaching two different values when they are approached from two sides.

    "g(x):=piecewise(0<x<2,1/(2 x-x^(2)),2 <x<=3,2 -x,3<x<4,x-4, 4<=x,Pi,undefined);"

    proc (x) options operator, arrow, function_assign; piecewise(0 < x and x < 2, 1/(2*x-x^2), 2 < x and x <= 3, 2-x, 3 < x and x < 4, x-4, 4 <= x, Pi, undefined) end proc

    (9)

    plot(g(x), x = -10 .. 10, y = -1 .. 10, discont = true, color = "Red")

     

    Suppose we want to approach 2 and see the function’s limit. This naturally leads to directions from which we can approach. Left-hand side and the right-hand side limits.

    The right-hand side limit is the value of the function that it takes while approaching it from the right-hand side of the desired point. Similarly, the left-hand side limit is the value of function while approaching it from the left-hand side.

    eval(g(x), x = 2)

    undefined

    (10)

    limit(g(x), x = 2, left)

    infinity

    (11)

    limit(g(x), x = 2, right)

    0

    (12)

    limit(g(x), x = 2)

    undefined

    (13)

    And the ordinary limit "does not exist".

    g(4)

    Pi

    (14)

    limit(g(x), x = 4, left)

    0

    (15)

    limit(g(x), x = 4, right)

    Pi

    (16)

    limit(g(x), x = 4)

    undefined

    (17)

    And the ordinary limit "does not exist".

    with(Student[Calculus1]); LimitTutor()

    Example 3

    Estimate the value of the following limit limit(h(x)*where, x = 2), h(x) = piecewise(x <> 2, x+12, x = 2, 4).

    "h(x):={[[x+12,x<>2],[4,x=2]];"

    proc (x) options operator, arrow, function_assign; piecewise(x <> 2, x+12, x = 2, 4) end proc

    (18)

    plot(h(x), x = -10 .. 10, discont = true, color = "#40e0d0")

     

    limit(h(x), x = 2)

    14

    (19)

    The limit is NOT 2025!Remember from the first example that limits do not care what the function is actually doing at the point in question. Limits are only concerned with what is going on around the point. Since the only thing about the function that we actually changed was its behavior at x = 2 this will not change the limit.

    Example 4

    " w(x):=piecewise( x<0,-x+5,x>=0,2 x);"

    proc (x) options operator, arrow, function_assign; piecewise(x < 0, -x+5, 0 <= x, 2*x) end proc

    (20)

    plot(w(x), x = -10 .. 10, y = -10 .. 10, discont = true, color = "Blue")

     

    limit(w(x), x = 5)

    10

    (21)

    limit(w(x), x = 6, left)

    12

    (22)

    limit(w(x), x = 1, right)

    2

    (23)

    Example 5

    " k(x):=piecewise( x<5,x+4,x>=5, x^(2)-2);"

    proc (x) options operator, arrow, function_assign; piecewise(x < 5, x+4, 5 <= x, x^2-2) end proc

    (24)

    plot(k(x), x = -10 .. 10, discont = true, color = orange)

     

    limit(k(x), x = 2)

    6

    (25)

    limit(k(x), x = 5, left)

    9

    (26)

    limit(k(x), x = 5, right)

    23

    (27)

    limit(k(x), x = 5)

    undefined

    (28)

    limit(k(x), x = 6)

    34

    (29)

    Example 6

    restart

    " l(x):=piecewise( x<=1,(x-8)/(x-3),x>=3, sqrt(x^(2)+x+2), undefined);"

    proc (x) options operator, arrow, function_assign; piecewise(x <= 1, (x-8)/(x-3), 3 <= x, sqrt(x^2+x+2), undefined) end proc

    (30)

    plot(l(x), x = -10 .. 10, discont = true, color = "Blue")

     

    limit(l(x), x = 0)

    8/3

    (31)

    limit(l(x), x = 1, left)

    7/2

    (32)

    limit(l(x), x = 1, right)

    undefined

    (33)

    limit(l(x), x = 2)

    undefined

    (34)

    Example 7

    Estimate the value of the following limit. limit(H(t), t = 0)where, H(t) = piecewise(t < 0, 0, t >= 0, 1)

    "  H(t):=piecewise( t<0,0,t>=0, 1);"

    proc (t) options operator, arrow, function_assign; piecewise(t < 0, 0, 0 <= t, 1) end proc

    (35)

    This function is often called either the Heaviside or step function. We could use a table of values to estimate the limit, but it’s probably just as quick in this case to use the graph so let’s do that. Below is the graph of this function.

    plot(H(t), t = -10 .. 10, discont = true, color = "Blue")

     

    limit(H(t), t = 0, left)

    0

    (36)

    limit(H(t), t = 0, right)

    1

    (37)

    We can see from the graph that if we approach t = 0from the right side the function is moving in towards a yvalue of 1. Well actually it’s just staying at 1, but in the terminology that we’ve been using in this section it’s moving in towards 1.

    Also, if we move in towards t = 0 from the left the function is moving in towards a yvalue of 0.

    According to our definition of the limit the function needs to move in towards a single value as we move in towards t = a (from both sides). This isn’t happening in this case and so in this example we will also say that the limit doesn’t exist.

     

    NULL

    Download limits.mw

    As some of you may have noticed, MaplePrimes was down for a few days, and some customers may have encountered a few problems with our web site, as well.

    As you can see, MaplePrimes is up and running again. We are still experiencing a few lingering issues on our main web site, but we expect these will be resolved shortly.

    Apologies to those who were inconvenienced by the outage. We appreciate your patience.

    fyi, there is new video showing Maple's 2025 new interface

    It seems oriented to document mode which I do not use. May be they also improved worksheet mode.

    I am still getting my Maple desktop getting shuffled few times each day where I have to close Maple and reopen it to clear it. I hope they fixed this in the new interface.,

    Hello everyone,

    I have created a Maple worksheet titled "ΕΜΒΑΔΟΝ ΕΠΙΠΕΔΟΥ ΧΩΡΙΟΥ", designed to help my students prepare for their final exams as they qualify for university. This worksheet focuses on area calculations in plane geometry, using Maple to visualize and solve problems efficiently.

    This worksheet is aimed at high school students preparing for university entrance exams, as well as teachers who want to integrate Maple into their teaching.

    I would love to hear your thoughts and feedback!

    Have you used Maple for similar exam preparation?
    εμβαδόν_χωρίου.mw

    MaplePrimes offers spell check and correction

    for the function "Contact Author"

     

    Suggestion:
    For the sake of message clarity (and to save time) it would be desireable to have spell checking and correction in other MaplePrimes message bodies of as well.

     

    Recently, @zenterix asked a question related to solving the quantum mechanics of the finite-wall-height particle-in-a-box problem. In the last two decades or so I've been teaching a quantum mechanics course every second year, in which I sketch the method of solution of this problem for the class, I but do not get into the mathematical details. A few times I've started to work on it in Maple, but abandoned it because of lack of time. So I decided to persist this time.

    This post is more about some of the challenges and tradeoffs in making it work, hence the title. The title is also a nod to the fact that this problem appears in Engel's textbook [1] in a chapter entitled "Applying the particle in a box model to real-world topics". You don't need to know much about quantum mechanics to understand it. From the mathematical point of view you are solving three ordinary differential equations (in three regions) and stitching the solutions together so that the solution, a wavefunction, is continuous with continuous derivative and tends to zero at plus infinity and minus infinity.

    The three regions and the wavefunctions for three eigenvalues (energies) are shown here in the figure, which is the final output of the worksheet. Next is the worksheet and I'll comment further on it below.

    (It's a quirk of quantum mechanics that we plot two things with different units (wavefunctions and energy) on the same axis, and worse, we offset one by the other.)

    Bound states for particle in finite height box of wall height
    V(x) = piecewise(x < -(1/2)*a, V__0, -(1/2)*a <= x and x <= (1/2)*a, 0, (1/2)*a < x, V__0)
    I'll call these regions left, box and right respectively.

    restart

    Schroedinger equation. Here as elsewhere "=0" is implied.

    Schroedinger := -`&hbar;`^2*(diff(psi(x), x, x))/(2*m)+V*psi(x)-E*psi(x)

    -(1/2)*`&hbar;`^2*(diff(diff(psi(x), x), x))/m+V*psi(x)-E*psi(x)

    (1)

    Nondimensionalize length as units of box length, X = x/a. Since psi(x)^2*dx is unitless probability, psi(x) has dimensions 1/length^(1/2)so we define a dimensionless Phi = a^(1/2)*psi

    tr := {x = a*X, psi(x) = Phi(X)/sqrt(a)}; NDSchroedinger := PDETools:-dchange(tr, Schroedinger, [X, Phi(X)], params = [`&hbar;`, m, a], simplify)

    {x = a*X, psi(x) = Phi(X)/a^(1/2)}

     

    -((1/2)*`&hbar;`^2*(diff(diff(Phi(X), X), X))+a^2*m*Phi(X)*(E-V))/(a^(5/2)*m)

    (2)

    This will also change the normalization integrals, e.g., over the box region:

    normint := Int(psi(x)^2, x = -(1/2)*a .. (1/2)*a); NDnormint := PDETools:-dchange(tr, normint, [X, Phi(X)], params = [`&hbar;`, m, a], simplify)

    Int(psi(x)^2, x = -(1/2)*a .. (1/2)*a)

     

    Int(Phi(X)^2, X = -1/2 .. 1/2)

    (3)

    Outside the box we have 0 <= E and E < V with V = V__0. Define a dimensionless positive constant kappa:

    `&kappa;__eqn` := kappa = a*sqrt(2*m*(V__0-E)/`&hbar;`^2); de_outside := simplify(sqrt(a)*kappa^2*(eval(NDSchroedinger, {isolate(`&kappa;__eqn`, m), V = V__0}))/(E-V__0)); outside := rhs(dsolve(de_outside))

    kappa = a*2^(1/2)*(m*(V__0-E)/`&hbar;`^2)^(1/2)

     

    diff(diff(Phi(X), X), X)-kappa^2*Phi(X)

     

    c__1*exp(-kappa*X)+c__2*exp(kappa*X)

    (4)

    Inside the box we have V = 0 and E > 0. Define a dimensionless positive constant k.

    k__eqn := k = a*sqrt(2*m*E/`&hbar;`^2); de_box := simplify(-sqrt(a)*k^2*(eval(NDSchroedinger, {isolate(k__eqn, m), V = 0}))/E); box := rhs(eval(dsolve(de_box), {c__1 = A, c__2 = B}))

    k = a*2^(1/2)*(m*E/`&hbar;`^2)^(1/2)

     

    diff(diff(Phi(X), X), X)+k^2*Phi(X)

     

    A*sin(k*X)+B*cos(k*X)

    (5)

    We have seven unknowns {A, B, E, c__1(L), c__1(R), c__2(L), c__2(R)}, where for example c__1(L) means Maple's c__1 for the left region. We need seven equations: goes to 0 at x = -infinityand at x = infinity(or Phi(X) wouldn't be square integrable), continuity at the two box boundaries, slopes continuous at the two box boundaries, and normalization. We deal with the ones at `&+-`(infinity)first.

    At X = -infinity the wavefunction will be unbounded unless one of the constants goes to zero. We'll rename the other one A__L or A__R. The code here is just to handle the fact that c__1 and c__2 can be the the opposite way around, depending on the session.

    `assuming`(['limit(outside, X = -infinity)' = limit(outside, X = -infinity)], [kappa > 0]); left0 := eval(outside, `~`[`=`](`minus`(indets(rhs(%), name), {infinity}), 0)); left1 := eval(left0, `~`[`=`](indets(left0, suffixed(c)), A__L)); `assuming`(['limit(outside, X = infinity)' = limit(outside, X = infinity)], [kappa > 0]); right0 := eval(outside, `~`[`=`](`minus`(indets(rhs(%), name), {infinity}), 0)); right1 := eval(right0, `~`[`=`](indets(right0, suffixed(c)), A__R))

    limit(c__1*exp(-kappa*X)+c__2*exp(kappa*X), X = -infinity) = signum(c__1)*infinity

     

    A__L*exp(kappa*X)

     

    limit(c__1*exp(-kappa*X)+c__2*exp(kappa*X), X = infinity) = signum(c__2)*infinity

     

    A__R*exp(-kappa*X)

    (6)

    Now we have five unknowns {A, A__L, A__R, B, E}. The strategy will be to eliminate A, A__L, A__R and then find the eigenvalues E. Then we will have expressions for all three regions containing the unknown B, which we will find by the normalization condition (The normalization condition fixes the scale of the wavefunction so that the probability of finding the particle somewhere is one.)

    Continuity of the wavefunction and its derivative at the left boundary are below
    We eliminate A__L and A using this boundary

    bcleft := {eval(left1-box, X = -1/2), eval(diff(left1, X)-(diff(box, X)), X = -1/2)}; sol := solve(bcleft, {A, A__L}); left2 := eval(left1, sol); box2 := eval(box, sol)

    {A__L*exp(-(1/2)*kappa)+A*sin((1/2)*k)-B*cos((1/2)*k), A__L*kappa*exp(-(1/2)*kappa)-A*k*cos((1/2)*k)-B*k*sin((1/2)*k)}

     

    {A = -B*(sin((1/2)*k)*k-cos((1/2)*k)*kappa)/(sin((1/2)*k)*kappa+cos((1/2)*k)*k), A__L = B*k*(sin((1/2)*k)^2+cos((1/2)*k)^2)/(exp(-(1/2)*kappa)*(sin((1/2)*k)*kappa+cos((1/2)*k)*k))}

     

    B*k*(sin((1/2)*k)^2+cos((1/2)*k)^2)*exp(kappa*X)/(exp(-(1/2)*kappa)*(sin((1/2)*k)*kappa+cos((1/2)*k)*k))

     

    -B*(sin((1/2)*k)*k-cos((1/2)*k)*kappa)*sin(k*X)/(sin((1/2)*k)*kappa+cos((1/2)*k)*k)+B*cos(k*X)

    (7)

    Eliminate A__R at the right box boundary

    bcright := {eval(box2-right1, X = 1/2), eval(diff(box2, X)-(diff(right1, X)), X = 1/2)}; elim := eliminate(bcright, A__R); right2 := eval(right1, elim[1]); energy_eqn := elim[2][]/B

    {-B*(sin((1/2)*k)*k-cos((1/2)*k)*kappa)*sin((1/2)*k)/(sin((1/2)*k)*kappa+cos((1/2)*k)*k)+B*cos((1/2)*k)-A__R*exp(-(1/2)*kappa), -B*(sin((1/2)*k)*k-cos((1/2)*k)*kappa)*k*cos((1/2)*k)/(sin((1/2)*k)*kappa+cos((1/2)*k)*k)-B*k*sin((1/2)*k)+A__R*kappa*exp(-(1/2)*kappa)}

     

    [{A__R = B*(2*sin((1/2)*k)*cos((1/2)*k)*kappa+2*cos((1/2)*k)^2*k-k)/(exp(-(1/2)*kappa)*(sin((1/2)*k)*kappa+cos((1/2)*k)*k))}, {-2*(sin((1/2)*k)*cos((1/2)*k)*k^2-sin((1/2)*k)*cos((1/2)*k)*kappa^2-2*cos((1/2)*k)^2*k*kappa+k*kappa)*B}]

     

    B*(2*sin((1/2)*k)*cos((1/2)*k)*kappa+2*cos((1/2)*k)^2*k-k)*exp(-kappa*X)/(exp(-(1/2)*kappa)*(sin((1/2)*k)*kappa+cos((1/2)*k)*k))

     

    -2*sin((1/2)*k)*cos((1/2)*k)*k^2+2*sin((1/2)*k)*cos((1/2)*k)*kappa^2+4*cos((1/2)*k)^2*k*kappa-2*k*kappa

    (8)

    Now we will get the quantized energies by finding which values will make energy_eqn zero. We recall the relationships of `&kappa;__` and k to energy:

    k__eqn; `&kappa;__eqn;`

    k = a*2^(1/2)*(m*E/`&hbar;`^2)^(1/2)

     

    kappa = a*2^(1/2)*(m*(V__0-E)/`&hbar;`^2)^(1/2)

    (9)

    The energy scale is set by V__0, so can get a non-dimensionalized energy e = E/V__0 and a non-dimensional parameter b = a*sqrt(2*m*V__0/`&hbar;`^2).

    b__eqn := b = a*sqrt(2*m*V__0/`&hbar;`^2); e__eqn := e = E/V__0; eqns := `assuming`([{simplify(eval(eval(k__eqn, isolate(e__eqn, E)), isolate(b__eqn, V__0))), simplify(eval(eval(`&kappa;__eqn`, isolate(e__eqn, E)), isolate(b__eqn, V__0)))}], [a > 0])

    b = a*2^(1/2)*(m*V__0/`&hbar;`^2)^(1/2)

     

    e = E/V__0

     

    {k = (e*b^2)^(1/2), kappa = (-b^2*(e-1))^(1/2)}

    (10)

    So now for any values of the box length and height, one finds the parameter b describing the problem and solves for the possible e values.

    energy_eqn2 := `assuming`([simplify((eval(energy_eqn, eqns))/b^2)], [b > 0, e > 0, e < 1])

    2*e^(1/2)*(-e+1)^(1/2)*cos(b*e^(1/2))-2*e*sin(b*e^(1/2))+sin(b*e^(1/2))

    (11)

    Maple cannot find an analytical solution, so we will find the eigenvalues numerically. The zero solution is unphysical.

    solve(energy_eqn2, e)

    0, RootOf(-2*(-(_Z^2-b^2)/b^2)^(1/2)*(_Z^2/b^2)^(1/2)*b^2+2*tan(_Z)*_Z^2-tan(_Z)*b^2)^2/b^2

    (12)

    For any b, read off the possible energies off the vertical line for that b, e.g. there are 3 bound states for b = 9. For other values, choose bval and nsols on the next line appropriately.

    bval := 9; nsols := 3; plots:-implicitplot([b = bval, energy_eqn2], b = 0 .. 15, e = 0 .. 1)

    9

     

    3

     

     

    evals := [fsolve(eval(energy_eqn2, b = bval), e = 0 .. 1, maxsols = nsols)]

    [0.8115199822e-1, .3189598393, .6884466500]

    (13)

    In terms of the parameters b and e and normalization constant B, the non-dimensionalized wavefunctions of the 3 regions are as below. The scale factor gives a simpler form and prevents 0/0 problems later in the calculation.

    scale := `assuming`([denom(simplify(eval(left2, eqns)))], [positive]); left3 := `assuming`([simplify(eval(scale*left2, eqns))], [positive]); box3 := `assuming`([simplify(eval(scale*box2, eqns))], [positive]); right3 := `assuming`([simplify(eval(scale*right2, eqns))], [positive])

    sin((1/2)*b*e^(1/2))*(-e+1)^(1/2)+cos((1/2)*b*e^(1/2))*e^(1/2)

     

    B*e^(1/2)*exp((1/2)*b*(-e+1)^(1/2)*(2*X+1))

     

    -(sin((1/2)*b*e^(1/2))*(e^(1/2)*sin(b*e^(1/2)*X)-(-e+1)^(1/2)*cos(b*e^(1/2)*X))-cos((1/2)*b*e^(1/2))*(cos(b*e^(1/2)*X)*e^(1/2)+sin(b*e^(1/2)*X)*(-e+1)^(1/2)))*B

     

    B*exp(-(1/2)*b*(-e+1)^(1/2)*(-1+2*X))*(e^(1/2)*cos(b*e^(1/2))+(-e+1)^(1/2)*sin(b*e^(1/2)))

    (14)

    It remains to find the normalization constant B. The three parts of the normalization integral are:

    P__L := `assuming`([int(left3^2, X = -infinity .. -1/2)], [positive]); P__box := `assuming`([int(box3^2, X = -1/2 .. 1/2)], [positive]); P__R := `assuming`([int(right3^2, X = 1/2 .. infinity)], [positive])

    (1/2)*B^2*e/(b*(-e+1)^(1/2))

     

    -(1/4)*B^2*(2*e^(1/2)*(-e+1)^(1/2)*cos(2*b*e^(1/2))-2*e^(1/2)*(-e+1)^(1/2)-2*b*e^(1/2)-2*e*sin(2*b*e^(1/2))+sin(2*b*e^(1/2)))/(b*e^(1/2))

     

    (1/2)*B^2*(e^(1/2)*cos(b*e^(1/2))+(-e+1)^(1/2)*sin(b*e^(1/2)))^2/(b*(-e+1)^(1/2))

    (15)

    Solve for B. There are two (messy) solutions for B of opposite sign; it is conventional to choose the sign such that the ground state has positive values.

    Bval := sort([solve(P__L+P__box+P__R = 1, B)])[2]

    Assemble it into a single piecewise function.

    `&Phi;__all` := eval(piecewise(X < -1/2, left3, `and`(X >= -1/2, X <= 1/2), box3, X > 1/2, right3), B = Bval)

    Plot. Plots offset vertically by the energy as usual

    Xmax := 1.5; psiscale := .15; colors := [red, blue, magenta]; wavefns := plot([seq(eval(psiscale*`&Phi;__all`, {b = bval, e = evals[i]})+evals[i], i = 1 .. nsols)], X = -Xmax .. Xmax, color = colors); walls := plot([-1/2, y, y = 0 .. 1], color = black), plot([1/2, y, y = 0 .. 1], color = black), plot(1, X = -Xmax .. -1/2, color = black), plot(0, X = -1/2 .. 1/2, color = black), plot(1, X = 1/2 .. Xmax, color = black); evalplot := plot(evals, X = -Xmax .. Xmax, linestyle = dash, color = colors); plots:-display(wavefns, walls, evalplot, axes = boxed, labels = [X, psiscale*Phi+E/V__0])

     

    NULL

    Download finite_box_non-dim2.mw

    Comments
    1. The first general comment is about the issue of how much Maple can do. In principle one should be able to give dsolve the three ODEs and the boundary conditions and it should output the result. We are still far away from that. Even for one region, the boundary conditions at infinity are not handled by dsolve. The other extreme is to solve the ODEs for their general solutions, and follow the algebraic manipulations for implementing the boundary conditions as one might do it on paper. Engel for example has a comment "At this point we notice that by dividing the equations in each pair, the coefficients can be eliminated to give [...]". Maple will not easily do that trick or the manipulations that led to the special form on which the trick is applied. Levine's textbook [2] gives a different non-intuitive set of steps.

    But Maple can solve complicated equations, so I wanted a more general strategy that avoids as much as possible the requirement to manipulate into special forms. On the other hand, it is tempting (as @zenterix tried) to just give the three general solutions with 6 arbitrary constants and the boundary conditions to Maple's solve, and solve for the six constants. One finds that all six are zero. This us a consequence of the fact that the ODEs are linear, and misses the crucial point that it is an eigenvalue problem. So there must be some sort of step-by-step guidance for Maple and there is a general but non-obvious strategy to solving such problems.

    2. The second general comment is that there is a trade-off between presenting the solution steps without any arcane Maple code, and getting to the solution efficiently and programmatically, without too many manual manipulations. One of my suggestions to @zenterix involved manually dividing both sides of an equations by a fairly complicated expression, to which @acer expressed a preference for programmatic solutions. At the time, I mainly did that to keep close the to existing solution, and was thinking that I usually set things up so that is not necessary. But I found here that I do it a lot, and it is a natural consequence of the interactive way I use Maple. Here's two examples:

    In the second line of label (4), I originally got a more complicated form of de_outside, then manually figured out the factor to put inside simplify to get the form you see. I do this frequently with Maple, especially when I use dchange.

    In the last line of label (8), I manually divided by B to remove it. Had B been in the denominator, I could have removed it programmatically with numer (an operation I use frequently, though it somewhat recklessly ignores denominator zeroes).

    Programmatically doing things can be relatively unreadable and disrupt the readability for the non-Maple reader. For example in label (6) there is "extraneous" code to find which coefficient to set to zero so that the solution goes to zero at +/- infinity. This was forced on me since after generating the nice figure I re-ran the worksheet only to have multiple errors. I originally used eval(outside, {c__1 = B__L, c__2 = A__L}); to change Maple's dsolve constants to read the way I wanted. But I realized that dsolve does not reproducibly assign c__1 and c__2 to the same term each time, raising the general issue of: 

    3. Session to session reproducibility. If I am only going to use a worksheet myself, I don't usually worry too much about this, since I can usually debug this fairly easily. On the other hand if the worksheet is for others, it needs to be foolproof. I have taken to always sorting the output of solve, e.g., the assignment to bval in the execution group after label (15). The code added for the dsolve constants works, but is non-trivial Maple (suffixed type testing). (Solving this issue for Eigenvectors is yet more difficult.)

    4. Some efforts to make the worksheet more readable were hard to do.

    - For some reason, Engel chose to call the two constants in the left region B and B'. Entering B*exp(-k*x) + B'*exp(k*x) in 2D leads to B(x)*exp(-k*x) + diff(B(x), x)*exp(k*x). No doubt this can be fixed, but I didn't pursue this.

    - I would normally use upper case Psi for the non-dimensionalized psi, but Psi is the name of a special function in Maple. But wait, now we can use "local Psi;" to solve this problem. However diff(Psi(x),x,x) displays as Psi(2,x), so I had to settle for Phi. (SCR submitted.)

    - I build up the left wavefunctions incrementally as expressions assigned to left1, left2, etc, which is rather ugly. Why not use arrow procedures so we can use psi(x), D(psi)(x) etc. Aside from some awkwardness in updating such procedures by redefining them as needed, I ran into the following problem:

    limit(A*exp(-k*x) + B*exp(k*x), x = infinity) assuming k > 0; gives as expected
    signum(B)*infinity;

    but

    psi := x -> A*exp(-k*x) + B*exp(k*x);
    limit(psi(x), x = infinity) assuming k > 0;
    returns unevaluated. What happened to full evaluation?

    - I'm used to entering hbar as hbar in 1D; it displays as h with the bar through. In 2-D entering hbar^2 is displayed nicely, but internally it is `&hbar;`. We can have the following puzzle:

    (hbar/m)^2-`&hbar;`^2/m^2

    hbar^2/m^2-`&hbar;`^2/m^2

    NULL

    - I originally wanted a vertical axis label Psi, E/V__0, but labels = [X, Psi, E/V__0] obviously won't work. It took some time to figure out the answer is to make "Psi, E/V__0" an atomic variable. I finally settled on the more honest, if cryptic, label shown.

    5. solve gives an "invalid" solution. Originally I worked with solutions left2 etc (see label (6)) after simplify(eval(left2, eqns)). Every second wavefunction was about 10^9 times higher than the others, which was difficult to debug. It turns out the denominator of those wavefunctions was exactly zero but numerically 10^(-10), meaning they plotted as large-scaled versions of the right shape, without any division-by-zero error. The fix is to remove the denominators by scaling (see label (14)).

    This is just a consequence of Maple giving generic solutions, and the only way around it is to be aware of it.

    Summary

    To solve such problems with Maple requires one to play around in an interactive way. After hacking to a solution, it can take some effort to make it presentable and usable to others. But the final result is pleasing, and is certainly much easier than working through the problem on paper. The ability to manipulate complicated expressions means circumventing tricks that might be needed in manual solutions. 

    References
    [1] T. Engel, Quantum Chemistry and spectroscopy. Pearson 2019, 4th ed. Sec. 5.1, Prob. 5.3. 
    [2] I. N. Levine, Quantum Chemistry. Pearson 2009, 6th ed. Sec. 2.4. 

    To Maplesoft,

    Please consider changing the name from Maple to something else.

    It is almost impossible to search for anything related to maple, since google keeps giving results about trees called Maple in Canada and about some maple syrup products which I have no interest in at all.

    One has to go through pages and pages of links looking for a real Maple software hit.

    I know the name Maple has been around for long time, but a new unique name will make searching easier and people will get used to the new name very quickly (may be in 1-2 years).

    Some examples:

     

    I know when Maple was created almost 40 years ago, the inernet itself was not even here (I forgot when VP. Gore created the internet but I think that was in early 90's), and search was not thought about then.

    But these days, the ability to search for something and to easily find it is very important for companies and having a unique name for Maple will make it also much more popular and easier to find things about it instead of finding  information about Maple syrup and Maple trees all the time.

    This is posted  under product suggestions.

    Hi Maple lovers, and others,

    I wrote up some simple Maple code that 
    checks if an expression (a^3+2) is a prime number.

    The output fits on one page of paper.
    n_cubed_plus_2.mw

    n_cubed_plus_2.pdf

    Also, I have an account at mersenneforum.org

    Kindest Regards,

    Matt

    PS Does this help?  The code works.

    1 2 3 4 5 6 7 Last Page 2 of 305