Software for approximation of zeros of polynomials

In this webpage, a set of software modules is described, which can be downloaded from this website. This webpage describes how to compile and run the software.

There are multiple modules, suitable for different use cases. All of the modules should be used from Java programs. The modules only do the actual arithmetic. Things like plotting, printing solutions and so on are not part of the modules. A helper module is provided, which can be used to make simple graphs and diagrams. Below follow descriptions of all different modules and how they can be compiled and used. Here, a complete download list is given, but links are provided further below as well in the sections about each particular piece of software:

• Helper class for generating PNG Images
• Helper classes for complex and 105 bit precision
• Java port of RPoly, CPoly and MPSolve 1996 algorithms
• Test code for RPoly, CPoly and MPSolve 1996
• Demo code using RPoly, CPoly and MPSolve 1996
  
• Source code for C-part of MPSolve 2012 algorithm
• Source code for Java-part of MPSolve 2012 algorithm
  

A little helper class for generating PNG images

Arithmetic helper classes: complex and 105-bit precision

• The class Complex provides basic complex arithmetic on standard 53 bit doubles.
• The class DoubleDouble provides basic arithmetic on 105 bit double double numbers. The underlying mathematics is quite advanced, but using this class is simple. If you have used a class like BigDecimal in Java, then the use of the DoubleDouble class should be easy for you.
• The class DoubleComplex provides basic arithmetic on 105 bit complex numbers. It uses the DoubleDouble class for the underlying 105 bit arithmetic.

Java ports of fixed precision polynomial solvers in a single archive

MPSolve 1996 version ported to Java

MPSolve 1996 is the first version of the MPSolve algorithm, developed by D.A. Bini. It is written in Fortran 90. It is one of the more robust algorithms for approximating the zeros of polynomial equations presently available. This algorithm was developed in the 1990's and a paper is published in 1996. Unfortunately, no port of this algorithm to another language than Fortran is available and hence I decided to make the port myself. The program is ported to Java and is available in standard double precision at 53 bits and in so-called double double precision at 105 bits. The port is a full Java implementation and it is completely self-contained, not depending on any special third-party library. The remarkable property of this algorithm is that it allows switching from a lower to a higher precision when computation at the lower precision does not help further improving the solutions or finding additional solutions. This is the basis for efficient high precision code.

The package net.oelen.polsolve.pzeros provides a single public class PZeros, which can be used to find the zeros of polynomials. The other classes are internal to the package. The class  PZeros provides a method for solving real or complex double precision polynomials, and a method for solving real or complex double double precision polynomials. This package uses the complex and double double arithmetic helper classes, described above.

A test program is written, which finds the zeros of 100000 polynomials with double double precision complex coefficients, having zeros in the square [-5, 5] + i∙[-5, 5]. The degree of these polynomials is taken randomly from 12 to 50. On a 3.4 GHz system, this took just over 80 seconds, meaning that appr. 1250 such polynomial equations can be solved per second in a single thread.

There are no special requirements for compiling this code. Using any IDE like Eclipse or NetBeans, you can compile the code without issue and make a JAR-file, containing the solver code. Be sure to have the complex and double double helper classes in the class path.

Jenkins Traub RPoly and CPoly ported to Java

In the 1960's, an algorithm was developed by Jenkins and Traub. This algorithm is implemented in Fortran and both a general version, named CPoly, was written and a specialized version for polynomials with real coefficients, named RPoly. RPoly and CPoly are published in the Netlib repository, published by the ACM.

From the original CPoly program there is a literal port from Fortran to Java, written by Joseph A. Huwaldt in 2000. I added some additional code to make the use of this port easier. The CPoly program works for general complex polynomials. I added a minor tweak to have it iterate once more once it finds a root, to make it somewhat more accurate. This works fine and is tested with millions of random polynomials with varying type of distribution of zeros.

The original RPoly program unfortunately was seriously flawed and results in inaccurate results, even for well-conditioned polynomial equations of low degree. An open source Java port of this algorithm is available, but it has the same flaw as the Fortran code. Commercial Java implementations are available, apparently these ports are better, but I have no personal experience with them. There is a C++ port as well. The C++ port has the same flaw as the original Fortran code. In my port from the Fortran code, the flaw is resolved. Accuracy still can be somewhat low (i.e. only 3 to 4 accurate digits for some zeros), but this only happens in rare ill-conditioned cases (e.g. multiple roots or very close roots) and I have not found any polynomial for which the algorithm fails to converge. Its speed is great, the RPoly software can solve well over 35000 polynomials with a degree randomly taken from the range 8 . . . 22 and zeros in the square [-5, 5] + i∙[-5, 5] per second on a 3.4 GHz Core i5 in a single thread. I also compared the speed of my port with the above mentioned Java-port and my port is appr. 5 times as fast for the same polynomial.

