# Solving Equations

Maxima can provide analytical solutions to a wide range of equations
and to many systems of equations.  Where analyical solutions cannot be
produced, Maxima can usually produce numerical approximations to
solutions for equations and systems of equations.

## Single Equations

Use the solve command to solve y=x^3+a*x^2+b*x+c. Assign the list 
of solutions the name soln. Which of the three solutions is likely 
to apply to economic models? Why? 

/*Input. Execute (ctrl-enter) for results */
y:x^3+a*x^2+b*x+c;
soln:solve(y,x); /*Only one solution is real*/ end;

Use the subst command to evaluate the solution when
a=-2, b=10, and c = 100. To ignore the two solutons
that involve imaginary numbers use a command like this
one: subst([list], soln[3]);. Convert the output to a 
floating-point representation. Set fpprintprec to 5.
Assign the name soln0 to this result.

/*Input. Execute (ctrl-enter) for results */
fpprintprec:5$
soln0:subst([a=-2, b=10, c = 100],soln[3])$ float(%);

To confirm that this value of x yields a value of y that
(very nearly) equals zero, enter 
subst([a=-2, b=10, c = 100, soln0],y);. The result is a 
long expression whose value is not obvious. Convert this
value to a floating-point representation.

Is the result very nearly zero? Explain.

/*Input. Execute (ctrl-enter) for results */
subst([a=-2, b=10, c = 100, soln0],y);      float(%);

