Solving Complex Mathematical Tasks With a Simple Recursive Solution

April 27, 2004
by Lee Cullens

The idea behind this article is to convey a very simple solution technique for an otherwise more involved mathematical task. Such is not going to hold your interest though, if you don't appreciate the elegance and versatility of the mathematical approach. So for the "un-indoctrinated" we'll begin with the fundamentals.

Going beyond pretty artwork and simple animation generally involves a little math, which in turn, not infrequently, involves problems with more than one unknown. For other than an arbitrary solution, such problems depend on multiple sets of like unknowns. These types of math problems are encountered in business, industry and the sciences. More specifically, a sampling includes graphics (e.g. intersecting, shaping/fitting, transforming, projecting), production management (e.g. what products to produce when and how best to assign production tasks), gaming theory (e.g. Markov chains involving permutations and combinations), municipal government (e.g. traffic flow networks), and on and on.

One can't escape the necessity of understanding a problem area in formulating solutions, but solving the equations for the unknowns is simple. No "high brow" math to it, just high school algebra like what I had in the early fifties.

In this approach when more than one set of like unknowns is considered, if the number of sets is less than the number of unknowns, there are, possibly, infinitely many solutions. If the number of sets equals the number of unknowns, though, then there is possibly a unique solution. Refining our terminology a little and substituting the word "equations" for "sets", herein we will discuss problems with an equal number of equations and unknowns.

This is part of what is called linear algebra and we are going to use matrix notation, which simplifies explanation. So, let's start with a likely familiar simple problem for our first example.

We want to find the intersection of two lines in the same plane. The first line passes through points (0, -3) and (6, 3), which in the form of a line equation is x - y = 3 (i.e. both points satisfy the equation). The second line passes through points (5, -4) and (0, 6) so the second equation is 2x + y = 6. To use matrix notation we construct an n by n (number of equations by number of unknowns) matrix of the equation coefficients (unknowns multipliers) and a n by 1 matrix of the constants (i.e. the equal to values). We will use rows for the equations and columns for the like equation coefficients, so each column in the coefficients matrix represents an unknown.

In abstract form our two equations would be a11x + a12y = k1 and a21x + a22y = k2, and the 2 by 2 matrix (with row column subscripts) would look like:

```        |a11 a12|       |k1|
|a21 a22|       |k2|
```
Thus, our sample problem in matrix form is:
```        |1 -1|  |3|
|2  1|  |6|
```

Back to our abstract form, we are going to solve this problem with something called determinates, which (in the interest of simplicity) is a shorthand for the multiplication and subtraction elimination of longhand algebra. We will calculate a determinate for the equation coefficients matrix and one each for the equation coefficients matrix with a column replaced with the constants matrix.

```        D =     |a11 a12|       Dx = |k1 a12|   Dy = |a11 k1|
|a21 a22|            |k2 a22|        |a21 k2|
```

The actual computation is D = a11* a22 - a21* a12 and likewise for Dx and Dy, and then the unknowns are x = Dx/D and y = Dy/D.

So, let's complete our sample problem, which you can verify graphically if you wish.

```        D =     |1 -1|                          Dx = |3 -1|     Dy = |1 3|
|2  1|                               |6  1|          |2 6|
D = (1* 1) - (2* -1)  = 3
Dx = (3* 1) - (6* -1)  = 9
Dy = (1* 6) - (2* 3)  = 0
x = 9/3  = 3   and   y = 0/3  = 0
```

Thus the intersection of the two lines is at point (3, 0).

Earlier in this article I said "there is possibly a unique solution." Now that we've gotten through the first sample problem, an explanation is due. Had the two lines been separate parallel lines, there is no intersection (i.e. no solution), and had the two lines been the same line (e.g. x - y = 3 and 2x - 2y = 6) there are infinitely many solutions. If the equations are inconsistent (i.e. no solution) or dependent (many solutions), the equation coefficients matrix will have a determinate of zero (i.e. D = 0).

If we get into 3D space, or any other problem that involves more than two unknowns, the determinate calculation is a little more involved. As defined, it is for a 2 by 2 matrix, so we need to work our way down through the matrix till we get to a 2 by 2 matrix, then back out with the products. This is called a "minors" approach (for the mathematicians among us: (-1)i+jaijAij). The idea is to work our way across the first row (or alternately down the first column) and multiply each element of such by its minor. The minor of an element is the determinate that remains when the row and column in which the element appears are omitted. So for a 3 by 3 matrix the determinate is:

```|a11 a12 a13 |
|a21 a22 a23 | = a11|a22 a23 | - a12|a21 a23 | + a13|a22 a23 |
|a31 a32 a33 |      |a32 a33 |      |a31 a33 |      |a32 a33 |
```

To make sure we understand, let's do one of the old "apples and oranges" problems. A child checks its piggy bank and counts 33 coins in nickels, dimes and quarters, which total to \$3.75. Interestingly (to the child) there are twice as many dimes as quarters. How many nickels, dimes and quarters does the child have?

Our equations for this problem are the total value (.05n + .10d + .25q = 3.75), the total quantity (n + d + q = 33) and the proportionality (d - 2q = 0), so our problem in matrix form is:

```| .05  .10   .25|   | 3.75|
|1.0  1.0   1.0 |   |33.0 |
|0.0  1.0  -2.0 |   | 0.0 |
```

The determinates then are:

