ecterrab

14540 Reputation

24 Badges

20 years, 24 days

MaplePrimes Activity


These are Posts that have been published by ecterrab

The computation of traces of products of Dirac matrices was implemented years ago - see Physics,Trace .

 

The simplification of products of Dirac matrices, however, was not. Now it is, and illustrating this new feature is the matter of this post. To reproduce the results below please update the Physics library with the one distributed at the Maplesoft R&D Physics webpage.

with(Physics)

 

First of all, when loading Physics, a frequent question is about the signature, the default is (- - - +)

Setup(signature)

[signature = `- - - +`]

(1)

This is important because the convention for the Algebra of Dirac Matrices depends on the signature. With the signatures (- - - +) as well as (+ - - -), the sign of the timelike component is 1

Library:-SignOfTimelikeComponent()

1

(2)

With the signatures (+ + + -) as well as (- + + +), the sign of the timelike component is of course -1

Library:-SignOfTimelikeComponent(`+ + + -`)

-1

(3)

The simplification of products of Dirac Matrices, illustrated below with the default signature, works fine with any of these signatures, and works without having to set a representation for the Dirac matrices -- all the results are representation-independent.

 

The examples below, however, also illustrate a new feature of Physics, for now implemented as a Library:-PerformMatrixOperations command (there is a related, also new, command, Library:-RewriteInMatrixForm, to just present the underlying matrix operations, without performing them). To illustrate this other new functionality , set a representation for the Dirac matrices, say the standard one

 

Setup(Dgamma = standard, math = true)

`* Partial match of  'Physics:-Dgamma' against keyword 'Dgammarepresentation'`

 

`* Partial match of  'math' against keyword 'mathematicalnotation'`

 

`Setting lowercaselatin letters to represent spinor indices `

 

`Defined Dirac gamma matrices (Dgamma) in standard representation`, gamma[1], gamma[2], gamma[3], gamma[4]

 

__________________________________________________

 

[Dgammarepresentation = standard, mathematicalnotation = true]

(4)

The four Dirac matrices are

TensorArray(Dgamma[`~mu`])

Array(%id = 18446744078360533342)

(5)

The definition of the Dirac matrices is implemented in Maple following the conventions of Landau books ([1] Quantum Electrodynamics, V4), and  does not depend on the signature, ie the form of these matrices is

"Library:-RewriteInMatrixForm(?)"

Array(%id = 18446744078360529726)

(6)

With the default signature, the space part components of  gamma[mu] change sign when compared with corresponding ones from gamma[`~mu`] while the timelike component remains unchanged

TensorArray(Dgamma[mu])

Array(%id = 18446744078565663678)

(7)

"Library:-RewriteInMatrixForm(?)"

Array(%id = 18446744078677131982)

(8)

For the default signature, the algebra of the Dirac Matrices, loaded by default when Physics is loaded, is (see page 80 of [1])

(%AntiCommutator = AntiCommutator)(Dgamma[`~mu`], Dgamma[`~nu`])

%AntiCommutator(Physics:-Dgamma[`~mu`], Physics:-Dgamma[`~nu`]) = 2*Physics:-g_[`~mu`, `~nu`]

(9)

When the sign of the timelike component of the signature is -1, we have a -1 factor on the right-hand side of (9).

 

Note as well that in (9) the right-hand side has no matrix elements. This is standard in particle physics where the computations are performed algebraically, without performing the matrix operations. For the purpose of actually performing the underlying matrix operations, however, one may want to rewrite this algebra including a 4x4 identity matrix. For that purpose, see Algebra of Dirac Matrices with an identity matrix on the right-hand side. For the purpose of this illustration, below we proceed with the algebra as shown in (9), interpreting right-hand sides as if they involve an identity matrix.

 

Verify the algebra rule by performing all the involved matrix operations

expand(%AntiCommutator(Physics[Dgamma][`~mu`], Physics[Dgamma][`~nu`]) = 2*Physics[g_][`~mu`, `~nu`])

Physics:-`*`(Physics:-Dgamma[`~mu`], Physics:-Dgamma[`~nu`])+Physics:-`*`(Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~mu`]) = 2*Physics:-g_[`~mu`, `~nu`]

(10)

Note that, regarding the spacetime indices, this is a 4x4 matrix, whose elements are in turn 4x4 matrices. Compute first the external 4x4 matrix related to mu and nu

TensorArray(Physics[`*`](Physics[Dgamma][`~mu`], Physics[Dgamma][`~nu`])+Physics[`*`](Physics[Dgamma][`~nu`], Physics[Dgamma][`~mu`]) = 2*Physics[g_][`~mu`, `~nu`])

Matrix(%id = 18446744078587020822)

(11)

Perform now all the matrix operations involved in each of the elements of this 4x4 matrix

"Library:-PerformMatrixOperations(?)"

Matrix(%id = 18446744078743243942)

(12)

By eye everything checks OK.NULL

 

Consider now the following five products of Dirac matrices

e0 := Dgamma[mu]^2

Physics:-`*`(Physics:-Dgamma[mu], Physics:-Dgamma[`~mu`])

(13)

e1 := Dgamma[mu]*Dgamma[`~nu`]*Dgamma[mu]

Physics:-`*`(Physics:-Dgamma[mu], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~mu`])

(14)

e2 := Dgamma[mu]*Dgamma[`~lambda`]*Dgamma[`~nu`]*Dgamma[mu]

Physics:-`*`(Physics:-Dgamma[mu], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~mu`])

(15)

e3 := Dgamma[mu]*Dgamma[`~lambda`]*Dgamma[`~nu`]*Dgamma[`~rho`]*Dgamma[mu]

Physics:-`*`(Physics:-Dgamma[mu], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~mu`])

(16)

e4 := Dgamma[mu]*Dgamma[`~lambda`]*Dgamma[`~nu`]*Dgamma[`~rho`]*Dgamma[`~sigma`]*Dgamma[mu]

Physics:-`*`(Physics:-Dgamma[mu], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~sigma`], Physics:-Dgamma[`~mu`])

(17)

New: the simplification of these products is now implemented

e0 = Simplify(e0)

Physics:-`*`(Physics:-Dgamma[mu], Physics:-Dgamma[`~mu`]) = 4

(18)

Verify this result performing the underlying matrix operations

T := SumOverRepeatedIndices(Physics[`*`](Physics[Dgamma][mu], Physics[Dgamma][`~mu`]) = 4)

Physics:-`*`(Physics:-Dgamma[1], Physics:-Dgamma[`~1`])+Physics:-`*`(Physics:-Dgamma[2], Physics:-Dgamma[`~2`])+Physics:-`*`(Physics:-Dgamma[3], Physics:-Dgamma[`~3`])+Physics:-`*`(Physics:-Dgamma[4], Physics:-Dgamma[`~4`]) = 4

(19)

Library:-PerformMatrixOperations(T)

Matrix(%id = 18446744078553169662) = 4

(20)

The same with the other expressions

e1 = Simplify(e1)

Physics:-`*`(Physics:-Dgamma[mu], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~mu`]) = -2*Physics:-Dgamma[`~nu`]

(21)

SumOverRepeatedIndices(Physics[`*`](Physics[Dgamma][mu], Physics[Dgamma][`~nu`], Physics[Dgamma][`~mu`]) = -2*Physics[Dgamma][`~nu`])

Physics:-`*`(Physics:-Dgamma[1], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~1`])+Physics:-`*`(Physics:-Dgamma[2], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~2`])+Physics:-`*`(Physics:-Dgamma[3], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~3`])+Physics:-`*`(Physics:-Dgamma[4], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~4`]) = -2*Physics:-Dgamma[`~nu`]

(22)

T := TensorArray(Physics[`*`](Physics[Dgamma][1], Physics[Dgamma][`~nu`], Physics[Dgamma][`~1`])+Physics[`*`](Physics[Dgamma][2], Physics[Dgamma][`~nu`], Physics[Dgamma][`~2`])+Physics[`*`](Physics[Dgamma][3], Physics[Dgamma][`~nu`], Physics[Dgamma][`~3`])+Physics[`*`](Physics[Dgamma][4], Physics[Dgamma][`~nu`], Physics[Dgamma][`~4`]) = -2*Physics[Dgamma][`~nu`])

Array(%id = 18446744078695012102)

(23)

Library:-PerformMatrixOperations(T)

Array(%id = 18446744078701714238)

(24)

For e2

e2 = Simplify(e2)

Physics:-`*`(Physics:-Dgamma[mu], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~mu`]) = 4*Physics:-g_[`~lambda`, `~nu`]

(25)

SumOverRepeatedIndices(Physics[`*`](Physics[Dgamma][mu], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~mu`]) = 4*Physics[g_][`~lambda`, `~nu`])

