Maxima Crash Course

1 Maxima as a programming environment

Computer algebra systems generally follow a functional programming paradigm. A number of helper and higher-order functions are available in Maxima.
Map and apply are typical of the functional paradigm.
--> map("=",[x,y],[a,b]);

\[\]\[\tag{%o1} \left[ x\mathop{=}a\mathop{,}y\mathop{=}b\right] \]

Example: Use "apply" to sum the elements of a list.
--> apply("+",[1,2,3]);

\[\]\[\tag{%o2} 6\]

Makelist creates a list using a general form by running through an index.
--> makelist(j^2,j,1,10);

\[\]\[\tag{%o3} \left[ 1\mathop{,}4\mathop{,}9\mathop{,}16\mathop{,}25\mathop{,}36\mathop{,}49\mathop{,}64\mathop{,}81\mathop{,}100\right] \]

Evaluation vs substitution.
--> f: sin(x^5)$
--> ev(f,x=5);

\[\]\[\tag{%o5} \sin{(3125)}\]

--> subst(5,x,f);

\[\]\[\tag{%o6} \sin{(3125)}\]

--> subst(x,5,f);

\[\]\[\tag{%o7} \sin{\left( {{x}^{x}}\right) }\]

--> subst(cos(t),5,f);

\[\]\[\tag{%o8} \sin{\left( {{x}^{\cos{(t)}}}\right) }\]

--> g: sqrt(x^2 + 2·x + 3);

\[\]\[\tag{} \sqrt{{{x}^{2}}\mathop{+}2 x\mathop{+}3}\]

--> ratsubst(14+y^2,x^2 + 2·x,g);

\[\]\[\tag{%o10} \sqrt{{{y}^{2}}\mathop{+}17}\]

--> subst(14+y^2,x^2+2·x,g);

\[\]\[\tag{%o11} \sqrt{{{x}^{2}}\mathop{+}2 x\mathop{+}3}\]

The last two show that subst only works on atoms or complete subexpressions.
--> atom(x);

\[\]\[\tag{%o12} true\]

--> atom(cot(x));

\[\]\[\tag{%o13} false\]

The quote operator prevents evaluation of an expression but not simplification.
--> 'sin(0);

\[\]\[\tag{%o1} 0\]

Maxima consider the above simplification rather than evaluation.
--> f(x) := 2·x;

\[\]\[\tag{%o17} \mathop{f}(x)\mathop{:=}2 x\]

--> 'f(2);

\[\]\[\tag{%o18} \mathop{f}(2)\]

Get documentation on a particular function.
--> ? map;