For real polynomials of up to degrees 30 or so, the program RPoly is one of the best available programs if only standard double precision is available. It is tested with hundreds of millions of random polynomials of all kinds (test program is included) and it has proven to be very reliable. Many images from the experiments web page were produced with the RPoly program. For many engineering problems (e.g. in electronics) this program also is extremely useful.

The CPoly program is approximately 2.5 times as slow as the RPoly program for the same polynomial, but it can be used for general polynomials with complex coefficients.

Compiling RPoly and CPoly does not require anything special. Any environment, capable of building Java programs, can be used to compile RPoly and CPoly into a single JAR file. They do not depend on any other classes, they have no imports at all.

Test code for MPSolve 1996 and for RPoly/CPoly

For RPoly, CPoly and MPSolve 1996 (PZeros), test programs are written, which generate lots of random polynomials and compute zeros of these. Some of the tests already are described above. The test programs contain a main() method, which can call other specific methods which perform certain tests. Besides test programs for RPoly, CPoly and MPSolve 1996 the package contains a test suite of some interesting polynomials.

The test programs can be downloaded here and it is recommended to have a look at them, modify them and run them in your own environment.

Demo code, using RPoly, CPoly and MPSolve 1996

MPSolve 2012 version

The previous programs are useful for computing zeros of low to moderate degree polynomials at high speed and with 105 bits precision one can go quite far with polynomials of degree 100 or so. The MPSolve 2012 algorithm, however, takes solving polynomial equations to the next level. There are virtually no limits in degree and no limits in accuracy. The previous programs are useful for all kinds of engineering problems and that also is their main intended use. The MPSolve 2012 package, however, is especially intended for doing research on all kinds of special polynomials, studying asymptotic behavior of special functions and other research topics from fundamental mathematics. This algorithm is used for solving special polynomial equations of degrees up to several millions, requiring precisions of thousands or even tens of thousands of bits.  There are two packages to be downloaded:

• The C-software, used for building a JNI library
• The Java software, which makes the C-software available from Java-programs.
  

Unfortunately, no native Java version is available of this, and there also never will be. The software is written in C and for its multiprecision computations, it uses the GMP package. Porting this to Java is not feasible and porting the MPSolve 2012 code itself to Java also is not feasible. The package is large and it is remarkably complicated.

The algorithm is developed by Bini and Robol and they provide a C-implementation. This package, however, contains quite a few serious issues. It has several memory leaks, sometimes simply crashes, and for higher precisions it becomes slow due to a serious bug in the refinement code. The MPSolve 2012 algorithm is a masterpiece of mathematics and logical thinking, but the software implementation is not that good. The code, provided here, is loosely based on the software implementation, provided by Bini and Robol, but it is largely refactored and some parts are rewritten. The original software can be used to solve equations with integer coefficients, rational coefficients and floating point coefficients. The refactored and rewritten version only supports integer and rational coefficients (both real and complex). The floating point support of the original software is of such low quality that I decided to leave that out completely and maybe in a future version I add my own support for this. There is no real need for that. For any given precision, a floating point number can be approximated by a quotient P/Q, with P and Q being integers.

Part of the software is written in C, and other parts are written in Java. My refactored version of the software is intended to be used from within Java code, but the software can also be used from within C. For all experiments, described on this website, however, all usage of the software is from Java code. The C-library is made available to Java runtime environments by means of JNI. A more formal specification of JNI can be found on the official Oracle page. For compiling and using the software there is no need to study the JNI-specification, this info is purely added for background purposes and for the interested reader/developer.

Building the MPSolve 2012 package requires the use of a build environment, which contains several programs. Below it is described how to setup a build environment for building the C-part of the software. The final result of this build-stage is a loadable module, which can be loaded into a JVM for usage in Java programs. Below, it is described how to setup the build-environment for Linux/Ubuntu on x86/x86-64, Linux/Ubuntu on ARM, and Windows 8.x/10 on x86/x86-64.

Preparing the system for building the C-part on Linux/Ubuntu x86

Before you can build the MPSolve 2012 software, you have to setup a Linux/Ubuntu system with the suitable build tools. Below, a procedure is given for setting up such a system on Ubuntu 16.01 server, but the same procedure also can be used on a desktop version of Ubuntu and also on previous versions of Ubuntu (e.g. 14.04) or other Debian-based systems.

Once you have a working installation and have access to a shell prompt, perform the following steps as super user (a.k.a. "root"). You can obtain root-privileges with the command:

• sudo su -
  

First update the system to the newest bugfixes and security fixes, by issuing the following commands while having access to the internet:

