Multiv Polyn

# Multivariate Polynomials  for the HP-41

Overview

1°)  Two Variables

1-a)  Evaluation of Polynomial Functions
1-b)  Degree of a Polynomial
1-c)  Product of 2 Polynomials

2°)  N Variables

2-a)  Evaluation of Polynomial Functions
2-b)  Product of 2 Polynomials
2-c)  Deleting the Zeros
2-d)  Simplifications
2-e)  Partial Derivatives

1°)  Two Variables

-In this paragraph, the coefficients  a , b , c , d , e , f , g , h , i , j , ....  of a polynomial p(x,y) are to be stored into consecutive registers as follows:

p(x,y) = a + b x + c y + c x2 + d x.y + e y2 + f x3 + g x2y + h x.y2 + i y3 + j x4 + k x3y + l x2y2 + m x.y3 + n y4 +......

-Thus, we start with the constant term, then the monomials of degree 1 , the monomials of degree 2 , ... etc ...
-The polynomial is identified by the control number of its block of registers:  bbb.eee  that begins with Rbb and ends with Ree.
-All the coefficients must be stored ( including the zeros )

-So, we use:

1 register for a polynomial of degree  0
3 registers for a polynomial of degree  1
6 registers for a polynomial of degree  2
10 registers for a polynomial of degree  3   and, more generally,  (p+1)(p+2)/2  registers for a polynomial of degree p

-Inversely,  the degree of the polynomial = (1/2)(-3+(8n+1)1/2)  where n is the number of coefficients in the block of registers ( cf §1-b) ).

a) Evaluation of Polynomial Functions

Data Registers:              R00 = bbb.eee .....  1+eee.eee                            ( Registers Rbb thru Ree are to be initialized before executing "PXY" )

•  Rbb  ......  •  Ree = coefficients of the polynomial
Flags: /
Subroutines: /

-Synthetic registers M N O P may be replaced by R01  R02  R03  R04
-In this case, delete line 42 and replace line 02 by   0   STO 04   RDN

 01  LBL "PXY" 02  CLA 03  STO M 04  RDN 05  STO N 06  X<>Y 07  STO 00 08  RCL IND X 09  STO O 10  ISG 00 11  FS? 30 12  GTO 03 13  LBL 01 14  RCL P 15  INT 16   E3 17  / 18  STO P 19  SIGN 20  RCL IND 00 21  LBL 02 22  ISG 00 23  RCL M 24  * 25  X<>Y 26  RCL N 27  * 28  STO Z 29  RCL IND 00 30  * 31  + 32  ISG P 33  GTO 02 34  ST+ O 35  ISG 00 36  GTO 01 37  LBL 03 38  RCL M         39  SIGN 40  X<> N 41  RCL O 42  CLA 43  END

( 75 bytes )

 STACK INPUTS OUTPUTS Z bbb.eee / Y y y X x p(x,y) L / x

with  bbb > 000

Example:      p(x,y) = 1 + 2 x + 3 y + 3 x2 + 2 x.y - 4 y2 + 2 x3 + 4 x2y - x.y2 + 3 y3      >>>>   Evaluate  p(3,4)

-Suppose we choose the first register = R11

-Store    1       2      3       3      2     -4      2       4     -1      3
into   R11  R12  R13  R14  R15  R16  R17  R18  R19  R20   respectively    ( control number = 11.020 )

11.020   ENTER^
4        ENTER^
3        XEQ "PXY"   >>>>   p(3;4)  =  348    ( in 4 seconds )

b)  Degree of a Polynomial

-This short routine calculates the degree of a polynomial from its control number and it produces DATA ERROR ( line 20 ) if this control number is invalid.
-It is used to compute the product of 2 polynomials in paragraph 1-c)

Data Registers: /
Flags: /
Subroutines: /

 01  LBL "DEGXY"  02  INT  03  LASTX  04  FRC  05   E3  06  *  07  1  08  +  09  X<>Y  10  -  11  8  12  *  13  1  14  +  15  SQRT  16  3  17  -  18  2  19  /  20  FACT  21  X<> L  22  END

( 34 bytes )

Examples:

11.020   XEQ "DEGXY"  yields  3
1.028           R/S              yields   6
1.029           R/S              yields   DATA ERROR  1.029 is not a valid control-number

c) Product of 2 Polynomials

-Storing the coefficients as explained above allows to write a relatively fast program to multiplicate 2 polynomials of 2 variables:

Data Registers:              R00 to R07: temp               ( Registers Rbb thru Ree & Rbb' thru Ree' are to be initialized before executing "PROXY" )