\[\]\[ -- Function: map (< f> , < expr\_ 1> , ..., < expr\_ n> ) \]\[ Returns an expression whose leading operator is the same as that of \]\[ the expressions < expr\_ 1> , ..., < expr\_ n> but whose subparts are the \]\[ results of applying < f> to the corresponding subparts of the \]\[ expressions. < f> is either the name of a function of n arguments \]\[ or is a ‘lambda’ form of n arguments. \]\[ ‘maperror’ - if ‘false’ will cause all of the mapping functions to \]\[ (1) stop when they finish going down the shortest < expr\_ i> if not \]\[ all of the < expr\_ i> are of the same length and (2) apply < f> to \]\[ [< expr\_ 1> , < expr\_ 2> , ...] if the < expr\_ i> are not all the same type \]\[ of object. If ‘maperror’ is ‘true’ then an error message will be \]\[ given in the above two instances. \]\[ One of the uses of this function is to ‘map’ a function (e.g. \]\[ ‘partfrac’) onto each term of a very large expression where it \]\[ ordinarily wouldn’t be possible to use the function on the entire \]\[ expression due to an exhaustion of list storage space in the course \]\[ of the computation. \]\[ See also ‘scanmap’, ‘maplist’, ‘outermap’, ‘matrixmap’ and ‘apply’. \]\[ (\% i1) map(f,x+a\ensuremath{\cdot}y+b\ensuremath{\cdot}z); \]\[ (\% o1) f(b z) + f(a y) + f(x) \]\[ (\% i2) map(lambda([u],partfrac(u,x)),x+1/(x\textasciicircum3+4\ensuremath{\cdot}x\textasciicircum2+5\ensuremath{\cdot}x+2)); \]\[ 1 1 1 \]\[ (\% o2) ----- - ----- + -------- + x \]\[ x + 2 x + 1 2 \]\[ (x + 1) \]\[ (\% i3) map(ratsimp, x/(x\textasciicircum2+x)+(y\textasciicircum2+y)/y); \]\[ 1 \]\[ (\% o3) y + ----- + 1 \]\[ x + 1 \]\[ (\% i4) map("=",[a,b],[-0.5,3]); \]\[ (\% o4) [a = - 0.5, b = 3] \]\[ There are also some inexact matches for `map'. \]\[ Try `?? map' to see them.\]

\[\]\[\tag{%o19} true\]

2 Calculus

Evaluate limits.
--> limit((1+1/n)^n,n,inf);

\[\]\[\tag{%o20} \% e\]

Derivatives.
--> f: sin(2·x) + exp(−x)$
--> diff(f,x);

\[\]\[\tag{%o22} 2 \cos{\left( 2 x\right) }\mathop{-}{{\% e}^{-x}}\]

Partial derivatives and the total differential.
--> g: exp(−x·y)·cos(sqrt(y))$
--> diff(g,y);

\[\]\[\tag{%o24} \mathop{-}\left( \frac{\sin{\left( \sqrt{y}\right) } {{\% e}^{-\left( x y\right) }}}{2 \sqrt{y}}\right) \mathop{-}x \cos{\left( \sqrt{y}\right) } {{\% e}^{-\left( x y\right) }}\]

--> diff(g);

\[\]\[\tag{%o25} \left( \mathop{-}\left( \frac{\sin{\left( \sqrt{y}\right) } {{\% e}^{-\left( x y\right) }}}{2 \sqrt{y}}\right) \mathop{-}x \cos{\left( \sqrt{y}\right) } {{\% e}^{-\left( x y\right) }}\right) \mathop{del}(y)\mathop{-}\cos{\left( \sqrt{y}\right) } y {{\% e}^{-\left( x y\right) }} \mathop{del}(x)\]

One could use the total differential and 'solve' function to do implicit differentiation.
--> expr: x^2 + x·y^2 = 9$
--> diff(expr);

\[\]\[\tag{%o27} 2 x y \mathop{del}(y)\mathop{+}\left( {{y}^{2}}\mathop{+}2 x\right) \mathop{del}(x)\mathop{=}0\]

--> solve(%,del(y));

\[\]\[\tag{%o28} \left[ \mathop{del}(y)\mathop{=}\mathop{-}\left( \frac{\left( {{y}^{2}}\mathop{+}2 x\right) \mathop{del}(x)}{2 x y}\right) \right] \]

--> %/del(x);

\[\]\[\tag{%o29} \left[ \frac{\mathop{del}(y)}{\mathop{del}(x)}\mathop{=}\mathop{-}\left( \frac{{{y}^{2}}\mathop{+}2 x}{2 x y}\right) \right] \]

Antiderivatives and definite integrals analytically.
--> f: cos(x)·sin(x)^2;

\[\]\[\tag{} \cos{(x)} {{\sin{(x)}}^{2}}\]

--> g: exp(−t)/t;

\[\]\[\tag{} \frac{{{\% e}^{-t}}}{t}\]

--> integrate(f,x);

\[\]\[\tag{%o32} \frac{{{\sin{(x)}}^{3}}}{3}\]

--> integrate(f,x,0,1);

\[\]\[\tag{%o33} \frac{{{\sin{(1)}}^{3}}}{3}\]

--> integrate(g,t);

\[\]\[\tag{%o34} \mathop{-}\mathop{gamma\_ incomplete}\left( 0\mathop{,}t\right) \]

This integral was recently reported as a bug case.
--> f: log(5 − 4·cos(x))$
--> integrate(f,x,0,%pi);

\[\]\[\tag{%o37} \mathop{-}\left( \frac{\pi \log{(9)}\mathop{+}\pi \log{\left( \frac{9}{4}\right) }\mathop{-}2 \% i {{\ensuremath{\mathrm{li}}}_2}\left( \mathop{-}\left( \frac{1}{2}\right) \right) \mathop{-}2 \% i {{\ensuremath{\mathrm{li}}}_2}\left( \mathop{-}2\right) \mathop{-}\% i {{\pi }^{2}}}{2}\right) \mathop{+}\pi \log{(9)}\mathop{+}\frac{\% i \left( 6 {{\log{(2)}}^{2}}\mathop{-}12 {{\ensuremath{\mathrm{li}}}_2}(2)\mathop{-}{{\pi }^{2}}\right) }{12}\mathop{-}\frac{\% i {{\pi }^{2}}}{3}\]

--> quad_qags(f,x,0,%pi);

\[\]\[\tag{%o38} \left[ 4.3551721806072035\mathop{,}3.895181855479842 {{10}^{-12}}\mathop{,}63\mathop{,}0\right] \]

--> first(%);

\[\]\[\tag{%o39} 4.3551721806072035\]

--> guess_exact_value(%);

\[\]\[\tag{%o40} 2 \pi \log{(2)}\]

The quote operator has a cousin which is frequently useful. Consider the following example of defining one function in terms of a derivative.
--> y(t) := t^2;

\[\]\[\tag{%o41} \mathop{y}(t)\mathop{:=}{{t}^{2}}\]

--> h(t) := diff(y(t),t);

\[\]\[\tag{%o42} \mathop{h}(t)\mathop{:=}\frac{d}{d t} \mathop{y}(t)\]

--> h(2);

\[\]\[diff: second argument must be a variable; found 2 \]\[\neq 0: h(t=2) \]\[\] \texttt{%error -- an error. To debug this try: debugmode(true);}\[\]

--> h(t) := ''(diff(y(t),t));

\[\]\[\tag{%o44} \mathop{h}(t)\mathop{:=}2 t\]

--> h(2);

\[\]\[\tag{%o45} 4\]

This works because the quote-quote operator causes evaluation when the function is defined, so the derivative is computed before substituting t=2.
How to calculate a path integral? Compute along γ from 0 to 1:
--> f: 4·x^3·y$
--> γ: [1−3·t,2−3·t]$
--> dγ: diff(γ,t);

\[\]\[\tag{\ensuremath{\gamma} \left[ \mathop{-}3\mathop{,}\mathop{-}3\right] \]

--> dS: sqrt(dγ . dγ);

\[\]\[\tag{} 3 \sqrt{2}\]

--> map("=",[x,y],γ);

\[\]\[\tag{%o50} \left[ x\mathop{=}1\mathop{-}3 t\mathop{,}y\mathop{=}2\mathop{-}3 t\right] \]

--> sublis(%,f);

\[\]\[\tag{%o51} 4 {{\left( 1\mathop{-}3 t\right) }^{3}} \left( 2\mathop{-}3 t\right) \]

--> integrate(%·dS,t,0,1);

\[\]\[\tag{%o52} \frac{57 \sqrt{2}}{5}\]

3 Matrices and linear algebra

--> load("linearalgebra")$
Define matrices.
--> A: matrix([2,4],[1,2]);

\[\]\[\tag{} \begin{pmatrix}2 & 4\\ 1 & 2\end{pmatrix}\]

--> B: matrix([6/7,4,29],[11,0,1],[0,3,0]);

\[\]\[\tag{} \begin{pmatrix}\frac{6}{7} & 4 & 29\\ 11 & 0 & 1\\ 0 & 3 & 0\end{pmatrix}\]

--> I: ident(5);

\[\]\[\tag{} \begin{pmatrix}1 & 0 & 0 & 0 & 0\\ 0 & 1 & 0 & 0 & 0\\ 0 & 0 & 1 & 0 & 0\\ 0 & 0 & 0 & 1 & 0\\ 0 & 0 & 0 & 0 & 1\end{pmatrix}\]

The dot operator " . " is used for matrix multiplication and the standard inner product.
--> A . A;

\[\]\[\tag{%o57} \begin{pmatrix}8 & 16\\ 4 & 8\end{pmatrix}\]

The ordinary multiplication operator performs the Hadamard product (element-wise).
--> D: identfor(B)·[x,1,1/y];

\[\]\[\tag{} \begin{pmatrix}x & 0 & 0\\ 0 & 1 & 0\\ 0 & 0 & \frac{1}{y}\end{pmatrix}\]

Difference between matrix^n and matrix^^n:
--> A^2;

\[\]\[\tag{%o59} \begin{pmatrix}4 & 16\\ 1 & 4\end{pmatrix}\]

--> A^^2;

\[\]\[\tag{%o60} \begin{pmatrix}8 & 16\\ 4 & 8\end{pmatrix}\]

There are two available functions to compute the determinant of a square matrix. "newdet" is supposed to be faster.
--> determinant(B);

\[\]\[\tag{%o61} \frac{6681}{7}\]

Matrix inversion.
--> invert(B);

\[\]\[\tag{%o62} \begin{pmatrix}\mathop{-}\left( \frac{7}{2227}\right) & \frac{203}{2227} & \frac{28}{6681}\\ 0 & 0 & \frac{1}{3}\\ \frac{77}{2227} & \mathop{-}\left( \frac{6}{2227}\right) & \mathop{-}\left( \frac{308}{6681}\right) \end{pmatrix}\]

The Moore-Penrose pseudoinverse is available via the linearalgebra share package.
--> moore_penrose_pseudoinverse(A);

\[\]\[\tag{%o64} \begin{pmatrix}\frac{2}{25} & \frac{1}{25}\\ \frac{4}{25} & \frac{2}{25}\end{pmatrix}\]

--> rank(A);

\[\]\[\tag{%o65} 1\]

The matrix exponential.
--> matrixexp(A);

\[\]\[\tag{%o66} \begin{pmatrix}\frac{{{\% e}^{4}}\mathop{+}1}{2} & {{\% e}^{4}}\mathop{-}1\\ \frac{{{\% e}^{4}}\mathop{-}1}{4} & \frac{{{\% e}^{4}}\mathop{+}1}{2}\end{pmatrix}\]

More generally, one may compute an analytic function of a square matrix, provided the corresponding power series of the matrix argument is convergent.
--> load("diag")$
--> mat_function(sin,A);

\[\]\[\tag{%o68} \begin{pmatrix}\frac{\sin{(4)}}{2} & \sin{(4)}\\ \frac{\sin{(4)}}{4} & \frac{\sin{(4)}}{2}\end{pmatrix}\]

Eigensystem of a matrix.
--> load("eigen")$
--> eigenvalues(A);

\[\]\[\tag{%o70} \left[ \left[ 0\mathop{,}4\right] \mathop{,}\left[ 1\mathop{,}1\right] \right] \]

--> eigenvectors(A);

\[\]\[\tag{%o71} \left[ \left[ \left[ 0\mathop{,}4\right] \mathop{,}\left[ 1\mathop{,}1\right] \right] \mathop{,}\left[ \left[ \left[ 1\mathop{,}\mathop{-}\left( \frac{1}{2}\right) \right] \right] \mathop{,}\left[ \left[ 1\mathop{,}\frac{1}{2}\right] \right] \right] \right] \]

Rather than have a separate function to diagonalize a matrix, Maxima just has the Jordan normal form via the diag package.
--> load("diag")$
--> jordan(B);

\[\]\[\tag{%o73} \]

--> dispJordan(%);

\[\]\[\tag{%o74} \begin{pmatrix}\frac{2315 \left( \frac{\sqrt{3} \% i}{2}\mathop{+}\frac{\mathop{-}1}{2}\right) }{147 {{\left( \frac{\sqrt{59719596463}}{98 {{3}^{\frac{3}{2}}}}\mathop{+}\frac{331991}{686}\right) }^{\frac{1}{3}}}}\mathop{+}{{\left( \frac{\sqrt{59719596463}}{98 {{3}^{\frac{3}{2}}}}\mathop{+}\frac{331991}{686}\right) }^{\frac{1}{3}}} \left( \frac{\mathop{-}1}{2}\mathop{-}\frac{\sqrt{3} \% i}{2}\right) \mathop{+}\frac{2}{7} & 0 & 0\\ 0 & {{\left( \frac{\sqrt{59719596463}}{98 {{3}^{\frac{3}{2}}}}\mathop{+}\frac{331991}{686}\right) }^{\frac{1}{3}}} \left( \frac{\sqrt{3} \% i}{2}\mathop{+}\frac{\mathop{-}1}{2}\right) \mathop{+}\frac{2315 \left( \frac{\mathop{-}1}{2}\mathop{-}\frac{\sqrt{3} \% i}{2}\right) }{147 {{\left( \frac{\sqrt{59719596463}}{98 {{3}^{\frac{3}{2}}}}\mathop{+}\frac{331991}{686}\right) }^{\frac{1}{3}}}}\mathop{+}\frac{2}{7} & 0\\ 0 & 0 & {{\left( \frac{\sqrt{59719596463}}{98 {{3}^{\frac{3}{2}}}}\mathop{+}\frac{331991}{686}\right) }^{\frac{1}{3}}}\mathop{+}\frac{2315}{147 {{\left( \frac{\sqrt{59719596463}}{98 {{3}^{\frac{3}{2}}}}\mathop{+}\frac{331991}{686}\right) }^{\frac{1}{3}}}}\mathop{+}\frac{2}{7}\end{pmatrix}\]

There are some high-performance numerical linear algebra functions via LAPACK (Linear Algebra PAckage), which was translated to Common Lisp from FORTRAN. (Matrices for these functions should have float or integer entries).
--> load("lapack")$
--> dgesvd(float(B));

\[\]\[\tag{%o76} \left[ \left[ 29.324140322483178\mathop{,}10.954531767377269\mathop{,}2.971148259913842\right] \mathop{,}false\mathop{,}false\right] \]

Compare to analytically calculating the singular values of B.
--> transpose(B) . B;

\[\]\[\tag{%o77} \begin{pmatrix}\frac{5965}{49} & \frac{24}{7} & \frac{251}{7}\\ \frac{24}{7} & 25 & 116\\ \frac{251}{7} & 116 & 842\end{pmatrix}\]

--> eigenvalues(%);

\[\]\[\tag{%o81} \]

Should there be complex eigenvalues? These are eigenvalues of a real symmetric matrix, so they should have zero imaginary part. The reason this occurs is that imaginary numbers may arise in using Cardano's formula for the cubic equation. We show that the imaginary part is indeed zero:
--> trigreduce(rectform(%));

\[\]\[\tag{%o82} \]

--> float(%);

\[\]\[\tag{%o83} \left[ \left[ 120.00176624247769\mathop{,}8.827721982389022\mathop{,}859.9052056526842\right] \mathop{,}\left[ 1.0\mathop{,}1.0\mathop{,}1.0\right] \right] \]

And the singular values are the square roots:
--> sqrt(%);

\[\]\[\tag{%o84} \left[ \left[ 10.954531767377265\mathop{,}2.971148259913837\mathop{,}29.324140322483185\right] \mathop{,}\left[ 1.0\mathop{,}1.0\mathop{,}1.0\right] \right] \]

4 Some fun: Lorenz attractor

--> bounds(x,y,z) := (
   makelist(
        [apply(min, l),apply(max, l)],l,[x,y,z]
    )
)$
--> round_bounds(L) := (
makelist(
    [floor(L[i][1]), ceiling(L[i][2])], i, 1, 3)
)$
--> kill(x,y,z)$
[sigma, rho, beta]: [10, 28, 8/3]$
eq: [sigma·(y−x), x·(rho−z)−y, x·y−beta·z]$
sol: rk(eq, [x, y, z], [1, 0, 0], [t, 0, 50, 1/80])$
len: length(sol)$
x: makelist(sol[k][2], k, len)$
y: makelist(sol[k][3], k, len)$
z: makelist(sol[k][4], k, len)$
L: round_bounds(bounds(x,y,z))$
wxanimate_framerate:20$
wxplot_size:[840,700]$
with_slider_draw3d(
    k,makelist(i,i,1,len,20),
    title = concat("Lorenz Attractor"),
    points_joined=true,
    color = midnight_blue,
    point_type = −1,
   xrange = [L[1][1],L[1][2]], yrange = [L[2][1],L[2][2]], zrange = [0,L[3][2]],
    points(firstn(x,k),firstn(y,k),firstn(z,k)),
    view = [75,35],
    axis_3d = true
)$

Created with wxMaxima.

The source of this Maxima session can be downloaded here.