c@427
|
1 /* ieeeck.f -- translated by f2c (version 20061008).
|
c@427
|
2 You must link the resulting object file with libf2c:
|
c@427
|
3 on Microsoft Windows system, link with libf2c.lib;
|
c@427
|
4 on Linux or Unix systems, link with .../path/to/libf2c.a -lm
|
c@427
|
5 or, if you install libf2c.a in a standard place, with -lf2c -lm
|
c@427
|
6 -- in that order, at the end of the command line, as in
|
c@427
|
7 cc *.o -lf2c -lm
|
c@427
|
8 Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
|
c@427
|
9
|
c@427
|
10 http://www.netlib.org/f2c/libf2c.zip
|
c@427
|
11 */
|
c@427
|
12
|
c@427
|
13 #include "f2c.h"
|
c@427
|
14 #include "blaswrap.h"
|
c@427
|
15
|
c@427
|
16 integer ieeeck_(integer *ispec, real *zero, real *one)
|
c@427
|
17 {
|
c@427
|
18 /* System generated locals */
|
c@427
|
19 integer ret_val;
|
c@427
|
20
|
c@427
|
21 /* Local variables */
|
c@427
|
22 real nan1, nan2, nan3, nan4, nan5, nan6, neginf, posinf, negzro, newzro;
|
c@427
|
23
|
c@427
|
24
|
c@427
|
25 /* -- LAPACK auxiliary routine (version 3.2) -- */
|
c@427
|
26 /* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
|
c@427
|
27 /* November 2006 */
|
c@427
|
28
|
c@427
|
29 /* .. Scalar Arguments .. */
|
c@427
|
30 /* .. */
|
c@427
|
31
|
c@427
|
32 /* Purpose */
|
c@427
|
33 /* ======= */
|
c@427
|
34
|
c@427
|
35 /* IEEECK is called from the ILAENV to verify that Infinity and */
|
c@427
|
36 /* possibly NaN arithmetic is safe (i.e. will not trap). */
|
c@427
|
37
|
c@427
|
38 /* Arguments */
|
c@427
|
39 /* ========= */
|
c@427
|
40
|
c@427
|
41 /* ISPEC (input) INTEGER */
|
c@427
|
42 /* Specifies whether to test just for inifinity arithmetic */
|
c@427
|
43 /* or whether to test for infinity and NaN arithmetic. */
|
c@427
|
44 /* = 0: Verify infinity arithmetic only. */
|
c@427
|
45 /* = 1: Verify infinity and NaN arithmetic. */
|
c@427
|
46
|
c@427
|
47 /* ZERO (input) REAL */
|
c@427
|
48 /* Must contain the value 0.0 */
|
c@427
|
49 /* This is passed to prevent the compiler from optimizing */
|
c@427
|
50 /* away this code. */
|
c@427
|
51
|
c@427
|
52 /* ONE (input) REAL */
|
c@427
|
53 /* Must contain the value 1.0 */
|
c@427
|
54 /* This is passed to prevent the compiler from optimizing */
|
c@427
|
55 /* away this code. */
|
c@427
|
56
|
c@427
|
57 /* RETURN VALUE: INTEGER */
|
c@427
|
58 /* = 0: Arithmetic failed to produce the correct answers */
|
c@427
|
59 /* = 1: Arithmetic produced the correct answers */
|
c@427
|
60
|
c@427
|
61 /* .. Local Scalars .. */
|
c@427
|
62 /* .. */
|
c@427
|
63 /* .. Executable Statements .. */
|
c@427
|
64 ret_val = 1;
|
c@427
|
65
|
c@427
|
66 posinf = *one / *zero;
|
c@427
|
67 if (posinf <= *one) {
|
c@427
|
68 ret_val = 0;
|
c@427
|
69 return ret_val;
|
c@427
|
70 }
|
c@427
|
71
|
c@427
|
72 neginf = -(*one) / *zero;
|
c@427
|
73 if (neginf >= *zero) {
|
c@427
|
74 ret_val = 0;
|
c@427
|
75 return ret_val;
|
c@427
|
76 }
|
c@427
|
77
|
c@427
|
78 negzro = *one / (neginf + *one);
|
c@427
|
79 if (negzro != *zero) {
|
c@427
|
80 ret_val = 0;
|
c@427
|
81 return ret_val;
|
c@427
|
82 }
|
c@427
|
83
|
c@427
|
84 neginf = *one / negzro;
|
c@427
|
85 if (neginf >= *zero) {
|
c@427
|
86 ret_val = 0;
|
c@427
|
87 return ret_val;
|
c@427
|
88 }
|
c@427
|
89
|
c@427
|
90 newzro = negzro + *zero;
|
c@427
|
91 if (newzro != *zero) {
|
c@427
|
92 ret_val = 0;
|
c@427
|
93 return ret_val;
|
c@427
|
94 }
|
c@427
|
95
|
c@427
|
96 posinf = *one / newzro;
|
c@427
|
97 if (posinf <= *one) {
|
c@427
|
98 ret_val = 0;
|
c@427
|
99 return ret_val;
|
c@427
|
100 }
|
c@427
|
101
|
c@427
|
102 neginf *= posinf;
|
c@427
|
103 if (neginf >= *zero) {
|
c@427
|
104 ret_val = 0;
|
c@427
|
105 return ret_val;
|
c@427
|
106 }
|
c@427
|
107
|
c@427
|
108 posinf *= posinf;
|
c@427
|
109 if (posinf <= *one) {
|
c@427
|
110 ret_val = 0;
|
c@427
|
111 return ret_val;
|
c@427
|
112 }
|
c@427
|
113
|
c@427
|
114
|
c@427
|
115
|
c@427
|
116
|
c@427
|
117 /* Return if we were only asked to check infinity arithmetic */
|
c@427
|
118
|
c@427
|
119 if (*ispec == 0) {
|
c@427
|
120 return ret_val;
|
c@427
|
121 }
|
c@427
|
122
|
c@427
|
123 nan1 = posinf + neginf;
|
c@427
|
124
|
c@427
|
125 nan2 = posinf / neginf;
|
c@427
|
126
|
c@427
|
127 nan3 = posinf / posinf;
|
c@427
|
128
|
c@427
|
129 nan4 = posinf * *zero;
|
c@427
|
130
|
c@427
|
131 nan5 = neginf * negzro;
|
c@427
|
132
|
c@427
|
133 nan6 = nan5 * 0.f;
|
c@427
|
134
|
c@427
|
135 if (nan1 == nan1) {
|
c@427
|
136 ret_val = 0;
|
c@427
|
137 return ret_val;
|
c@427
|
138 }
|
c@427
|
139
|
c@427
|
140 if (nan2 == nan2) {
|
c@427
|
141 ret_val = 0;
|
c@427
|
142 return ret_val;
|
c@427
|
143 }
|
c@427
|
144
|
c@427
|
145 if (nan3 == nan3) {
|
c@427
|
146 ret_val = 0;
|
c@427
|
147 return ret_val;
|
c@427
|
148 }
|
c@427
|
149
|
c@427
|
150 if (nan4 == nan4) {
|
c@427
|
151 ret_val = 0;
|
c@427
|
152 return ret_val;
|
c@427
|
153 }
|
c@427
|
154
|
c@427
|
155 if (nan5 == nan5) {
|
c@427
|
156 ret_val = 0;
|
c@427
|
157 return ret_val;
|
c@427
|
158 }
|
c@427
|
159
|
c@427
|
160 if (nan6 == nan6) {
|
c@427
|
161 ret_val = 0;
|
c@427
|
162 return ret_val;
|
c@427
|
163 }
|
c@427
|
164
|
c@427
|
165 return ret_val;
|
c@427
|
166 } /* ieeeck_ */
|