3
$\begingroup$

I have a list of polynomials polys.

I want a set of $f(x) \in$ polys which have roots strictly in $\{S^1,0\}$, with $S^1$ being the unit sphere.

They are all degree $n$, so I have crudely picked out those $f(x)$ such that the sum of the norms of the roots of $f$ is $\leq n$.

rts = Table[Values[{ToRules[Roots[Part[polys, i] == 0, x]]}], {i, 1, Length[polys]}] rtsAbsSum = Flatten[Table[Sum[Abs[Part[Part[rts, i], j]], {j, 1, Length[Part[rts, i]]}], {i, 1, Length[polys]}]] parts = Select[Range[Length[polys]], Part[rtsAbsSum, #] <= n &] polysIWant = Table[Part[polys, i], {i, parts}]

Now the list polysIWant contains the polynomials I want. However, it also contains many false positives.

This code is crude. What its doing:

  1. Takes each polynomial of polys and replaces it with a list of its roots (multiplicity not counted.)
  2. Adds together the Abs of the elements of each list.
  3. Finds the polynomials whose sum of root norms is $\leq n$.

Is there an easier way? ie. Instead of summing roots, checking the inequality, etc., how can I write a program that tosses out the polynomials in polys with a root whose norm is neither $0$ nor $1$, and keeps the rest?

Improve this question
$\endgroup$
5
  • $\begingroup$ If you have a function roots that gets the roots of a polynomial into a list form of consistent dimension, then why not something like: Select[polys, Not[0!=Norm[roots[#]]!=1] &] $\endgroup$ Commented Aug 6, 2020 at 18:08
  • $\begingroup$ I assume that by $S^1$ you mean that unit circle in the complex plane? So you are allowing complex roots? Just for clarification. $\endgroup$ Commented Aug 6, 2020 at 20:34
  • $\begingroup$ Apparently there is something called the Jury stability criterion which you can use to determine if the roots of a polynomial are all inside the unit disk. $\endgroup$ Commented Aug 6, 2020 at 20:42
  • 2
    $\begingroup$ Are "roots of undesirable magnitude" related to "rodents of unusual size"? $\endgroup$ Commented Aug 8, 2020 at 16:33
  • 2
    $\begingroup$ R.O.U.M.s? I do not think they exist. $\endgroup$ Commented Aug 11, 2020 at 2:58

1 Answer 1

Reset to default
6
$\begingroup$

Here is an implementation of the Jury stability criterion (also known as the Jury-Marden criterion; there is also the Schur-Cohn criterion), which determines if the roots of a real-coefficient polynomial are all in the open unit disk. The function returns True if all of the roots of the polynomial are in the open unit disk and returns False otherwise.

juryStabilityCheck[poly_, var_] := juryStabilityCheck[Reverse@CoefficientList[poly, var]] juryStabilityCheck[coefs_List] := Module[ {lst = coefs}, If[First@lst < 0, lst = -lst]; While[lst != {} && First@lst > 0, lst = lst - Last@lst/First@lst Reverse@lst // Most]; lst == {} ]

Here is a polynomial whose roots are all in the open unit disk:

poly = (z - 0.5) (z - I/2) (z + I/2) (z - 0.3 I) (z + 0.3 I);

If we run the check on it, we get:

juryStabilityCheck[poly, z] (* True *)

It works for any scalar multiple of that list:

juryStabilityCheck[-10 poly, z] (* True *)

If some roots are outside the unit disk:

poly = (z - 0.5) (z - 3 I/2) (z + 3 I/2) (z - 0.3 I) (z + 0.3 I); juryStabilityCheck[poly, z] (* False *)

If a root is on the unit circle (i.e. the boundary of the unit disk), then the check finds that they are outside the open unit disk:

poly = (z - 1) (z - I/2) (z + I/2) (z - 0.3 I) (z + 0.3 I); juryStabilityCheck[poly, z] (* False *)

Now, if you want to include polynomials that have a root on the unit circle, the quick fix I've come up with is to "enlarge" the circle slightly by "shrinking" the z's. You must choose a tolerance, and be okay with polynomials being included that have roots within this tolerance outside the unit circle.

Add to the code the definition

juryStabilityCheck[poly_, var_, tolerance_] := juryStabilityCheck[poly /. var -> var (1 + tolerance), var]

As an example:

poly = (z - 1) (z - I/2) (z + I/2) (z - 0.3 I) (z + 0.3 I) // Expand; juryStabilityCheck[poly, z] juryStabilityCheck[poly, z, 10^-9] (* False *) (* True *)

Of course, this will also include a polynomial whose roots are within the tolerance but outside the unit circle:

poly = (z - 1 + 10^-10) (z - I/2) (z + I/2) (z - 0.3 I) (z + 0.3 I); juryStabilityCheck[poly, z, 10^-9] (* True *)
Improve this answer
$\endgroup$
2
  • $\begingroup$ Thank you, this is excellent. I've made a one-character change to your code that I believe was a typo (at least, everything worked AFTER I made it.) $\endgroup$ Commented Aug 6, 2020 at 23:02
  • $\begingroup$ For searching purposes: this is also sometimes referred to as the "Jury-Marden" criterion. There is also the more classical Schur-Cohn criterion for doing the same test. $\endgroup$ Commented Aug 7, 2020 at 9:12

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.