Physics:-`*`(Physics:-Dgamma[1], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~1`])+Physics:-`*`(Physics:-Dgamma[2], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~2`])+Physics:-`*`(Physics:-Dgamma[3], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~3`])+Physics:-`*`(Physics:-Dgamma[4], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~4`]) = 4*Physics:-g_[`~lambda`, `~nu`]

(26)

T := TensorArray(Physics[`*`](Physics[Dgamma][1], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~1`])+Physics[`*`](Physics[Dgamma][2], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~2`])+Physics[`*`](Physics[Dgamma][3], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~3`])+Physics[`*`](Physics[Dgamma][4], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~4`]) = 4*Physics[g_][`~lambda`, `~nu`])

Matrix(%id = 18446744078470204942)

(27)

Library:-PerformMatrixOperations(T)

Matrix(%id = 18446744078550068870)

(28)

For e3 we have

e3 = Simplify(e3)

Physics:-`*`(Physics:-Dgamma[mu], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~mu`]) = -2*Physics:-`*`(Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~lambda`])

(29)

Verify this result,

SumOverRepeatedIndices(Physics[`*`](Physics[Dgamma][mu], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~rho`], Physics[Dgamma][`~mu`]) = -2*Physics[`*`](Physics[Dgamma][`~rho`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~lambda`]))

Physics:-`*`(Physics:-Dgamma[1], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~1`])+Physics:-`*`(Physics:-Dgamma[2], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~2`])+Physics:-`*`(Physics:-Dgamma[3], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~3`])+Physics:-`*`(Physics:-Dgamma[4], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~4`]) = -2*Physics:-`*`(Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~lambda`])

(30)

In this case, with three free spacetime indices lambda, nu, rho, the spacetime components form an array 4x4x4 of 64 components, each of which is a matrix equation

T := TensorArray(Physics[`*`](Physics[Dgamma][1], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~rho`], Physics[Dgamma][`~1`])+Physics[`*`](Physics[Dgamma][2], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~rho`], Physics[Dgamma][`~2`])+Physics[`*`](Physics[Dgamma][3], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~rho`], Physics[Dgamma][`~3`])+Physics[`*`](Physics[Dgamma][4], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~rho`], Physics[Dgamma][`~4`]) = -2*Physics[`*`](Physics[Dgamma][`~rho`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~lambda`]))

Array(%id = 18446744078647326830)

(31)

For instance, the first element is

T[1, 1, 1]

Physics:-`*`(Physics:-Dgamma[1], Physics:-`^`(Physics:-Dgamma[`~1`], 4))+Physics:-`*`(Physics:-Dgamma[2], Physics:-`^`(Physics:-Dgamma[`~1`], 3), Physics:-Dgamma[`~2`])+Physics:-`*`(Physics:-Dgamma[3], Physics:-`^`(Physics:-Dgamma[`~1`], 3), Physics:-Dgamma[`~3`])+Physics:-`*`(Physics:-Dgamma[4], Physics:-`^`(Physics:-Dgamma[`~1`], 3), Physics:-Dgamma[`~4`]) = -2*Physics:-`^`(Physics:-Dgamma[`~1`], 3)

(32)

and it checks OK:

Library:-PerformMatrixOperations(T[1, 1, 1])

Matrix(%id = 18446744078647302614) = Matrix(%id = 18446744078647302974)

(33)

How can you test the 64 components of T all at once?

1. Compute the matrices, without displaying the whole thing, take the elements of the array and remove the indices (ie take the right-hand side); call it M

 

M := map(rhs, ArrayElems(Library:-PerformMatrixOperations(T)))

 

For instance,

M[1]

Matrix(%id = 18446744078629635726) = Matrix(%id = 18446744078629636206)

(34)

Now verify all these matrix equations at once: take the elements of the arrays on each side of the equations and verify that the are the same: we expect for output just {true}

 

map(proc (u) options operator, arrow; evalb(map(ArrayElems, u)) end proc, M)

{true}

(35)

The same for e4

e4 = Simplify(e4)

Physics:-`*`(Physics:-Dgamma[mu], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~sigma`], Physics:-Dgamma[`~mu`]) = 2*Physics:-`*`(Physics:-Dgamma[`~sigma`], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`])+2*Physics:-`*`(Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~sigma`])

(36)

SumOverRepeatedIndices(Physics[`*`](Physics[Dgamma][mu], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~rho`], Physics[Dgamma][`~sigma`], Physics[Dgamma][`~mu`]) = 2*Physics[`*`](Physics[Dgamma][`~sigma`], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~rho`])+2*Physics[`*`](Physics[Dgamma][`~rho`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~sigma`]))

Physics:-`*`(Physics:-Dgamma[1], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~sigma`], Physics:-Dgamma[`~1`])+Physics:-`*`(Physics:-Dgamma[2], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~sigma`], Physics:-Dgamma[`~2`])+Physics:-`*`(Physics:-Dgamma[3], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~sigma`], Physics:-Dgamma[`~3`])+Physics:-`*`(Physics:-Dgamma[4], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~sigma`], Physics:-Dgamma[`~4`]) = 2*Physics:-`*`(Physics:-Dgamma[`~sigma`], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~rho`])+2*Physics:-`*`(Physics:-Dgamma[`~rho`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~sigma`])

(37)

Regarding the spacetime indices this is now an array 4x4x4x4

T := TensorArray(Physics[`*`](Physics[Dgamma][1], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~rho`], Physics[Dgamma][`~sigma`], Physics[Dgamma][`~1`])+Physics[`*`](Physics[Dgamma][2], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~rho`], Physics[Dgamma][`~sigma`], Physics[Dgamma][`~2`])+Physics[`*`](Physics[Dgamma][3], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~rho`], Physics[Dgamma][`~sigma`], Physics[Dgamma][`~3`])+Physics[`*`](Physics[Dgamma][4], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~rho`], Physics[Dgamma][`~sigma`], Physics[Dgamma][`~4`]) = 2*Physics[`*`](Physics[Dgamma][`~sigma`], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~rho`])+2*Physics[`*`](Physics[Dgamma][`~rho`], Physics[Dgamma][`~nu`], Physics[Dgamma][`~lambda`], Physics[Dgamma][`~sigma`]))

Array(%id = 18446744078730196382)

(38)

For instance the first of these 256 matrix equations

T[1, 1, 1, 1]

Physics:-`*`(Physics:-Dgamma[1], Physics:-`^`(Physics:-Dgamma[`~1`], 5))+Physics:-`*`(Physics:-Dgamma[2], Physics:-`^`(Physics:-Dgamma[`~1`], 4), Physics:-Dgamma[`~2`])+Physics:-`*`(Physics:-Dgamma[3], Physics:-`^`(Physics:-Dgamma[`~1`], 4), Physics:-Dgamma[`~3`])+Physics:-`*`(Physics:-Dgamma[4], Physics:-`^`(Physics:-Dgamma[`~1`], 4), Physics:-Dgamma[`~4`]) = 4*Physics:-`^`(Physics:-Dgamma[`~1`], 4)

(39)

verifies OK:

Library:-PerformMatrixOperations(Physics[`*`](Physics[Dgamma][1], Physics[`^`](Physics[Dgamma][`~1`], 5))+Physics[`*`](Physics[Dgamma][2], Physics[`^`](Physics[Dgamma][`~1`], 4), Physics[Dgamma][`~2`])+Physics[`*`](Physics[Dgamma][3], Physics[`^`](Physics[Dgamma][`~1`], 4), Physics[Dgamma][`~3`])+Physics[`*`](Physics[Dgamma][4], Physics[`^`](Physics[Dgamma][`~1`], 4), Physics[Dgamma][`~4`]) = 4*Physics[`^`](Physics[Dgamma][`~1`], 4))

Matrix(%id = 18446744078727227382) = Matrix(%id = 18446744078727227862)

(40)

Now all the 256 matrix equations verified at once as done for e3

 

M := map(rhs, ArrayElems(Library:-PerformMatrixOperations(T)))

map(proc (u) options operator, arrow; evalb(map(ArrayElems, u)) end proc, M)

{true}

(41)

Finally, although there is more work to be done here, let's define some tensors and contract their product with these expressions involving products of Dirac matrices.

 

For example,

Define(A, B)

`Defined as tensors`

 

{A[nu], B[lambda], Physics:-Dgamma[mu], Physics:-Psigma[mu], Physics:-d_[mu], Physics:-g_[mu, nu], Physics:-KroneckerDelta[mu, nu], Physics:-LeviCivita[alpha, beta, mu, nu]}

(42)

Contract with e1 and e2 and simplify

A[nu]*e1; % = Simplify(%)

A[nu]*Physics:-`*`(Physics:-Dgamma[mu], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~mu`]) = -2*A[`~nu`]*Physics:-Dgamma[nu]

(43)

