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
  • Yahoo Finance recently discontinued their (largely undocumented) historical stock quote API.

    Previously, you get simply send a HTTP:-Get request like this…


    …and get historical OHLCV (open, high, low, close, trading volume) data in your worksheet (in this case for AAPL between 1 January 2016 and 1 January 2017).

    This no longer works! Yahoo shut the door on this easy-to-use and widely disseminated API.

    You can still download historical stock quotes from Yahoo Finance into Maple, but the process is now somewhat more involved. My complete code in this worksheet but I'll step through the process below.

    If you visit the updated Yahoo Finance website and download historical data for a ticker, you see a URL like this in the status bar of your browser

    Let's examine how ths URL is constructed.

    • period1 and period2 are Unix time stamps for your start and end date
    • interval is the data retrieval interval (this can be either 1d, 1w or 1m)
    • crumb is an alphanumeric code that’s periodically regenerated every time you download new historical data from from the Yahoo Finance website using your browser. Moreover, crumb is paired with a cookie that’s stored by your browser.

    Here’s how to extract and supply the cookie-crumb pair to Yahoo Finance so you can still use Maple to retrieve historical stock quotes

    Send a dummy request to get a cookie-crumb pair


    Grab the crumb from the response

    crumbValue := res[2][i+22..i+32]
                      crumbValue := "btW01FWTBn3"

    Store the cookie from the response

        cookieHeader := "B=702eqhdcmq7cl&b=3&s=0t; expires=Mon,17-Jul-2018 20:27:01 GMT; path=/;

    Construct the URL

    • Your desired start and end dates have to be defined as Unix time stamps. Converting a human readable date (like 1st January 2017) to a Unix timestamp is simple, so I won't cover it here.
    • The previously retrieved crumb has to be added to the URL.
    p1 := 1497709183:
    p2 := 1500301183:
    url:=cat("",ticker,"?period1=",p1,"&period2=",p2,"&interval=1d&events=history&crumb=", crumbValue):

    Send the request to Yahoo Finance, including the cookie in the header

    data:=HTTP:-Get(url,headers = ["Cookie" = cookieHeader])

    Your historical data is now returned

    The historical data is now easily parsed into a matrix.

    Please note that any use of Yahoo Finance has to be consistent with their terms of service.

    The Lattice package to investigate particle accelerator magnet lattices (original post) has been updated to V1.1. This is a significant update, addressing a number of inaccuracies and bugs of V1.0 as well as introducing new elements: Octupole, Fringe effect in dipoles, a MatchedSection allowing to insert a piece of beamline when the details are irrelevant, and a few experimental elements like WireQuad. New functions include the 6th synchrotron-radiation integral I6x, momentum compaction alphap and TaylorMap, which allows to compute the Taylor expansion of  the non-linear map to any degree.

    The code and documentation are available in the Application Center.

    U. Wienands, aka Mac Dude


    A couple of weeks ago, I recorded a short video that discussed various applications for the Statistics:-Fit command. One of the more interesting examples examined how manually adjusting the number of parameters used for a regression model affected the resulting adjusted r-squared value.

    I won’t go into detail about r-squared here, but to briefly summarize: In a linear regression model, r-squared measures the proportion of the variation in a model's dependent variable explained by the independent variables. Basically, r-squared gives a statistical measure of how well the regression line approximates the data. R-squared values usually range from 0 to 1 and the closer it gets to 1, the better it is said that the model performs as it accounts for a greater proportion of the variance (an r-squared value of 1 means a perfect fit of the data). When more variables are added, r-squared values typically increase. They can never decrease when adding a variable; and if the fit is not 100% perfect, then adding a variable that represents random data will increase the r-squared value with probability 1. The adjusted r-squared attempts to account for this phenomenon by adjusting the r-squared value based on the number of independent variables in the model.

    The formula for the adjusted r-squared is:


    n is the number of points in the data sample

    k is the number of independent variables in the model excluding the constant

    By taking the number of independent variables into consideration, the adjusted r-squared behaves different than r-squared; adding more variables doesn’t necessarily produce better fitting models. In many cases, more variables can often lead to lower adjusted r-squared values. In particular, if you add a variable representing random data, the expected change in the adjusted r-squared is 0.

    As such, the adjusted r-squared has a slightly different interpretation than the r-squared. While r-squared is perceived to give an indication of the measure of fit for a chosen regression model, the adjusted r-squared is perceived more as a comparative tool that can be useful for picking variables and designing models that may require less predictors than other models. The science of “gaming” models is a broad topic, so I won’t go into any more detail here, but there’s lots of great information out there if you are looking to learn more (here’s a good place to start).

    The following example adjusts a fitted model by adding or removing variables in order to find better adjusted r-squared values.


    The Import command reads a datafile into a new DataFrame.

    ExperimentalData := Import(FileTools:-JoinPath(["Excel", "ExperimentalData.xls"], base = datadir));

    The dataset has seven variables: time and experimental readings for 6 various concentrations. Removing “time” from our variable set, the convert command converts the values in the DataFrame to a Matrix of values.

    ExMat := convert( ExperimentalData, Matrix )[..,2..7];

    We start by fitting a model that includes predicting variables for each of the columns of data. We mark “Concentration A” as our dependent variable.

    Fit( C + C2*v + C3*w + C4*x + C5*y + C6*z, ExMat[..,2..6], ExMat[..,1], [v,w,x,y,z], summarize=embed ):

    From the above, we can observe that both the r-squared and adjusted r-squared are reasonably high, however only one of the coefficient values has a significant p-value, C3.

    Note: Maple shows all p-values less than 0.05 in bold.

    Let's try to fit the data again, this time keeping the two coefficients with the lowest p-values and the intercept.

    Fit( C + C3*v + C5*w, ExMat[..,[3,5]], ExMat[..,1], [v,w], summarize=embed ):

    From the above, we can see that the r-squared value does go down, however the adjusted r-squared goes up! Let's fit the model one last time to see if removing C5 increases or decreases the adjusted r-squared.

    Fit( C + C3*v, ExMat[..,3], ExMat[..,1], [v], summarize=embed ):

    We can see that the final adjusted r-squared value is lower than the previous two, so we are probably better to keep the additional C5 coefficient value.

    You can see this example as well as a couple of other examples of using the Fit command in the following video:

    You can download the worksheet here:

    We have just released the 3rd edition of the Mathematics Survival Kit – Maple Edition.

    The Math Survival Kit helps students get unstuck when they are stuck. Sometimes students are prevented from solving a problem, not because they haven’t understood the new concept, but because they forget how to do one of the steps, like completely the square, or dealing with log properties.  That’s where this interactive e- book comes in. It gives students the opportunity to review exactly the concept or technique they are stuck on, work through an example, practice as much (or as little) as they want using randomly generated, automatically graded questions on that exact topic, and then continue with their homework.

    This book covers over 150 topics known to cause students grief, from dividing fractions to integration by parts. This 3rd edition contains 31 additional topics, deepening the coverage of mathematical topics at every level, from pre-high school to university.

    See the Mathematics Survival Kit for more information about this updated e-book, including the complete list of topics.


    This post is the answer to this question.

    The procedure named  IntOverDomain  finds a double integral over an arbitrary domain bounded by a non-selfintersecting piecewise smooth curve. The code of the procedure uses the well-known Green's theorem.

    Each section in the border should be specified by a list in the following formats :    
    1. If a section is given parametrically, then  [[f(t), g(t)], t=t1..t2]    
    2. If several consecutive sections of the border or the entire border is a broken line, then it is sufficient to set vertices of this broken line  [ [x1,y1], [x2,y2], .., [xn,yn] ] (for the entire border should be  [xn,yn]=[x1,y1] ).

    Required parameters of the procedure:  f  is an expression in variables  x  and  y , L  is the list of all the sections. The sublists of the list  L  must follow in the positive direction (counterclockwise).

    The code of the procedure:

    IntOverDomain := proc(f, L) 
    local n, i, j, m, yk, yb, xk, xb, Q, p, P, var;
    for i from 1 to n do 
    if type(L[i], listlist(algebraic)) then
    for j from 1 to m-1 do
    yk:=L[i,j+1,2]-L[i,j,2]; yb:=L[i,j,2];
    xk:=L[i,j+1,1]-L[i,j,1]; xb:=L[i,j,1];
    P[i]:=add(p[j],j=1..m-1) else
    var := lhs(L[i, 2]);
    P[i]:=int(eval(Q*diff(L[i,1,2],var),[x=L[i,1,1],y=L[i,1,2]]),L[i,2]) fi;
    add(P[i], i = 1 .. n); 
    end proc:


    Examples of use.

    1. In the first example, we integrate over a quadrilateral:

    with(plottools): with(plots):
    display(polygon([[0,0],[3,0],[0,3],[1,1]], color="LightBlue"));  
    # Visualization of the domain of integration
    IntOverDomain(x^2+y^2, [[[0,0],[3,0],[0,3],[1,1],[0,0]]]);  # The value of integral


    2. In the second example, some sections of the boundary of the domain are curved lines:

    display(inequal({{y<=sqrt(x),y>=sin(Pi*x/3)/2,y<=3-x}, {y>=-2*x+3,y>=sqrt(x),y<=3-x}}, x=0..3,y=0..3, color="LightGreen", nolines), plot([[t,sqrt(t),t=0..1],[t,-2*t+3,t=0..1],[t,3-t,t=0..3],[t,sin(Pi*t/3)/2,t=0..3]], color=black, thickness=2));
    f:=x^2+y^2: L:=[[[t,sin(Pi*t/3)/2],t=0..3],[[3,0],[0,3],[1,1]], [[t,sqrt(t)],t=1..0]]:
    IntOverDomain(f, L);


    3. If  f=1  then the procedure returns the area of the domain:

    IntOverDomain(1, L);  # The area of the above domain



    The help page for  ?procedure has the following example.

    addList := proc(a::list,b::integer)::integer;
        local x,i,s;
        description "add a list of numbers and multiply by a constant";
        for i in a do
        end do;
    end proc:




    Of course it's a bug (actually a typo) here:  s:=s+a[i];   should be   s:=s+i;

    A strange/funny fact is that the result is correct but, for almost any other list the call will produce an error or a wrong answer.

    This app shows the calculation of the final speed of a body after it made contact with a variable force; Taking as reference the initial velocity, mass and graph of the variation of F as a function of time. That is, given the variable forces represented by the lines in the time intervals, we will show the equation of momentum and momentum; With their respective values, followed by their response.

    Lenin Araujo C.

    Ambassador of Maple







    When we first started trying to use Maple to create a maple leaf like the one in the Canada 150 logo, we couldn’t find any references online to the exact geometry, so we went back to basics. With our trusty ruler and protractor, we mapped out the geometry of the maple leaf logo by hand.

    Our first observation was that the maple leaf could be viewed as being comprised of 9 kites. You can read more about the meaning of these shapes on the Canada 150 site (where they refer to the shapes as diamonds).

    We also observed that the individual kites had slightly different scales from one another. The largest kites were numbers 3, 5 and 7; we represented their length as 1 unit of length. Also, each of the kites seemed centred at the origin, but was rotated about the y-axis at a certain angle.

    As such, we found the kites to have the following scales and rotations from the vertical axis:


    1, 9: 0.81 at +/- Pi/2

    2, 8: 0.77 at +/- 2*Pi/5

    3, 5, 7: 1 at +/-Pi/4, 0

    4, 6: 0.93 at +/- Pi/8

    This can be visualized as follows:

    To draw this in Maple we put together a simple procedure to draw each of the kites:

    # Make a kite shape centred at the origin.
    opts := thickness=4, color="#DC2828":
    MakeKite := proc({scale := 1, rotation := 0})
        local t, p, pts, x;

        t := 0.267*scale;
        pts := [[0, 0], [t, t], [0, scale], [-t, t], [0, 0]]:
        p := plot(pts, opts);
        if rotation<>0.0 then
            p := plottools:-rotate(p, rotation);
        end if;
        return p;
    end proc:


    The main idea of this procedure is that we draw a kite using a standard list of points, which are scaled and rotated. Then to generate the sequence of plots:

    shapes := MakeKite(rotation=-Pi/4),
              MakeKite(scale=0.77, rotation=-2*Pi/5),

              MakeKite(scale=0.81, rotation=-Pi/2),
              MakeKite(scale=0.93, rotation=-Pi/8),
              MakeKite(scale=0.93, rotation=Pi/8),
              MakeKite(scale=0.81, rotation=Pi/2),
              MakeKite(scale=0.77, rotation=2*Pi/5),
              plot([[0,-0.5], [0,0]], opts): #Add in a section for the maple leaf stem
    plots:-display(shapes, scaling=constrained, view=[-1..1, -0.75..1.25], axes=box, size=[800,800]);

    This looked pretty similar to the original logo, however the kites 2, 4, 6, and 8 all needed to be moved behind the other kites. This proved somewhat tricky, so we just simply turned on the point probe in Maple and drew in the connected lines to form these points.

    shapes := MakeKite(rotation=-Pi/4),

              MakeKite(scale=0.81, rotation=-Pi/2),
              MakeKite(scale=0.81, rotation=Pi/2),
              plot([[0,-0.5], [0,0]], opts):
    plots:-display(shapes, scaling=constrained, view=[-1..1, -0.75..1.25], axes=box, size=[800,800]);

    Happy Canada Day!

    With Maple, you can create amazing visualizations that go far beyond the standard mathematical plots that you might typically expect (I wince every time I see yet another sine curve).

    At your fingertips, you have

    • plotting primitives that can be assembled in new and novel ways
    • precise control over coloring (yay for ColorTools) and placement
    • an interactive coding environment with inline plots, giving you quick visual feedback over aesthetic changes
    • and a comprehensive mathematical programming language to glue everything together

    Here, I thought I'd share a few of the visualizations I've really enjoyed creating over the last few years (and I'd like to emphasize 'enjoy' - doing this stuff is fun!)

    Let me know if you want any of the worksheets.


    Psychrometric chart with historical weather data for Waterloo, Ontario.


    Ternary plot of the color of gold-silver-copper alloys


    Spectrogram of a violin note played with vibrato


    Colored zoom of the Mandelbrot set


    Reporting dashboard for an Organic Rankine Cycle


    Temperature-entropy plot of an ideal Rankine Cycle


    Quaternion fractal


    Historical sunpot data


    Earthquake data


    African literacy rates

    Do you have Maple content that you want to protect from editing and viewing, while still allowing others to execute the code within and obtain results? In Maple, worksheets can be password protected so the users of your Maple application can benefit from the specialized routines you've created while the details remain hidden.

    The password protection feature can be useful for a variety of situations, such as:

    • • Providing a Maple-based solution while protecting the intellectual property embodied in your algorithms
    • • Ensuring the users of your application can not accidentally make changes that break your code


    To learn more about this feature in Maple, you can download the free Tips & Techniques from the Application Center.

    We have just released an update to Maple, Maple 2017.1.  It includes improvements to the display on high resolution monitors for the debugger, MapleCloud, and help system table of contents. It also contains a variety of small improvements to the math engine, including in limit, series, Physics, typesetting, and PackageTools. This update is available through Tools>Check for Updates in Maple, and is also available from our website on the Maple 2017.1 download page.


    This Saturday is Canada’s 150th birthday. As you can imagine, the country has been paying a lot more attention to this year’s anniversary than our usual low key approach, and as a Canadian company, we at Maplesoft decided to join in the fun.

    And what better way for Maplesoft to celebrate Canada’s birthday than to create a maple leaf in Maple! 

    So here is a maple leaf inspired by the Canada 150 logo, which was created by Ariana Cuvin, a student at the University of Waterloo and former co-op here at Maplesoft:

    Here’s the code to reproduce this plot (more details can be found in this follow up post):





    Know other ways to plot a maple leaf in Maple?  If so, please share them below - we’d love to see them!





    Maple provides a state-of-the-art environment for algebraic and tensorial computations in Physics, with emphasis on ensuring that the computational experience is as natural as possible.


    The theme of the Physics project for Maple 2017 has been the consolidation of the functionality introduced in previous releases, together with significant enhancements and new functionality in General Relativity, in connection with classification of solutions to Einstein's equations and tensor representations to work in an embedded 3D curved space - a new ThreePlusOne  package. This package is relevant in numerical relativity and a Hamiltonian formulation of gravity. The developments also include first steps in connection with computational representations for all the objects entering the Standard Model in particle physics.

    Classification of solutions to Einstein's equations and the Tetrads package


    In Maple 2016, the digitizing of the database of solutions to Einstein's equations  was finished, added to the standard Maple library, with all the metrics from "Stephani, H.; Kramer, D.; MacCallum, M.; Hoenselaers, C.; and Herlt, E., Exact Solutions to Einstein's Field Equations". These metrics can be loaded to work with them, or change them, or searched using g_  (the Physics command representing the spacetime metric that also sets the metric to your choice in one go) or using the command DifferentialGeometry:-Library:-MetricSearch .

    In Maple 2017, the Physics:-Tetrads  package has been vastly improved and extended, now including new commands like PetrovType  and SegreType  to classify these metrics, and the TransformTetrad  now has an option canonicalform to automatically derive a transformation and put the tetrad in canonical form (reorientation of the axis of the local system of references), a relevant step in resolving the equivalence between two metrics.



    Petrov and Segre types, tetrads in canonical form


    Equivalence for Schwarzschild metric (spherical and Kruskal coordinates)


    Formulation of the problem (remove mixed coordinates)


    Solving the Equivalence


    The ThreePlusOne (3 + 1) new Maple 2017 Physics package


    ThreePlusOne , is a package to cast Einstein's equations in a 3+1 form, that is, representing spacetime as a stack of nonintersecting 3-hypersurfaces Σ. This 3+1 description is key in the Hamiltonian formulation of gravity as well as in the study of gravitational waves, black holes, neutron stars, and in general to study the evolution of physical system in general relativity by running numerical simulations as traditional initial value (Cauchy) problems. ThreePlusOne includes computational representations for the spatial metric gamma[i, j] that is induced by g[mu, nu] on the 3-dimensional hypersurfaces, and the related covariant derivative, Christoffel symbols and Ricci and Riemann tensors, the Lapse, Shift, Unit normal and Time vectors and Extrinsic curvature related to the ADM equations.


    The following is a list of the available commands:















    The other four related new Physics  commands:



    Decompose , to decompose 4D tensorial expressions (free and/or contracted indices) into the space and time parts.


    gamma_ , representing the three-dimensional metric tensor, with which the element of spatial distance is defined as  `#mrow(msup(mi("dl"),mrow(mo("&InvisibleTimes;"),mn("2"))),mo("&equals;"),msub(mi("&gamma;",fontstyle = "normal"),mrow(mi("i"),mo("&comma;"),mi("j"))),mo("&InvisibleTimes;"),msup(mi("dx"),mi("i")),mo("&InvisibleTimes;"),msup(mi("dx"),mi("j")))`.


    Redefine , to redefine the coordinates and the spacetime metric according to changes in the signature from any of the four possible signatures(− + + +), (+ − − −), (+ + + −) and ((− + + +) to any of the other ones.


    EnergyMomentum , is a computational representation for the energy-momentum tensor entering Einstein's equations as well as their 3+1 form, the ADMEquations .




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

    `Default differentiation variables for d_, D_ and dAlembertian are: `*{X = (x, y, z, t)}


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


    [coordinatesystems = {X}]



    `Setting lowercaselatin_is letters to represent space indices `


    0, "%1 is not a command in the %2 package", ThreePlusOne, Physics


    `Changing the signature of spacetime from `(`- - - +`)*` to `(`+ + + -`)*` in order to match the signature customarily used in the ADM formalism`


    [ADMEquations, Christoffel3, D3_, ExtrinsicCurvature, Lapse, Ricci3, Riemann3, Shift, TimeVector, UnitNormalVector, gamma3_]


    Note the different color for gamma[mu, nu], now a 4D tensor representing the metric of a generic 3-dimensional hypersurface induced by the 4D spacetime metric g[mu, nu]. All the ThreePlusOne tensors are displayed in black to distinguish them of the corresponding 4D or 3D tensors. The particular hypersurface gamma[mu, nu] operates is parameterized by the Lapse  alpha and the Shift  beta[mu].

    The induced metric gamma[mu, nu]is defined in terms of the UnitNormalVector  n[mu] and the 4D metric g[mu, nu] as


    Physics:-ThreePlusOne:-gamma3_[mu, nu] = Physics:-ThreePlusOne:-UnitNormalVector[mu]*Physics:-ThreePlusOne:-UnitNormalVector[nu]+Physics:-g_[mu, nu]


    where n[mu] is defined in terms of the Lapse  alpha and the derivative of a scalar function t that can be interpreted as a global time function


    Physics:-ThreePlusOne:-UnitNormalVector[mu] = -Physics:-ThreePlusOne:-Lapse*Physics:-D_[mu](t)


    The TimeVector  is defined in terms of the Lapse  alpha and the Shift  beta[mu] and this vector  n[mu] as


    Physics:-ThreePlusOne:-TimeVector[mu] = Physics:-ThreePlusOne:-Lapse*Physics:-ThreePlusOne:-UnitNormalVector[mu]+Physics:-ThreePlusOne:-Shift[mu]


    The ExtrinsicCurvature  is defined in terms of the LieDerivative  of  gamma[mu, nu]


    Physics:-ThreePlusOne:-ExtrinsicCurvature[mu, nu] = -(1/2)*Physics:-LieDerivative[Physics:-ThreePlusOne:-UnitNormalVector](Physics:-ThreePlusOne:-gamma3_[mu, nu])


    The metric gamma[mu, nu]is also a projection tensor in that it projects 4D tensors into the 3D hypersurface Σ. The definition for any 4D tensor that is also a 3D tensor in Σ, can thus be written directly by contracting their indices with gamma[mu, nu]. In the case of Christoffel3 , Ricci3  and Riemann3,  these tensors can be defined by replacing the 4D metric g[mu, nu] by gamma[mu, nu] and the 4D Christoffel symbols GAMMA[mu, nu, alpha] by the ThreePlusOne GAMMA[mu, nu, alpha] in the definitions of the corresponding 4D tensors. So, for instance


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



    Physics:-ThreePlusOne:-Ricci3[mu, nu] = Physics:-d_[alpha](Physics:-ThreePlusOne:-Christoffel3[`~alpha`, mu, nu], [X])-Physics:-d_[nu](Physics:-ThreePlusOne:-Christoffel3[`~alpha`, mu, alpha], [X])+Physics:-ThreePlusOne:-Christoffel3[`~beta`, mu, nu]*Physics:-ThreePlusOne:-Christoffel3[`~alpha`, beta, alpha]-Physics:-ThreePlusOne:-Christoffel3[`~beta`, mu, alpha]*Physics:-ThreePlusOne:-Christoffel3[`~alpha`, nu, beta]



    Physics:-ThreePlusOne:-Riemann3[mu, nu, alpha, beta] = Physics:-g_[mu, lambda]*(Physics:-d_[alpha](Physics:-ThreePlusOne:-Christoffel3[`~lambda`, nu, beta], [X])-Physics:-d_[beta](Physics:-ThreePlusOne:-Christoffel3[`~lambda`, nu, alpha], [X])+Physics:-ThreePlusOne:-Christoffel3[`~lambda`, upsilon, alpha]*Physics:-ThreePlusOne:-Christoffel3[`~upsilon`, nu, beta]-Physics:-ThreePlusOne:-Christoffel3[`~lambda`, upsilon, beta]*Physics:-ThreePlusOne:-Christoffel3[`~upsilon`, nu, alpha])


    When working with the ADM formalism, the line element of an arbitrary spacetime metric can be expressed in terms of the differentials of the coordinates dx^mu, the Lapse , the Shift  and the spatial components of the 3D metric gamma3_ . From this line element one can derive the relation between the Lapse , the spatial part of the Shift , the spatial part of the gamma3_  metric and the g[0, j] components of the 4D spacetime metric.

    For this purpose, define a tensor representing the differentials of the coordinates and an alias  dt = `#msup(mi("dx"),mn("0"))`


    `Defined objects with tensor properties`


    {Physics:-ThreePlusOne:-D3_[mu], Physics:-Dgamma[mu], Physics:-Psigma[mu], Physics:-ThreePlusOne:-Ricci3[mu, nu], Physics:-ThreePlusOne:-Shift[mu], Physics:-d_[mu], dx[mu], Physics:-g_[mu, nu], Physics:-ThreePlusOne:-gamma3_[mu, nu], Physics:-gamma_[i, j], Physics:-ThreePlusOne:-Christoffel3[mu, nu, alpha], Physics:-KroneckerDelta[mu, nu], Physics:-LeviCivita[alpha, beta, mu, nu], Physics:-ThreePlusOne:-Riemann3[mu, nu, alpha, beta], Physics:-ThreePlusOne:-TimeVector[mu], Physics:-ThreePlusOne:-ExtrinsicCurvature[mu, nu], Physics:-ThreePlusOne:-UnitNormalVector[mu], Physics:-SpaceTimeVector[mu](X)}


    "alias(dt = dx[~0]):"

    The expression for the line element in terms of the Lapse  and Shift   is (see [2], eq.(2.123))

    ds^2 = (-Lapse^2+Shift[i]^2)*dt^2+2*Shift[i]*dt*dx[`~i`]+gamma_[i, j]*dx[`~i`]*dx[`~j`]

    ds^2 = (-Physics:-ThreePlusOne:-Lapse^2+Physics:-ThreePlusOne:-Shift[i]*Physics:-ThreePlusOne:-Shift[`~i`])*dt^2+2*Physics:-ThreePlusOne:-Shift[i]*dt*dx[`~i`]+Physics:-gamma_[i, j]*dx[`~i`]*dx[`~j`]


    Compare this expression with the 3+1 decomposition of the line element in an arbitrary system. To avoid the automatic evaluation of the metric components, work with the inert form of the metric %g_

    ds^2 = %g_[mu, nu]*dx[`~mu`]*dx[`~nu`]

    ds^2 = %g_[mu, nu]*dx[`~mu`]*dx[`~nu`]


    Decompose(ds^2 = %g_[mu, nu]*dx[`~mu`]*dx[`~nu`])

    ds^2 = %g_[0, 0]*dt^2+%g_[0, j]*dt*dx[`~j`]+%g_[i, 0]*dt*dx[`~i`]+%g_[i, j]*dx[`~i`]*dx[`~j`]


    The second and third terms on the right-hand side are equal

    op(2, rhs(ds^2 = dt^2*%g_[0, 0]+dt*%g_[0, j]*dx[`~j`]+dt*%g_[i, 0]*dx[`~i`]+%g_[i, j]*dx[`~i`]*dx[`~j`])) = op(3, rhs(ds^2 = dt^2*%g_[0, 0]+dt*%g_[0, j]*dx[`~j`]+dt*%g_[i, 0]*dx[`~i`]+%g_[i, j]*dx[`~i`]*dx[`~j`]))

    %g_[0, j]*dt*dx[`~j`] = %g_[i, 0]*dt*dx[`~i`]


    subs(%g_[0, j]*dt*dx[`~j`] = %g_[i, 0]*dt*dx[`~i`], ds^2 = dt^2*%g_[0, 0]+dt*%g_[0, j]*dx[`~j`]+dt*%g_[i, 0]*dx[`~i`]+%g_[i, j]*dx[`~i`]*dx[`~j`])

    ds^2 = %g_[0, 0]*dt^2+2*%g_[i, 0]*dt*dx[`~i`]+%g_[i, j]*dx[`~i`]*dx[`~j`]


    Taking the difference between this expression and the one in terms of the Lapse  and Shift  we get

    simplify((ds^2 = dt^2*%g_[0, 0]+2*dt*%g_[i, 0]*dx[`~i`]+%g_[i, j]*dx[`~i`]*dx[`~j`])-(ds^2 = (-Physics:-ThreePlusOne:-Lapse^2+Physics:-ThreePlusOne:-Shift[i]*Physics:-ThreePlusOne:-Shift[`~i`])*dt^2+2*Physics:-ThreePlusOne:-Shift[i]*dt*dx[`~i`]+Physics:-gamma_[i, j]*dx[`~i`]*dx[`~j`]))

    0 = (Physics:-ThreePlusOne:-Lapse^2-Physics:-ThreePlusOne:-Shift[i]*Physics:-ThreePlusOne:-Shift[`~i`]+%g_[0, 0])*dt^2+2*dx[`~i`]*(%g_[i, 0]-Physics:-ThreePlusOne:-Shift[i])*dt-dx[`~i`]*dx[`~j`]*(Physics:-gamma_[i, j]-%g_[i, j])


    Taking coefficients, we get equations for the Shift , the Lapse  and the spatial components of the metric gamma3_

    eq[1] := coeff(coeff(rhs(0 = (Physics:-ThreePlusOne:-Lapse^2-Physics:-ThreePlusOne:-Shift[i]*Physics:-ThreePlusOne:-Shift[`~i`]+%g_[0, 0])*dt^2+2*dx[`~i`]*(%g_[i, 0]-Physics:-ThreePlusOne:-Shift[i])*dt-dx[`~i`]*dx[`~j`]*(Physics:-gamma_[i, j]-%g_[i, j])), dt), dx[`~i`]) = 0

    2*%g_[i, 0]-2*Physics:-ThreePlusOne:-Shift[i] = 0


    eq[2] := coeff(rhs(0 = (Physics:-ThreePlusOne:-Lapse^2-Physics:-ThreePlusOne:-Shift[i]*Physics:-ThreePlusOne:-Shift[`~i`]+%g_[0, 0])*dt^2+2*dx[`~i`]*(%g_[i, 0]-Physics:-ThreePlusOne:-Shift[i])*dt-dx[`~i`]*dx[`~j`]*(Physics:-gamma_[i, j]-%g_[i, j])), dt^2)

    Physics:-ThreePlusOne:-Lapse^2-Physics:-ThreePlusOne:-Shift[i]*Physics:-ThreePlusOne:-Shift[`~i`]+%g_[0, 0]


    eq[3] := coeff(coeff(rhs(0 = (Physics:-ThreePlusOne:-Lapse^2-Physics:-ThreePlusOne:-Shift[i]*Physics:-ThreePlusOne:-Shift[`~i`]+%g_[0, 0])*dt^2+2*dx[`~i`]*(%g_[i, 0]-Physics:-ThreePlusOne:-Shift[i])*dt-dx[`~i`]*dx[`~j`]*(Physics:-gamma_[i, j]-%g_[i, j])), dx[`~i`]), dx[`~j`]) = 0

    -Physics:-gamma_[i, j]+%g_[i, j] = 0


    Using these equations, these quantities can all be expressed in terms of the time and space components of the 4D metric g[0, 0] and g[i, j]

    isolate(eq[1], Shift[i])

    Physics:-ThreePlusOne:-Shift[i] = %g_[i, 0]


    isolate(eq[2], Lapse^2)

    Physics:-ThreePlusOne:-Lapse^2 = Physics:-ThreePlusOne:-Shift[i]*Physics:-ThreePlusOne:-Shift[`~i`]-%g_[0, 0]


    isolate(eq[3], gamma_[i, j])

    Physics:-gamma_[i, j] = %g_[i, j]




    [1] Landau, L.D., and Lifshitz, E.M. The Classical Theory of Fields, Course of Theoretical Physics Volume 2, fourth revised English edition. Elsevier, 1975.


    [2] Alcubierre, M., Introduction to 3+1 Numerical Relativity, International Series of Monographs on Physics 140, Oxford University Press, 2008.


    [3] Baumgarte, T.W., Shapiro, S.L., Numerical Relativity, Solving Einstein's Equations on a Computer, Cambridge University Press, 2010.


    [4] Gourgoulhon, E., 3+1 Formalism and Bases of Numerical Relativity, Lecture notes, 2007,


    [5] Arnowitt, R., Dese, S., Misner, C.W., The Dynamics of General Relativity, Chapter 7 in Gravitation: an introduction to current research (Wiley, 1962),



    Examples: Decompose, gamma_



    Setup(mathematicalnotation = true)

    [mathematicalnotation = true]


    Define  now an arbitrary tensor A


    `Defined objects with tensor properties`


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


    So A^mu is a 4D tensor with only one free index, where the position of the time-like component is the position of the different sign in the signature, that you can query about via


    [signature = `- - - +`]


    To perform a decomposition into space and time, set - for instance - the lowercase latin letters from i to s to represent spaceindices and

    Setup(spaceindices = lowercase_is)

    [spaceindices = lowercaselatin_is]


    Accordingly, the 3+1 decomposition of A^mu is


    Array(%id = 18446744078724512334)


    The 3+1 decomposition of the inert representation %g_[mu,nu] of the 4D spacetime metric; use the inert representation when you do not want the actual components of the metric appearing in the output

    Decompose(%g_[mu, nu])

    Matrix(%id = 18446744078724507998)


    Note the position of the component %g_[0, 0], related to the trailing position of the time-like component in the signature "(- - - +)".

    Compare the decomposition of the 4D inert with the decomposition of the 4D active spacetime metric


    g[mu, nu] = (Matrix(4, 4, {(1, 1) = -1, (1, 2) = 0, (1, 3) = 0, (1, 4) = 0, (2, 1) = 0, (2, 2) = -1, (2, 3) = 0, (2, 4) = 0, (3, 1) = 0, (3, 2) = 0, (3, 3) = -1, (3, 4) = 0, (4, 1) = 0, (4, 2) = 0, (4, 3) = 0, (4, 4) = 1}))


    Decompose(g_[mu, nu])

    Matrix(%id = 18446744078724494270)


    Note that in general the 3D space part of g[mu, nu] is not equal to the 3D metric gamma[i, j] whose definition includes another term (see [1] Landau & Lifshitz, eq.(84.7)).


    Physics:-gamma_[i, j] = -Physics:-g_[i, j]+Physics:-g_[0, i]*Physics:-g_[0, j]/Physics:-g_[0, 0]


    The 3D space part of -g[`~mu`, `~nu`] is actually equal to the 3D metric "gamma[]^(i,j)"


    Physics:-gamma_[`~i`, `~j`] = -Physics:-g_[`~i`, `~j`]


    To derive the formula  for the covariant components of the 3D metric, Decompose into 3+1 the identity

    %g_[`~alpha`, `~mu`]*%g_[mu, beta] = KroneckerDelta[`~alpha`, beta]

    %g_[`~alpha`, `~mu`]*%g_[mu, beta] = Physics:-KroneckerDelta[beta, `~alpha`]


    To the side, for illustration purposes, these are the 3 + 1 decompositions, first excluding the repeated indices, then excluding the free indices

    Eq := Decompose(%g_[`~alpha`, `~mu`]*%g_[mu, beta] = Physics[KroneckerDelta][beta, `~alpha`], repeatedindices = false)

    Matrix(%id = 18446744078132963318)


    Eq := Decompose(%g_[`~alpha`, `~mu`]*%g_[mu, beta] = Physics[KroneckerDelta][beta, `~alpha`], freeindices = false)

    %g_[0, beta]*%g_[`~alpha`, `~0`]+%g_[i, beta]*%g_[`~alpha`, `~i`] = Physics:-KroneckerDelta[beta, `~alpha`]


    Compare with a full decomposition

    Eq := Decompose(%g_[`~alpha`, `~mu`]*%g_[mu, beta] = Physics[KroneckerDelta][beta, `~alpha`])

    Matrix(%id = 18446744078724489454)


    Eq is a symmetric matrix of equations involving non-contracted occurrences of `#msup(mi("g"),mrow(mn("0"),mo("&comma;"),mn("0")))`, `#msup(mi("g"),mrow(mi("j"),mo("&comma;"),mo("0")))` and `#msup(mi("g"),mrow(mi("j"),mo("&comma;"),mi("i")))`. Isolate, in Eq[1, 2], `#msup(mi("g"),mrow(mi("j"),mo("&comma;"),mo("0")))`, that you input as %g_[~j, ~0], and substitute into Eq[1, 1]

    "isolate(Eq[1, 2], `%g_`[~j, ~0]);"

    %g_[`~j`, `~0`] = -%g_[i, 0]*%g_[`~j`, `~i`]/%g_[0, 0]


    subs(%g_[`~j`, `~0`] = -%g_[i, 0]*%g_[`~j`, `~i`]/%g_[0, 0], Eq[1, 1])

    -%g_[0, k]*%g_[i, 0]*%g_[`~j`, `~i`]/%g_[0, 0]+%g_[i, k]*%g_[`~j`, `~i`] = Physics:-KroneckerDelta[k, `~j`]


    Collect `#msup(mi("g"),mrow(mi("j"),mo("&comma;"),mi("i")))`, that you input as %g_[~j, ~i]

    collect(-%g_[0, k]*%g_[i, 0]*%g_[`~j`, `~i`]/%g_[0, 0]+%g_[i, k]*%g_[`~j`, `~i`] = Physics[KroneckerDelta][k, `~j`], %g_[`~j`, `~i`])

    (-%g_[0, k]*%g_[i, 0]/%g_[0, 0]+%g_[i, k])*%g_[`~j`, `~i`] = Physics:-KroneckerDelta[k, `~j`]


    Since the right-hand side is the identity matrix and, from , `#msup(mi("g"),mrow(mi("i"),mo("&comma;"),mi("j")))` = -`#msup(mi("&gamma;",fontstyle = "normal"),mrow(mi("i"),mo("&comma;"),mi("j")))`, the expression between parenthesis, multiplied by -1, is the reciprocal of the contravariant 3D metric `#msup(mi("&gamma;",fontstyle = "normal"),mrow(mi("i"),mo("&comma;"),mi("j")))`, that is the covariant 3D metric gamma[i, j], in accordance to its definition for the signature `- - - +`


    Physics:-gamma_[i, j] = -Physics:-g_[i, j]+Physics:-g_[0, i]*Physics:-g_[0, j]/Physics:-g_[0, 0]





    [1] Landau, L.D., and Lifshitz, E.M. The Classical Theory of Fields, Course of Theoretical Physics Volume 2, fourth revised English edition. Elsevier, 1975.

    Example: Redefine


    Tensors in Special and General Relativity


    A number of relevant changes happened in the tensor routines of the Physics package, towards making the routines pack more functionality both for special and general relativity, as well as working more efficiently and naturally, based on Maple's Physics users' feedback collected during 2016.

    New functionality


    Implement conversions to most of the tensors of general relativity (relevant in connection with functional differentiation)


    New setting in the Physics Setup  allows for specifying the cosmologicalconstant and a default tensorsimplifier



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

    In this file you will be able to observe and analyze how the exercises and problems of Kinematics and Dynamics are solved using the commands and operators through a very well-structured syntax; Allowing me to save time and use it in interpretation. I hope you can share and spread to break the traditional and unnecessary myths. Only for Engineering and Science. Share if you like.

    In Spanish.

    Lenin Araujo Castillo

    Ambassador of Maple

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