~ubuntu-branches/ubuntu/karmic/python-scipy/karmic

« back to all changes in this revision

Viewing changes to Lib/integrate/mach/d1mach.f

  • Committer: Bazaar Package Importer
  • Author(s): Daniel T. Chen (new)
  • Date: 2005-03-16 02:15:29 UTC
  • Revision ID: james.westby@ubuntu.com-20050316021529-xrjlowsejs0cijig
Tags: upstream-0.3.2
ImportĀ upstreamĀ versionĀ 0.3.2

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
      DOUBLE PRECISION FUNCTION D1MACH(I)
 
2
      INTEGER I
 
3
C
 
4
C  DOUBLE-PRECISION MACHINE CONSTANTS
 
5
C  D1MACH( 1) = B**(EMIN-1), THE SMALLEST POSITIVE MAGNITUDE.
 
6
C  D1MACH( 2) = B**EMAX*(1 - B**(-T)), THE LARGEST MAGNITUDE.
 
7
C  D1MACH( 3) = B**(-T), THE SMALLEST RELATIVE SPACING.
 
8
C  D1MACH( 4) = B**(1-T), THE LARGEST RELATIVE SPACING.
 
9
C  D1MACH( 5) = LOG10(B)
 
10
C
 
11
      INTEGER SMALL(2)
 
12
      INTEGER LARGE(2)
 
13
      INTEGER RIGHT(2)
 
14
      INTEGER DIVER(2)
 
15
      INTEGER LOG10(2)
 
16
      INTEGER SC, CRAY1(38), J
 
17
      COMMON /D9MACH/ CRAY1
 
18
      SAVE SMALL, LARGE, RIGHT, DIVER, LOG10, SC
 
19
      DOUBLE PRECISION DMACH(5)
 
20
      EQUIVALENCE (DMACH(1),SMALL(1))
 
21
      EQUIVALENCE (DMACH(2),LARGE(1))
 
22
      EQUIVALENCE (DMACH(3),RIGHT(1))
 
23
      EQUIVALENCE (DMACH(4),DIVER(1))
 
24
      EQUIVALENCE (DMACH(5),LOG10(1))
 
25
C  THIS VERSION ADAPTS AUTOMATICALLY TO MOST CURRENT MACHINES.
 
26
C  R1MACH CAN HANDLE AUTO-DOUBLE COMPILING, BUT THIS VERSION OF
 
27
C  D1MACH DOES NOT, BECAUSE WE DO NOT HAVE QUAD CONSTANTS FOR
 
28
C  MANY MACHINES YET.
 
29
C  TO COMPILE ON OLDER MACHINES, ADD A C IN COLUMN 1
 
30
C  ON THE NEXT LINE
 
31
      DATA SC/0/
 
32
C  AND REMOVE THE C FROM COLUMN 1 IN ONE OF THE SECTIONS BELOW.
 
33
C  CONSTANTS FOR EVEN OLDER MACHINES CAN BE OBTAINED BY
 
34
C          mail netlib@research.bell-labs.com
 
35
C          send old1mach from blas
 
36
C  PLEASE SEND CORRECTIONS TO dmg OR ehg@bell-labs.com.
 
37
C
 
38
C     MACHINE CONSTANTS FOR THE HONEYWELL DPS 8/70 SERIES.
 
39
C      DATA SMALL(1),SMALL(2) / O402400000000, O000000000000 /
 
40
C      DATA LARGE(1),LARGE(2) / O376777777777, O777777777777 /
 
41
C      DATA RIGHT(1),RIGHT(2) / O604400000000, O000000000000 /
 
42
C      DATA DIVER(1),DIVER(2) / O606400000000, O000000000000 /
 
43
C      DATA LOG10(1),LOG10(2) / O776464202324, O117571775714 /, SC/987/
 
44
C
 
45
C     MACHINE CONSTANTS FOR PDP-11 FORTRANS SUPPORTING
 
46
C     32-BIT INTEGERS.
 
47
C      DATA SMALL(1),SMALL(2) /    8388608,           0 /
 
48
C      DATA LARGE(1),LARGE(2) / 2147483647,          -1 /
 
49
C      DATA RIGHT(1),RIGHT(2) /  612368384,           0 /
 
50
C      DATA DIVER(1),DIVER(2) /  620756992,           0 /
 