A[nu]*B[lambda]*e2; % = Simplify(%)

A[nu]*B[lambda]*Physics:-`*`(Physics:-Dgamma[mu], Physics:-Dgamma[`~lambda`], Physics:-Dgamma[`~nu`], Physics:-Dgamma[`~mu`]) = 4*B[`~nu`]*A[nu]

(44)

 


 

Download DiracMatricesAndPerformMatrixOperation.mw

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

In a previous Mapleprimes question related to Dirac Matrices, I was asked how to represent the algebra of Dirac matrices with an identity matrix on the right-hand side of  %AntiCommutator(Physics:-Dgamma[j], Physics:-Dgamma[k]) = 2*g[j, k]. Since this is a hot-topic in general, in that, making it work, involves easy and useful functionality however somewhat hidden, not known in general, it passed through my mind that this may be of interest in general. (To reproduce the computations below you need to update your Physics library with the one distributed at the Maplesoft R&D Physics webpage.)

 

restart

with(Physics)

 

First of all, this shows the default algebra rules loaded when you load the Physics package, for the Pauli  and Dirac  matrices

Library:-DefaultAlgebraRules()

%Commutator(Physics:-Psigma[j], Physics:-Psigma[k]) = (2*I)*(Physics:-Psigma[1]*Physics:-LeviCivita[4, j, k, `~1`]+Physics:-Psigma[2]*Physics:-LeviCivita[4, j, k, `~2`]+Physics:-Psigma[3]*Physics:-LeviCivita[4, j, k, `~3`]), %AntiCommutator(Physics:-Psigma[j], Physics:-Psigma[k]) = 2*Physics:-KroneckerDelta[j, k], %AntiCommutator(Physics:-Dgamma[j], Physics:-Dgamma[k]) = 2*Physics:-g_[j, k]

(1)

Now, you can always overwrite these algebra rules.

 

For instance, to represent the algebra of Dirac matrices with an identity matrix on the right-hand side, one can proceed as follows.

First create the identity matrix. To emulate what we do with paper and pencil, where we write I to represent an identity matrix without having to see the actual table 2x2 with the number 1 in the diagonal and a bunch of 0, I will use the old matrix command, not the new Matrix (see more comments on this at the end). One way of entering this identity matrix is

`𝕀` := matrix(4, 4, proc (i, j) options operator, arrow; KroneckerDelta[i, j] end proc)

array( 1 .. 4, 1 .. 4, [( 4, 1 ) = (0), ( 1, 2 ) = (0), ( 2, 3 ) = (0), ( 1, 3 ) = (0), ( 2, 2 ) = (1), ( 4, 2 ) = (0), ( 3, 4 ) = (0), ( 1, 4 ) = (0), ( 3, 1 ) = (0), ( 4, 4 ) = (1), ( 3, 2 ) = (0), ( 1, 1 ) = (1), ( 2, 1 ) = (0), ( 4, 3 ) = (0), ( 3, 3 ) = (1), ( 2, 4 ) = (0)  ] )

(2)

The most important advantage of the old matrix command is that I is of type algebraic and, consequently, this is the important thing, one can operate with it algebraically and its contents are not displayed:

type(`𝕀`, algebraic)

true

(3)

`𝕀`

`𝕀`

(4)

And so, most commands of the Maple library, that only work with objects of type algebraic, will handle the task. The contents are displayed only on demand, for instance using eval

eval(`𝕀`)

array( 1 .. 4, 1 .. 4, [( 4, 1 ) = (0), ( 1, 2 ) = (0), ( 2, 3 ) = (0), ( 1, 3 ) = (0), ( 2, 2 ) = (1), ( 4, 2 ) = (0), ( 3, 4 ) = (0), ( 1, 4 ) = (0), ( 3, 1 ) = (0), ( 4, 4 ) = (1), ( 3, 2 ) = (0), ( 1, 1 ) = (1), ( 2, 1 ) = (0), ( 4, 3 ) = (0), ( 3, 3 ) = (1), ( 2, 4 ) = (0)  ] )

(5)

Returning to the topic at hands: set now the algebra the way you want, with an I matrix on the right-hand side, and without seeing a bunch of 0 and 1

%AntiCommutator(Dgamma[mu], Dgamma[nu]) = 2*g_[mu, nu]*`𝕀`

%AntiCommutator(Physics:-Dgamma[mu], Physics:-Dgamma[nu]) = 2*Physics:-g_[mu, nu]*`𝕀`

(6)

Setup(algebrarules = (%AntiCommutator(Physics[Dgamma][mu], Physics[Dgamma][nu]) = 2*Physics[g_][mu, nu]*`𝕀`))

[algebrarules = {%AntiCommutator(Physics:-Dgamma[mu], Physics:-Dgamma[nu]) = 2*Physics:-g_[mu, nu]*`𝕀`}]

(7)

And that is all.

 

Check it out

(%AntiCommutator = AntiCommutator)(Dgamma[mu], Dgamma[nu])

%AntiCommutator(Physics:-Dgamma[mu], Physics:-Dgamma[nu]) = 2*Physics:-g_[mu, nu]*`𝕀`

(8)

Set now a Dirac spinor; this is how you could do that, step-by-step.

 

Again you can use {vector, matrix, array} or {Vector, Matrix, Array}, and again, if you use the Upper case commands, you always have the components visible, and cannot compute with the object normally since they are not of type algebraic. So I use matrix, not Matrix, and matrix instead of vector so that the Dirac spinor that is both algebraic and matrix, is also displayed in the usual display as a "column vector"

 

_local(Psi)

Setup(anticommutativeprefix = {Psi, psi})

[anticommutativeprefix = {_lambda, psi, :-Psi}]

(9)

In addition, following your question, in this example I explicitly specify the components of the spinor, in any preferred way, for example here I use psi[j]

Psi := matrix(4, 1, [psi[1], psi[2], psi[3], psi[4]])

array( 1 .. 4, 1 .. 1, [( 4, 1 ) = (psi[4]), ( 3, 1 ) = (psi[3]), ( 1, 1 ) = (psi[1]), ( 2, 1 ) = (psi[2])  ] )

(10)

Check it out:

Psi

Psi

(11)

type(Psi, algebraic)

true

(12)

Let's see all this working together by multiplying the anticommutator equation by Psi

(%AntiCommutator(Physics[Dgamma][mu], Physics[Dgamma][nu]) = 2*Physics[g_][mu, nu]*`𝕀`)*Psi

Physics:-`*`(%AntiCommutator(Physics:-Dgamma[mu], Physics:-Dgamma[nu]), Psi) = 2*Physics:-g_[mu, nu]*Physics:-`*`(`𝕀`, Psi)

(13)

Suppose now that you want to see the matrix form of this equation

Library:-RewriteInMatrixForm(Physics[`*`](%AntiCommutator(Physics[Dgamma][mu], Physics[Dgamma][nu]), Psi) = 2*Physics[g_][mu, nu]*Physics[`*`](`𝕀`, Psi))

Physics:-`.`(%AntiCommutator(Physics:-Dgamma[mu], Physics:-Dgamma[nu]), array( 1 .. 4, 1 .. 1, [( 4, 1 ) = (psi[4]), ( 3, 1 ) = (psi[3]), ( 1, 1 ) = (psi[1]), ( 2, 1 ) = (psi[2])  ] )) = 2*Physics:-g_[mu, nu]*Physics:-`.`(array( 1 .. 4, 1 .. 4, [( 4, 1 ) = (0), ( 1, 2 ) = (0), ( 2, 3 ) = (0), ( 1, 3 ) = (0), ( 2, 2 ) = (1), ( 4, 2 ) = (0), ( 3, 4 ) = (0), ( 1, 4 ) = (0), ( 3, 1 ) = (0), ( 4, 4 ) = (1), ( 3, 2 ) = (0), ( 1, 1 ) = (1), ( 2, 1 ) = (0), ( 4, 3 ) = (0), ( 3, 3 ) = (1), ( 2, 4 ) = (0)  ] ), array( 1 .. 4, 1 .. 1, [( 4, 1 ) = (psi[4]), ( 3, 1 ) = (psi[3]), ( 1, 1 ) = (psi[1]), ( 2, 1 ) = (psi[2])  ] ))

(14)

The above has the matricial operations delayed; unleash them

%

Physics:-`.`(%AntiCommutator(Physics:-Dgamma[mu], Physics:-Dgamma[nu]), array( 1 .. 4, 1 .. 1, [( 4, 1 ) = (psi[4]), ( 3, 1 ) = (psi[3]), ( 1, 1 ) = (psi[1]), ( 2, 1 ) = (psi[2])  ] )) = 2*Physics:-g_[mu, nu]*(array( 1 .. 4, 1 .. 1, [( 4, 1 ) = (psi[4]), ( 3, 1 ) = (psi[3]), ( 1, 1 ) = (psi[1]), ( 2, 1 ) = (psi[2])  ] ))

