Polynomials

Introduction

Polynomials are well treated in J as can be seen by here .

Multiplication and division of polynomials however is not as efficient as could be.

Multiplication

Normally, +//.@:(*/) is the usual way of multiplying polynomials (look for 'product').

   1 _2 1 +//.@:(*/) 2 _6 _6 2
2 _10 8 8 _10 2

'X Y'=: 3000 4000 <@(?@$"0)1000

   ts' X +//.@:(*/) Y'
0.29188701 67143552

   # X +//.@:(*/) Y
6999

However, look at this

   rank' X +//.@:(*/) Y';'X pm Y'
+----+------+----+------+
|rank|tm*sz |time|size  |
+----+------+----+------+
| 1  |646.79|1.63|398.00|
+----+------+----+------+
| 0  |  1.00|1.00|  1.00|
+----+------+----+------+

   X (+//.@:(*/) -: pm) Y
1

As is obvious here, pm is a lot leaner and more than a bit faster. Here is the code

pm=: 4 :0
 if. x -: y do. (+/@:(*|.)\y) , }.+/@:(*|.)\.y
 else.
  z=. i.2 0
  'x y'=. x (,|.&.>)/@(; \:,&#) y
  t=. # y
  for_j. i.&.<: t do. z=. z,. ((j {. x) +/@:* (-j) {. y) , ((-j) {. x) +/@:* j {. y end.
  z ({.@[ , ] , |.@{:@[) t +/@:(*&y)\ x 
 end.
)

For squaring a polynomial, ie. multiplying it with itself, the performance improvement is even more impressive.

   rank'+//.@:(*/)~ Y';'pm~ Y'
+----+-------+----+------+
|rank|tm*sz  |time|size  |
+----+-------+----+------+
| 1  |1675.79|3.38|496.51|
+----+-------+----+------+
| 0  |   1.00|1.00|  1.00|
+----+-------+----+------+

   (+//.@:(*/)~ -: pm~) Y
1

Sparse multiplication

For polynomials with much coefficients equal to zero one can do better ie. faster multiplication than pm.

spm=: 4 : '(] /:"1 {.) x ( ,@:(+/)&{. (~.@[ ,: +//.) ,@:(*/)&{: ) y'
NB. x, y and the result are in form indexes ,: values like
NB. 8 15 21 39 45 50 64 73 89 96
NB. 0 13 12 46 23 43 30  5 20 44
NB. or represented as 4&(,@:$.),:5&(,@:$.)

   'X Y'=:50 100 <@(>:@?@#~@[ ,:~ ?@# + (*i.)~)"(0) 10

   $&.>X;Y
+----+-----+
|2 50|2 100|
+----+-----+

   expnd=: ({: #^:_1~ (e.~ i.@>:@(>./))@{.) :. ([:(] /:"1 {.) ](I.@] ,: #~)|@*)

   'Xe Ye'=: expnd&.> X;Y

   $&.>Xe;Ye
+---+---+
|494|991|
+---+---+

   (X;Y)-:]&.(expnd&.>) X;Y
1

   rank'Xe pm Ye';'X spm Y'
+----+-----+-----+----+
|rank|tm*sz|time |size|
+----+-----+-----+----+
| 1  |5.06 |22.60|1.00|
+----+-----+-----+----+
| 0  |1.00 | 1.00|4.46|
+----+-----+-----+----+

spm is the sparse polynomial multiplication, X and Y are the test polynomials and expnd expands a sparse polynomial in a normal one, and back. Some other results are.

   'X Y'=:300 500 <@(>:@?@#~@[ ,:~ ?@# + (*i.)~)"(0) 10
   $&.>' Xe Ye'=:expnd &.> X;Y
+----+----+
|2993|4994|
+----+----+
   rank'Xe pm Ye';'X spm Y'
+----+-----+-----+-----+
|rank|tm*sz|time |size |
+----+-----+-----+-----+
| 1  |1.15 |21.55| 1.00|
+----+-----+-----+-----+
| 0  |1.00 | 1.00|18.70|
+----+-----+-----+-----+

   'X Y'=:300 500 <@(?@#~@[ ,:~ ?@# + (*i.)~)"(0) 20
   $&.>' Xe Ye'=:expnd &.> X;Y
+----+----+
|5995|9981|
+----+----+
   rank'Xe pm Ye';'X spm Y'
++----+-----+-----+----+
|rank|tm*sz|time |size|
+----+-----+-----+----+
| 1  |8.67 |82.24|1.00|
+----+-----+-----+----+
| 0  |1.00 | 1.00|9.49|
+----+-----+-----+----+

   'X Y'=:300 500 <@(?@#~@[ ,:~ ?@# + (*i.)~)"(0) 30
   $&.>' Xe Ye'=:expnd &.> X;Y
+----+-----+
|9000|14991|
+----+-----+
   rank'Xe pm Ye';'X spm Y'
+----+-----+------+----+
|rank|tm*sz|time  |size|
+----+-----+------+----+
| 1  |24.41|116.64|1.00|
+----+-----+------+----+
| 0  | 1.00|  1.00|4.78|
+----+-----+------+----+

   'X Y'=:100 150 <@(?@#~@[ ,:~ ?@# + (*i.)~)"(0) 50
   $&.>' Xe Ye'=:expnd &.> X;Y
+----+----+
|4982|7500|
+----+----+
   rank'Xe pm Ye';'X spm Y'
+----+------+------+----+
|rank|tm*sz |time  |size|
+----+------+------+----+
| 1  |144.04|193.60|1.00|
+----+------+------+----+
| 0  |  1.00|  1.00|1.34|
+----+------+------+----+

So we can conclude that is is worthwile to use sparse multiplication if 90% or more of the coefficients are zero.

Extended precision

For extended precision pm appears to be only much leaner, but a bit slower.

   'X Y'=:900 1200 <@(?@$"0)1000x

   rank' X +//.@:(*/) Y';'X pm Y'
+----+------+----+------+
|rank|tm*sz |time|size  |
+----+------+----+------+
| 1  |279.75|1.00|329.00|
+----+------+----+------+
| 0  |  1.00|1.18|  1.00|
+----+------+----+------+

Contrary to the sparse case, where the perfromance ratio is as with normal precision.

   'X Y'=:30 50 <@(?@#~@[ ,:~ ?@# + (*i.)~)"(0) 30x

   $&.>' Xe Ye'=:expnd &.> X;Y
+---+----+
|880|1479|
+---+----+

   rank'Xe pm Ye';'X spm Y'
+----+------+------+----+
|rank|tm*sz |time  |size|
+----+------+------+----+
| 1  |189.70|200.20|1.00|
+----+------+------+----+
| 0  |  1.00|  1.00|1.06|
+----+------+------+----+

Divison

On the PolynomialsRational (look for 'rational conjunction')division of polynomials is defined a bit hidden. But also this can be done much better.