Math Development Program

Before we get into the “main course” for this issue, I would like to present a brief discussion of what “Simultaneous Linear Equations” are, for the benefit of those who may not be familiar with these concepts, or need to brush up a hazy memory from high-school or college days. It’s not intended as a course in S.L.E., as we just don’t have the space for such an undertaking. I do hope to give enough information to spur your interest to do as Mr. Bent suggests and bone up on the subject using a good text or even an outline book if all you need is a refresher.

Let’s start out with a simple mathematical relation expressed in words. Suppose someone told you he bought 20 apples and 5 mangos for $5.75. You could not deduce much about the price of apples and mangos as there is an infinite number of possibilities (assuming that fractional cent values are allowed). But suppose he also told you that had he bought 6 apples and 10 mangos it would have cost $6.40. Ah, now you can sit down with pencil and paper and figure out how much each costs. Before you read on, try to work out the solution to this little problem.

As it turns out, the word problem above can very easily be converted to algebraic equations (how convenient!).

20*X+5*Y = 5.75 and 6*X+10*Y = 6.40

where X is the price of apples and Y is the price of mangos.

This pair of equations is termed a “system” of equations. In order for both equations to hold true at the same time (hence “simultaneous”), the values for X and Y must “satisfy” both; i.e., when we plug them into either equation we should have both sides of the equal sign the same. Put the values you found into the equations as a check on your work; both equations should balance. (If you’ve forgotten how to solve this type of system, here’s a hint: multiply each term in the first equation by 2, then subtract the second equation term by term from the first. This results in an equation from which Y has been eliminated, and you should have no trouble finding X=.15 and then put that value back into either of the original equations and get Y=.55.)

OK so far, but why “linear?” Linear means having to do with straight lines, and linear equations are ones that produce straight lines when graphed. As it turns out, this is only true of equations in which the variables are not raised to powers or have other functional operations (e.g. SIN, ATN, LN, etc.) performed on them. In this installment we will only deal with this type of system. If we graph our fruity problem above, we come up with a plot like Fig. 1. Line 1 represents all the combinations of x and y which satisfy the first equation, and line two is all the points which Satisfy the second. The point where they intersect is the solution – the only point which simultaneously satisfies both equations

At this point, let’s assume that the second equation were 21*X + 5*Y = 5.90 instead. The system still gives the same solution, as you can verify. A graph of this modified system is shown in Fig, 2. Again there’s one discrete solution, but in this case the lines are very nearly parallel. Such a system is called “ill-conditioned,” and is less likely to give accurate results when using a computer. Just as it’s harder to discern the intersection point on the graph, so also the computer has difficulty resolving the exact solution. The program includes a routine that checks for ill-conditioning and yields a number between 0 (parallel lines or the same line and possessing no discrete solution) and 1 (perpendicular lines, and most accurate) to rate the system.

So big deal. We don’t need a computer to solve simple systems like this. The rub comes when we try to solve systems with more than two variables – try solving a six-variable system sometime! Note that for a solution we must have at least as many equations as we have unknowns, so you can see that things rapidly get out of hand. What we need is a systematic approach the computer can just crunch away at. Such systematic approaches were developed long before computers were, so don’t get the idea that all this is a new development. With a computer, though, even the huge task of solving 10th order S.E. becomes routine.

The first step is to convert our system of equations into an array called a matrix. Matrices are mathematical entities that can be operated on according to certain rules. Study Mr. Bent’s multiplication example to see why we can express our system as follows:

We thus have two arrays of numbers, one on either side of the equal sign, that specify the system. The matrix on the left side contains the multipliers (coefficients) of the unknowns, and the one on the right contains the constants each equation adds up to. For convenience, we can combine the two matrices into one 2×3 matrix, as this makes it easier to input the numbers. Just remember that the last column is the constants, the others are the coefficients.

Now that we have our matrix, we need some way to boil it down into answers. One way is by breaking it down into “determinants,” which is another form of array. Unlike a matrix, though, a determinant boils down to a single number, when evaluated according to certain simple rules, The bigger it is, the more calculations it takes, but you’ll always end up with a number that can be used to obtain a solution. To keep this distinction in view, we’ll indicate matrices with brackets, and determinants with straight lines. Arrays that are essentially a “butt-splice” of two matrices — Le, the way you input the numbers to the computer, are shown with curved brackets as they are neither matrices nor determinants but a shortcut for convenience. These will usually be N*N+1 in dimension.