(15)

Or directly perform in one go the matrix operations behind (13)

Library:-PerformMatrixOperations(Physics[`*`](%AntiCommutator(Physics[Dgamma][mu], Physics[Dgamma][nu]), Psi) = 2*Physics[g_][mu, nu]*Physics[`*`](`𝕀`, Psi))

Physics:-`.`(%AntiCommutator(Physics:-Dgamma[mu], Physics:-Dgamma[nu]), array( 1 .. 4, 1 .. 1, [( 4, 1 ) = (psi[4]), ( 3, 1 ) = (psi[3]), ( 1, 1 ) = (psi[1]), ( 2, 1 ) = (psi[2])  ] )) = 2*Physics:-g_[mu, nu]*(array( 1 .. 4, 1 .. 1, [( 4, 1 ) = (psi[4]), ( 3, 1 ) = (psi[3]), ( 1, 1 ) = (psi[1]), ( 2, 1 ) = (psi[2])  ] ))

(16)

REMARK: As shown above, in general, the representation using lowercase commands allows you to use `*` or `.` depending on whether you want to represent the operation or perform the operation. For example this represents the operation, as an exact mimicry of what we do with paper and pencil, both regarding input and output

`𝕀`*Psi

Physics:-`*`(`𝕀`, Psi)

(17)

And this performs the operation

`𝕀`.Psi

array( 1 .. 4, 1 .. 1, [( 4, 1 ) = (psi[4]), ( 3, 1 ) = (psi[3]), ( 1, 1 ) = (psi[1]), ( 2, 1 ) = (psi[2])  ] )

(18)

Or to only displaying the operation

Library:-RewriteInMatrixForm(Physics[`*`](`𝕀`, Psi))

Physics:-`.`(array( 1 .. 4, 1 .. 4, [( 4, 1 ) = (0), ( 1, 2 ) = (0), ( 2, 3 ) = (0), ( 1, 3 ) = (0), ( 2, 2 ) = (1), ( 4, 2 ) = (0), ( 3, 4 ) = (0), ( 1, 4 ) = (0), ( 3, 1 ) = (0), ( 4, 4 ) = (1), ( 3, 2 ) = (0), ( 1, 1 ) = (1), ( 2, 1 ) = (0), ( 4, 3 ) = (0), ( 3, 3 ) = (1), ( 2, 4 ) = (0)  ] ), array( 1 .. 4, 1 .. 1, [( 4, 1 ) = (psi[4]), ( 3, 1 ) = (psi[3]), ( 1, 1 ) = (psi[1]), ( 2, 1 ) = (psi[2])  ] ))

(19)

And to perform all these matricial operations within an algebraic expression,

Library:-PerformMatrixOperations(Physics[`*`](`𝕀`, Psi))

Matrix(%id = 18446744079185513758)

(20)

``

 


 

Download DiracAlgebraWithIdentityMatrix.mw

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


 

Quantum Commutation Rules Basics

 

Pascal Szriftgiser1 and Edgardo S. Cheb-Terrab2 

(1) Laboratoire PhLAM, UMR CNRS 8523, Université Lille 1, F-59655, France

(2) Maplesoft

NULL

NULL

In Quantum Mechanics, in the coordinates representation, the component of the momentum operator along the x axis is given by the differential operator


 "`p__x`=-i `ℏ`(∂)/(∂x)  "

 

The purpose of the exercises below is thus to derive the commutation rules, in the coordinates representation, between an arbitrary function of the coordinates and the related momentum, departing from the differential representation

 

p[n] = -i*`ℏ`*`∂`[n]

These two exercises illustrate how to have full control of the computational process by using different elements of the Maple language, including inert representations of abstract vectorial differential operators, Hermitian operators, algebra rules, etc.

 

These exercises also illustrate a new feature of the Physics package, introduced in Maple 2017, that is getting refined (the computation below requires the Maplesoft updates of the Physics package) which is the ability to perform computations algebraically, using the product operator, but with differential operators, and transform the products into the application of the operators only when we want that, as we do with paper and pencil.

 

%Commutator(g(x, y, z), p_) = I*`ℏ`*Nabla(F(X))

 

restart; with(Physics); with(Physics[Vectors]); interface(imaginaryunit = i)

 

Start setting the problem:

– 

 all ofx, y, z, p__x, p__y, p__z are Hermitian operators

– 

 all of x, y, z commute between each other

– 

 tell the system only that the operators x, y, z are the differentiation variables of the corresponding (differential) operators p__x, p__y, p__z but do not tell what is the form of the operators

 

Setup(mathematicalnotation = true, differentialoperators = {[p_, [x, y, z]]}, hermitianoperators = {p, x, y, z}, algebrarules = {%Commutator(x, y) = 0, %Commutator(x, z) = 0, %Commutator(y, z) = 0}, quiet)

[algebrarules = {%Commutator(x, y) = 0, %Commutator(x, z) = 0, %Commutator(y, z) = 0}, differentialoperators = {[p_, [x, y, z]]}, hermitianoperators = {p, x, y, z}, mathematicalnotation = true]

(1.1)

Assuming F(X) is a smooth function, the idea is to apply the commutator %Commutator(F(X), p_) to an arbitrary ket of the Hilbert space Ket(psi, x, y, z), perform the operation explicitly after setting a differential operator representation for `#mover(mi("p",mathcolor = "olive"),mo("→"))`, and from there get the commutation rule between F(X) and `#mover(mi("p",mathcolor = "olive"),mo("→"))`.

 

Start introducing the commutator, to proceed with full control of the operations we use the inert form %Commutator

alias(X = (x, y, z))

CompactDisplay(F(X))

` F`(X)*`will now be displayed as`*F

(1.2)

%Commutator(F(X), p_)*Ket(psi, X)

Physics:-`*`(%Commutator(F(X), p_), Physics:-Ket(psi, x, y, z))

(1.3)

For illustration purposes only (not necessary), expand this commutator

Physics[`*`](%Commutator(F(X), p_), Physics[Ket](psi, x, y, z)) = expand(Physics[`*`](%Commutator(F(X), p_), Physics[Ket](psi, x, y, z)))

Physics:-`*`(%Commutator(F(X), p_), Physics:-Ket(psi, x, y, z)) = Physics:-`*`(F(X), p_, Physics:-Ket(psi, x, y, z))-Physics:-`*`(p_, F(X), Physics:-Ket(psi, x, y, z))

(1.4)

Note that  `#mover(mi("p",mathcolor = "olive"),mo("→"))`, F(X) and the ket Ket(psi, x, y, z) are operands in the products above and that they do not commute: we indicated that the coordinates x, y, z are the differentiation variables of `#mover(mi("p",mathcolor = "olive"),mo("→"))`. This emulates what we do when computing with these operators with paper and pencil, where we represent the application of a differential operator as a product operation.

 

This representation can be transformed into the (traditional in computer algebra) application of the differential operator when desired, as follows:

Physics[`*`](%Commutator(F(X), p_), Physics[Ket](psi, x, y, z)) = Library:-ApplyProductsOfDifferentialOperators(Physics[`*`](%Commutator(F(X), p_), Physics[Ket](psi, x, y, z)))

Physics:-`*`(%Commutator(F(X), p_), Physics:-Ket(psi, x, y, z)) = Physics:-`*`(F(X), p_(Physics:-Ket(psi, x, y, z)))-p_(Physics:-`*`(F(X), Physics:-Ket(psi, x, y, z)))

(1.5)

Note that, in `#mover(mi("p",mathcolor = "olive"),mo("→"))`(F(X)*Ket(psi, x, y, z)), the application of `#mover(mi("p",mathcolor = "olive"),mo("→"))` is not expanded: at this point nothing is known about  `#mover(mi("p",mathcolor = "olive"),mo("→"))` , it is not necessarily a linear operator. In the Quantum Mechanics problem at hands, however, it is. So give now the operator  `#mover(mi("p",mathcolor = "olive"),mo("→"))` an explicit representation as a linear vectorial differential operator (we use the inert form %Nabla, %Nabla, to be able to proceed with full control one step at a time)

p_ := proc (f) options operator, arrow; -I*`ℏ`*%Nabla(f) end proc

proc (f) options operator, arrow; -Physics:-`*`(Physics:-`*`(I, `ℏ`), %Nabla(f)) end proc

(1.6)

The expression (1.5) becomes

Physics[`*`](%Commutator(F(X), p_), Physics[Ket](psi, x, y, z)) = Physics[`*`](F(X), p_(Physics[Ket](psi, x, y, z)))-p_(Physics[`*`](F(X), Physics[Ket](psi, x, y, z)))