•  Rbb  ......  •  Ree  = coefficients of the 1st polynomial     bb > 07
•  Rbb'  .....  •  Ree' = coefficients of the 2nd polynomial    bb' > 07
Flags: /
Subroutine:  "DEGXY"

-If you don't have an HP-41CX, replace lines 27-28-29  by

SIGN
CLX
STO 00
LBL 00
STO IND L
ISG L
GTO 00

 01  LBL "PROXY" 02  STO 03 03  RDN 04  STO 02 05  X<>Y 06  STO 01 07  XEQ "DEGXY" 08  1 09  + 10  X<>Y 11  XEQ "DEGXY" 12  ST+ Z 13  + 14  1 15  + 16  * 17  2 18  / 19  RCL 03 20  + 21  DSE X              22   E3 23  / 24  RCL 03 25  + 26  STO 03 27  CLRGX    28  CLX       29  STO 00   30  LBL 01 31  RCL 03 32  RCL 00 33  +       34  STO 04 35  LASTX             36  8 37  * 38  1 39  ST+ 00 40  + 41  SQRT 42  1 43  - 44  2 45  / 46  INT 47  STO 05 48  CLX 49  STO 06 50  RCL 02 51  STO 07 52  LBL 02 53  RCL IND 01    54  RCL IND 07 55  * 56  ST+ IND 04 57  ISG 06 58  GTO 03 59  RCL 05 60  ST+ 04 61  RCL 06 62  INT 63   E3 64  / 65  STO 06 66  LBL 03             67  ISG 04 68  CLX 69  ISG 07 70  GTO 02 71  ISG 01 72  GTO 01 73  RCL 03 74  END

( 113 bytes )

 STACK INPUTS OUTPUTS Z bbb.eee / Y bbb.eee' / X bbb'' bbb.eee''

with  bbb and bbb' > 007

Example:

p(x,y) = 1 + 2 x + 3 y + 3 x2 + 2 x.y - 4 y2 + 2 x3 + 4 x2y - x.y2 + 3 y3
q(x,y) = 2 + 3 x + 4 y + 5 x2 + 6 x.y + 7 y2

-If we choose Rbb = R11 for p(x,y)   and   Rbb' = R21 for q(x,y)

-Store    1       2      3       3      2     -4      2       4     -1      3
into   R11  R12  R13  R14  R15  R16  R17  R18  R19  R20   respectively    ( control number = 11.020 )

and     2       3       4       5        6       7
into   R21  R22   R23   R24   R25   R26   respectively   ( control number = 21.026 )

-We can take for instance  Rbb" = R31

11.020  ENTER^
21.026  ENTER^
31     XEQ "PROXY"  >>>>   bbb.eee'' =  31.051   ( in 47 seconds )   and we have the coefficients of the product in registers R31 thru R51, namely:

p(x,y) q(x,y) =  2 + 7 x + 10 y + 17 x2 + 27 x.y + 11 y2 + 23 x3 + 53 x2y + 26 x.y2 + 11 y3 + 21 x4 + 48 x3y + 26 x2y2 - 5 x.y3 -16 y4
+ 10 x5 + 32 x4y + 33 x3y2 + 37 x2y3 + 11 x.y4 + 21 y5

2°)  N Variables

-In the following programs, each monomial  µ xa yb zc td ue  is stored into 2 consecutive registers
-The coefficient µ is stored in one data register  and the exponents are coded in the next register under the form  aa.bbccddee

-For instance, to store  12 x2 y3 z4 t15 u6  in registers  R14-R15 ,

12          STO 14
2.03041506  STO 15

-Thus, we assume that all the ( partial ) degrees are smaller than 99 and that the number of variables does not exceed 5.
-However, if we are sure that all the degrees will be smaller than 10, we could deal with polynomials of up to 10 variables ( x2 y3 z4 t5 u6 v7 is coded 2.34567 )

-On the other hand, if you want to handle degrees up to 999, the maximum number of variables will be reduced to 3 and each exponent will be coded with 3 digits.
( x21 y34 z471 will be coded 21.034471 )

a) Evaluation of Polynomial Functions

Data Registers:              R00: temp              ( Registers R01 thru Rnn  &  Rbb thru Ree are to be initialized before executing "PEVAL" )

•  R01 = x    •  R02 = y   ......
•  Rbb  ......  •  Ree = coefficients of the polynomial
Flags: /
Subroutines: /

-If you don't want to use the synthetic registers  M N O , replace them by R06 R07 R08
-In this case, delete line 32 and replace line 01 by  0  STO 06  RDN