I hope this has whetted your appetite for the article. A few short notes on the program: Enter RUN to start, push F or S for fast or slow mode. The main menu appears, options 1-3 allow matrix operations, 4 evaluates a determinant, and 5 goes on to the main S.L.E. portion. Selecting this gives an SE menu that is a veritable smorgasbord of possibilities. I particularly like combination plate #7, I think you will too. After completing analysis by Gauss Elimination and/or Gauss Seidel any size, you can repeat the problem using Cramers Rule by selecting 5, then 2, and pressing “N” on “NEW DATA” prompt. This is the only program I’ve seen allowing all three methods, even though CR any size “cheats” by using GE to do the dirty work. For this reason it usually gives the same results as GE (actually results ar more accurate, but the time scale is considerably larger), but I added it to demonstrate row-swapping of the last row with the others to get the determinants for each answer.

By Thomas Bent

First, a hearty hello to all you ZX/TSers out there! I would like to say a few words about myself before going too far. My computer background consists of a couple years of Basic and Fortran with Pascal as a new addition. For those of you familiar with DEC equipment, some Focal language also. My job doesn’t consist of a lot of programming, but the need to use a computer to solve problems for myself (and my boss) sent me in search of this programming tool. Since I do quite a bit of materials testing, I do a fair amount of set-up and program design on an assortment of computers. I am also a professional pilot holding a commercial, multi-engine rating.

Many of those who initially purchase the ZX/TS computer buy it with little knowledge of how computers work. Most of the rest, who know how they work, can’t afford anything more expensive. Then there are those who hold out for a more expensive computer simply because they can’t believe that a $30.00 computer can do what they want.

I have been pleasantly surprised with my Sinclair in the past two and a half years, and I look forward to a computerized world in which I know I will be able to hold up my end of the load. What I hope to do in a series of upcoming articles is to bring newcomers “up to speed” in some basic techniques of programming math routines on the ZX/TS, so that you won’t have to shy away from a project simply for a lack of knowledge on how to program these routines. Those who are familiar with this material can incorporate the routines into their own programs conveniently.

What I plan to cover will include but not be limited to more integration, numerical differentiation (polynomial and real-time,) curvefitting of various types, finding roots of equations, and simultaneous equations. Some of these techniques are hard to understand at first, but very easy to program. Please try to follow the development of the first subject, Simultaneous Linear Equations, because these routines will be used as subroutines in many programs that lie ahead.

Most high-school and college students know how to solve linear simultaneous equations (sim. eq.) and I’m sure they will agree that they have better things to do with their time. This first part will provide some useful programs for solving these types of problems, and will greatly speed up that engineering and/or math homework.

Before we start, let’s get one thing straight; all computers are notorious liars, Apple, Atari and VIC computers (6502 chip) do not handle floating point mathematics well. One reason that the ZX/TS is much slower when grinding numbers is because it has much better resolution (5-byte floating point as opposed to 4-byte for the others, or two extra decimal places accuracy). However if you want speed then you can always get a compiler (e.g. Bob Berch, 19 Jacques St., Rochester, NY 14620) and speed up those “integer” calculations anyway. REMEMBER!! There are many large holes in the computer number line and many numbers are stored inaccurately (e.g. Pi, .2, 1.3). This is quantizing error. Addition, subtraction, etc. are also done inaccurately (try taking the cube root [or square root] of eight, cube it [square it], and subtract 8 on your calculator). This is roundoff error, Heaven forbid you should add two inaccurate numbers many times and say, “the answer is….”!! One last point, if your data is accurate to two decimal places then your answer after several calculations certainly is no better. Throw away all those insignificant figures, they do you more harm than good.

Alright now, on with the show.

The good news is, this series will provide a complete computer math program and hopefully some understanding behind how the routines work. The bad news is that the material tends to be rather dry (and may require some WORK). The laughs will come later when you see how easy the homework becomes. The future chapters of this report will be integrated with this first section. I think if you should be very easy.

The Sim. Eq. techniques (S.E.) we will delve into this time around are Gauss-Seidel, Cramer’s Rule, and Gauss elimination, This is a tough subject to explain and to grasp, fortunately it is easier to program. It may require some extra homework, but if there are any questions don’t hesitate to write (this still holds-ed.) If S.E. are understood then many upcoming topics will be a snap.

