(Copy of MPFUN page before significant editing.) |
|||
| Line 5: | Line 5: | ||
}} | }} | ||
===Introduction=== | ===Introduction=== | ||
| + | |||
MPFUN90 is a Fortran-90 arbitrary precision package (i.e. user specifies the desired precision). The library is written by David H. Bailey et. al from the Lawrence Berkeley National Laboratory. If 64 digits of precision are sufficient for an application then users should use QD instead of MPFUN which would run much faster. | MPFUN90 is a Fortran-90 arbitrary precision package (i.e. user specifies the desired precision). The library is written by David H. Bailey et. al from the Lawrence Berkeley National Laboratory. If 64 digits of precision are sufficient for an application then users should use QD instead of MPFUN which would run much faster. | ||
=== Version Selection === | === Version Selection === | ||
| − | + | * This page details how to use MPFUN90 (version 20100825) that is available as a SHARCNET module. | |
| + | ** '''NOTE:''' The only SHARCNET module currently available supports the Intel v15.0.3 Fortran compiler. | ||
| + | * There is an older MPFUN library version 6.1.23 for the intel compiler available on SHARCNET. See the [[MPFUN]] page for those details. | ||
===Job Submission=== | ===Job Submission=== | ||
| − | To compile an | + | To compile an MPFUN90 job you need to have the intel/15.0.3 module loaded. After loading it, use the following |
| − | command to compile an | + | command to compile an MPFUN90 source file, my_mpfun.f90: |
| + | |||
| + | ifort my_mpfun.f90 | ||
| − | |||
| − | |||
To submit the job to the serial queue use: | To submit the job to the serial queue use: | ||
sqsub -r 10m --mpp=2G -n 1 -o OUTPUT_MY_MPFUN ./a.out | sqsub -r 10m --mpp=2G -n 1 -o OUTPUT_MY_MPFUN ./a.out | ||
| Line 101: | Line 104: | ||
===OpenMP with MPFUN and QD=== | ===OpenMP with MPFUN and QD=== | ||
| − | |||
| − | |||
The Multi Precision packages MPFUN90 and QD use "user defined data type" for their variables. | The Multi Precision packages MPFUN90 and QD use "user defined data type" for their variables. | ||
Revision as of 14:24, 28 September 2015
| MPFUN90 |
|---|
| Description: A Fortran-90 arbitrary precision package |
| SHARCNET Package information: see MPFUN90 software page in web portal |
| Full list of SHARCNET supported software |
Contents
Introduction
MPFUN90 is a Fortran-90 arbitrary precision package (i.e. user specifies the desired precision). The library is written by David H. Bailey et. al from the Lawrence Berkeley National Laboratory. If 64 digits of precision are sufficient for an application then users should use QD instead of MPFUN which would run much faster.
Version Selection
- This page details how to use MPFUN90 (version 20100825) that is available as a SHARCNET module.
- NOTE: The only SHARCNET module currently available supports the Intel v15.0.3 Fortran compiler.
- There is an older MPFUN library version 6.1.23 for the intel compiler available on SHARCNET. See the MPFUN page for those details.
Job Submission
To compile an MPFUN90 job you need to have the intel/15.0.3 module loaded. After loading it, use the following command to compile an MPFUN90 source file, my_mpfun.f90:
ifort my_mpfun.f90
To submit the job to the serial queue use:
sqsub -r 10m --mpp=2G -n 1 -o OUTPUT_MY_MPFUN ./a.out
Testing the MPFUN90 functions
This example shows how to compile the test program testmp90.f90 using the SHARCNET mpfun build. This program is used to verify DHB's multiprecision computation package MPFUN. It exercises most routines and verifies that they are working properly. Note that MPFUN was compiled with the intel compiler on each cluster shown in the availability table therefore its necessary to also compile the below example with the appropriate compiler. The environment variable $FC handles this automatically as follows:
Following commands copy the testmp90.f90 file into a user subdirectory and run the test.
[roberpj@hnd50:~/samples/mpfun] cp /opt/sharcnet/mpfun/6.1.23/f90/src/testmp90.f90 .
[roberpj@hnd50:~/samples/mpfun] $FC testmp90.f90 -I/opt/sharcnet/mpfun/6.1.23/f90/include \
-L/opt/sharcnet/mpfun/6.1.23/f90/lib -lmp
[roberpj@hnd50:~/samples/mpfun] ./a.out
Completed Test 1
Completed Test 2
etc
Completed Test 60
Completed Test 62
Simple MPFUN example
Here is a simple MPFUN90 example and the commands required to compile and execute it:
! file = factorial_100.f90
program factorial_100_pgm
use mpmodule
implicit none
integer, parameter :: Nmx = 100
integer :: AllocateStatus
integer :: J
type(mp_real) zero, one
type(mp_real), dimension (:), allocatable :: fac
call mpinit(800)
call mpsetprec(800)
call mpsetoutputprec(700)
zero='0.'
one ='1.'
allocate(fac(0:Nmx),STAT=AllocateStatus)
IF (AllocateStatus /= 0) STOP "Not enough memory for fac "
fac(0)=one
do J=1,Nmx
fac(J)=fac(J-1)*J
end do
J=4
write(6,1001) J
1001 format(/"fac(",i3.3,") =",$)
call mpwrite(6,fac(J))
J=100
write(6,1001) J
call mpwrite(6,fac(J))
write(6,1005)
1005 format(/"DONE"//)
end
Compile (copy/paste the three lines)
$FC factorial_100.f90 \ -I/opt/sharcnet/mpfun/6.1.23/f90/include \ -L/opt/sharcnet/mpfun/6.1.23/f90/lib -lmp >& AOUT1
Execute code
./a.out fac(004) =10 ^ 1 x 2.4, fac(100) =10 ^ 157 x 9.3326215443944152681699238856266700490715968264381621468592 963895217599993229915608941463976156518286253697920827223758251185210916864, DONE
OpenMP with MPFUN and QD
The Multi Precision packages MPFUN90 and QD use "user defined data type" for their variables. Since the OpenMP reduction clause cannot be used for "user defined data types" to run MPFUN90 or QD with OpenMP the reduction must be done manually in a PARALLEL OpenMP region.
Detailed documentation on how to use OpenMP with MPFUN90 and QD are presented in the following document:
https://www.sharcnet.ca/help/index.php/OpenMP_Reduction_for_User_Defined_Data_Types
Handling large arrays in an MPFUN program
Declaring large arrays in an MPFUN program often leads to compile errors like:
"relocation truncated to fit r_x86_64_pc32 against .bss"
Even with the flag '-mcmodel large' in the compile command the error will persist, because there is a limit to the amount of data you can have allocated in static arrays for x86_64 architecture. This data segment is limited to 2GB. The data in COMMON plus any other statically declared arrays PLUS your code must fit in 2GB.
In order to overcome this problem use allocatable arrays declared in a module together with a short subroutine which allocates the arrays.
We illustrate these concepts by starting with a program 'betax.f90' which work fine for values of parameter 'Nsz' less than 37. Thus, if you compile it with
integer, parameter :: Nsz=36
it works, but for any higher values it fails with error message:
"relocation truncated to fit: R_X86_64_PC32 against symbol '...'"
To overcome this problem, we first copy program betax.f90 into a new file, say, betay.f90, and write a module called 'mp_data' and a subroutine called 'allocate_arrays' at the top of the new file. We move all the declarations of large arrays into the module, and declare them as 'allocatable'.
In the subroutine 'allocate_arrays' we do the allocations and in the main program we add the statements:
USE mp_data
...
call allocate_arrays
...
Below you will find the listings for the files 'betax.f90', 'betay.f90' and the 'HIST' file containing the commands to compile these files and submit jobs:
! file betax.f90
program BETAX
use mpmodule
implicit none
! integer, parameter :: Nsz=12 ! compilation fails for 37
integer, parameter :: Nsz=37
INTEGER :: I,J,K,L
type(mp_real) zero, one
type(mp_real) fac(0:Nsz),fci(0:Nsz)
type(mp_real) BETA2(Nsz,Nsz)
type(mp_real) BETA3(Nsz,Nsz,Nsz)
type(mp_real) BETA4(Nsz,Nsz,Nsz,Nsz)
call mpinit(800)
call mpsetprec(800)
call mpsetoutputprec(700)
zero='0.'
one ='1.'
fac(0)=one
fci(0)=one
do j=1,Nsz
fac(j)=fac(j-1)*j
fci(j)=one/fac(j)
end do
I=Nsz/2
J=I-1
K=J-1
L=K-1
write(6,1000) I,J,K,L
1000 format(/"I,J,K,L = ",4i3)
write(6,1001) I
1001 format(/"fac(",i3.3,") =",$)
call mpwrite(6,fac(I))
BETA2(I,J) =fac(I)*fac(J)*fci(I+J)
BETA3(I,J,K) =fac(I)*fac(J)*fac(K)*fci(I+J)*fci(J+K)
BETA4(I,J,K,L) =fac(I)*fac(J)*fac(K)*fac(L)* &
& fci(I+J)*fci(J+K)*fci(K+L)
write(6,1002) I,J
1002 format(/"BETA2(",i3.3,",",i3.3,") =")
call mpwrite(6,BETA2(I,J))
write(6,1003) I,J,K
1003 format(/"BETA3(",i3.3,",",i3.3,",",i3.3,") =")
call mpwrite(6,BETA3(I,J,K))
write(6,1004) I,J,K,L
1004 format(/"BETA4(",i3.3,",",i3.3,",",i3.3,",",i3.3,") =")
call mpwrite(6,BETA4(I,J,K,L))
write(6,1005)
1005 format(/"DONE"//)
stop
end
! file betay.f90
module mp_data
use mpmodule
use mpfunmod
implicit none
integer, parameter :: Nsz =48
integer :: AllocateStatus
! type(mp_real) zero, one
!
! type(mp_real) fac(0:Nsz),fci(0:Nsz)
! type(mp_real) BETA2(Nsz,Nsz)
! type(mp_real) BETA3(Nsz,Nsz,Nsz)
! type(mp_real) BETA4(Nsz,Nsz,Nsz,Nsz)
type(mp_real) zero, one
type(mp_real), dimension (:), allocatable :: fac,fci
type(mp_real), dimension (:,:), allocatable :: BETA2
type(mp_real), dimension (:,:,:), allocatable :: BETA3
type(mp_real), dimension (:,:,:,:), allocatable :: BETA4
end module mp_data
subroutine allocate_arrays
USE mp_data
implicit none
allocate(fac(0:Nsz),STAT=AllocateStatus)
IF (AllocateStatus /= 0) STOP "Not enough memory for fac "
allocate(fci(0:Nsz),STAT=AllocateStatus)
IF (AllocateStatus /= 0) STOP "Not enough memory for fci "
allocate(BETA2(Nsz,Nsz),STAT=AllocateStatus)
IF (AllocateStatus /= 0) STOP "Not enough memory for BETA2"
allocate(BETA3(Nsz,Nsz,Nsz),STAT=AllocateStatus)
IF (AllocateStatus /= 0) STOP "Not enough memory for BETA3"
allocate(BETA4(Nsz,Nsz,Nsz,Nsz),STAT=AllocateStatus)
IF (AllocateStatus /= 0) STOP "Not enough memory for BETA4"
print *,"Done with allocations"
return
END
! ----------------------------------------------------------------------
program BETAY
USE mp_data
implicit none
INTEGER :: I,J,K,L
call mpinit(800)
call mpsetprec(800)
call mpsetoutputprec(700)
zero='0.'
one ='1.'
call allocate_arrays
fac(0)=one
fci(0)=one
do j=1,Nsz
fac(j)=fac(j-1)*j
fci(j)=one/fac(j)
end do
I=Nsz/2
J=I-1
K=J-1
L=K-1
write(6,1000) I,J,K,L
1000 format(/"I,J,K,L = ",4i3)
write(6,1001) I
1001 format(/"fac(",i3.3,") =",$)
call mpwrite(6,fac(I))
BETA2(I,J) =fac(I)*fac(J)*fci(I+J)
BETA3(I,J,K) =fac(I)*fac(J)*fac(K)*fci(I+J)*fci(J+K)
BETA4(I,J,K,L) =fac(I)*fac(J)*fac(K)*fac(L)* &
& fci(I+J)*fci(J+K)*fci(K+L)
write(6,1002) I,J
1002 format(/"BETA2(",i3.3,",",i3.3,") =")
call mpwrite(6,BETA2(I,J))
write(6,1003) I,J,K
1003 format(/"BETA3(",i3.3,",",i3.3,",",i3.3,") =")
call mpwrite(6,BETA3(I,J,K))
write(6,1004) I,J,K,L
1004 format(/"BETA4(",i3.3,",",i3.3,",",i3.3,",",i3.3,") =")
call mpwrite(6,BETA4(I,J,K,L))
write(6,1005)
1005 format(/"DONE"//)
stop
end
# file HIST
#
# Set Nsz=12
vi betax.f90
and set Nsz=12 in file betax.f90
# Compile (copy/paste the three lines)
ifort betax.f90 \
-I/opt/sharcnet/mpfun/6.1.23/f90/include \
-L/opt/sharcnet/mpfun/6.1.23/f90/lib -lmp >& AOUT1
# Execute code (it should work)
./a.out
# Set Nsz=40
vi betax.f90
and set Nsz=40 in file betax.f90
# Re-compile and look at file AOUT1
ifort betax.f90 \
-I/opt/sharcnet/mpfun/6.1.23/f90/include \
-L/opt/sharcnet/mpfun/6.1.23/f90/lib -lmp >& AOUT1
============================================================================
# File betay.f90 was copied from betax.f90 and modified as described earlier.
# Compile betay.f90
ifort betay.f90 -mcmodel medium -shared-intel \
-I/opt/sharcnet/mpfun/6.1.23/f90/include \
-L/opt/sharcnet/mpfun/6.1.23/f90/lib -lmp >& AOUT2
============================================================================
# The executable might require more than 2.0G of memory:
# This will pruduce an error:
sqsub -r 5m --mpp=4.0G -o OUTPUTFILE4 ./a.out
# This should work:
sqsub -r 5m --mpp=8.0G -o OUTPUTFILE8 ./a.out
============================================================================
MPFUN's TOOLKIT
This toolkit provides and extensive example of how the mpfun package can be used. Read the header of files
mathinit.f90
and
mathtool.f90
found under
/opt/sharcnet/mpfun/current/toolkit
for an explanation what the package does. The Makefile in this directory shows how the mpfun libraries were linked to compile this sample program pair. The following steps shows the first initialization step followed by actually running the mathtool program to perform some examples such as performing an integration and displaying the value of pie. Note that the displayed zero cputimes are due to a bug in the mpfun current build which will be fixed in future updates on the clusters.
Initializing Mathtool
To run Mathtool you must first generate following two input files
const.dat and quadts.dat
which can be accomplished by running following command:
[roberpj@bal34:~/samples/mpfun] /opt/sharcnet/mpfun/current/toolkit/mathinit mathinit: start const complete file const.dat written cpu time = 0.E+0 initqts: Tanh-sinh quadrature initialization 0 18632 1000 18632 2000 18632 3000 18632 4000 18632 5000 18632 6000 18632 7000 18632 8000 18632 initqts: Table spaced used = 8177 initqts complete file quadts.dat written cpu time = 0.E+0 total cpu time = 0.E+0 [roberpj@bal34:~/samples/mpfun] ls const.dat quadts.dat
Running Mathtool
Use the files created in the previous section to run mathtool as follows:
[roberpj@bal34:~/samples/mpfun] /opt/sharcnet/mpfun/current/toolkit/mathtool
Welcome to the Experimental Mathematician's Toolkit
Initializing...
Current settings:
Debug level = 2
Primary precision level = 100 digits
Secondary precision level = 200 digits
Primary epsilon = 10^ -100
Secondary epsilon = 10^ -200
PSLQ bound = 100
PSLQ level = 1
Quadrature level = 6
Quadrature type = 3 (Tanh-Sinh)
-- snip --
TYPE: integrate[1/(1+x^2), {x, 0, 1}] (hit enter)
quadts: Iteration 1 of 6; est error = 10^ 0; approx value =
7.853157298876894309110693190063938396203990510665917998384243e-1
quadts: Iteration 2 of 6; est error = 10^ 0; approx value =
7.853981605125528448099840689980111720633592558858121242731082e-1
quadts: Iteration 3 of 6; est error = 10^ -17; approx value =
7.853981633974483082028533923930157568609007576776270726821597e-1
quadts: Iteration 4 of 6; est error = 10^ -36; approx value =
7.853981633974483096156608458198757190586422063114153956005278e-1
quadts: Iteration 5 of 6; est error = 10^ -71; approx value =
7.853981633974483096156608458198757210492923498437764552437361e-1
quadts: Iteration 6 of 6; est error = 10^ -101; approx value =
7.853981633974483096156608458198757210492923498437764552437361e-1
Result[ 1] =
7.85398163397448309615660845819875721049292349843776455243736148076954e-1
CPU time = 0.0000
TYPE: functions (hit enter)
Abs Arccos Arcsin Arctan
Arctan2 Bessel Besselexp Binomial
Cos Erf Exp Factorial
Gamma Integrate Log Max
Min Polyroot Pslq Result
Sin Sqrt Sum Table
Tan Zeta Zetap Zetaz
TYPE: variables (hit enter)
E Log2 Log10 Pi
Catalan Eulergamma Infinity Arg1
Arg2 Arg3 Arg4 Arg5
Arg6 Arg7 Arg8 Arg9
Debug Digits Digits2 Epsilon
Epsilon2 Pslqbound Pslqlevel Quadlevel
Quadtype x
TYPE: Pi (hit enter)
Result[ 2] =
3.14159265358979323846264338327950288419716939937510582097494459230781e0
CPU time = 0.0000
More Mathtool commands
To start mathtool (assuming you have issued the command /opt/sharcnet/mpfun/current/toolkit/mathinit and have these two files: const.dat HIST quadts.dat - you type:
/opt/sharcnet/mpfun/current/toolkit/mathtool
Then try following commands:
integrate[2*x,{x,0,1}]
sqrt[1600]
help Binomial
Binomial[3,2]
hypo[x,y] = sqrt[x^2+y^2]
z=hypo[3,4]
z
exit
This is what you should get from above commands:
integrate[2*x,{x,0,1}]
quadts: Iteration 1 of 6; est error = 10^ 0; approx value =
1.000003359570811224742398159504692881030640064061042259417294e0
quadts: Iteration 2 of 6; est error = 10^ 0; approx value =
1.000000000000036584185557611216071731517571134591404403031033e0
quadts: Iteration 3 of 6; est error = 10^ -27; approx value =
1.000000000000000000000000000001046972220803162709182580101627e0
quadts: Iteration 4 of 6; est error = 10^ -60; approx value =
1.000000000000000000000000000000000000000000000000000000000000e0
quadts: Iteration 5 of 6; est error = 10^ -101; approx value =
9.999999999999999999999999999999999999999999999999999999999999e-1
Result[ 1] =
9.99999999999999999999999999999999999999999999999999999999999999999999e-1
CPU time = 0.0000
sqrt[1600]
Result[ 2] =
4.00000000000000000000000000000000000000000000000000000000000000000000e1
CPU time = 0.0000
help Binomial Binomial[m,n] computes the binomial coefficient of integers (m,n).
Binomial[3,2]
Result[ 3] =
3.00000000000000000000000000000000000000000000000000000000000000000000e0
CPU time = 0.0000
hypo[x,y] = sqrt[x^2+y^2] parse: new function name = hypo CPU time = 0.0000
z=hypo[3,4] parse: new variable name = z CPU time = 0.0000
z
Result[ 4] =
5.00000000000000000000000000000000000000000000000000000000000000000000e0
CPU time = 0.0000
exit
References
o High-Precision Software Directory
http://crd.lbl.gov/~dhbailey/mpdist/index.html