Essays/PolynomialsIntersection

From J Wiki
Jump to: navigation, search

Find the intersection of polynomials

Two polynomials

Given two polynomials, for example:

f(x)= 3 + x (a)
f(x)= 1 + 2x + x^2 (b)

how can we find their intersection?

The coefficients of the above polynomials from lowest to highest order are:

a=: 3 1
b=: 1 2 1

The dyadic form of the J primitive p. (Polynomial) evaluates a polynomial of order #x with coefficients x for the argument y. So f(-2) for polynomial b is

   1 2 1 p. _2
1
Abintersect.png
plot _3 3;'a p. y ` b p. y'

We can plot these polynomials as follows:

   load 'plot'
   plot _3 3;'3 1 p. y ` 1 2 1 p. y'
NB. or
   plot _3 3;'a p. y ` b p. y'

From the plot we can see that the two polynomials intersect when x is _2 or 1. How can we get that result using J?

This is equivalent to finding the roots (the values for x where a polynomial is zero) of a polynomial formed by subtracting polynomial a from polynomial b.

Firstly subtract one polynomial from another:

   a (-/@,:) b
2 _1 _1

i.e.    f(x)= 2 - x - x^2      (c)

Then find the roots of the resulting polynomial c using the monadic form of the J primitive p. (Roots):

   p. 2 _1 _1
┌──┬────┐
│_1│_2 1│
└──┴────┘

The roots are contained in the 2nd box (the first contains the multiplier) so we can put these ideas all together to give a vector of the x values where two polynomials intersect:

   findIntersect=: 1 {:: [: p. -/@,:

   a findIntersect b
_2 1

Where the roots of the polynomial formed by subtraction are complex, findIntersect will return the complex roots.

   6 2 1 findIntersect 1 2
0j2.236068 0j_2.236068

If we only wished to return real roots we could extend findIntersect as follows:

   6 2 1 (#~ (= +))@findIntersect 1 2

Two or more polynomials

The examples given above are probably the best way to find the intersection of 2 polynomials. If you want the common intersection of several, here is a less-succinct method.

You can test whether f_0, f_1, ... are all equal by forming the sum of (f_i - f_j)^2 and setting it to zero.

ppr   =: +//.@(*/)  NB. polynomial product
pdiff =: -/@,:      NB. polynomial difference
pps   =: ppr~      NB. polynomial square

comb=: 4 : 0
 k=. i.>:d=.y-x
 z=. (d$<i.0 0),<i.1 0
 for. i.x do. z=. k ,.&.> ,&.>/\. >:&.> z end.
 ; z
)

NB. findIntersectM <list of boxes of coefficients>
NB. returns x-coordinates of common intersection points
findIntersectM=:3 : 0
a=. >y
c=. 2 comb #a
~. (#~ (= +)) 1{:: p. +/ pps@pdiff /"2  c { a
)

For example to find the intersection of the polynomials:

f(x)= x^2 (d)
f(x)= 2 - x^2 (e)
f(x)= 1 (f)
d=: 0 0 1
e=: 2 0 _1
f=: 1
Cdeintersect.png
plot _3 3;'d p. y`e p. y`f p. y'

We can plot these as follows.

   require 'plot'
   plot _3 3;'d p. y ` e p. y ` f p. y'

By inspection these polynomials intersect when x is 1 or _1.

Given a list of boxed coefficients for each polynomial we can find the x coordinates where they all intersect.

   ]r=: findIntersectM c;d;e
1 _1

We can evaluate each polynomial at those x values to show that they all have the same f(x).

   d p. r
1 1
   e p. r
1 1
   f p. r
1 1

Contributions

This essay was compiled by Ric Sherlock from posts in this forum thread by John Randall, RaulMiller and Henry Rich.

See Also

  • The Polynomials lab available from the J session (Studio | Labs... | Polynomials).