Curve Fit for HP-41cv/cx (PPC)

This program will determine a curve of best fit to a set of data points.  It handles the four standard curve types:\(\)

Curve types
name formula restriction
linear \(y=bx + a\)
exponential \(y=ae^{bx}\) for \(a>0\)
logarithmic \(y=b\ln{x}+1\)
power \(y=ax^{b}\) for \(a>0\)

The program will compute the coefficients a and b for the curve types as well as compute the coefficient of determination (r2) which is a measure of the goodness of fit.  Once a set of data has been fit to a given curve type, a prediction may be made for the y-value given a new x-value, or a prediction may be made for the x-value given a new y-value.

History

The code originates from the PPC ROM Module released by the international group People Programming Computers on December 1981 [1,2]. That in turn was based on Gary Tenzer’s curve fit program published in PPC Journal.

Instructions

The functions map to the top row on the USER keyboard.  These same functions are referenced in the examples by enclosing them in a blue box.

Keyboard overlay for CV

The key mappings are shown in the table below.

Key mappings
USER keyboard description usage
∑+ ∑+ Add a data point x ENTER  y ∑+
∑- ∑- Delete a data point x ENTER  y ∑+
SOLVEj 1/x Solve type j j SOLVEj
y^ √x Predict y value x y^
x^ LOG Predict x value y x^
SOLVE LN Solve best type
INIT ex Initialize

Operation

  1. GTO “CV”
  2. SIZE 027
  3. USER mode
  4. INIT, the display will show 1.
  5. key in the data pairs as x ENTER y and push ∑+.  The display will stop with the count of the number of the next data pair.  This feature makes it possible to enter only the y-values when the x-values are consecutive integers starting at 1.  In this case the display provides the x-values which need to be entered.  If an improper data pair has just been input with the ∑+ key, then immediately press R/S will delete the pair.  Otherwise an improper or undesired data pair can be deleted by re-entering both x and y and pressing ∑-.
  6. repeat step until all data pairs are entered.
  7. If any x-values are negative or zero, then only types 1 and 2 are feasible curves.  If any y-values are negative or zero, then only types 1 and 3 are feasible curves.  If in any data pair both and x and y are negative or zero, the type 1 is the only feasible curve.  The a coefficient must be positive for curve type 2 and 4.
  8. Select the desired curve type.  Once the curve type is selected, the HP-41 should not be interrupted.
    • To fit a particular curve type, key in the number 1-4 for that type and press SOLVEj.  The stack returns with:
      Stack usage
      stack register description
      Z R10 R2
      Y R09 a coefficient
      X R08 b coefficient

      The closer R2 is to one of the extremes (-1 or 1), the better the fit.  The sign indicates whether the data is positively or negatively skewed.
      Note that this step may be repeated any number of times of any of the four curve types.

    • If all data data input is positive, then pressing SOLVE will automatically choose the curve of best fit according to the curve type with largest absolute value of r2.  In this case the stack returns with:
      Stack usage
      stack register description
      Z R10 R2
      Z R09 a coefficient
      Y R08 b coefficient
      X R07 j = best curve type
  9. Predictions for new x or y values may be made only after the previous step has been completed.  Predictions for new values are based on curve type (encoded in flags F08 and F09).
    • To predict y given x, key in x and press y^
    • To predict x given y, key in y and press x^
  10. New data may be added or deleted at any time via the ∑+ and ∑- keys.  However the “solve” must be preformed after updating the data before any new predictions can be made.  The parameters a and b automatically destroyed after input of new data.

Examples

We reuse the examples from the PPC ROM User Manual. [hp41.com]

Straight line

Find the straight line which best fits the data set shown below, and predict y when \(x=20\) and predict \(x\) when \(y=25\).

$$(1.1, 5.2), (4.5, 12.6), (8, 20), (10, 23), (15.6, 34)$$

Minimum size is SIZE 25.

Straight line
key strokes display
GTO “CV”
SIZE 027
USER
INIT 1.0000
1.1 ENTER  5.2 ∑+ 2.0000
4.5 ENTER  12.6 ∑+ 3.0000
8 ENTER  20 ∑+ 4.0000
10 ENTER  23 ∑+ 5.0000
15.6 ENTER  34 6.0000
1 SOLVEj 1.9720

After solving for linear (type 1), the values for \(a\), \(b\) and \(R^2\) are on the stack.  We find:

\(y=1.972047542x + 3.499147270\)

To determine the points along that line:

Straight line
key strokes display
20 y^ 42.94009811
25 x^ 10.90280649

Linear or exponential curve

The data shown below fits either a linear or exponential curve.  Determine which is more appropriate.

$$(2,12), (-1,2), (3,17), (5,23)$$