Physics:-`*`(%Commutator(F(X), p_), Physics:-Ket(psi, x, y, z)) = -I*`ℏ`*Physics:-`*`(F(X), %Nabla(Physics:-Ket(psi, x, y, z)))+I*`ℏ`*%Nabla(Physics:-`*`(F(X), Physics:-Ket(psi, x, y, z)))

(1.7)

Activate now the inert operator VectorCalculus[Nabla] and simplify taking into account the algebra rules for the coordinate operators {%Commutator(x, y) = 0, %Commutator(x, z) = 0, %Commutator(y, z) = 0}

Simplify(value(Physics[`*`](%Commutator(F(X), p_), Physics[Ket](psi, x, y, z)) = -I*`ℏ`*Physics[`*`](F(X), %Nabla(Physics[Ket](psi, x, y, z)))+I*`ℏ`*%Nabla(Physics[`*`](F(X), Physics[Ket](psi, x, y, z)))))

Physics:-`*`(Physics:-Commutator(F(X), p_), Physics:-Ket(psi, x, y, z)) = I*`ℏ`*_i*Physics:-`*`(diff(F(X), x), Physics:-Ket(psi, x, y, z))+I*`ℏ`*_j*Physics:-`*`(diff(F(X), y), Physics:-Ket(psi, x, y, z))+I*`ℏ`*_k*Physics:-`*`(diff(F(X), z), Physics:-Ket(psi, x, y, z))

(1.8)

To make explicit the gradient in disguise on the right-hand side, factor out the arbitrary ket Ket(psi, x, y, z)

Factor(Physics[`*`](Physics[Commutator](F(X), p_), Physics[Ket](psi, x, y, z)) = I*`ℏ`*_i*Physics[`*`](diff(F(X), x), Physics[Ket](psi, x, y, z))+I*`ℏ`*_j*Physics[`*`](diff(F(X), y), Physics[Ket](psi, x, y, z))+I*`ℏ`*_k*Physics[`*`](diff(F(X), z), Physics[Ket](psi, x, y, z)))

Physics:-`*`(Physics:-Commutator(F(X), p_), Physics:-Ket(psi, x, y, z)) = I*`ℏ`*Physics:-`*`((diff(F(X), y))*_j+(diff(F(X), z))*_k+(diff(F(X), x))*_i, Physics:-Ket(psi, x, y, z))

(1.9)

Combine now the expanded gradient into its inert (not-expanded) form

subs((Gradient = %Gradient)(F(X)), Physics[`*`](Physics[Commutator](F(X), p_), Physics[Ket](psi, x, y, z)) = I*`ℏ`*Physics[`*`]((diff(F(X), y))*_j+(diff(F(X), z))*_k+(diff(F(X), x))*_i, Physics[Ket](psi, x, y, z)))

Physics:-`*`(Physics:-Commutator(F(X), p_), Physics:-Ket(psi, x, y, z)) = I*`ℏ`*Physics:-`*`(%Gradient(F(X)), Physics:-Ket(psi, x, y, z))

(1.10)

Since (1.10) is true for allKet(psi, x, y, z), this ket can be removed from both sides of the equation. One can do that either taking coefficients (see Coefficients ) or multiplying by the "formal inverse" of this ket, arriving at the (expected) form of the commutation rule between F(X) and `#mover(mi("p",mathcolor = "olive"),mo("→"))`

(Physics[`*`](Physics[Commutator](F(X), p_), Ket(psi, x, y, z)) = I*`ℏ`*Physics[`*`](%Gradient(F(X)), Ket(psi, x, y, z)))*Inverse(Ket(psi, x, y, z))

Physics:-Commutator(F(X), p_) = I*`ℏ`*%Gradient(F(X))

(1.11)

Tensor notation, "[`X__m`,P[n]][-]=i `ℏ` g[m,n]"

 

The computation rule for position and momentum, this time in tensor notation, is performed in the same way, just that, additionally, specify that the space indices to be used are lowercase latin letters, and set the relationship between the differential operators and the coordinates directly using tensor notation.

You can also specify that the metric is Euclidean, but that is not necessary: the default metric of the Physics package, a Minkowski spacetime, includes a 3D subspace that is Euclidean, and the default signature, (- - - +), is not a problem regarding this computation.

 

restart; with(Physics); interface(imaginaryunit = i)

Setup(mathematicalnotation = true, coordinates = cartesian, spaceindices = lowercaselatin, algebrarules = {%Commutator(x, y) = 0, %Commutator(x, z) = 0, %Commutator(y, z) = 0}, hermitianoperators = {P, X, p}, differentialoperators = {[P[m], [x, y, z]]}, quiet)

[algebrarules = {%Commutator(x, y) = 0, %Commutator(x, z) = 0, %Commutator(y, z) = 0}, coordinatesystems = {X}, differentialoperators = {[P[m], [x, y, z]]}, hermitianoperators = {P, p, t, x, y, z}, mathematicalnotation = true, spaceindices = lowercaselatin]

(2.1)

Define now the tensor P[m]

Define(P[m], quiet)

{Physics:-Dgamma[mu], P[m], Physics:-Psigma[mu], Physics:-d_[mu], Physics:-g_[mu, nu], Physics:-gamma_[a, b], Physics:-KroneckerDelta[mu, nu], Physics:-LeviCivita[alpha, beta, mu, nu], Physics:-SpaceTimeVector[mu](X)}

(2.2)

Introduce now the Commutator, this time in active form, to show how to reobtain the non-expanded form at the end by resorting the operands in products

Commutator(X[m], P[n])*Ket(psi, x, y, z)

Physics:-`*`(Physics:-Commutator(Physics:-SpaceTimeVector[m](X), P[n]), Physics:-Ket(psi, x, y, z))

(2.3)

Expand first (not necessary) to see how the operator P[n] is going to be applied

Physics[`*`](Physics[Commutator](Physics[SpaceTimeVector][m](X), P[n]), Ket(psi, x, y, z)) = expand(Physics[`*`](Physics[Commutator](Physics[SpaceTimeVector][m](X), P[n]), Ket(psi, x, y, z)))

Physics:-`*`(Physics:-Commutator(Physics:-SpaceTimeVector[m](X), P[n]), Physics:-Ket(psi, x, y, z)) = Physics:-`*`(Physics:-SpaceTimeVector[m](X), P[n], Physics:-Ket(psi, x, y, z))-Physics:-`*`(P[n], Physics:-SpaceTimeVector[m](X), Physics:-Ket(psi, x, y, z))

(2.4)

Now expand and directly apply in one ago the differential operator P[n]

Physics[`*`](Physics[Commutator](Physics[SpaceTimeVector][m](X), P[n]), Ket(psi, x, y, z)) = Library:-ApplyProductsOfDifferentialOperators(Physics[`*`](Physics[Commutator](Physics[SpaceTimeVector][m](X), P[n]), Ket(psi, x, y, z)))

Physics:-`*`(Physics:-Commutator(Physics:-SpaceTimeVector[m](X), P[n]), Physics:-Ket(psi, x, y, z)) = Physics:-`*`(Physics:-SpaceTimeVector[m](X), P[n](Physics:-Ket(psi, x, y, z)))-P[n](Physics:-`*`(Physics:-SpaceTimeVector[m](X), Physics:-Ket(psi, x, y, z)))

(2.5)

Introducing the explicit differential operator representation for P[n], here again using the inert %d_[n] to keep control of the computations step by step

P[n] := proc (f) options operator, arrow; -I*`ℏ`*%d_[n](f) end proc

proc (f) options operator, arrow; -Physics:-`*`(Physics:-`*`(I, `ℏ`), %d_[n](f)) end proc

(2.6)

The expanded and applied commutator (2.5) becomes

Physics[`*`](Physics[Commutator](Physics[SpaceTimeVector][m](X), P[n]), Ket(psi, x, y, z)) = Physics[`*`](Physics[SpaceTimeVector][m](X), P[n](Ket(psi, x, y, z)))-P[n](Physics[`*`](Physics[SpaceTimeVector][m](X), Ket(psi, x, y, z)))

Physics:-`*`(Physics:-Commutator(Physics:-SpaceTimeVector[m](X), P[n]), Physics:-Ket(psi, x, y, z)) = -I*`ℏ`*Physics:-`*`(Physics:-SpaceTimeVector[m](X), %d_[n](Physics:-Ket(psi, x, y, z)))+I*`ℏ`*%d_[n](Physics:-`*`(Physics:-SpaceTimeVector[m](X), Physics:-Ket(psi, x, y, z)))

(2.7)

Activate now the inert operators %d_[n] and simplify taking into account Einstein's rule for repeated indices

Simplify(value(Physics[`*`](Physics[Commutator](Physics[SpaceTimeVector][m](X), P[n]), Ket(psi, x, y, z)) = -I*`ℏ`*Physics[`*`](Physics[SpaceTimeVector][m](X), %d_[n](Ket(psi, x, y, z)))+I*`ℏ`*%d_[n](Physics[`*`](Physics[SpaceTimeVector][m](X), Ket(psi, x, y, z)))))

