You can make all sorts of mathematical and graphical computations with this web interface to Maxima.

Write in the text area de code you want to be executed and then press the button below.

Examples included in this page can be copied and pasted (CTRL-C and CTRL-V) in the form before pressing the button. You can make your own variations, since the form is editable.

The documentation in this page is very basic. For a deeper introduction on Maxima's syntax:

- Official reference: Maxima manual
- From official web site: Tutorials

The output obtained is like a dialog. Our questions are preceded by label **%i** followed by numbers and the results are labeled with **%o** followed by the corresponding number.

Integers

You can use integer numbers as large as you need. Here is the factorial of 100,

100! ;

As in any other programming language, you can add (**+**), subtract (**-**), multiply (*****), divide (**/**), exponentiate (**^**), and calculate square roots (**sqrt**).

3^4 * (5-91+(-24)) + 75 / sqrt(25) ;

Let's know the list of divisors of a given integer and let's factor it in primes,

/* we first set variable a to the number */ a: 169996284 $ divisors(a) ; factor(a) ;

Ask Maxima whether 1502933 is a prime number.

primep(1502933);

We use **gcd** to calculate the *greatest common divisor* of two numbers, and **lcm** for the *least common multiple* of an arbitrary set of numbers, but take into account that for **lcm** we have to load first package **functs**.

gcd(16479936, 43375500000); load("functs")$ lcm(16479936, 43375500000, 111375);

Function **divide** returns a list with two elements, the quotient and the rest (modulus) of an integer division.

divide(34, 5);

Rationals

Maxima likes exact results. Here we make rational operations. Sentences must end with semicolons (**;**).

1 + 5 * (5 - 7/9)^4 ;

But it is also possible to get floting point approximations. Here we store in variable **x** (with a colon **:**) a rational expression and then we call function **float**.

x : 1 + 5 * (5 - 7/9)^4 ; float(x);

The same result can be obtained with only one sentence.

float (1 + 5 * (5 - 7/9)^4) ;

Another example with variables. Since comments are written between symbols **/*** and ***/**, this is a self explained example.

/* Bind variable 'a' to a value */ a : 8 / 24; /* Bind variable 'b' to another value */ b: 1/5 ; /* Squaring the sum of 'a' and 'b' */ (a + b)^2;

Function **float** returns numbers in double precission. Let's calculate number π with five thousand decimal digits. First, we give global variable
**fpprec** the desired precission and then we call function **bfloat**. Remember that the Maxima symbol for π is **%pi**. We also calculate the square root of 2 with 5000 decimal digits.

fpprec : 5000; bfloat(%pi); /* The square root of two */ bfloat( sqrt(2) );

Complex

The imaginary unit is represented as **%i**. With this symbol, we can write any complex number and make some operations with them.

z1: 3 + 5*%i; z2: 2/3 - %i; z1 + z2; z1 - z2;

Sometimes, Maxima returns complex expressions with parentheses. If we want to simplify these expressions, we can call function **expand**.

z1: 3 + 5*%i $ z2: 2/3 - %i $ 3*z1 - 7*z2; expand(3*z1 - 7*z2); z1 * z2; expand(z1 * z2);

We can transform complex numbers in binomial form into polar form and vice-versa.

z: 3 + 5*%i ; p: polarform(z); b: rectform(p);

Let's divide two complex numbers and give the result in binomial form,

d: (3 + 5*%i) / (7/8 - 1/4*%i) ; rectform(d);

Given a complex number, we can obtain its real and imaginary parts, as well as its modulus and angle.

z: 7/8 - 1/4*%i; realpart(z); imagpart(z); cabs(z); carg(z);

Conversion factors

Package **ezunits** provides all we need for unit conversions. Let's see some examples.

load("ezunits") $ /* convert m/s to km/hour */ 5 ` m/s `` km/hour ; /* convert N/cm^2 to Pa */ 32 ` N/cm^2 `` Pa ; /* convert m to inch in two steps */ d : 10 ` m $ d `` inch; /* list all conversion factors */ for k in known_unit_conversions do print(k)$

You can declare your own conversion factors.

load("ezunits") $ declare_unit_conversion ( Pa = 9.869*10^-6 * atm, bar = 100000 * Pa ) $ 5.3 ` atm `` bar ;

We can convert any physical units to fundamental units.

load("ezunits") $ f: 7 ` kg*km/minute^2; /* tell me the dimensions */ dimensions (f); /* tell me the fundamental units */ fundamental_units (f); /* make conversion */ f `` fundamental_units (f);

Polynomials

We can make operations with variables. Look how we transform algebraic expressions.

1 + x + 5 * y - 2/3 * x + y + 2 ; (5*x^2 * u * r) / (u^2 * t) ;

Expanding products.

expand( (x+1)^2 * (x+5) * (x-a) ) ;

Factoring polynomials.

factor(u^6 - 14*u^5 - 23*u^4 + 808*u^3 - 320*u^2 - 12800*u) ;

Make use of function **subst** if you want to calculate the numerical value of an algebraic expression. Square brackets in Maxima are very important, since we can build lists with them.

f : 6.67 * 10^-11 * m1 * m2 / r^2 ; subst([m1 = 4, m2 = 5, r = 7], f) ; subst([m1 = mass, m2 = mass], f) ;

Equations

In general, we use function **solve**. Two arguments are needed, the list of equations followed by the list of unknowns.

A simple linear equation.

solve([5 * (x + 1/2) = 1 + x /3], [x]) ;

Linear equation with parameter.

solve([5 * (x + 1/2) = k + x /3], [x]) ;

A quadratic equation. Remember that **sqrt** is the square root and **%i** is the Maxima symbol for the imaginary unit.

solve([2*x^2 - 3*x + 5 = 0], [x]) ;

A system of equations. We can edit more than one line for entering a sentence. We also call for the floating point result.

sol : solve([5*x + y/8 = 4, x^2 = y], [x,y]) ; float(sol);

When algebraic methods are impossible to apply, we can make use of numerical routines. Library **mnewton** can solve equations and systems of equations by Newton's method. We call function **mnewton** with three arguments: the list of equations, the list of unknowns, and the coordinates of the seed point.

load("mnewton") $ /* equation with sigle unknown */ mnewton([2*a^a = 5],[a],[1]); /* multiple unknowns */ mnewton([2*3^u-v/u-5 = 0, u+2^v = 4], [u, v], [2, 2]);

Inequalities

Solving inequalities require package **fourier_elim**. We first solve a linear inequality with only one unknown. Function **fourier_elim** needs two arguments, the list of inequalities, and the list of unknowns.

load("fourier_elim") $ ine : (5*x-16)/6 + (x+8)/12 < (x+1)/3; fourier_elim([ine], [x]);

An inequality of degree four.

load("fourier_elim") $ fourier_elim([x^4+5*x^3+5*x^2-5*x-6 > 0],[x]);

A system of inequalities with one unknown.

load("fourier_elim") $ fourier_elim( [(x+2)/4 <= x/2-3, (8-x)/3 < (1+x)/2-1], [x]);

A system of inequalities with two unknowns.

load("fourier_elim") $ des1: 3*x - 5*y < 2; des2: x+y> (1+x)/3; des3: y < 2; fourier_elim([des1, des2, des3],[x,y]);

Vectors

Vectors are represented as list with square brackets. We can make with them several common operations.

v1: [1, 2, 3] $ v2: [a, b, c] $ /* addition and subtraction */ v1 + v2; v1 - v2; /* multiplication by a number */ k * v1; /* dot multiplication */ v1 . v2;

Cross multiplication is defined in package **vector3d**.

load(vector3d) $ cross([a, b, c] , [x, y, z]);

Here is an example of how to define a new function in Maxima. We want to write a function to calculate the modulus of a vector; see that functions are defined by means of the **:=** operator. Ok, don't worry if you don't understand the inner part of the function.

modulus(v):= sqrt(apply("+", v^2)) $ modulus([a, 5, c]);

Matrices

Function **matrix** builds matrices; its arguments are the rows of the matrix in list form. We define two matrices and then we add and subtract them.

a : matrix([1, 2, 3], [4, 5, 6]) ; b : matrix([7, 1/9, 2], [z, x+y, -5]) ; a + b ; a - b ;

The operator for the matrix product is the dot (**.**)

a : matrix([1, 2, 3], [4, 5, 6]) ; b : matrix([7, 1/9], [z, -5], [1, 0]) ; a . b ;

The operator for matrix exponentiation is **^^**

a : matrix([1, 2, 3], [4, 5, 6], [7, 8, 9]) ; a^^3 ;

Other operations with matrices.

a : matrix([1, 4, 3], [7, 5, 6], [3, 8, 9]) ; /* Inverse */ invert(a) ; /* Transpose */ transpose(a) ; /* Determinant */ determinant(a) ;

There are some special matrices, like the identity and the zero matrices.

ident(4); zeromatrix(4,2);

Given a matrix, we calculate its characteristic polynomial, eigenvalues, and eigenvectors.

A: matrix([1,3], [5,7]); /* the second argument is the name of the variable */ charpoly(A, x); /* the eigenvalues and their multiplicities */ eigenvalues(A); /* the eigenvalues and their corresponding vectors */ eigenvectors(A);

Limits

Function **limit** lets you calculate limits.

Infinity is represented by **inf** and minus infinity by **minf**.

f: (x-1) / (2*x + 3) ; /* x approaches 1 */ limit(f, x, 1); /* x approaches positive infinity */ limit(f, x, inf); /* x approaches negative infinity */ limit(f, x, minf);

You can approximate the limit from the right or from the left.

g: 1 / (x-1); limit(g, x, 1); /* approaching 1 from the right */ limit(g, x, 1, plus); /* approaching 1 from the left */ limit(g, x, 1, minus);

Derivatives

We have to call function **diff**.

fun: sin(x) * x^4 ; diff(f, x);

Calculating derivatives of order greater than one.

/* order three */ diff(log(x)*x^x, x, 3);

If we want to calculate implicit derivatives, we have to tell Maxima which variables depend on the independent variable.

/* y is a constant */ i: x^4 - x*y = cos(y^2+x); diff(i, x); /* now, y depends on x */ depends(y,x) $ diff(i, x);

Integrals

We call **integrate** for indefinite and definite integrals.

/* indefinite integral */ integrate(sin (2 + 3 * x), x) ; /* definite integral */ integrate( x * exp(x), x, 1, 2) ;

When symbolic methods can't be applied, we can make use of some numerical routines. One of them is **quad_qag**; it needs five arguments, the integrand, the independent variable, the two limits of integration, and an integer between 1 and 6, which is a parameter of the procedure (3 is a good choice). What we obtain is a list of numbers, the first one is the value of the definite integral, and the second is an estimation of the error.

quad_qag(exp(sin(x)),x,2,7,3);

Differential equations

When working with differential equations, we need to tell Maxima which variable depends on the independent variable. We have to call function **ode2** in order to solve ordinary differential equations.

depends(y, x) $ /* an equation with separate variables */ e: (x-1)*y^3+(y-1)*x^3*diff(y,x)=0; r: ode2(e,y,x); /* apply initial conditions */ ic1(r,x=2,y=-3);

Solving a differential equation of Bernoulli type.

depends(y, x) $ r: ode2(diff(y,x)-y+sqrt(y),y,x); /* isolate dependent variable with solve */ solve(r, y);

Together with function **ode2**, package **contrib_ode** adds more capabilities for solving differential equations, please consult available documentation.

Solving a second order differential equation, namely the *simple harmonic motion* equation.

depends(x, t) $ e: diff(x,t,2) + 3 * x = 0; r: ode2(e, x, t); /* apply initial conditions */ ic2(r, t=0, x=1, diff(x,t)=1/2); /* or apply boundary conditions */ s: bc2(r, t=0, x=1, x=3, t=-1); /* let's plot it! */ draw2d(explicit(rhs(s), t, 0, 10))$

With function **desolve** we can solve systems of differential equations, but in this case we have to explicity indicate function dependencies.

desolve( [diff(f(x),x)=3*f(x)-2*g(x), diff(g(x),x)=2*f(x)-2*g(x)], [f(x),g(x)]);

If we want to apply initial conditions, we have to declare them with **atvalue** before calling **desolve**.

atvalue(f(x),x=0,-1)$ atvalue(g(x),x=0,3)$ atvalue(h(x),x=0,1)$ desolve( [diff(f(x),x)=f(x)+g(x)+3*h(x), diff(g(x),x)=g(x)-2*h(x), diff(h(x),x)=f(x)+h(x)], [f(x),g(x),h(x)]);

When symbolic methods are not allowed, we can make use of the Runge-Kutta (**rk**) method. Here, we also plot the resulting function.

/* solving dy/dt = 2 y^2-exp(-3 t) */ sol: rk(2*y^2-exp(-3*t), y, 0, [t, 0, 10, 0.2]) $ draw2d( points_joined = true, points(sol))$

The best reference for package **draw** is the html document A Maxima-Gnuplot interface, which contains lots of examples. You can find more examples in the Maxima manual.

2D graphics

2D drawings are made with function **draw2d**.

Options are written in the form **option = value** and graphic objects as **object(arg _{1}, arg_{2}, ...)**.

Here is an explicit function in red. Explicit functions need an expression, the name of the independent variable, and the extrems of the domine to be plotted.

draw2d ( color = red, explicit(x^3-2*x^2+x-1, x, -3, 3) ) $

Options are read sequentially, and they affect the objects declared after them. Her is a an explicit function together with a parametric one. A grid is added in the background.

draw2d ( explicit(x^3-2*x^2+x-1, x, -3, 3), grid = true, color = red, line_width = 4, parametric(sin(t), t^2, t, -2*%pi, 2*%pi) ) $

We can plot isolated points and polygonal lines. Coordinate points can be passed to object **points** as a list of pairs, writing the list of abscissas followed by the list of ordinates, or giving only the list of ordinates. See also that colors can be declared in hexadecimal form, and how point styles, except the color, affect the two last objects.

draw2d ( /* default style */ points([[2,3], [7,4], [0,5], [5,2/3], [7,-1], [8,3], [-2,10]]), point_type = circle, point_size = 3, color = brown, points_joined = true, points([3, 4, 6, 8], [-2, -5, 3, 1]), color = "#ff670d", points([2, 6, 2, 8, 9]) ) $

An implicit function, a function defined in polar coordinates and a polygon.

draw2d( /* global options */ title = "Write here your title", xlabel = "x-axis", ylabel = "y-axis", grid = true, dimensions = [500,500], /* implicit function */ key = "Implicit", implicit(y^2=x^3-2*x+1, x, -4,4, y, -4,6), /* polar coordinates */ key = "Polar", color = red, polar(10/theta,theta,1,3*%pi), /* polygon */ key = "Polygon", fill_color = sea_green, color = aquamarine, line_width = 6, polygon([[1,1],[3,0],[4,5],[0,3]]) ) $

3D graphics

3D drawings are made with function **draw3d**.

Three styles for drawing the same explicit surface.

/* wired surface */ draw3d( explicit(20*exp(-x^2-y^2)-10,x,-3,3,y,-3,3)) $ /* hidden surface with red lines */ draw3d( surface_hide = true, color = red, explicit(20*exp(-x^2-y^2)-10,x,-3,3,y,-3,3)) $ /* surface colored as a function of point coordinates */ draw3d( enhanced3d = [sin(x*y), x, y, z], explicit(20*exp(-x^2-y^2)-10,x,-3,3,y,-3,3)) $

We can also plot parametric surfaces.

draw3d( /* fix the number of nodes to be calculated */ xu_grid = 50, yv_grid = 15, surface_hide = true, parametric_surface( cos(a)*(3+b*cos(a/2)), sin(a)*(3+b*cos(a/2)), b*sin(a/2), a,-%pi,%pi,b,-1,1) )$

A parametric curve. See how the color depends on parameter *u*.

draw3d( /* fix the number of points to be calculated */ nticks = 60, title = "Colored parametric curve", line_width = 2, enhanced3d = (u-1)^2, dimensions = [500, 500], /* image size */ parametric(cos(5*u)^2,sin(7*u),u-2,u,0,2)) $

An implicit surface.

draw3d( enhanced3d = [x-z,x,y,z], palette = gray, dimensions = [500, 600], implicit(min(x^2+y^2+z^2,2*x^2+2*y^2)=1, x,-1.5,1.5,y,-1.5,1.5,z,-1.5,1.5)) $

Isolated points and polygonal lines.

xx: [1,2,3,4,5]$ yy: [4,5,7,1,2]$ zz: [3,2,1,9,4]$ xyz: [[10,15,12],[11,13,16],[16,10,11],[15,15,15]]$ draw3d( point_type = filled_square, color = blue, points(xx,yy,zz), point_type = circle, color = green, points_joined = true, points(xyz)) $

A surface defined in spherical coordinates.

draw3d( color = green, surface_hide = true, axis_3d = false, xtics = none, ytics = none, ztics = none, title = "Nautilus", spherical(a+z,a,0,3*%pi,z,0,%pi))$

Family of surfaces defined in cylindrical coordinates.

draw3d( proportional_axes = xyz, dimensions = [500,500], surface_hide = true, axis_3d = false, xtics = none, ytics = none, ztics = none, color = blue, cylindrical(z,z,-2,2,a,0,2*%pi), /*cone*/ color = brown, cylindrical(3,z,-2,2,az,0,%pi), /*cylinder*/ color = green, cylindrical(sqrt(25-z^2),z,-5,5,a,0,%pi) /*sphere*/ )$

Multiplots

In order to draw several independent scenes, you have to call function **draw**, passing as many 2D and 3D graphics as you need, with objects **gr2d** and **gr3d**, respectively. Here is a simple example, where we fix the dimensions to our needs.

draw( dimensions = [400, 600], gr2d( color = red, explicit(sin(x), x, -5, 5)), gr3d( color = black, surface_hide = true, explicit(sin(x)*cos(y), x, -5, 5, y, -5, 5)) ) $

By default, multiple plots are arranged in only one column, but you can change this behaviour with option **columns**.

draw( columns = 2, dimensions = [600, 300], gr2d( color = red, explicit(sin(x), x, -5, 5)), gr3d( color = black, surface_hide = true, explicit(sin(x)*cos(y), x, -5, 5, y, -5, 5)) ) $

This last example shows how to plot a scene inside another one. Option **allocation** needs a list of two pairs of numbers: the coordinates of the lower-left corner of the scene, and its width and height; all values must be relative, from 0 to 1.

draw( dimensions = [500, 400], gr2d( background_color = cyan, color = red, explicit(x^2,x,-1,1)), gr2d( allocation = [[0.30, 0.50],[0.50, 0.40]], background_color = green, explicit(x^3,x,-1,1), grid = true) ) $

Animations

For animations, we have to use function **draw** and set **terminal = animated_gif**.

draw(terminal = animated_gif, delay = 40, dimensions = [300,300], makelist(gr2d(explicit(x^(k/10),x,0,1)), k, 1, 10) ) $

Statistical functions are defined in package **descriptive** and we have to load them before performing our calculations.

Statiscal graphics defined in package **descriptive** must be written with first letter in upper case: **Scatterplot**, **Histogram**, **Barsplot**, **Piechart**, **Boxplot**, and **Starplot**.

Our records here are stored in a list.

load(descriptive) $ sample : [4, 7, 6, 1, 5, 10, 3, 6, 6, 6, 9, 9, 5, 2, 2, 7, 7, 4, 6, 7, 8, 4, 10, 10, 4] $ mean(sample) ; median(sample) ; /* standard deviation */ std(sample) ; /* piechart */ Piechart(sample) $

In case of bivariate records, pairs are inroduced as a two columns matrix.

load(descriptive) $ sample: matrix( [125.1,140.7], [130.6,155.1], [135.1,160.3], [140.0,167.2], [145.4,169.8], [142.7,168.5])$ /* mean vector */ mean(sample) ; /* covariance vector */ cov(sample) ; /* correlations matrix */ cor(sample) ; /* dispersion cloud */ draw2d( point_type = circle, point_size = 3, color = navy, points(sample)) $

Mario RodrÃguez Riotorto