Again, to solve this, we use the CV program.

Linear or exponential curve
key strokes display description
INIT 1.0000
2 ENTER  12 ∑+ 2.0000
1 CHS ENTER  2 ∑+ 3.0000
3 ENTER  17 ∑+ 4.0000
5 ENTER  23 ∑+ 5.0000
1 SOLVEj 3.5467 linear, b coefficient
R  5.5200 linear, a coefficient
R  0.9976 linear, R2
2 SOLVEj 3.5467 exponential, b coefficient
R 3.8262 exponential, a coefficient
R 0.9586 exponential, R2

Choosing the best \(R^2\), we find a linear fit is more appropriate.

\(y=3.5467x+5.5200\)

Predict \(y\) when \(x=-10\).  Since we just finished the exponential fit, the exponential parameters are still in the machine, and hence we must go back and solve for linear to predict \(y\) for \(x=-10\).

Predict \(y\) for \(x=-10\)
key strokes display description
1 SOLVEj 3.5467
10 CHS y^ -29.9467

Add the additional data points shown below and solve the same problem.

$$(-4, 0.713), (2.5, 10.93), (6, 47.53), (10, 254.95)$$

Note that the display should show 6 after entering the first new data pair.  We try a linear and exponential fit.

Add points and resolve
key strokes display description
4 CHS ENTER 0.713 ∑+ 6.0000
2.5 ENTER 10.93 ∑+ 7.0000
6 ENTER 47.53 ∑+ 8.0000
10 ENTER 254.95 ∑+ 9.0000
1 SOLVEj 15.3315 linear, b coefficient
R 0.9790 linear, a coefficient
R 0.7657 linear, R2
2 SOLVEj 0.4199 exponential, b coefficient
R 3.8256 exponential, a coefficient
R 0.9936 exponential, R2

Now choosing the best \(R^2\), we see that the new data reflects a change in the curve type.

\(y=3.8256e^{0.4199}\)

Since the exponential parameters are still in the machine, we can predict \(y\) for \(x=-10\) as

Predict y for x=-10
key strokes display
10 CHS y^ 0.0574

Best curve

Fit the best curve to the data points shown below.

$$(1,2), (2,2.828), (3, 3.464), (4, 4), (5, 4.472), (6, 4.899), (7, 5.292), (8, 5.657) and (9, 6)$$

In this example the \(x\)-coordinates start counting from \(1\) and are consecutive integers.  So we need only input the y-coordinates, but they must be in the proper order.  The count in the display will serve as the x-coordinates.

Find the best curve
key strokes display
INIT 1.0000
2 ∑+ 2.0000
2.828 ∑+ 3.0000
3.464 ∑+ 4.0000
4 ∑+ 5.0000
4.472 ∑+ 6.0000
4.899 ∑+ 7.0000
5.292 ∑+ 8.0000
5.657 ∑+ 9.0000
6 ∑+ 10.0000
SOLVE 4.0000

We could use the best fit function because all the values were positive.  After solving, the stack represents j, b, a and r2.  The value j=4 indicates a power curve.  Rounded to 2 decimals, this is

\(y = 2.00 x^{0.50}\)
using zumzum.com

Source code

  • Requires X-Functions module on the HP-41cv, and a minimum SIZE 25.
  • Available as source code, raw code for the V41 emulator and bar code for the HP Wand (HP82153A)
  • Size 321 Bytes, 3 magnetic tracks (321/112=2.866)

Listing

;  /---------------------------------------------------------------------\
;  |                       C u r v e   F i t                             |
;  |                                                                     |
;  |                         for the HP-41                               |
;  \---------------------------------------------------------------------/
;
;                               1.00
;                             PPC ROM
;
;         http://www.coertvonk.com/technology/hp41/curve-fit-4581

; BLOCK CLEAR ROUTINE, stores zeroes in a block of registers
;
; uses the complete form of the general block control word bbb.eeeii an
; can thus be used to clear blocks of consecutive registers or can be used
; to sip over registers within a block.

01      LBL "BC"
02      LBL 16
03      SIGN
04      CLX
05      LBL 20
06      STO IND L
07      ISG L
08      GTO 20
09      RTN

10      LBL 17
11      LBL "CV"
12      GTO IND 06      ; provides access to all numeric labels within CV

; INPUT A DATA POINT
;
; All summations are updated when this routine is called.  These summations
; include sums of x, x^2, y, y^2, xy, ln(x), ln(x)^2, ln(y)^2, ln(x).ln(y),
; x.ln(y), y.ln(x).