51
C      DATA LOG10(1),LOG10(2) / 1067065498, -2063872008 /, SC/987/
 
52
C
 
53
C     MACHINE CONSTANTS FOR THE UNIVAC 1100 SERIES.
 
54
C      DATA SMALL(1),SMALL(2) / O000040000000, O000000000000 /
 
55
C      DATA LARGE(1),LARGE(2) / O377777777777, O777777777777 /
 
56
C      DATA RIGHT(1),RIGHT(2) / O170540000000, O000000000000 /
 
57
C      DATA DIVER(1),DIVER(2) / O170640000000, O000000000000 /
 
58
C      DATA LOG10(1),LOG10(2) / O177746420232, O411757177572 /, SC/987/
 
59
C
 
60
C     ON FIRST CALL, IF NO DATA UNCOMMENTED, TEST MACHINE TYPES.
 
61
      IF (SC .NE. 987) THEN
 
62
         DMACH(1) = 1.D13
 
63
         IF (      SMALL(1) .EQ. 1117925532
 
64
     *       .AND. SMALL(2) .EQ. -448790528) THEN
 
65
*           *** IEEE BIG ENDIAN ***
 
66
            SMALL(1) = 1048576
 
67
            SMALL(2) = 0
 
68
            LARGE(1) = 2146435071
 
69
            LARGE(2) = -1
 
70
            RIGHT(1) = 1017118720
 
71
            RIGHT(2) = 0
 
72
            DIVER(1) = 1018167296
 
73
            DIVER(2) = 0
 
74
            LOG10(1) = 1070810131
 
75
            LOG10(2) = 1352628735
 
76
         ELSE IF ( SMALL(2) .EQ. 1117925532
 
77
     *       .AND. SMALL(1) .EQ. -448790528) THEN
 
78
*           *** IEEE LITTLE ENDIAN ***
 
79
            SMALL(2) = 1048576
 
80
            SMALL(1) = 0
 
81
            LARGE(2) = 2146435071
 
82
            LARGE(1) = -1
 
83
            RIGHT(2) = 1017118720
 
84
            RIGHT(1) = 0
 
85
            DIVER(2) = 1018167296
 
86
            DIVER(1) = 0
 
87
            LOG10(2) = 1070810131
 
88
            LOG10(1) = 1352628735
 
89
         ELSE IF ( SMALL(1) .EQ. -2065213935
 
90
     *       .AND. SMALL(2) .EQ. 10752) THEN
 
91
*               *** VAX WITH D_FLOATING ***
 
92
            SMALL(1) = 128
 
93
            SMALL(2) = 0
 
94
            LARGE(1) = -32769
 
95
            LARGE(2) = -1
 
96
            RIGHT(1) = 9344
 
97
            RIGHT(2) = 0
 
98
            DIVER(1) = 9472
 
99
            DIVER(2) = 0
 
100
            LOG10(1) = 546979738
 
101
            LOG10(2) = -805796613
 
102
         ELSE IF ( SMALL(1) .EQ. 1267827943
 
103
     *       .AND. SMALL(2) .EQ. 704643072) THEN
 
104
*               *** IBM MAINFRAME ***
 
105
            SMALL(1) = 1048576
 
106
            SMALL(2) = 0
 
107
            LARGE(1) = 2147483647
 
108
            LARGE(2) = -1
 
109
            RIGHT(1) = 856686592
 
110
            RIGHT(2) = 0
 
111
            DIVER(1) = 873463808
 
112
            DIVER(2) = 0
 
113
            LOG10(1) = 1091781651
 
114
            LOG10(2) = 1352628735
 
115
         ELSE IF ( SMALL(1) .EQ. 1120022684
 
116
     *       .AND. SMALL(2) .EQ. -448790528) THEN
 
117
*           *** CONVEX C-1 ***
 
118
            SMALL(1) = 1048576
 
119
            SMALL(2) = 0
 
120
            LARGE(1) = 2147483647
 
121
            LARGE(2) = -1
 
122
            RIGHT(1) = 1019215872
 
123
            RIGHT(2) = 0
 
124
            DIVER(1) = 1020264448
 
125
            DIVER(2) = 0
 
126
            LOG10(1) = 1072907283
 
127
            LOG10(2) = 1352628735
 