-Replace line 21 by  10  if you code the exponents with one digit - provided they are smaller than 10
or  by  E3  if you code the exponents with  3  digits - provided they are smaller than 1000

 01  LBL "PEVAL" 02  CLA 03  STO 00 04  STO O 05  LBL 01 06  CLX 07  STO N 08  RCL IND 00 09  ISG 00 10  RCL IND 00 11  LBL 02 12  ISG N 13  CLX 14  FRC 15  RCL IND N 16  LASTX 17  INT 18  Y^X 19  ST* Z 20  CLX 21   E2 22  *          23  X#0? 24  GTO 02         25  X<>Y 26  ST+ M 27  ISG 00 28  GTO 01 29  RCL O          30  SIGN 31  X<> M 32  CLA 33  END

( 58 bytes )

 STACK INPUTS OUTPUTS X bbb.eee p(x,y,z,...) L / bbb.eee

with bbb > nnn  ( n = number of variables )

Example:      p(x,y,z) = 3 x y4 z7 + 4 x2 y2 z3 - 7 y5 z + 2      Evaluate  p(2,3,4)

-If you choose to store the first coefficient into R11

3   STO 11    1.0407  STO 12
4   STO 13    2.0203  STO 14
-7  STO 15    0.0501  STO 16
2   STO 17        0       STO 18            ( control number = 11.018 )

2  STO 01
3  STO 02
4  STO 03
11.018   XEQ "PEVAL"  >>>>   p(2,3,4) = 7965038   ( in 8 seconds )

Note:

-There will be a DATA ERROR line 18 if a variable and a corresponding exponent equal zero because the HP-41 refuses to compute 00
-So, you can add after line 17    X#0?  GTO 02  SIGN  STO Y  LBL 02     and you'll get for instance the correct result  p(0,3,4) = -6802

-You can also replace line 18 by a M-Code routine  Y^X  that does produce 00 = 1 for example:

098   "X"
01E   "^"
019   "Y"
0F8   READ3(X)
10E   A=C  ALL
0B8   READ2(Y)
355   ?NCXQ          this subroutine checks A and C for alphadata
050   14D5              executes A<>C  and  SETDEC
128   WRIT4(L)      saves X in L-register
2EE   ?C#0?  ALL
02F   JC+ 05
04E   C=0  ALL      C
35C   PT= 12          =
050    LD@PT-1     1
023   JNC+04
044   CLRF 4          C
045   ?NCXQ          =
06C   1B11            A^C
0E8   WRIT3(X)
078   READ1(Z)
0A8   WRIT2(Y)
046   C=0 S&X
270   RAMSLCT
038   READDATA
068   WRIT1(Z)
3E0   RTN

XEQ "Y^X"  will work like the built-in function Y^X but it yields 00 = 1

-However, this M-Code routine does not check for overflows and underflows.

b) Product of 2 Polynomials

Data Registers:              R00 to R05: temp              ( Registers Rbb thru Ree  &  Rbb' thru Ree' are to be initialized before executing "PMUL" )

•  Rbb  ......  •  Ree = coefficients of the 1st polynomial           bb > 05
•  Rbb'  .....  •  Ree' = coefficients of the 2nd polynomial        bb' > 05
Flag:  F07 is cleared
Subroutine:  "DEL0"  ( cf next paragraph )

 01  LBL "PMUL" 02  STO 03 03  DSE X 04   E-3 05  STO 05 06  ST+ 05 07  * 08  ST+ 03 09  RDN 10  STO 02 11  STO 04 12  X<>Y 13  STO 01 14  SIGN 15  2 E-5 16  ST+ 01 17  + 18  STO 00 19  SF 07 20  LBL 01 21  RCL IND 01 22  RCL IND 04 23  * 24  RCL 00 25  RCL 01 26  + 27  ISG 04 28  RDN 29  RCL IND T 30  RCL IND 04 31  + 32  RCL 00 33  RCL 03 34  + 35  X<>Y 36  FS?C 07 37  GTO 03 38  LBL 02 39  RCL IND Y 40  X=Y? 41  GTO 04 42  RDN 43  ISG Y 44  GTO 02 45  LBL 03 46  STO IND Y 47  SIGN 48  - 49  X<>Y 50  STO IND Y  51  RCL 05 52  ST+ 03 53  GTO 05 54  LBL 04 55  RDN 56  SIGN 57  - 58  X<>Y 59  ST+ IND Y 60  LBL 05 61  ISG 04 62  GTO 01 63  RCL 02 64  STO 04 65  ISG 01 66  GTO 01 67  RCL 03 68  XEQ "DEL0" 69  END

( 113 bytes )

 STACK INPUTS OUTPUTS Z bbb.eee / Y bbb.eee' / X bbb'' bbb.eee''

