Multiple Precision Computation

Web search keyword: dsmithfmlibrary

August, 2007: An error was corrected in the polygamma function (routine FMPGAM). Polygamma(N,X) returned the wrong result when N = 0 and X was very close to a root. The only positive root of Polygamma(0,X) is near 1.46, and when using 50-digit precision the failure happened only when X was within 1.0E-10 of a root. Polygamma(0,X) is the same as Psi(X), and Psi did not have this error. An error was corrected in the integer division with remainder routine IMDVIR. This is a specialized routine for dividing a multiple precision integer by a machine precision integer, giving a multiple precision quotient and a machine precision remainder. The error occurs rarely, but is more likely when a small base like 2 is used and the multiple precision input is about the same size as the machine precision input. The corrected versions of FMPGAM and IMDVIR have been put into the 8-21-01 version of FM.f90 available below.

February, 2005: An error was corrected in routine FMADDI. When the sign bit for FM numbers was moved from element (2) of the array to element (-1) in version 1.2, one of the cases was not updated properly. When the FM number in a call to FMADDI was negative and between 1 and the base being used for the arithmetic, a wrong result could be returned. The corrected version of FMADDI has been put into the 8-21-01 version of FM.f90 available below.

August, 2001: Version 1.2 of the FMLIB multiple-precision package is available. The new version has functions for real Gamma and related functions. The code has been updated to Fortran-90 free source form, although the older Fortran-77 version 1.1 is still available below.

FMLIB is a package of Fortran routines for real and complex arithmetic and elementary functions. The precision and base for the arithmetic can be set by the user, and routines are available for floating-point arithmetic, conversion and input/output operations, trigonometric, exponential, logarithmic, and hyperbolic functions.

The package also has routines for integer multiple-precision arithmetic and functions, including GCD, modular products and powers, and a random number generator based on 49-digit primes.

FMZM90 is a Fortran-90 module that defines three derived types: multiple-precision real, multiple-precision integer, and multiple-precision complex. Interfaces are provided so that a program can declare variables to be multiple-precision types and then write the code for multiple-precision operations using normal Fortran syntax.

Example using FMLIB without derived types:

...

USE FMVALS

DOUBLE PRECISION A(-1:LUNPCK),B(-1:LUNPCK),C(-1:LUNPCK)

DOUBLE PRECISION D(-1:LUNPCK),T1(-1:LUNPCK),T2(-1:LUNPCK)

...

CALL FMSQR(A,T1)

CALL FMSQR(B,D)

CALL FMADD(T1,D,T2)

CALL FMSQR(C,D)

CALL FMADD(T2,D,T1)

CALL FMSQRT(T1,D)

...

The same example using derived types:

...

USE FMZM

TYPE ( FM ) A,B,C,D

...

D = SQRT( A**2 + B**2 + C**2 )

...

List of files for FMLIB 1.2: (8-21-01 version)

1-4-06: The server is not able to find files with "f90" or other uncommon filetypes. Until this is fixed, here is a workaround:

Click a given file. If you directly download it, there may be an extra ".txt" filetype at the end of the file name. Delete it. (You might also have to change "f90" to "f95" or something else, depending on your compiler.)

FM.f90
Subroutine library for multiple-precision operations. 36,894 lines of code

FMZM90.f90
Module for derived type interfaces. 7,396 lines of code

FMSAVE.f90
Module for FM internal global variables. 544 lines of code

TestFM.f90
Test program for most of the FM routines. 10,944 lines of code

SampleFM.f90
Small sample program using FM. 639 lines of code

SAMPLE.CHK
Expected output file from SampleFM.f90. 107 lines

ReadMe
User's Guide for the package, along with a list of the files, some troubleshooting advice, and an example set of compiler/linker commands for building the programs. 1,626 lines

FMLIB 1.1: (5-19-97 Fortran-77 version)

FMLIB 1.1

PAPERS ONLINE (in pdf format)

All of these are technical papers explaining the algorithms used by the multiple-precision package, except the last one.

"Using Multiple-Precision Arithmetic" gives some samples and discusses the package from a user's point of view.

Efficient Multiple-Precision Evaluation of Elementary Functions
Mathematics of Computation 52 (1989) 131 -- 134

A Fortran Package For Floating-Point Multiple-Precision Arithmetic
Transactions on Mathematical Software 17 (1991) 273 -- 283

A Multiple-Precision Division Algorithm
Mathematics of Computation 66 (1996) 157 -- 163

Multiple Precision Complex Arithmetic and Functions
Transactions on Mathematical Software 24 (1998) 359 -- 367

Multiple-Precision Gamma Function and Related Functions
Transactions on Mathematical Software 27 (2001) 377 -- 387

Using Multiple-Precision Arithmetic
Computing in Science and Engineering 5 (July, 2003) 88 -- 93

PAPERS ONLINE (individual pages can be viewed with a browser)

Efficient Multiple-Precision Evaluation of Elementary Functions
Mathematics of Computation 52 (1989) 131 -- 134
Page 1, 2, 3, 4

A Fortran Package For Floating-Point Multiple-Precision Arithmetic
Transactions on Mathematical Software 17 (1991) 273 -- 283
Page 1, 2, 3, 4, 5, 6, 7, 8, 9, 10

A Multiple-Precision Division Algorithm
Mathematics of Computation 66 (1996) 157 -- 163
Page 1, 2, 3, 4, 5, 6, 7, 8

Multiple Precision Complex Arithmetic and Functions
Transactions on Mathematical Software 24 (1998) 359 -- 367
Page 1, 2, 3, 4, 5, 6, 7, 8, 9, 10

Multiple-Precision Gamma Function and Related Functions
Transactions on Mathematical Software 27 (2001) 377 -- 387
Page 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11

Using Multiple-Precision Arithmetic
Computing in Science and Engineering 5 (July, 2003) 88 -- 93
Page 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11