• apt-get update
• apt-get dist-upgrade
  

Then reboot the system by issuing the command reboot, still logged in as root. Once the system is up again, login and become super user again. Issue the following two commands for installation of the compiler tool chain and the m4 preprocessor:

• apt-get install build-essential
• apt-get install m4

• mkdir /opt/jdk
• cd /opt/jdk
• cp ..../jdk-8uxxx-linux-x64.tar.gz .
• tar zxvf jdk-8uxxx-linux-x64.tar.gz

Here .... is the path where you have the jdk-8uxxx-linux-x64.tar.gz file.

Finally, make a symbolic link to the created directory in the /opt/jdk directory. Give the following command while you are still in /opt/jdk:

• ln -s jdk1.8.0-xxx current

The use of this current symbolic link has the advantage, that if you at a later point in time update your JDK to a newer version, you only have to unpack the newly downloaded archive file and have to create a new symbolic link with the name current and all other parts of the system can be left unchanged.

Now you have the system ready for building JNI connectors for access from Java code to C/C++ code. For convenience it is best to do some path and variable settings in the system, so that each time you login, you have the environment set for immediate use of Java.

• Edit file /etc/environment: Add path "/opt/jdk/current/bin" to start of PATH variable.
• Add line "JAVA_HOME=/opt/jdk/current" to the same /etc/environment file.
  

Reboot the system again and you are ready to build the C-part of MPSolve, including the JNI library for access from Java packages. Further instructions are in the downloaded source code.

Preparing the system for building the C-part on Linux/Ubuntu ARM

• jdk-8uxxx-linux-arm32-vfp-hflt.tar.gz for 32 bit processors (e.g. Odroid U3, Odroid C1+, Raspberry Pi 2)
• jdk-8uxxx-linux-arm64-vfp-hflt.tar.gz for 64 bit processors (e.g. Odroid C2, Raspberry Pi 3)

All steps are the same as for Linux/Ubuntu on an x86 system.

Preparing the system for building the C-part on Windows

A suitable and very good build environment for Windows is MSYS2, which is based on the Cygwin and MinGW build environments. If you go to that page, then you can download a 32 bit installation executable, or a 64-bit installation executable. The name of these executables is of the form msys2-i686-yyyymmdd.exe for the 32 bit version and msys2-x86_64-yyyymmdd.exe for the 64 bit version (here, yyyymmdd stands for year, month, day, e.g. 20161025 for October 25 in the year 2016). For all modern versions of Windows you can take the 64 bit version, 32 bits Windows is becoming rare these days. Download the file and execute the installer by double clicking it in an explorer. Below, only the installation for 64 bits MinGW is given. If you really need 32 bits, then replace "x86_64" by "i686" in the names of files or packages.

When the installable executable is started, then it prompts for a path, where to install the software. The default location is c:\msys64, which is perfectly fine in most situations. When installation is done, then the program allows you to run msys2 straight away. When this is done, then a command-line window appears with a $ prompt. If by accident you skipped the option to start msys2, then open a CMD.EXE command line window and go to c:\msys64 and give the msys2_shell command while being in the msys64 directory. Next, the following steps must be performed:

• In the msys2 command line window give the command pacman -Syu, which updates the package database. When this command is finished, you cannot use the msys2 window anymore.
• Close the msys2 command line window and open a new one, as described above.
  
• In the new window give the command pacman -Su, which updates several parts of MSYS2.
  
• Finally, in the same window, give the following commands:
• pacman -S git
• pacman -S tar
• pacman -S m4
• pacman -S mingw-w64-x86_64-toolchain
• pacman -S make
• Edit the PATH environment variable under Windows. At the end of the path, add the following directory: <basedir>\mingw64\bin. Here, the directory <basedir> must be set to the location of msys2, which was set in the install wizard. If the default value was accepted then <basedir> is equal to c:\msys64 and then the path must be set to c:\msys64\mingw64\bin.
  

After setting the PATH variable, the 64 bit environment can be used by going to the <basedir> directory and starting the mingw64 program. This starts a shell in which all 64 bit development tools are available, needed for building the JNI connectors for Java programs. Package management can be done in the msys2 shell, but when you want to build software, then you need to start the program mingw64 in the <basedir> directory (for a default install this is the directory c:\msys64).

For the Windows platform, the following JDK version should be downloaded from the Oracle site:

• jdk-8uxxx-windows-i586.exe for 32 bit windows (this is becoming increasingly rare and one usually needs the version from the next bullet)
• jdk-8uxxx-windows-x64.exe for 64 bit windows (Windows 7, 8, 8.1, 10)

Building the C-part of MPSolve 2012

Using MPSolve 2012 from Java programs

back to top of page

back to main polynomial page

back to main experiments page