Manipulating matrices is the main scientific use of computers, You can add and subtract matrices of the same size, or multiply them by a constant. You can also multiply square or rectangular matrices where the number of columns. of the first matrix equals the number of rows of the second. Multiplying matrices is only tricky at first, but as you can see in the program, this is done with three loops and one calculation line, Please take the time to work through the given multiplication example so that the following S.E. routines will be more familiar and easier to debug.

You multiply rectangular matrices as follows:

A 4x2 * 2x4 = 4x4
A 2x4 * 4x2 = 2x2

I have worked out an example of the first case to show all the calculations involved. In the program, you input the size of matrix A and the computer only accepts input for the number of columns of matrix B. The product, matrix C has as many rows as A and as many columns as B.

It is much easier for the ZX/TS to do this kind of grunt work, and I think you’ll agree.

Determinants (det) are SQUARE arrays of any size, e.g. 2×2, 10×10, with some unique properties, A det is a single number in the form of an array of numbers, An Nth order det generates a predetermined number of calculations to solve, namely N factorial (N!) times (N-1). A 3rd order (3×3) det generates 12 calculations (3*2*1)*(3-1). A 10th order det generates 10!9 or in excess of 32,600,000 calculations! There are identities (det = 1) and inverse dets, among others. Please refer to a good math book or even flip through a Shaum’s outline book at your nearest book store for the properties of dets. Two important properties we must program for are:

1. If you interchange two rows of a det, the absolute value does not change, but the sign is reversed,
2. Adding one row of a det to another does not change the value of the determinant.

Evaluating a small det is easy and is done as follows:

In array form it looks like this:

Larger ones just get more cumbersome and tougher to keep straight.

The large det just gets broken down into smaller units and manipulated. By the way, if one row is a multiple of another row, the determinant = 0.

Taking what has just been discussed, we will transform our det (since we know some key properties) to reduce the astronomical number of calculations to something more manageable. In order to reduce round-off error, we should insure that the diagonal terms [A(1,1), A(2,2), etc.] contain the largest absolute value in their respective columns, for the terms below them. By this I mean that after evaluating col 1 row 1, leave that row alone [the col values below A(1,1) will all be zero.] Of course the last row gets whatever values are left.

What happens in the program is, the det is INPUT and col 1 is checked for the largest term. A zero term is checked for as this would mean an infinite number of solutions. The row with the largest value in col 1 is swapped with row 1. The new row 1 is multiplied by a different constant for each remaining row, so that every value in col 1 below the first term is zero when row 1 is subtracted from the row in question. Row 1 is now complete, as is col 1. Our “Pivot” now moves diagonally A(2,2), and the same trick is done again on the smaller array and so on until the last row (where there is only one term left). You must keep track of the number of times you swap rows (sign change) and simply multiply all of the diagonal terms together with the value of the sign (+ or -1) and voila! A 10×10 array takes less than 1000 calculations. A 6×6 array would not be uncommon for an engineer or engineering student. Most of us usually only see a 3×3 (because of the complexity and grunt work).

Simultaneous Equations (SE) will typically be 3×4 [or generally, Nx(N+1)] arrays and are handled much the same as determinants. Now we no longer have one number to deal with and must “back-calculate” the values that will solve the system of equations in question. The matrix is reduced the same way as the det is, except now we rearrange (put everything on the other side of the equality) our array to produce the values for X and Y [X(1) and X(2).]

The general form is:

AX + BY = E 
CX + DY = F

A, B, C, D, E and F are constants that are entered. X and Y we must solve for.

After reducing the matrix by multiplying row 1 by C/A and subtracting from the second row, the CX term cancels out and the second row becomes DY – CB/AY =F (the first row doesn’t change value.) The computer can quickly solve for Y and put that value in equation 1 and back out the value of X. There are no sign changes in SE to keep track of.

The computer will look for equations like this to solve for:

A(1,1) + A(1,2) + A(1,3) = A(1,4)
A(2,1) + A(2,2) + A(2,3) = A(2,4)
A(3,1) + A(3,2) + A(3,3) = A(3,4)

In the programs we dimension A(M,N) where M is the number of rows and N is the number of columns. This will give generality for any size array you have to solve for. X is dimensioned by M also, as linear SE are only valid for one unknown per equation. The computer will respond with one X for each column, that is, the answer for column 1 is X(1), etc. and that X is supposed to satisfy all M equations. (Hope and pray!)

Cramer’s Rule is very easy to program for small arrays (2×3), but gets tougher quickly thereafter. For 2 equations & 2 unknowns, in fact, it consists of only 3 program lines. Cramer’s Rule (CR) breaks the matrix into determinants (only one number) and solves for each unknown independently. This method is very close to how you solve linear SE by hand. Look at the following equations:

