Mercurial > hg > octave-nkf
annotate liboctave/lo-ieee.cc @ 9825:7483fe200fab
narrow complex values with negative zero imaginary parts
| author | Jaroslav Hajek <highegg@gmail.com> |
|---|---|
| date | Fri, 13 Nov 2009 12:34:46 +0100 |
| parents | bb70d16cca3b |
| children | 4c0cdbe0acca |
| rev | line source |
|---|---|
| 1967 | 1 /* |
| 2 | |
| 7017 | 3 Copyright (C) 1996, 1997, 1998, 1999, 2000, 2002, 2003, 2004, 2005, |
|
9592
5828d64ca004
lo-ieee.cc: include cstdlib, update copyright date
John W. Eaton <jwe@octave.org>
parents:
9591
diff
changeset
|
4 2007, 2008, 2009 John W. Eaton |
| 1967 | 5 |
| 6 This file is part of Octave. | |
| 7 | |
| 8 Octave is free software; you can redistribute it and/or modify it | |
| 9 under the terms of the GNU General Public License as published by the | |
| 7016 | 10 Free Software Foundation; either version 3 of the License, or (at your |
| 11 option) any later version. | |
| 1967 | 12 |
| 13 Octave is distributed in the hope that it will be useful, but WITHOUT | |
| 14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or | |
| 15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License | |
| 16 for more details. | |
| 17 | |
| 18 You should have received a copy of the GNU General Public License | |
| 7016 | 19 along with Octave; see the file COPYING. If not, see |
| 20 <http://www.gnu.org/licenses/>. | |
| 1967 | 21 |
| 22 */ | |
| 23 | |
| 24 #ifdef HAVE_CONFIG_H | |
| 25 #include <config.h> | |
| 26 #endif | |
| 27 | |
| 28 #include <cfloat> | |
|
9592
5828d64ca004
lo-ieee.cc: include cstdlib, update copyright date
John W. Eaton <jwe@octave.org>
parents:
9591
diff
changeset
|
29 #include <cstdlib> |
| 1967 | 30 |
| 31 #if defined (HAVE_FLOATINGPOINT_H) | |
| 32 #include <floatingpoint.h> | |
| 33 #endif | |
| 34 | |
| 3268 | 35 #if defined (HAVE_IEEEFP_H) |
| 36 #include <ieeefp.h> | |
| 37 #endif | |
| 38 | |
| 2598 | 39 #if defined (HAVE_NAN_H) |
| 40 #if defined (SCO) | |
| 2560 | 41 #define _IEEE 1 |
| 2598 | 42 #endif |
| 2508 | 43 #include <nan.h> |
| 2598 | 44 #if defined (SCO) |
| 2560 | 45 #undef _IEEE |
| 2508 | 46 #endif |
| 2598 | 47 #endif |
| 2508 | 48 |
| 4698 | 49 #include "lo-error.h" |
| 1967 | 50 #include "lo-ieee.h" |
| 7231 | 51 #include "lo-math.h" |
| 4025 | 52 #include "mach-info.h" |
| 1967 | 53 |
| 54 void | |
| 55 octave_ieee_init (void) | |
| 56 { | |
| 4698 | 57 // Default values. DBL_MAX is not right for NaN and NA, but do you |
| 58 // have a better suggestion? If you don't have IEEE floating point | |
| 59 // values, there are many parts of Octave that will not work | |
| 60 // correctly. | |
| 61 | |
| 62 octave_Inf = octave_NaN = octave_NA = DBL_MAX; | |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
63 octave_Float_Inf = octave_Float_NaN = octave_Float_NA = FLT_MAX; |
| 4698 | 64 |
| 4601 | 65 oct_mach_info::float_format ff = oct_mach_info::native_float_format (); |
| 66 | |
| 4698 | 67 switch (ff) |
| 4601 | 68 { |
| 4698 | 69 case oct_mach_info::flt_fmt_ieee_big_endian: |
| 70 case oct_mach_info::flt_fmt_ieee_little_endian: | |
| 71 { | |
| 72 // Don't optimize away tmp_inf / tmp_inf to generate octave_NaN. | |
| 73 | |
| 74 volatile double tmp_inf; | |
| 1967 | 75 |
| 2713 | 76 #if defined (SCO) |
| 4698 | 77 volatile double tmp = 1.0; |
| 78 tmp_inf = 1.0 / (tmp - tmp); | |
| 4164 | 79 #elif defined (__alpha__) && defined (__osf__) |
| 4698 | 80 extern unsigned int DINFINITY[2]; |
| 81 tmp_inf = (*(X_CAST(double *, DINFINITY))); | |
| 1967 | 82 #else |
| 4698 | 83 double tmp = 1e+10; |
| 84 tmp_inf = tmp; | |
| 85 for (;;) | |
| 86 { | |
| 87 tmp_inf *= 1e+10; | |
| 88 if (tmp_inf == tmp) | |
| 89 break; | |
| 90 tmp = tmp_inf; | |
| 91 } | |
| 1967 | 92 #endif |
| 93 | |
| 4164 | 94 #if defined (__alpha__) && defined (__osf__) |
| 4698 | 95 extern unsigned int DQNAN[2]; |
| 96 octave_NaN = (*(X_CAST(double *, DQNAN))); | |
|
9441
160c564d5d25
initialize floating point values properly for NetBSD systems
Aleksej Saushev <asau@inbox.ru>
parents:
8920
diff
changeset
|
97 #elif defined (__NetBSD__) |
|
160c564d5d25
initialize floating point values properly for NetBSD systems
Aleksej Saushev <asau@inbox.ru>
parents:
8920
diff
changeset
|
98 octave_NaN = nan (""); |
| 1967 | 99 #else |
| 4698 | 100 octave_NaN = tmp_inf / tmp_inf; |
|
8029
090001c04619
initialization check for correct NaN sign
Jaroslav Hajek <highegg@gmail.com>
parents:
7814
diff
changeset
|
101 // try to ensure that lo_ieee_sign gives false for a NaN. |
|
090001c04619
initialization check for correct NaN sign
Jaroslav Hajek <highegg@gmail.com>
parents:
7814
diff
changeset
|
102 if (lo_ieee_signbit (octave_NaN)) |
|
090001c04619
initialization check for correct NaN sign
Jaroslav Hajek <highegg@gmail.com>
parents:
7814
diff
changeset
|
103 octave_NaN = -octave_NaN; |
|
090001c04619
initialization check for correct NaN sign
Jaroslav Hajek <highegg@gmail.com>
parents:
7814
diff
changeset
|
104 |
| 1967 | 105 #endif |
| 106 | |
| 4698 | 107 octave_Inf = tmp_inf; |
| 108 | |
| 109 // This is patterned after code in R. | |
| 110 | |
| 111 if (ff == oct_mach_info::flt_fmt_ieee_big_endian) | |
| 112 { | |
| 113 lo_ieee_hw = 0; | |
| 114 lo_ieee_lw = 1; | |
| 115 } | |
| 116 else | |
| 117 { | |
| 118 lo_ieee_hw = 1; | |
| 119 lo_ieee_lw = 0; | |
| 120 } | |
| 4025 | 121 |
| 4698 | 122 lo_ieee_double t; |
| 123 t.word[lo_ieee_hw] = LO_IEEE_NA_HW; | |
| 124 t.word[lo_ieee_lw] = LO_IEEE_NA_LW; | |
| 125 | |
| 126 octave_NA = t.value; | |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
127 |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
128 volatile float float_tmp_inf; |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
129 |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
130 #if defined (SCO) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
131 volatile float float_tmp = 1.0; |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
132 float_tmp_inf = 1.0 / (float_tmp - float_tmp); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
133 #else |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
134 float float_tmp = 1e+10; |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
135 float_tmp_inf = float_tmp; |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
136 for (;;) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
137 { |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
138 float_tmp_inf *= 1e+10; |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
139 if (float_tmp_inf == float_tmp) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
140 break; |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
141 float_tmp = float_tmp_inf; |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
142 } |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
143 #endif |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
144 |
|
9441
160c564d5d25
initialize floating point values properly for NetBSD systems
Aleksej Saushev <asau@inbox.ru>
parents:
8920
diff
changeset
|
145 #if defined (__NetBSD__) |
|
160c564d5d25
initialize floating point values properly for NetBSD systems
Aleksej Saushev <asau@inbox.ru>
parents:
8920
diff
changeset
|
146 octave_Float_NaN = nanf (""); |
|
160c564d5d25
initialize floating point values properly for NetBSD systems
Aleksej Saushev <asau@inbox.ru>
parents:
8920
diff
changeset
|
147 #else |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
148 octave_Float_NaN = float_tmp_inf / float_tmp_inf; |
|
9441
160c564d5d25
initialize floating point values properly for NetBSD systems
Aleksej Saushev <asau@inbox.ru>
parents:
8920
diff
changeset
|
149 #endif |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7231
diff
changeset
|
150 octave_Float_Inf = float_tmp_inf; |
|
7814
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
151 |
|
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
152 lo_ieee_float tf; |
|
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
153 tf.word = LO_IEEE_NA_FLOAT; |
|
87865ed7405f
Second set of single precision test code and fix of resulting bugs
David Bateman <dbateman@free.fr>
parents:
7789
diff
changeset
|
154 octave_Float_NA = tf.value; |
| 4698 | 155 } |
| 156 break; | |
| 4025 | 157 |
| 4698 | 158 case oct_mach_info::flt_fmt_cray: |
| 159 case oct_mach_info::flt_fmt_vax_d: | |
| 160 case oct_mach_info::flt_fmt_vax_g: | |
| 161 default: | |
| 162 // If the format is unknown, then you will probably not have a | |
|
9805
bb70d16cca3b
fail at configure time if IEEE floating point format is not detected
John W. Eaton <jwe@octave.org>
parents:
9592
diff
changeset
|
163 // useful system, so we will abort here. Anyone wishing to |
|
bb70d16cca3b
fail at configure time if IEEE floating point format is not detected
John W. Eaton <jwe@octave.org>
parents:
9592
diff
changeset
|
164 // experiment with building Octave on a system without IEEE |
|
bb70d16cca3b
fail at configure time if IEEE floating point format is not detected
John W. Eaton <jwe@octave.org>
parents:
9592
diff
changeset
|
165 // floating point should be capable of removing this check and |
|
bb70d16cca3b
fail at configure time if IEEE floating point format is not detected
John W. Eaton <jwe@octave.org>
parents:
9592
diff
changeset
|
166 // the configure test. |
|
9591
264fb5520973
abort if floating point format is not recognized as IEEE
John W. Eaton <jwe@octave.org>
parents:
9441
diff
changeset
|
167 (*current_liboctave_error_handler) |
|
264fb5520973
abort if floating point format is not recognized as IEEE
John W. Eaton <jwe@octave.org>
parents:
9441
diff
changeset
|
168 ("lo_ieee_init: floating point format is not IEEE! Maybe DLAMCH is miscompiled, or you are using some strange system without IEEE floating point math?"); |
|
264fb5520973
abort if floating point format is not recognized as IEEE
John W. Eaton <jwe@octave.org>
parents:
9441
diff
changeset
|
169 abort (); |
| 4601 | 170 } |
| 1967 | 171 } |
| 172 | |
| 173 /* | |
| 174 ;;; Local Variables: *** | |
| 175 ;;; mode: C++ *** | |
| 176 ;;; End: *** | |
| 177 */ |