Physics:-`*`(Physics:-Commutator(Physics:-SpaceTimeVector[m](X), P[n]), Physics:-Ket(psi, x, y, z)) = I*`ℏ`*Physics:-g_[m, n]*Physics:-Ket(psi, x, y, z)

(2.8)

Since the ket Ket(psi, x, y, z) is arbitrary, we can take coefficients (or multiply by the formal Inverse  of this ket as done in the previous section). For illustration purposes, we use   Coefficients  and note hwo it automatically expands the commutator

Coefficients(Physics[`*`](Physics[Commutator](Physics[SpaceTimeVector][m](X), P[n]), Ket(psi, x, y, z)) = I*`ℏ`*Physics[g_][m, n]*Ket(psi, x, y, z), Ket(psi, x, y, z))

Physics:-`*`(Physics:-SpaceTimeVector[m](X), P[n])-Physics:-`*`(P[n], Physics:-SpaceTimeVector[m](X)) = I*`ℏ`*Physics:-g_[m, n]

(2.9)

One can undo this (frequently undesired) expansion of the commutator by sorting the products on the left-hand side using the commutator between X[m] and P[n]

Library:-SortProducts(Physics[`*`](Physics[SpaceTimeVector][m](X), P[n])-Physics[`*`](P[n], Physics[SpaceTimeVector][m](X)) = I*`ℏ`*Physics[g_][m, n], [P[n], X[m]], usecommutator)

Physics:-Commutator(Physics:-SpaceTimeVector[m](X), P[n]) = I*`ℏ`*Physics:-g_[m, n]

(2.10)

And that is the result we wanted to compute.

 

Additionally, to see this rule in matrix form,

TensorArray(-(Physics[Commutator](Physics[SpaceTimeVector][m](X), P[n]) = I*`ℏ`*Physics[g_][m, n]))

Matrix(%id = 18446744078261558678)

(2.11)

In the above, we use equation (2.10) multiplied by -1 to avoid a minus sign in all the elements of (2.11), due to having worked with the default signature (- - - +); this minus sign is not necessary if in the Setup at the beginning one also sets  signature = `+ + + -`

 

For display purposes, to see this matrix expressed in terms of the geometrical components of the momentum `#mover(mi("p",mathcolor = "olive"),mo("→"))` , redefine the tensor P[n] explicitly indicating its Cartesian components

Define(P[m] = [p__x, p__y, p__z], quiet)

{Physics:-Dgamma[mu], P[m], Physics:-Psigma[mu], Physics:-d_[mu], Physics:-g_[mu, nu], Physics:-gamma_[a, b], Physics:-KroneckerDelta[mu, nu], Physics:-LeviCivita[alpha, beta, mu, nu], Physics:-SpaceTimeVector[mu](X)}

(2.12)

TensorArray(-(Physics[Commutator](Physics[SpaceTimeVector][m](X), P[n]) = I*`ℏ`*Physics[g_][m, n]))

Matrix(%id = 18446744078575996430)

(2.13)

Finally, in a typical situation, these commutation rules are to be taken into account in further computations, and for that purpose they can be added to the setup via

"Setup(?)"

[algebrarules = {%Commutator(x, p__x) = I*`ℏ`, %Commutator(x, p__y) = 0, %Commutator(x, p__z) = 0, %Commutator(x, y) = 0, %Commutator(x, z) = 0, %Commutator(y, p__x) = 0, %Commutator(y, p__y) = I*`ℏ`, %Commutator(y, p__z) = 0, %Commutator(y, z) = 0, %Commutator(z, p__x) = 0, %Commutator(z, p__y) = 0, %Commutator(z, p__z) = I*`ℏ`}]

(2.14)

For example, from herein computations are performed taking into account that

(%Commutator = Commutator)(x, p__x)

%Commutator(x, p__x) = I*`ℏ`

(2.15)

NULL

NULL


 

Download DifferentialOperatorCommutatorRules.mw

 

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

This presentation is about magnetic traps for neutral particles, first achieved for cold neutrons and nowadays widely used in cold-atom physics. The level is that of undergraduate electrodynamics and tensor calculus courses. Tackling this topic within a computer algebra worksheet as shown below illustrates well the kind of advanced computations that can be done today with the Physics package. A new feature minimizetensorcomponents and related functionality is used along the presentation, that requires the updated Physics library distributed at the Maplesoft R&D Physics webpage.
 

 

Magnetic traps in cold-atom physics

 

Pascal Szriftgiser1 and Edgardo S. Cheb-Terrab2 

(1) Laboratoire PhLAM, UMR CNRS 8523, Université Lille 1, F-59655, France

(2) Maplesoft

 

We consider a device constructed with a set of electrical wires fed with constant electrical currents. Those wires can have an arbitrary complex shape. The device is operated in a regime such that, in some region of interest, the moving particles experience a magnetic field that varies slowly compared to the Larmor spin precession frequency. In this region, the effective potential is proportional to the modulus of the field: LinearAlgebra[Norm](`#mover(mi("B"),mo("→"))`(x, y, z)), this potential has a minimum and, close to this minimum, the device behaves as a magnetic trap.

 

 

 

Figure 1: Schematic representation of a Ioffe-Pritchard magnetic trap. It is made of four infinite rods and two coils.

_________________________________________

 

Following [1], we show that:

 

  

a) For a time-independent magnetic field  `#mover(mi("B"),mo("→"))`(x, y, z) in vacuum, up to order two in the relative coordinates X__i = [x, y, z] around some point of interest, the coefficients of orders 1 and 2 in this expansion, `v__i,j` and `c__i,j,k` , respectively the gradient and curvature, contain only 5 and 7 independent components.

  

b) All stationary points of LinearAlgebra[Norm](`#mover(mi("B"),mo("→"))`(x, y, z))^2 (nonzero minima and saddle points) are confined to a curved surface defined by det(`∂`[j](B[i])) = 0.

  

c) The effective potential, proportional to LinearAlgebra[Norm](`#mover(mi("B"),mo("→"))`(x, y, z)), has no maximum, only a minimum.

 

Finally, we draw the stationary condition surface for the case of the widely used Ioffe-Pritchard magnetic trap.

  

 

  

Reference

  

[1] R. Gerritsma and R. J. C. Spreeuw, Topological constraints on magnetostatic traps,  Phys. Rev. A 74, 043405 (2006)

  

 

The independent components of `v__i,j` and `c__i,j,k` entering B[i] = u[i]+v[i, j]*X[j]+(1/2)*c[i, j, k]*X[j]*X[k]

   

The stationary points are within the surface det(`∂`[j](B[i])) = 0

   

U = LinearAlgebra[Norm](`#mover(mi("B",fontweight = "bold"),mo("→",fontweight = "bold"))`)^2 has only minima, no maxima

   

Drawing the Ioffe-Pritchard Magnetic Trap

   


 

MagneticTraps.mw or in pdf format with the sections open: MagneticTraps.pdf

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


 

Maple 2017.1 and 2017.2 introduced several improvements in the solution of PDE & Boundary conditions problems (exact solutions).  Maple 2017.3 includes more improvements in this same area.

 

The following is a set of 25 examples of different PDE & Boundary Conditions problems that are solvable in Maple 2017.3 but not in Maple 2017.2 or previous releases. In the examples that follow, in some cases the PDE is different, in other cases the boundary conditions are of a different kind or solving the problem involves different computational strategies.

As usual, at the end there is a link pointing to the corresponding worksheet.

 

pde[1] := diff(u(x, t), t) = k*(diff(u(x, t), x, x)); bc[1] := u(x, 0) = 6+4*cos(3*Pi*x/L), (D[1](u))(0, t) = 0, (D[1](u))(L, t) = 0

diff(u(x, t), t) = k*(diff(diff(u(x, t), x), x))

 

u(x, 0) = 6+4*cos(3*Pi*x/L), (D[1](u))(0, t) = 0, (D[1](u))(L, t) = 0

(1)

`assuming`([pdsolve([pde[1], bc[1]], u(x, t))], [L > 0, k > 0])

u(x, t) = 6+4*cos(3*Pi*x/L)*exp(-9*k*Pi^2*t/L^2)

(2)

pde[2] := diff(g(t, x), t) = diff(g(t, x), x, x)+a*g(t, x); bc[2] := g(0, x) = 1

diff(g(t, x), t) = diff(diff(g(t, x), x), x)+a*g(t, x)

 

g(0, x) = 1

(3)

pdsolve([pde[2], bc[2]])

g(t, x) = exp(a*t)

(4)

pde[3] := diff(u(x, t), t) = diff(u(x, t), x, x)