AX + BY = C
DX + EY = F

Multiply the top equation by E, the bottom equation by B and subtract the second equation from the first.

+(EAX + EBY = EC)
-(BDX + BEY = BF)
-----------------
EAX - BDX = EC - BF

Then consolidate X and solve.

X * (EA-BD) = EC-BF or
X = (EC-BF) / (EA - BD)

If this does not look familiar to you at this point, have no fear. We shall repeat the same thing again below the way the computer sees it.

AX + BY = C
DX + EY = F
X = (C*E - B*F) / (A*E - B*D) and 
Y = (A*F - C*D) / (A*E - B*D) 
or,
X = (DET X) / (DET) 
Y = (DET Y) / (DET)

where

If you look back at the section on determinants, you should see the connection. The denominator is the same for both X and Y, and the numerator swaps terms. Watch that the-denominator is not zero as the solution is undefined. How this works is easy to see in a 3×4, but is tougher to program directly, as in the case of the 2×3.

Al + Bl + CI = Y1 
A2 + B2 + C2 = Y2 
A3 + B3 + C3 = Y3

The three answers are as follows:

I prefer Gauss Elimination when solving for any system of higher order than 2×3, as CR requires a lot more calculations (takes longer and accuracy suffers) and also requires more memory space.

This method is designed for a computer although it is an old technique. What it does is start with a guess at what the value is, plugs and chugs and goes back to the beginning with a new guess (the value it just solved for). Sounds great! Only one drawback; it doesn’t always work. However, we can stack the odds in our favor. Gauss-Seidel (GS) is guaranteed to work if the value of the diagonal terms are greater than the sum of the terms in the row plus the the sum of the terms in the column. Another way to help make this work will be to begin with a different starting value if necessary. GS is very memory efficient and is not susceptible to round-off error like the other techniques. Each iteration is just like starting over, if the present answer is not good enough then go back and try again. You can also watch the computer approach the answer rather than sitting and waiting.

Should GS diverge or jump around, simply BREAK in and assign X values directly (LET X(1) = -10 or some other value), then CONT to see if that helps. If that doesn’t work, then BREAK and RETURN to try another way. Initially all X values will be assigned zero. The equations are swapped using the GE subroutine (but not reduced), unless you input the equations in the proper order, in which case the subroutine can be bypassed (saving memory), The computer will manipulate the equations to simulate the following:

AX(1) + BX(2) = C 
DX(1) + EX(2) = F
Then,
X(1) = (C - BX(2)) / A 
X(2) = (F - DX(1)) / E 
and so on.

Actually, the X term being evaluated is on both sides of the equation, but is assigned to zero every time through the loop. There must be an outside loop controlling the maximum number of iterations in case of divergence. Normally, 30 iterations will be more than enough for 6-place accuracy, but 100 may be necessary depending on how ill-conditioned the system of equations is (see the graphs).

To test for ill-conditioning, you “normalize” the determinant of the system of equations. This means that you square each term in the row, add them up and divide each term in the row by the sum of the squares (more on this next time), GOSUB GE to evaluate the determinant and you have a number less than one. If it is less than .2 then you may have trouble trying to solve your SE by any other method than by hand (but not necessarily). If you are concerned about the accuracy of your answers then try changing your matrix slightly and see if your answers vary largely, If so then get out your pencil and paper.

10 REM SE VERSN 3.2 4/1/84
34 GOTO 1000
38 REM INPUT DATA POINTS
42 PRINT "HOW MANY POINTS?"
46 INPUT M
50 LET DP=M
54 DIM P(M)
58 DIM U(M)
62 DIM V(M)
66 DIM Y(M)
70 PRINT " INPUT X THEN Y VAL"
74 FAST
78 FOR I=A TO M
82 PRINT TAB TR;I;" ";
86 INPUT P(I)
90 INPUT Y(I)
94 PRINT P(I),Y(I)
98 IF I>18 THEN GOSUB 910
102 NEXT I
106 GOTO 666
110 PRINT ,," INPUT START / STOP FOR X THEN Y"
114 DIM E(D)
118 DIM E$(D,D)
122 FOR I=A TO D
126 IF I=A OR I=TR THEN PRINT ("X  " AND I=A)+("Y  " AND I=TR)
129 INPUT E$(I)
134 PRINT E$(I),
138 NEXT I
142 RETURN