128
         ELSE IF ( SMALL(1) .EQ. 815547074
 
129
     *       .AND. SMALL(2) .EQ. 58688) THEN
 
130
*           *** VAX G-FLOATING ***
 
131
            SMALL(1) = 16
 
132
            SMALL(2) = 0
 
133
            LARGE(1) = -32769
 
134
            LARGE(2) = -1
 
135
            RIGHT(1) = 15552
 
136
            RIGHT(2) = 0
 
137
            DIVER(1) = 15568
 
138
            DIVER(2) = 0
 
139
            LOG10(1) = 1142112243
 
140
            LOG10(2) = 2046775455
 
141
         ELSE
 
142
            DMACH(2) = 1.D27 + 1
 
143
            DMACH(3) = 1.D27
 
144
            LARGE(2) = LARGE(2) - RIGHT(2)
 
145
            IF (LARGE(2) .EQ. 64 .AND. SMALL(2) .EQ. 0) THEN
 
146
               CRAY1(1) = 67291416
 
147
               DO 10 J = 1, 20
 
148
                  CRAY1(J+1) = CRAY1(J) + CRAY1(J)
 
149
 10               CONTINUE
 
150
               CRAY1(22) = CRAY1(21) + 321322
 
151
               DO 20 J = 22, 37
 
152
                  CRAY1(J+1) = CRAY1(J) + CRAY1(J)
 
153
 20               CONTINUE
 
154
               IF (CRAY1(38) .EQ. SMALL(1)) THEN
 
155
*                  *** CRAY ***
 
156
                  CALL I1MCRY(SMALL(1), J, 8285, 8388608, 0)
 
157
                  SMALL(2) = 0
 
158
                  CALL I1MCRY(LARGE(1), J, 24574, 16777215, 16777215)
 
159
                  CALL I1MCRY(LARGE(2), J, 0, 16777215, 16777214)
 
160
                  CALL I1MCRY(RIGHT(1), J, 16291, 8388608, 0)
 
161
                  RIGHT(2) = 0
 
162
                  CALL I1MCRY(DIVER(1), J, 16292, 8388608, 0)
 
163
                  DIVER(2) = 0
 
164
                  CALL I1MCRY(LOG10(1), J, 16383, 10100890, 8715215)
 
165
                  CALL I1MCRY(LOG10(2), J, 0, 16226447, 9001388)
 
166
               ELSE
 
167
                  WRITE(*,9000)
 
168
                  STOP 779
 
169
                  END IF
 
170
            ELSE
 
171
               WRITE(*,9000)
 
172
               STOP 779
 
173
               END IF
 
174
            END IF
 
175
         SC = 987
 
176
         END IF
 
177
*    SANITY CHECK
 
178
      IF (DMACH(4) .GE. 1.0D0) STOP 778
 
179
      IF (I .LT. 1 .OR. I .GT. 5) THEN
 
180
         WRITE(*,*) 'D1MACH(I): I =',I,' is out of bounds.'
 
181
         STOP
 
182
         END IF
 
183
      D1MACH = DMACH(I)
 
184
      RETURN
 
185
 9000 FORMAT(/' Adjust D1MACH by uncommenting data statements'/
 
186
     *' appropriate for your machine.')
 
187
* /* Standard C source for D1MACH -- remove the * in column 1 */
 
188
*#include <stdio.h>
 
189
*#include <float.h>
 
190
*#include <math.h>
 
191
*double d1mach_(long *i)
 
192
*{
 
193
*       switch(*i){
 
194
*         case 1: return DBL_MIN;
 
195
*         case 2: return DBL_MAX;
 
196
*         case 3: return DBL_EPSILON/FLT_RADIX;
 
197
*         case 4: return DBL_EPSILON;
 
198
*         case 5: return log10(FLT_RADIX);
 
199
*         }
 
200
*       fprintf(stderr, "invalid argument: d1mach(%ld)\n", *i);
 
201
*       exit(1); return 0; /* some compilers demand return values */
 
202
*}
 
203
      END
 
204
      SUBROUTINE I1MCRY(A, A1, B, C, D)
 
205
**** SPECIAL COMPUTATION FOR OLD CRAY MACHINES ****
 
206
      INTEGER A, A1, B, C, D
 
207
      A1 = 16777216*B + C
 
208
      A = 16777216*A1 + D
 
209
      END