mpblas.f 2.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778
  1. *DECK MPBLAS
  2. SUBROUTINE MPBLAS (I1)
  3. C***BEGIN PROLOGUE MPBLAS
  4. C***SUBSIDIARY
  5. C***PURPOSE Subsidiary to DQDOTA and DQDOTI
  6. C***LIBRARY SLATEC
  7. C***TYPE ALL (MPBLAS-A)
  8. C***AUTHOR (UNKNOWN)
  9. C***DESCRIPTION
  10. C
  11. C This subroutine is called to set up Brent's 'mp' package
  12. C for use by the extended precision inner products from the BLAS.
  13. C
  14. C In the SLATEC library we require the Extended Precision MP number
  15. C to have a mantissa twice as long as Double Precision numbers.
  16. C The calculation of MPT (and MPMXR which is the actual array size)
  17. C in this routine will give 2x (or slightly more) on the machine
  18. C that we are running on. The INTEGER array size of 30 was chosen
  19. C to be slightly longer than the longest INTEGER array needed on
  20. C any machine that we are currently aware of.
  21. C
  22. C***SEE ALSO DQDOTA, DQDOTI
  23. C***REFERENCES R. P. Brent, A Fortran multiple-precision arithmetic
  24. C package, ACM Transactions on Mathematical Software 4,
  25. C 1 (March 1978), pp. 57-70.
  26. C R. P. Brent, MP, a Fortran multiple-precision arithmetic
  27. C package, Algorithm 524, ACM Transactions on Mathema-
  28. C tical Software 4, 1 (March 1978), pp. 71-81.
  29. C***ROUTINES CALLED I1MACH, XERMSG
  30. C***COMMON BLOCKS MPCOM
  31. C***REVISION HISTORY (YYMMDD)
  32. C 791001 DATE WRITTEN
  33. C 890531 Changed all specific intrinsics to generic. (WRB)
  34. C 890531 REVISION DATE from Version 3.2
  35. C 891214 Prologue converted to Version 4.0 format. (BAB)
  36. C 900402 Added TYPE section. (WRB)
  37. C 920501 Reformatted the REFERENCES section. (WRB)
  38. C 930124 Increased Array size in MPCON for SUN -r8, and calculate
  39. C size for Quad Precision for 2x DP. (RWC)
  40. C***END PROLOGUE MPBLAS
  41. COMMON /MPCOM/ MPB, MPT, MPM, MPLUN, MPMXR, MPR(30)
  42. C***FIRST EXECUTABLE STATEMENT MPBLAS
  43. I1 = 1
  44. C
  45. C For full extended precision accuracy, MPB should be as large as
  46. C possible, subject to the restrictions in Brent's paper.
  47. C
  48. C Statements below are for an integer wordlength of 48, 36, 32,
  49. C 24, 18, and 16. Pick one, or generate a new one.
  50. C 48 MPB = 4194304
  51. C 36 MPB = 65536
  52. C 32 MPB = 16384
  53. C 24 MPB = 1024
  54. C 18 MPB = 128
  55. C 16 MPB = 64
  56. C
  57. MPBEXP = I1MACH(8)/2-2
  58. MPB = 2**MPBEXP
  59. C
  60. C Set up remaining parameters
  61. C UNIT FOR ERROR MESSAGES
  62. MPLUN = I1MACH(4)
  63. C NUMBER OF MP DIGITS
  64. MPT = (2*I1MACH(14)+MPBEXP-1)/MPBEXP
  65. C DIMENSION OF R
  66. MPMXR = MPT+4
  67. C
  68. if (MPMXR.GT.30) THEN
  69. CALL XERMSG('SLATEC', 'MPBLAS',
  70. * 'Array space not sufficient for Quad Precision 2x ' //
  71. * 'Double Precision, Proceeding.', 1, 1)
  72. MPT = 26
  73. MPMXR = 30
  74. ENDIF
  75. C EXPONENT RANGE
  76. MPM = MIN(32767,I1MACH(9)/4-1)
  77. RETURN
  78. END