with  bbb and bbb' > 005

Example:            p(x,y,z) = 2 x + 3 x.y2z + 4 x3y.z4                                          q(x,y,z) = 4 - 6 y2z  + x2y.z4

-With  bb = 07 & bb' = 13

2  STO 07         3      STO 09          4      STO 11     control number             4  STO 13         -6     STO 15          1      STO 17          control number
1  STO 08    1.0201  STO 10     3.0104  STO 12          = 7.012                   0  STO 14     0.0201  STO 16    2.0104  STO 18             = 13.018

-If you choose  bb'' = 21

7.012   ENTER^
13.018  ENTER^
21      XEQ "PMUL"  >>>>   21.030   ( execution time = 15 seconds )    and we have:

R21 = 8      R23 = 4              R25 = 18             R27 = -18           R29 = -21
R22 = 1      R24 = 5.0208     R26 = 3.0104      R28 =  1.0402     R30 =  3.0305

-So,    p(x,y,z) q(x,y,z) =  8 x + 4 x5y2z8 + 18 x3y.z4 - 18 x.y4z2 - 21 x3y3z5

Notes:

-Each new monomial is compared to the other ones to collect like terms as soon as possible ( lines 38 to 44 & 54 to 59 )
-It reduces the number of registers during the calculations.
-A shorter program may be written if the HP-41 simplifies the polynomial at the end
but sometimes, ( too? ) many registers are used, especially if there are many like-terms.

-"SIM" is listed in paragraph 2-d)

 01  LBL "PMUL2" 02  STO 00 03  STO 03 04  RDN 05  STO 02 06  STO 04 07  X<>Y 08  STO 01 09  LBL 01 10  RCL IND 01 11  RCL IND 04 12  * 13  STO IND 03 14  1 15  ST+ 03 16  ST+ 04 17  RCL 01 18  + 19  RCL IND X 20  RCL IND 04 21  + 22  STO IND 03 23  ISG 03 24  CLX 25  ISG 04 26  GTO 01 27  RCL 02 28  STO 04         29  ISG 01 30  ISG 01 31  GTO 01 32  RCL 03 33  DSE X 34   E3 35  / 36  RCL 00 37  + 38  XEQ "SIM"   39  END

( 69 bytes )

 STACK INPUTS OUTPUTS Z bbb.eee / Y bbb.eee' / X bbb'' bbb.eee''

with  bbb and bbb' > 004

-The same example is solved in 21 seconds ( instead of 15 )
-Note that the monomials are stored in a different order.
-But if you are sure that there is no like-term in the product, delete line 38 and "PMUL2" will be much faster than "PMUL"

c) Deleting the Zeros

-This subroutine is called by "PMUL" and "SIM" to eliminate the monomials that equal zero.

Data Registers:              R00 to R02: temp              ( Registers Rbb thru Ree  are to be initialized before executing "DEL0" )

•  Rbb  ......  •  Ree = coefficients of the polynomial           bb > 03
Flags:  /
Subroutines:  /

 01  LBL "DEL0"  02  STO 00 03  STO 01 04  FRC 05   E3 06  * 07  STO 02 08  2 E-5 09  ST+ 00 10  LBL 01 11  ISG 00 12  FS? 30 13  GTO 02 14  RCL IND 00 15  X#0? 16  GTO 01 17  RCL 02 18  DSE X 19  RCL IND X 20  STO IND 00 21  RCL 00 22  1 23  + 24  RCL IND 02 25  STO IND Y  26  2 27  ST- 00 28  ST- 02 29   E3 30  / 31  ST- 00 32  ST- 01 33  GTO 01 34  LBL 02 35  RCL IND 01 36  X=0? 37  STO IND 02 38  RCL 01 39  END

( 71 bytes )

 STACK INPUTS OUTPUTS X bbb.eee bbb.eee'

with bbb > 003

-  bbb.eee = bbb.eee'  if all the coefficients are different from zero
-If all the coefficients equal zero  eee' = bbb+1 &  Rbb = Ree' = 0

d) Simplifications

-This subroutine collects like terms to simplify a polynomial expression.

Data Registers:              R00 to R04: temp              ( Registers Rbb thru Ree  are to be initialized before executing "SIM" )

•  Rbb  ......  •  Ree = coefficients of the polynomial           bb > 04
Flags:  /
Subroutine:  "DEL0"  ( cf § 2-c) )