```        D = .05|1.0  1.0| - .10|1.0  1.0| + .25|1.0 1.0|
|1.0 -2.0|      |0.0 -2.0|      |0.0 1.0|
= .05((1.0*-2.0) - (1.0*1.0))
- .10((1.0*-2.0) - (0.0*1.0))
+ .25((1.0*1.0) - (0.0*1.0))
= 0.3
Dx = 3.75|1.0  1.0| - .10|33.0  1.0| + .25|33.0 1.0|
|1.0 -2.0|      |0.0  -2.0|      | 0.0 1.0|
= 3.75((1.0*-2.0) - (1.0*1.0))
- .10((33.0*-2.0) - (0.0*1.0))
+ .25((33.0*1.0) - (0.0*1.0))
= 3.6
Dy = .05|33.0  1.0| - 3.75|1.0  1.0| + .25|1.0 33.0|
| 0.0 -2.0|       |0.0 -2.0|      |0.0  0.0|
= .05((33.0*-2.0) - (0.0*1.0))
- 3.75((1.0*-2.0) - (0.0*1.0))
+ .25((1.0*0.0) - (0.0*33.0))
= 4.2
Dz = .05|1.0 33.0| - .10|1.0 33.0| + 3.75|1.0 1.0|
|1.0  0.0|      |0.0  0.0|       |0.0 1.0|
= .05((1.0*0.0) - (1.0*33.0))
- .10((1.0*0.0) - (0.0*33.0))
+ 3.75((1.0*1.0) - (0.0*1.0))
= 2.1
```

So

```
x = 3.6 / 0.3  = 12 nickels
y = 4.2 / 0.3  = 14 dimes
z = 2.1 / 0.3  = 7 quarters
```

As you can easily see, the calculations expand rather quickly to more than you want to be coding out longhand. For having persevered this far, the reward is a very simple recursive solution, so that you can concentrate on setting up the problem without worrying about the mundane math of solving it.

If you copy out the following script, you can test it with the two sample problems and this additional system of equations:

```         x + 3y + 5z -  w =  4
3x -  y - 2z + 4w =  3
2x + 6y + 3z -  w = -1
5x + 7y + 2z + 2w =  1
```

In the message window enter:

```        TList = [[1.0, 3.0, 5.0, -1.0, 4.0],
[3.0, -1.0, -2.0, 4.0, 3.0],
[2.0, 6.0, 3.0, -1.0, -1.0],
[5.0, 7.0, 2.0, 2.0, 1.0]]
TV = SolveLinearSystem(TList)
put TV
-- [-4.0000, 1.0000, 2.0000, 5.0000]
```

----------------------------------------------------------------------
----------------------------------------------------------------------
-- Common use script to --
-- Solve an n by n Linear Equations system --
----------------------------------------------------------------------
--
-- Pass a list structure of the equations in the form
--
-- EqList = [[a11, a12, a13, ... a1n, k1],
-- [a21, a22, a23, ... a2n, k2],
-- ...
-- [an1, an2, an3, ... ann, kn]]
--
-- where n is the number of unknowns and the number of equations
--
-- A list of the solved unknowns is returned in the form
--
-- [x1, x2, x3, ... xn]
--
-- If the equation coefficients determinate is found to be
-- zero (no solution or many solutions) then the return value
-- will be VOID.
--
-- A simple determinate minors approach is employed
-- in a recursive process.
----------------------------------------------------------------------
----------------------------------------------------------------------
on SolveLinearSystem EqList

-- As an initial step in this handler you might want to check the
-- passed list structure to insure that it is a two dimensional
-- list (ilk(EqList, #list) and (ilk(EqList[n], #list)),
-- that it contains at least two rows (EqList.count >= 2),
-- that each of the rows has an equal and appropriate number of
-- columns ((EqList.count + 1) = EqList[n].count)), and that all
-- of the elements are real (float) numbers.

-- This script copies out the determinates from the passed list,
-- so the passed list is unaltered

-- Initialize results vector and working matrix
RV = []
WM = []

-- Define each top level determinate in turn, building results
-- after the first determinate
repeat with n = 1 to (EqList.count + 1)

if n = 1 then

-- on first loop build determinate of coefficients
repeat with nr = 1 to EqList.count
WM[nr] = []
repeat with nc = 1 to EqList.count
WM[nr][nc] = EqList[nr][nc]
end repeat
end repeat

-- calculate determinate of coefficients
DOC = DetCalc(WM)
if DOC = 0 then
RV = VOID
exit repeat
end if

else

-- on subsequent loops swap columns for unknown determinates
repeat with nr = 1 to EqList.count
WM[nr][n - 1] = EqList[nr][EqList[nr].count]
if n > 2 then
-- restore previous swap
WM[nr][n - 2] = EqList[nr][n - 2]
end if
end repeat

-- calculate determinate of variable
DOV = DetCalc(WM)

-- calculate variable
RV[n - 1] = DOV / DOC

end if

end repeat

return RV

end SolveLinearSystem

----------------------------------------------------------------------
----------------------------------------------------------------------
-- Recursive determinate calculation handler/function
----------------------------------------------------------------------
on DetCalc matrix

if matrix.count = 2 then

Dminor = (matrix * matrix) - \
(matrix * matrix)
return Dminor

else

Dminor = 0.0
repeat with i = 1 to matrix.count

minor = matrix.duplicate()
minor.deleteAt(i)
repeat with j = 1 to minor.count
minor[j].deleteAt(1)
end repeat
Dminor = Dminor + (power(-1, (i + 1)) * \
matrix[i] * DetCalc(minor))

end repeat

return Dminor

end if

end

----------------------------------------------------------------------
----------------------------------------------------------------------

Lee Cullens retired early from a software engineering career to pursue artistic interests. Then, creating distributed presentations of artwork evolved into a marriage of skills and he now provides authoring assistance where Macromedia Director is the prominent tool.