13      LBL A           ; [sigma+]
14      LBL 01
15      CF 10
16      LBL 06
17      STO 09
18      X<>Y
19      STO 08
20      sREG 13
21      FC? 10
22      s+
23      FS? 10
24      s-
25      RDN
26      RCL 08
27      ENTER^
28      X>0?
29      LN
30      ST* Z
31      RCL 09
32      X>0?
33      LN
34      ST* Z
35      X<>Y
36      sREG 19
37      FC? 10
38      s+
39      FS? 10
40      s-
41      R^
42      FS? 10
43      CHS
44      ST+ 12
45      R^
46      FS? 10
47      CHS
48      ST+ 11
49      X<> Z
50      SIGN
51      ST+ L
52      RCL 08
53      RCL 09
54      X<> L
55      TONE 9
56      RTN

57      RCL 08
58      RCL 09

; REMOVE A DATA POINT
;

59      LBL a
60      SF 10           ; signal "remove data point"
61      GTO 06          ; jump to "add data point"

; SOLVE SPECIFIED CURVE TYPE J
;
; b is stored in R08, a is stored in R09, r2 is stored in R10

62      LBL B           ; [SOLVEj]
63      LBL 02
64      CF 08
65      CF 09
66      STO 07
67      2
68      X<Y?
69      SF 09
70      /
71      FRC
72      X=0?
73      SF 08
74      8
75      ST+ 07
76      XEQ IND 07
77      RCL 17
78      RCL 13
79      RCL 15
80      STO 09
81      *
82      RCL 18
83      /
84      -
85      STO 10
86      RCL 14
87      RCL 13
88      X^2
89      RCL 18
90      /
91      -
92      STO Z
93      /
94      STO 08
95      RCL 13
96      *
97      ST- 09
98      X<>Y
99      RCL 16
100     RCL 15
101     X^2
102     RCL 18
103     ST/ 09
104     /
105     -
106     *
107     SQRT
108     ST/ 10
109     XEQ IND 07
110     8
111     ST- 07
112     RCL 10
113     RCL 09
114     FS? 08
115     E^X
116     STO 09
117     RCL 08
118     TONE 5
119     RTN

; A SERIES OF INTERTWINED SUBROUTINES, called in the curve fitting process.
;
; These routine simply perform a series of register exchanges which place
; the proper sums in the sigma registers for the calculation of the
; parameters a, b and r2 depending on the curve type selected.  Since the
; exchange is performed twice (line 076, line 109) all registers are
; returned to their original state.

120     LBL 10
121     RCL 11
122     X<> 17
123     STO 11
124     LBL 13
125     RCL 21
126     X<> 15
127     STO 21
128     RCL 22
129     X<> 16
130     STO 22
131     LBL 09
132     RTN

133     LBL 11
134     RCL 12
135     X<> 17
136     STO 12
137     LBL 14
138     RCL 19
139     X<> 13
140     STO 19
141     RCL 20
142     X<> 14
143     STO 20
144     RTN

145     LBL 12
146     RCL 23
147     X<> 17
148     STO 23
149     XEQ 14
150     GTO 13

; PREDICT THE Y VALUE

151     LBL C           ; [Y^]
152     LBL 03
153     FS? 09
154     LN
155     RCL 08
156     *
157     RCL 09
158     FS? 08
159     LN
160     +
161     FS? 08
162     E^X
163     RTN

; PREDICT THE X VALUE

164     LBL D           ; [X^]
165     LBL 04
166     FS? 08
167     LN
168     RCL 09
169     FS? 08
170     LN
171     -
172     RCL 08
173     /
174     FS? 09
175     E^X
176     RTN

; INITIALIZE, clear the data registers used to accumulate the sums

177     LBL e           ; [INIT]
178     LBL 00
179     11.024
180     XEQ 16
181      E
182     RTN

; SOLVE BEST CURVE TYPE

183     LBL E           ; [SOLVE]
184     LBL 05
185     .
186     STO 25
187     4
188     STO 07
189     LBL 07
190     RCL 07
191     XEQ B
192     RCL 25
193     RCL 10
194     ABS
195     X<=Y?
196     GTO 15
197     STO 25
198     RCL 07
199     STO 26
200     LBL 15
201     DSE 07
202     GTO 07
203     RCL 26
204     XEQ 02
205     RCL 26
206     TONE 5

207     END

References

[1] CV – Curve Fit
Gary Tenzer, Keith Jarret and John Kennedy, August 1981
PPC User Manual, page 110-111
[2] BC – Block Clear
John Kennedy, August 1981
PPC User Manual, page 50-51
Embedded software developer
Passionately curious and stubbornly persistent. Enjoys to inspire and consult with others to exchange the poetry of logical ideas.

Leave a Reply

Your email address will not be published. Required fields are marked *

 

This site uses Akismet to reduce spam. Learn how your comment data is processed.