-Line 51 may be replaced by  GTO "DEL0"

 01  LBL "SIM"    02  STO 01 03  1.00002 04  + 05  STO 02 06  RCL 01 07  FRC 08   E3 09  * 10  STO 00 11  LBL 01 12  RCL 02 13  STO 03 14  RCL IND X 15  LBL 02 16  ISG 03 17  FS? 30 18  GTO 03 19  RCL IND 03 20  X<>Y 21  X#Y? 22  GTO 02 23  STO 04 24  RCL IND 00 25  STO IND 03 26  RCL 02 27  1 28  ST- 00 29  ST- 03 30  - 31  RCL IND 00 32  X<> IND 03 33  ST+ IND Y 34  2 E-3 35  ST- 02 36  ST- 03 37  DSE 00         38  SIGN 39  ST- 03 40  RCL 04 41  GTO 02 42  LBL 03 43  ISG 02 44  GTO 01 45  RCL 00 46   E3 47  / 48  RCL 01 49  INT 50  + 51  XEQ "DEL0" 52  END

( 96 bytes )

 STACK INPUTS OUTPUTS X bbb.eee bbb.eee'

with bbb > 004

Example:      p(x,y) = 2 x + 3 x2y + 4 x.y + 5 x2y2 + 8 x.y + 7 x2y

-If we store this polynomial into R05 R06 .....

2  STO 05       3    STO 07       4    STO 09       5    STO 11       8     STO 13       7    STO 15
1  STO 06    2.01  STO 08    1.01  STO 10    2.02  STO 12     1.01  STO 14    2.01  STO 16

5.016  XEQ "SIM"  >>>>   5.012      ( in 8 seconds )      and we find

R05 = 2    R07 = 10      R09 = 12      R11 = 5
R06 = 1    R08 = 2.01   R10 = 1.01   R12 = 2.02

-So   p(x,y) = 2 x + 10 x2y + 12 x.y + 5 x2y2

Note:

-This program may be used to add 2 polynomials:
-Store the 2 polynomials into consecutive registers so that their control numbers (bbb.eee)1 and (bbb.eee)2 verify  bbb2 = 1+eee1
-Then,  key in  bbb1.eee2  XEQ "SIM"

e) Partial Derivatives

-This program calculates  a+b+c+...p / xaybzc....       [  d a+b+c+...p / dxadybdzc....  if the symbol  doesn't appear correctly ]

Data Registers:              R00 to R04: temp              ( Registers Rbb thru Ree  are to be initialized before executing "PDRV" )

•  Rbb  ......  •  Ree = coefficients of the polynomial           bb > 04
Flags:  /
Subroutine:  "SIM"  ( cf § 2-d) )

-Replace line 31 by  10  if you code the exponents with one digit - provided they are smaller than 10
or  by  E3  if you code the exponents with  3  digits - provided they are smaller than 1000

 01  LBL "PDRV" 02  STO 00 03  X<>Y 04  STO 01 05  STO 02 06  ISG 02 07  LBL 01 08  RCL IND 02 09  RCL 00 10  ST- IND 02 11  DSE 02 12  LBL 02 13  RCL Y 14  INT 15  ST- Z 16  RCL Y 17  INT 18  ST- Z 19  X=0? 20  GTO 04 21  X<>Y 22  LBL 03 23  ST* IND 02 24  DSE X 25  ""       26  DSE Y 27  GTO 03 28  LBL 04 29  RDN 30  CLX 31   E2       32  ST* Z  33  * 34  X#0? 35  GTO 02       36  2 37  ST+ 02 38  ISG 02 39  GTO 01 40  RCL 01 41  XEQ "SIM"  42  END

( 74 bytes )

 STACK INPUTS OUTPUTS Y bbb.eee / X aa.bbccddee bbb.eee'

with bbb > 004

Example:      p(x,y,z) = 2 x3 y4 z5 +  x2 y z3 + x y z3 + 3 x12 y7 z10 + 4       Calculate   6p / x2y z3   ( aa.bbcc = 2.0103 )

-If p(x,y,z) is stored in R05 R06 ....

2      STO 05         1        STO 07         1       STO 09          3        STO 11     4   STO 13
3.0405  STO 06     2.0103   STO 08    1.0103   STO 10     12.0710  STO 12     0   STO 14      control number = 5.014

5.014     ENTER^
2.0103   XEQ "PDRV"  >>>>   5.010   ( in 24 seconds )  and we have

R05 = 2880       R07 = 12    R09 = 1995840
R06 = 1.0302    R08 =  0     R10 = 10.0607

-Thus,    6p / x2y z3 = 2880 x y3 z2 + 12 + 1995840 x10 y6 z7

Note:

-This routine replaces the polynomial by its derivative.
-So you'll have to store it in a second block of registers if you want to save it.