Execute the cell below to see a graph of this expression.
Graphs are treated in a later exercise set, so just examine
the result and not the commands used to produce it (unless
you're interested in doing so).

wxdraw2d(xaxis=true,
explicit(subst([a=-2,b=10,c=100],y),x,-5,5))$

The solution above is analytical. Sometimes no analytical
solution can be found. In such cases, numerical solutions
can be useful. To solve this problem numerically, create
an expression named ynumer that is y with after the 
substitution of a=-2, b=10, and c=100. Use the command
find_root(ynumer,x,-5,2) to approximate the solution for 
the expression.

Confirm that the result agrees with the result of the 
analytical solution, to the available number of decimal 
points.

Maxima offers other numerical methods, but we limit this exercise
to find_root.

By the way, if one's interest is in finding accurate numerical roots,
then it is generally better to use a numerical root finder like 
find_root rather than substituting into the analytical solution.
See http://www.math.utexas.edu/pipermail/maxima/2013/031938.html
and the references therein.

/*Input. Execute (ctrl-enter) for results */
ynumer: subst([a=-2,b=10,c=100],y);
find_root(ynumer, x, -5, 2);
/*
These commands could have been nested, and the name 
need not be assigned
*/
find_root(subst([a=-2,b=10,c=100],y),x,-5,2);
/*
Or, we could have used the % operator to get a solution
without naming the expression
*/
subst([a=-2,b=10,c=100],y)$    find_root(%, x, -5, 2);

Execute this command: solve(A*x^a = B + b*x, x);. Confirm that
no unique value of x.

/*Input. Execute (ctrl-enter) for results */
solve(A*x^a = B + b*x, x);

Substitute A=1200, a=-2, B=5,and b=0.5 into the expression
that you tried to solve. Name this expression expr1.

/*Input. Execute (ctrl-enter) for results */
expr1: subst([A=1200,a=-2,B=5,b=0.5],A*x^a = B + b*x);

The graph below can guide our use of find_root to generate
a numerical solution to this equation. Again, just use the
graph; don't worry about the command. The graph indicates
a solution somewhere around x = 10.

wxdraw2d(yrange=[-10,10],xaxis=true,
explicit(1000*x^-2 - (5 + 0.5*x),x,.001,50));

Use the find_root command to find the value of x that satisfies
this expression. Use 5 and 15 as the endpoints for the command.
Any positives would work (x = 0 is not allowed; why?). Generally,
it is a good idea to use information from graphs when executing
find_root, for reasons that we see below.

/*Input. Execute (ctrl-enter) for results */
find_root(expr1,x,5, 15);

Suppose that the expression that we wish to examine is 
the cubic 1000 + 50*x - 3*x^2 + 0.03*x^3. The graph below
shows that it has three roots, and it provides some
guidance regarding the ranges over which one should
search to find numerical solutions.

wxdraw2d(yrange=[-500,1250],
explicit(1000 + 50*x - 3*x^2 + 0.03*x^3, x, -50,100));

Create an expression expr2:1000 + 50*x - 3*x^2 + 0.03*x^3.
Use find_root to find the three solutions. Name these 
solutions x1, x2, and x3.

/*Input. Execute (ctrl-enter) for results */
kill(all)$
expr1:1000 + 50*x - 3*x^2 + 0.03*x^3;
/*
These endpoint choices are somewhat arbitrary.
Use your reading of the graph to select values.
You need not place the solutions in a list.
*/
[x1: find_root(expr1,x,-20,15), 
x2: find_root(expr1,x,20,50),
x3: find_root(expr1,x,60,100)];

Substitute each of these values into expr1 and confirm that 
the result is (very nearly) zero.

/*Input. Execute (ctrl-enter) for results */
subst(x1,x,expr1); subst(x2,x,expr1); subst(x3,x,expr1);

## Systems of Equations

Create a linear equation expr1: y-x=3 and an equation for 
a circle expr2:x^2+y^2-5=0. (The =0 can be omitted; Maxima
assumes this is so unless specified otherwise.)

/*Input. Execute (ctrl-enter) for results */
kill(all)$
[expr1,expr2]:[y-x=3,x^2+y^2-5=0]$
soln: solve([expr1,expr2],[x,y]);

Substitute the list soln1 into the two expressions and confirm that
a solution has been found.

/*Input. Execute (ctrl-enter) for results */
subst(soln[1],[expr1,expr2]);

Do the same for soln[2].

/*Input. Execute (ctrl-enter) for results */
subst(soln[2],[expr1,expr2]);

The graph below shows the two solutions. As before, don't
worry about the graph's creation at this point.

Note that the circle does not look circular. Keep in mind that
x and y need not be in the same units. If they were, and if
the x and y scales were the same, then the graph would 
show a circle.

/*Input. Execute (ctrl-enter) for results */
wxdraw2d(
implicit(y-x-3,x,-3,3,y,-3,3),
implicit(x^2+y^2-5, x,-3,3,y,-3,3))$

Leave expr2 as is, but replace expr1 a follows: expr1:y=3*x^3-5.
Use solve to determine the solutions. Which of the solutions might
be pertinent in the context of economic analysis?

/*Input. Execute (ctrl-enter) for results */
kill(all)$ fpprintprec:5$
[expr1,expr2]:[y=3*x^3-5,x^2+y^2-5]$
solve([expr1,expr2],[x,y]);

Here's the graph.

wxdraw2d(
implicit(y-3*x^3+5,x,-3,3,y,-10,3),
implicit(x^2+y^2-5, x,-3,3,y,-10,3))$

Load the program mnewton. Then apply mnewton to the system
[expr1,expr2]. The syntax is mnewton([expr1,expr2],[x,y],[xguess,yguess]).
Apply mnewton twice. First enter [1,-2] to determine the bottom intersection
in the graph above. Then enter [1.5,2] to determine the top intersection.
Confirm that the values are the same as with solve, at least to the 
specified number of decimal places.

/*Input. Execute (ctrl-enter) for results */
load(mnewton)$ 
mnewton([expr1,expr2],[x,y],[1,-2]);
mnewton([expr1,expr2],[x,y],[1.5,2]);

As before, leave expr2 as is, but replace expr1 with y-log(x)=-1.
Now try applying solve( ). Confirm that solve( ) cannot solve this
pair of expressioins.

/*Input. Execute (ctrl-enter) for results */
kill(all)$
[expr1,expr2]:[y-log(x)=-1,x^2+y^2-5];
solve([expr1,expr2],[x,y]);

Here's the graph of this pair of expressions over the relevant range.
Intersections appear to occur at approximately x=-.025 and x=2.25.

wxdraw2d(
implicit(y-log(x)+1,x,.001,3,y,-3,3),
implicit(x^2+y^2-5, x,.001,3,y,-3,3))$

Use the same set of mnewton commands as above to find the 
solutions to this set of equations. Use the graph above to guide
your selection of initial guesses. You need not enter load(mnewton);
again.

/*Input. Execute (ctrl-enter) for results */
load(mnewton)$
soln1: mnewton([expr1,expr2],[x,y],[0.025,-2.25])[1];
soln2: mnewton([expr1,expr2],[x,y],[2.25,-0.5])[1];

Finally, make the appropriate substiutions to confirm that the two 
expressions are satisfied.

/*Input. Execute (ctrl-enter) for results */
subst(soln1,[expr1,expr2]);
subst(soln2,[expr1,expr2]);