bc[3] := u(x, 0) = f(x), u(-1, t) = 0, u(1, t) = 0

diff(u(x, t), t) = diff(diff(u(x, t), x), x)

 

u(x, 0) = f(x), u(-1, t) = 0, u(1, t) = 0

(5)

pdsolve([pde[3], bc[3]], u(x, t))

u(x, t) = Sum((Int(f(x)*sin(n*Pi*x), x = -1 .. 1))*sin(n*Pi*x)*exp(-Pi^2*n^2*t), n = 1 .. infinity)

(6)

pde[4] := diff(u(x, t), t) = k*(diff(u(x, t), x, x)); bc[4] := u(0, t) = 0, u(L, t) = 0, u(x, 0) = piecewise(0 < x and x <= (1/2)*L, 1, (1/2)*L < x and x < L, 2)

pde[4] := diff(u(x, t), t) = k*(diff(u(x, t), x, x))

 

u(0, t) = 0, u(L, t) = 0, u(x, 0) = piecewise(0 < x and x <= (1/2)*L, 1, (1/2)*L < x and x < L, 2)

(7)

`assuming`([pdsolve([pde[4], bc[4]], u(x, t))], [L > 0])

u(x, t) = Sum((2*cos((1/2)*Pi*n)+2+4*(-1)^(1+n))*sin(n*Pi*x/L)*exp(-k*Pi^2*n^2*t/L^2)/(Pi*n), n = 1 .. infinity)

(8)

pde[5] := diff(u(x, t), t) = k*(diff(u(x, t), x, x))

bc[5] := (D[1](u))(0, t) = 0, (D[1](u))(L, t) = 0, u(x, 0) = -3*cos(8*Pi*x/L)

diff(u(x, t), t) = k*(diff(diff(u(x, t), x), x))

 

(D[1](u))(0, t) = 0, (D[1](u))(L, t) = 0, u(x, 0) = -3*cos(8*Pi*x/L)

(9)

`assuming`([pdsolve([pde[5], bc[5]], u(x, t))], [0 < L, 0 < k])

u(x, t) = -3*cos(8*Pi*x/L)*exp(-64*k*Pi^2*t/L^2)

(10)

pde[6] := diff(u(x, t), t) = k*(diff(u(x, t), x, x))+f(x, t); bc[6] := u(0, t) = 0, u(l, t) = 0, u(x, 0) = g(x)

diff(u(x, t), t) = k*(diff(diff(u(x, t), x), x))+f(x, t)

 

u(0, t) = 0, u(l, t) = 0, u(x, 0) = g(x)

(11)

pdsolve([pde[6], bc[6]], u(x, t))

u(x, t) = Sum(2*(Int(g(tau1)*sin(Pi*n1*tau1/l), tau1 = 0 .. l))*sin(Pi*n1*x/l)*exp(-k*Pi^2*n1^2*t/l^2)/l, n1 = 1 .. infinity)+Int(Sum(2*(Int(f(x, tau1)*sin(Pi*n*x/l), x = 0 .. l))*sin(Pi*n*x/l)*exp(-k*Pi^2*n^2*(t-tau1)/l^2)/l, n = 1 .. infinity), tau1 = 0 .. t)

(12)

pde[7] := diff(u(x, t), t) = diff(u(x, t), x, x); bc[7] := u(x, 0) = f(x), u(-1, t) = 0, u(1, t) = 0

diff(u(x, t), t) = diff(diff(u(x, t), x), x)

 

u(x, 0) = f(x), u(-1, t) = 0, u(1, t) = 0

(13)

pdsolve([pde[7], bc[7]], u(x, t))

u(x, t) = Sum((Int(f(x)*sin(n*Pi*x), x = -1 .. 1))*sin(n*Pi*x)*exp(-Pi^2*n^2*t), n = 1 .. infinity)

(14)

pde[8] := diff(u(x, t), t) = k*(diff(u(x, t), x, x))-h*u(x, t); bc[8] := u(x, 0) = sin(x), u(-Pi, t) = u(Pi, t), (D[1](u))(-Pi, t) = (D[1](u))(Pi, t)

diff(u(x, t), t) = k*(diff(diff(u(x, t), x), x))-h*u(x, t)

 

u(x, 0) = sin(x), u(-Pi, t) = u(Pi, t), (D[1](u))(-Pi, t) = (D[1](u))(Pi, t)

(15)

pdsolve([pde[8], bc[8]], u(x, t))

u(x, t) = sin(x)*exp(-t*(k+h))

(16)

pde[9] := diff(u(x, t), t) = diff(u(x, t), x, x)

bc[9] := u(0, t) = 20, u(1, t) = 50, u(x, 0) = 0

diff(u(x, t), t) = diff(diff(u(x, t), x), x)

 

u(0, t) = 20, u(1, t) = 50, u(x, 0) = 0

(17)

pdsolve([pde[9], bc[9]], u(x, t))

u(x, t) = 20+Sum((-40+100*(-1)^n)*sin(n*Pi*x)*exp(-Pi^2*n^2*t)/(Pi*n), n = 1 .. infinity)+30*x

(18)

pde[10] := diff(u(x, y), x, x)+diff(u(x, y), y, y) = 0; bc[10] := u(x, 0) = 0, u(x, 1) = 0, u(0, y) = y^2-y, u(1, y) = 0

diff(diff(u(x, y), x), x)+diff(diff(u(x, y), y), y) = 0

 

u(x, 0) = 0, u(x, 1) = 0, u(0, y) = y^2-y, u(1, y) = 0

(19)

pdsolve([pde[10], bc[10]], u(x, y))

u(x, y) = Sum(-4*((-1)^n-1)*sin(Pi*y*n)*(exp(Pi*n*(2*x-1))-exp(Pi*n))*exp(-Pi*n*(x-1))/((exp(2*Pi*n)-1)*Pi^3*n^3), n = 1 .. infinity)

(20)

pde[11] := diff(u(x, y), x, x)+diff(u(x, y), y, y) = 0

bc[11] := (D[1](u))(0, y) = 0, (D[1](u))(L, y) = 0, u(x, H) = f(x), u(x, 0) = 0

diff(diff(u(x, y), x), x)+diff(diff(u(x, y), y), y) = 0

 

(D[1](u))(0, y) = 0, (D[1](u))(L, y) = 0, u(x, H) = f(x), u(x, 0) = 0

(21)

`assuming`([pdsolve([pde[11], bc[11]], u(x, y))], [0 < L, 0 < H])

u(x, y) = Sum(2*(Int(cos(Pi*x*n/L)*f(x), x = 0 .. L))*cos(Pi*x*n/L)*exp(Pi*n*(H-y)/L)*(exp(2*Pi*y*n/L)-1)/(L*(exp(2*Pi*H*n/L)-1)), n = 1 .. infinity)

(22)

pde[12] := diff(u(x, y), x, x)+diff(u(x, y), y, y) = 0

bc[12] := (D[1](u))(L, y) = 0, u(x, H) = 0, u(x, 0) = 0, (D[1](u))(0, y) = g(y)

diff(diff(u(x, y), x), x)+diff(diff(u(x, y), y), y) = 0

 

(D[1](u))(L, y) = 0, u(x, H) = 0, u(x, 0) = 0, (D[1](u))(0, y) = g(y)

(23)

`assuming`([pdsolve([pde[12], bc[12]], u(x, y))], [0 < x, x <= L, 0 < y, y <= H])

u(x, y) = Sum(-2*(Int(sin(Pi*y*n/H)*g(y), y = 0 .. H))*sin(Pi*y*n/H)*(exp(-Pi*n*(L-2*x)/H)+exp(Pi*L*n/H))*exp(Pi*n*(L-x)/H)/(Pi*n*(exp(2*Pi*L*n/H)-1)), n = 1 .. infinity)

(24)

pde[13] := diff(u(x, y), x, x)+diff(u(x, y), y, y) = 0

bc[13] := (D[1](u))(0, y) = 0, u(x, 0) = 0, u(x, H) = 0, u(L, y) = g(y)

diff(diff(u(x, y), x), x)+diff(diff(u(x, y), y), y) = 0

 

(D[1](u))(0, y) = 0, u(x, 0) = 0, u(x, H) = 0, u(L, y) = g(y)

(25)

`assuming`([pdsolve([pde[13], bc[13]], u(x, y))], [0 < L, 0 < H])

u(x, y) = Sum(2*(Int(sin(Pi*y*n/H)*g(y), y = 0 .. H))*sin(Pi*y*n/H)*exp(Pi*n*(L-x)/H)*(exp(2*Pi*x*n/H)+1)/(H*(exp(2*Pi*L*n/H)+1)), n = 1 .. infinity)

(26)

pde[14] := diff(u(x, y), x, x)+diff(u(x, y), y, y) = 0

bc[14] := u(0, y) = g(y), u(L, y) = 0, (D[2](u))(x, 0) = 0, u(x, H) = 0

diff(diff(u(x, y), x), x)+diff(diff(u(x, y), y), y) = 0

 

u(0, y) = g(y), u(L, y) = 0, (D[2](u))(x, 0) = 0, u(x, H) = 0

(27)

`assuming`([pdsolve([pde[14], bc[14]], u(x, y))], [0 < x, x <= L, 0 < y, y <= H])

u(x, y) = Sum(-2*(exp(-(L-2*x)*(1/2+n)*Pi/H)-exp((1/2)*Pi*(1+2*n)*L/H))*cos((1/2)*Pi*(1+2*n)*y/H)*(Int(cos((1/2)*Pi*(1+2*n)*y/H)*g(y), y = 0 .. H))*exp((1/2)*Pi*(1+2*n)*(L-x)/H)/(H*(exp(Pi*(1+2*n)*L/H)-1)), n = 0 .. infinity)

(28)

pde[15] := diff(u(x, y), x, x)+diff(u(x, y), y, y) = 0

bc[15] := u(0, y) = 0, u(L, y) = 0, u(x, 0) = (D[2](u))(x, 0), u(x, H) = f(x)

diff(diff(u(x, y), x), x)+diff(diff(u(x, y), y), y) = 0

 

u(0, y) = 0, u(L, y) = 0, u(x, 0) = (D[2](u))(x, 0), u(x, H) = f(x)

(29)

`assuming`([pdsolve([pde[15], bc[15]], u(x, y))], [0 < x, x <= L, 0 < y, y <= H])

u(x, y) = Sum(2*exp(Pi*n*(H-y)/L)*sin(Pi*x*n/L)*((Pi*n+L)*exp(2*Pi*y*n/L)+Pi*n-L)*(Int(sin(Pi*x*n/L)*f(x), x = 0 .. L))/(L*((Pi*n+L)*exp(2*Pi*H*n/L)+Pi*n-L)), n = 1 .. infinity)

(30)

pde[16] := diff(u(x, t), t, t) = c^2*(diff(u(x, t), x, x))

bc[16] := u(0, t) = 0, (D[1](u))(L, t) = 0, (D[2](u))(x, 0) = 0, u(x, 0) = f(x)

diff(diff(u(x, t), t), t) = c^2*(diff(diff(u(x, t), x), x))

 

u(0, t) = 0, (D[1](u))(L, t) = 0, (D[2](u))(x, 0) = 0, u(x, 0) = f(x)

(31)

`assuming`([pdsolve([pde[16], bc[16]], u(x, t))], [0 < x, x <= L])

u(x, t) = Sum(2*cos((1/2)*c*Pi*(1+2*n)*t/L)*(Int(sin((1/2)*Pi*(1+2*n)*x/L)*f(x), x = 0 .. L))*sin((1/2)*Pi*(1+2*n)*x/L)/L, n = 0 .. infinity)

(32)

pde[17] := diff(w(x1, x2, x3, t), t) = diff(w(x1, x2, x3, t), x1, x1)+diff(w(x1, x2, x3, t), x2, x2)+diff(w(x1, x2, x3, t), x3, x3)

bc[17] := w(x1, x2, x3, 0) = exp(x1)+x2*x3^5

diff(w(x1, x2, x3, t), t) = diff(diff(w(x1, x2, x3, t), x1), x1)+diff(diff(w(x1, x2, x3, t), x2), x2)+diff(diff(w(x1, x2, x3, t), x3), x3)

 

w(x1, x2, x3, 0) = exp(x1)+x2*x3^5

(33)

pdsolve([pde[17], bc[17]])

w(x1, x2, x3, t) = Sum(t^n*((proc (U) options operator, arrow; diff(diff(U, x1), x1)+diff(diff(U, x2), x2)+diff(diff(U, x3), x3) end proc)@@n)(exp(x1)+x2*x3^5)/factorial(n), n = 0 .. infinity)

(34)

pde[18] := diff(u(x, t), t)+u(x, t)*(diff(u(x, t), x)) = -x

bc[18] := u(x, 0) = x

diff(u(x, t), t)+u(x, t)*(diff(u(x, t), x)) = -x

 

u(x, 0) = x

(35)

pdsolve([pde[18], bc[18]], u(x, t))

u(x, t) = x/tan(t+(1/4)*Pi)

(36)

pde[19] := diff(u(x, t), t)-u(x, t)^2*(diff(u(x, t), x)) = 3*u(x, t)

bc[19] := u(x, 0) = x

diff(u(x, t), t)-u(x, t)^2*(diff(u(x, t), x)) = 3*u(x, t)

 

u(x, 0) = x

(37)

pdsolve([pde[19], bc[19]], u(x, t))

u(x, t) = -exp(3*t)*((-6*exp(6*t)*x+6*x+9)^(1/2)-3)/(exp(6*t)-1)

(38)

pde[20] := diff(u(x, t), t)-x*u(x, t)*(diff(u(x, t), x)) = 0

bc[20] := u(x, 0) = x

diff(u(x, t), t)-x*u(x, t)*(diff(u(x, t), x)) = 0

 

u(x, 0) = x

(39)

pdsolve([pde[20], bc[20]], u(x, t))

u(x, t) = -LambertW(-t*x)/t

(40)

pde[21] := diff(u(x, t), t)+u(x, t)*(diff(u(x, t), x)) = 0

bc[21] := u(x, 0) = x

diff(u(x, t), t)+u(x, t)*(diff(u(x, t), x)) = 0

 

u(x, 0) = x

(41)

pdsolve([pde[21], bc[21]], u(x, t))

u(x, t) = x/(t+1)

(42)

pde[22] := (t+u(x, t))*(diff(u(x, t), x))+t*(diff(u(x, t), t)) = 0; bc[22] := u(x, 1) = x

(t+u(x, t))*(diff(u(x, t), x))+t*(diff(u(x, t), t)) = 0

 

u(x, 1) = x

(43)

pdsolve([pde[22], bc[22]], u(x, t))

u(x, t) = (t-x-1)/(ln(1/t)-1)

(44)

pde[23] := (diff(r*(diff(u(r, theta), r)), r))/r+(diff(u(r, theta), theta, theta))/r^2 = 0; bc[23] := u(a, theta) = f(theta), u(r, -Pi) = u(r, Pi), (D[2](u))(r, -Pi) = (D[2](u))(r, Pi)

(diff(u(r, theta), r)+r*(diff(diff(u(r, theta), r), r)))/r+(diff(diff(u(r, theta), theta), theta))/r^2 = 0

 

u(a, theta) = f(theta), u(r, -Pi) = u(r, Pi), (D[2](u))(r, -Pi) = (D[2](u))(r, Pi)

(45)

`assuming`([pdsolve([pde[23], bc[23]], u(r, theta), HINT = boundedseries)], [a > 0])

u(r, theta) = (1/2)*(2*(Sum(r^n*(sin(n*theta)*(Int(f(theta)*sin(n*theta), theta = -Pi .. Pi))+cos(n*theta)*(Int(f(theta)*cos(n*theta), theta = -Pi .. Pi)))*a^(-n)/Pi, n = 1 .. infinity))*Pi+Int(f(theta), theta = -Pi .. Pi))/Pi

(46)

pde[24] := diff(g(t, x), t) = diff(g(t, x), x, x)+a*g(t, x); bc[24] := g(0, x) = f(x)

diff(g(t, x), t) = diff(diff(g(t, x), x), x)+a*g(t, x)

 

g(0, x) = f(x)

(47)

pdsolve([pde[24], bc[24]])

g(t, x) = exp(a*t)*invfourier(fourier(f(x), x, s1)*exp(-t*s1^2), s1, x)

(48)

pde[25] := diff(u(x, y), x, x)+diff(u(x, y), y, y) = 0; bc[25] := u(0, y) = y*(-y+1), u(1, y) = 0, (D[2](u))(x, 0) = 0, (D[2](u))(x, 1) = 0

diff(diff(u(x, y), x), x)+diff(diff(u(x, y), y), y) = 0

 

u(0, y) = y*(-y+1), u(1, y) = 0, (D[2](u))(x, 0) = 0, (D[2](u))(x, 1) = 0

(49)

pdsolve([pde[25], bc[25]], u(x, y))

u(x, y) = Sum(2*((-1)^n+1)*cos(Pi*y*n)*(exp(Pi*n*(2*x-1))-exp(Pi*n))*exp(-Pi*n*(x-1))/((exp(2*Pi*n)-1)*Pi^2*n^2), n = 1 .. infinity)

(50)

``


 

Download ImprovementsInPdsolve.mw

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

First 7 8 9 10 11 12 13 Last Page 9 of 18