From 66bee2b0ba122d4951d833caeeb0fb9ed7a9b310 Mon Sep 17 00:00:00 2001 From: jfrijters Date: Fri, 11 May 2007 08:55:35 +0000 Subject: [PATCH] Imported Sun's GPLed Double/Float toString and parsing code to fix the long standing incompatibilities. Update license and copyright texts to note the inclusion of Sun code. Note that FloatingDecimal has a couple of IKVM specific changes to work around a bug in the x64 CLR JIT. These changes are marked with /*IKVM*/ --- LICENSE | 364 ++- classpath/AssemblyInfo.java | 5 +- classpath/allsources.lst | 5 +- classpath/java/lang/DoubleToString.java | 663 ------ classpath/java/lang/VMDouble.java | 60 +- classpath/java/lang/VMFloat.java | 14 +- classpath/sun/misc/DoubleConsts.java | 118 + classpath/sun/misc/FloatConsts.java | 113 + classpath/sun/misc/FloatingDecimal.java | 2885 +++++++++++++++++++++++ classpath/sun/misc/FpUtils.java | 1230 ++++++++++ runtime/classpath.cs | 103 - 11 files changed, 4734 insertions(+), 826 deletions(-) delete mode 100644 classpath/java/lang/DoubleToString.java create mode 100644 classpath/sun/misc/DoubleConsts.java create mode 100644 classpath/sun/misc/FloatConsts.java create mode 100644 classpath/sun/misc/FloatingDecimal.java create mode 100644 classpath/sun/misc/FpUtils.java diff --git a/LICENSE b/LICENSE index d7ea9d02..9ff084a4 100644 --- a/LICENSE +++ b/LICENSE @@ -1,16 +1,24 @@ IMPORTANT NOTICE - Some files in this distribution are part of GNU Classpath and have a - different license. This applies in particular to: + Copyright (C) 1998-2007 Free Software Foundation, Inc. + Copyright (C) 1996-2007 Sun Microsystems, Inc. + Copyright (C) 2002-2007 Jeroen Frijters + + Some files in this distribution are part of GNU Classpath or OpenJDK and + are licensed under the GNU General Public License (GPL) version 2 + with "Classpath" exception. This applies in particular to: - IKVM.GNU.Classpath.dll - some of the *.java files (see each file header for license) See http://www.gnu.org/software/classpath/ for information on the - GNU Classpath license. + GNU Classpath license and "Classpath" exception. + + See below for a full copy of the GPL license and the Sun version of the + "Classpath" exception. ----------------------------------------------------------------------------- - Copyright (C) 2002, 2003, 2004, 2005 Jeroen Frijters + Copyright (C) 2002-2007 Jeroen Frijters This software is provided 'as-is', without any express or implied warranty. In no event will the authors be held liable for any damages @@ -32,3 +40,351 @@ IMPORTANT NOTICE jeroen@frijters.net ----------------------------------------------------------------------------- + +The GNU General Public License (GPL) + +Version 2, June 1991 + +Copyright (C) 1989, 1991 Free Software Foundation, Inc. +59 Temple Place, Suite 330, Boston, MA 02111-1307 USA + +Everyone is permitted to copy and distribute verbatim copies of this license +document, but changing it is not allowed. + +Preamble + +The licenses for most software are designed to take away your freedom to share +and change it. By contrast, the GNU General Public License is intended to +guarantee your freedom to share and change free software--to make sure the +software is free for all its users. This General Public License applies to +most of the Free Software Foundation's software and to any other program whose +authors commit to using it. (Some other Free Software Foundation software is +covered by the GNU Library General Public License instead.) You can apply it to +your programs, too. + +When we speak of free software, we are referring to freedom, not price. Our +General Public Licenses are designed to make sure that you have the freedom to +distribute copies of free software (and charge for this service if you wish), +that you receive source code or can get it if you want it, that you can change +the software or use pieces of it in new free programs; and that you know you +can do these things. + +To protect your rights, we need to make restrictions that forbid anyone to deny +you these rights or to ask you to surrender the rights. These restrictions +translate to certain responsibilities for you if you distribute copies of the +software, or if you modify it. + +For example, if you distribute copies of such a program, whether gratis or for +a fee, you must give the recipients all the rights that you have. You must +make sure that they, too, receive or can get the source code. And you must +show them these terms so they know their rights. + +We protect your rights with two steps: (1) copyright the software, and (2) +offer you this license which gives you legal permission to copy, distribute +and/or modify the software. + +Also, for each author's protection and ours, we want to make certain that +everyone understands that there is no warranty for this free software. If the +software is modified by someone else and passed on, we want its recipients to +know that what they have is not the original, so that any problems introduced +by others will not reflect on the original authors' reputations. + +Finally, any free program is threatened constantly by software patents. We +wish to avoid the danger that redistributors of a free program will +individually obtain patent licenses, in effect making the program proprietary. +To prevent this, we have made it clear that any patent must be licensed for +everyone's free use or not licensed at all. + +The precise terms and conditions for copying, distribution and modification +follow. + +TERMS AND CONDITIONS FOR COPYING, DISTRIBUTION AND MODIFICATION + +0. This License applies to any program or other work which contains a notice +placed by the copyright holder saying it may be distributed under the terms of +this General Public License. The "Program", below, refers to any such program +or work, and a "work based on the Program" means either the Program or any +derivative work under copyright law: that is to say, a work containing the +Program or a portion of it, either verbatim or with modifications and/or +translated into another language. (Hereinafter, translation is included +without limitation in the term "modification".) Each licensee is addressed as +"you". + +Activities other than copying, distribution and modification are not covered by +this License; they are outside its scope. The act of running the Program is +not restricted, and the output from the Program is covered only if its contents +constitute a work based on the Program (independent of having been made by +running the Program). Whether that is true depends on what the Program does. + +1. You may copy and distribute verbatim copies of the Program's source code as +you receive it, in any medium, provided that you conspicuously and +appropriately publish on each copy an appropriate copyright notice and +disclaimer of warranty; keep intact all the notices that refer to this License +and to the absence of any warranty; and give any other recipients of the +Program a copy of this License along with the Program. + +You may charge a fee for the physical act of transferring a copy, and you may +at your option offer warranty protection in exchange for a fee. + +2. You may modify your copy or copies of the Program or any portion of it, thus +forming a work based on the Program, and copy and distribute such modifications +or work under the terms of Section 1 above, provided that you also meet all of +these conditions: + + a) You must cause the modified files to carry prominent notices stating + that you changed the files and the date of any change. + + b) You must cause any work that you distribute or publish, that in whole or + in part contains or is derived from the Program or any part thereof, to be + licensed as a whole at no charge to all third parties under the terms of + this License. + + c) If the modified program normally reads commands interactively when run, + you must cause it, when started running for such interactive use in the + most ordinary way, to print or display an announcement including an + appropriate copyright notice and a notice that there is no warranty (or + else, saying that you provide a warranty) and that users may redistribute + the program under these conditions, and telling the user how to view a copy + of this License. (Exception: if the Program itself is interactive but does + not normally print such an announcement, your work based on the Program is + not required to print an announcement.) + +These requirements apply to the modified work as a whole. If identifiable +sections of that work are not derived from the Program, and can be reasonably +considered independent and separate works in themselves, then this License, and +its terms, do not apply to those sections when you distribute them as separate +works. But when you distribute the same sections as part of a whole which is a +work based on the Program, the distribution of the whole must be on the terms +of this License, whose permissions for other licensees extend to the entire +whole, and thus to each and every part regardless of who wrote it. + +Thus, it is not the intent of this section to claim rights or contest your +rights to work written entirely by you; rather, the intent is to exercise the +right to control the distribution of derivative or collective works based on +the Program. + +In addition, mere aggregation of another work not based on the Program with the +Program (or with a work based on the Program) on a volume of a storage or +distribution medium does not bring the other work under the scope of this +License. + +3. You may copy and distribute the Program (or a work based on it, under +Section 2) in object code or executable form under the terms of Sections 1 and +2 above provided that you also do one of the following: + + a) Accompany it with the complete corresponding machine-readable source + code, which must be distributed under the terms of Sections 1 and 2 above + on a medium customarily used for software interchange; or, + + b) Accompany it with a written offer, valid for at least three years, to + give any third party, for a charge no more than your cost of physically + performing source distribution, a complete machine-readable copy of the + corresponding source code, to be distributed under the terms of Sections 1 + and 2 above on a medium customarily used for software interchange; or, + + c) Accompany it with the information you received as to the offer to + distribute corresponding source code. (This alternative is allowed only + for noncommercial distribution and only if you received the program in + object code or executable form with such an offer, in accord with + Subsection b above.) + +The source code for a work means the preferred form of the work for making +modifications to it. For an executable work, complete source code means all +the source code for all modules it contains, plus any associated interface +definition files, plus the scripts used to control compilation and installation +of the executable. However, as a special exception, the source code +distributed need not include anything that is normally distributed (in either +source or binary form) with the major components (compiler, kernel, and so on) +of the operating system on which the executable runs, unless that component +itself accompanies the executable. + +If distribution of executable or object code is made by offering access to copy +from a designated place, then offering equivalent access to copy the source +code from the same place counts as distribution of the source code, even though +third parties are not compelled to copy the source along with the object code. + +4. You may not copy, modify, sublicense, or distribute the Program except as +expressly provided under this License. Any attempt otherwise to copy, modify, +sublicense or distribute the Program is void, and will automatically terminate +your rights under this License. However, parties who have received copies, or +rights, from you under this License will not have their licenses terminated so +long as such parties remain in full compliance. + +5. You are not required to accept this License, since you have not signed it. +However, nothing else grants you permission to modify or distribute the Program +or its derivative works. These actions are prohibited by law if you do not +accept this License. Therefore, by modifying or distributing the Program (or +any work based on the Program), you indicate your acceptance of this License to +do so, and all its terms and conditions for copying, distributing or modifying +the Program or works based on it. + +6. Each time you redistribute the Program (or any work based on the Program), +the recipient automatically receives a license from the original licensor to +copy, distribute or modify the Program subject to these terms and conditions. +You may not impose any further restrictions on the recipients' exercise of the +rights granted herein. You are not responsible for enforcing compliance by +third parties to this License. + +7. If, as a consequence of a court judgment or allegation of patent +infringement or for any other reason (not limited to patent issues), conditions +are imposed on you (whether by court order, agreement or otherwise) that +contradict the conditions of this License, they do not excuse you from the +conditions of this License. If you cannot distribute so as to satisfy +simultaneously your obligations under this License and any other pertinent +obligations, then as a consequence you may not distribute the Program at all. +For example, if a patent license would not permit royalty-free redistribution +of the Program by all those who receive copies directly or indirectly through +you, then the only way you could satisfy both it and this License would be to +refrain entirely from distribution of the Program. + +If any portion of this section is held invalid or unenforceable under any +particular circumstance, the balance of the section is intended to apply and +the section as a whole is intended to apply in other circumstances. + +It is not the purpose of this section to induce you to infringe any patents or +other property right claims or to contest validity of any such claims; this +section has the sole purpose of protecting the integrity of the free software +distribution system, which is implemented by public license practices. Many +people have made generous contributions to the wide range of software +distributed through that system in reliance on consistent application of that +system; it is up to the author/donor to decide if he or she is willing to +distribute software through any other system and a licensee cannot impose that +choice. + +This section is intended to make thoroughly clear what is believed to be a +consequence of the rest of this License. + +8. If the distribution and/or use of the Program is restricted in certain +countries either by patents or by copyrighted interfaces, the original +copyright holder who places the Program under this License may add an explicit +geographical distribution limitation excluding those countries, so that +distribution is permitted only in or among countries not thus excluded. In +such case, this License incorporates the limitation as if written in the body +of this License. + +9. The Free Software Foundation may publish revised and/or new versions of the +General Public License from time to time. Such new versions will be similar in +spirit to the present version, but may differ in detail to address new problems +or concerns. + +Each version is given a distinguishing version number. If the Program +specifies a version number of this License which applies to it and "any later +version", you have the option of following the terms and conditions either of +that version or of any later version published by the Free Software Foundation. +If the Program does not specify a version number of this License, you may +choose any version ever published by the Free Software Foundation. + +10. If you wish to incorporate parts of the Program into other free programs +whose distribution conditions are different, write to the author to ask for +permission. For software which is copyrighted by the Free Software Foundation, +write to the Free Software Foundation; we sometimes make exceptions for this. +Our decision will be guided by the two goals of preserving the free status of +all derivatives of our free software and of promoting the sharing and reuse of +software generally. + +NO WARRANTY + +11. BECAUSE THE PROGRAM IS LICENSED FREE OF CHARGE, THERE IS NO WARRANTY FOR +THE PROGRAM, TO THE EXTENT PERMITTED BY APPLICABLE LAW. EXCEPT WHEN OTHERWISE +STATED IN WRITING THE COPYRIGHT HOLDERS AND/OR OTHER PARTIES PROVIDE THE +PROGRAM "AS IS" WITHOUT WARRANTY OF ANY KIND, EITHER EXPRESSED OR IMPLIED, +INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND +FITNESS FOR A PARTICULAR PURPOSE. THE ENTIRE RISK AS TO THE QUALITY AND +PERFORMANCE OF THE PROGRAM IS WITH YOU. SHOULD THE PROGRAM PROVE DEFECTIVE, +YOU ASSUME THE COST OF ALL NECESSARY SERVICING, REPAIR OR CORRECTION. + +12. IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING WILL +ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MAY MODIFY AND/OR REDISTRIBUTE THE +PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES, INCLUDING ANY +GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING OUT OF THE USE OR +INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED TO LOSS OF DATA OR DATA +BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY YOU OR THIRD PARTIES OR A +FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER PROGRAMS), EVEN IF SUCH HOLDER +OR OTHER PARTY HAS BEEN ADVISED OF THE POSSIBILITY OF SUCH DAMAGES. + +END OF TERMS AND CONDITIONS + +How to Apply These Terms to Your New Programs + +If you develop a new program, and you want it to be of the greatest possible +use to the public, the best way to achieve this is to make it free software +which everyone can redistribute and change under these terms. + +To do so, attach the following notices to the program. It is safest to attach +them to the start of each source file to most effectively convey the exclusion +of warranty; and each file should have at least the "copyright" line and a +pointer to where the full notice is found. + + One line to give the program's name and a brief idea of what it does. + + Copyright (C) + + This program is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the Free + Software Foundation; either version 2 of the License, or (at your option) + any later version. + + This program is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for + more details. + + You should have received a copy of the GNU General Public License along + with this program; if not, write to the Free Software Foundation, Inc., 59 + Temple Place, Suite 330, Boston, MA 02111-1307 USA + +Also add information on how to contact you by electronic and paper mail. + +If the program is interactive, make it output a short notice like this when it +starts in an interactive mode: + + Gnomovision version 69, Copyright (C) year name of author Gnomovision comes + with ABSOLUTELY NO WARRANTY; for details type 'show w'. This is free + software, and you are welcome to redistribute it under certain conditions; + type 'show c' for details. + +The hypothetical commands 'show w' and 'show c' should show the appropriate +parts of the General Public License. Of course, the commands you use may be +called something other than 'show w' and 'show c'; they could even be +mouse-clicks or menu items--whatever suits your program. + +You should also get your employer (if you work as a programmer) or your school, +if any, to sign a "copyright disclaimer" for the program, if necessary. Here +is a sample; alter the names: + + Yoyodyne, Inc., hereby disclaims all copyright interest in the program + 'Gnomovision' (which makes passes at compilers) written by James Hacker. + + signature of Ty Coon, 1 April 1989 + + Ty Coon, President of Vice + +This General Public License does not permit incorporating your program into +proprietary programs. If your program is a subroutine library, you may +consider it more useful to permit linking proprietary applications with the +library. If this is what you want to do, use the GNU Library General Public +License instead of this License. + + +"CLASSPATH" EXCEPTION TO THE GPL + +Certain source files distributed by Sun Microsystems, Inc. are subject to +the following clarification and special exception to the GPL, but only where +Sun has expressly included in the particular source file's header the words +"Sun designates this particular file as subject to the "Classpath" exception +as provided by Sun in the LICENSE file that accompanied this code." + + Linking this library statically or dynamically with other modules is making + a combined work based on this library. Thus, the terms and conditions of + the GNU General Public License cover the whole combination. + + As a special exception, the copyright holders of this library give you + permission to link this library with independent modules to produce an + executable, regardless of the license terms of these independent modules, + and to copy and distribute the resulting executable under terms of your + choice, provided that you also meet, for each linked independent module, + the terms and conditions of the license of that module. An independent + module is a module which is not derived from or based on this library. If + you modify this library, you may extend this exception to your version of + the library, but you are not obligated to do so. If you do not wish to do + so, delete this exception statement from your version. diff --git a/classpath/AssemblyInfo.java b/classpath/AssemblyInfo.java index 02937d1c..d20abd15 100644 --- a/classpath/AssemblyInfo.java +++ b/classpath/AssemblyInfo.java @@ -1,5 +1,5 @@ /* - Copyright (C) 2006 Jeroen Frijters + Copyright (C) 2006, 2007 Jeroen Frijters This software is provided 'as-is', without any express or implied warranty. In no event will the authors be held liable for any damages @@ -23,9 +23,10 @@ */ @cli.System.Reflection.AssemblyCopyrightAttribute.Annotation( - "This software is licensed under the GNU General Public License + GNU Classpath exception.\r\n" + + "This software is licensed under the GNU General Public License version 2 + \"Classpath\" exception.\r\n" + "See http://www.gnu.org/software/classpath/license.html for details.\r\n" + "Copyright (C) 1998-2007 Free Software Foundation, Inc.\r\n" + + "Copyright (C) 1996-2007 Sun Microsystems, Inc.\r\n" + "Copyright (C) 2002-2007 Jeroen Frijters" ) diff --git a/classpath/allsources.lst b/classpath/allsources.lst index dc8cbe06..b6b35a8d 100644 --- a/classpath/allsources.lst +++ b/classpath/allsources.lst @@ -5002,7 +5002,6 @@ ikvm/runtime/Startup.java ikvm/runtime/Util.java java/io/VMFile.java java/io/VMObjectInputStream.java -java/lang/DoubleToString.java java/lang/ExceptionHelper.java java/lang/LibraryVMInterfaceImpl.java java/lang/management/VMManagementFactory.java @@ -5041,6 +5040,10 @@ java/util/concurrent/atomic/AtomicReference.java java/util/concurrent/atomic/AtomicReferenceArray.java java/util/concurrent/locks/LockSupport.java java/util/VMTimeZone.java +sun/misc/DoubleConsts.java +sun/misc/FloatConsts.java +sun/misc/FloatingDecimal.java +sun/misc/FpUtils.java sun/misc/Ref.java sun/misc/Unsafe.java sun/reflect/Reflection.java diff --git a/classpath/java/lang/DoubleToString.java b/classpath/java/lang/DoubleToString.java deleted file mode 100644 index 8b96db3c..00000000 --- a/classpath/java/lang/DoubleToString.java +++ /dev/null @@ -1,663 +0,0 @@ -// NOTE this code was adapted from source code accompanying the article -// http://www.onjava.com/pub/a/onjava/2000/12/15/formatting_doubles.html?page=2 -// by Jack Shirazi - -package java.lang; - -class DoubleToString -{ - //Hardcode some byte arrays to make them quickly available - public static final char[] INFINITY = {'I','n','f','i','n','i','t','y'}; - public static final char[] NaN = {'N','a','N'}; - public static final char[][] ZEROS = { - {}, - {'0'}, - {'0','0'}, - {'0','0','0'}, - {'0','0','0','0'}, - {'0','0','0','0','0'}, - {'0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0','0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0'}, - {'0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0','0'}, - }; - - private static final char[] charForDigit = { - '0','1','2','3','4','5','6','7','8','9','a','b','c','d','e','f','g','h', - 'i','j','k','l','m','n','o','p','q','r','s','t','u','v','w','x','y','z' - }; - - //And required double related constants. - private static final long DoubleSignMask = 0x8000000000000000L; - private static final long DoubleExpMask = 0x7ff0000000000000L; - private static final long DoubleFractMask= ~(DoubleSignMask|DoubleExpMask); - private static final int DoubleExpShift = 52; - private static final int DoubleExpBias = 1023; - - private static final double[] d_tenthPowers = { -1e-323D, 1e-322D, 1e-321D, 1e-320D, 1e-319D, 1e-318D, 1e-317D, 1e-316D, 1e-315D, 1e-314D, -1e-313D, 1e-312D, 1e-311D, 1e-310D, 1e-309D, 1e-308D, 1e-307D, 1e-306D, 1e-305D, 1e-304D, -1e-303D, 1e-302D, 1e-301D, 1e-300D, 1e-299D, 1e-298D, 1e-297D, 1e-296D, 1e-295D, 1e-294D, -1e-293D, 1e-292D, 1e-291D, 1e-290D, 1e-289D, 1e-288D, 1e-287D, 1e-286D, 1e-285D, 1e-284D, -1e-283D, 1e-282D, 1e-281D, 1e-280D, 1e-279D, 1e-278D, 1e-277D, 1e-276D, 1e-275D, 1e-274D, -1e-273D, 1e-272D, 1e-271D, 1e-270D, 1e-269D, 1e-268D, 1e-267D, 1e-266D, 1e-265D, 1e-264D, -1e-263D, 1e-262D, 1e-261D, 1e-260D, 1e-259D, 1e-258D, 1e-257D, 1e-256D, 1e-255D, 1e-254D, -1e-253D, 1e-252D, 1e-251D, 1e-250D, 1e-249D, 1e-248D, 1e-247D, 1e-246D, 1e-245D, 1e-244D, -1e-243D, 1e-242D, 1e-241D, 1e-240D, 1e-239D, 1e-238D, 1e-237D, 1e-236D, 1e-235D, 1e-234D, -1e-233D, 1e-232D, 1e-231D, 1e-230D, 1e-229D, 1e-228D, 1e-227D, 1e-226D, 1e-225D, 1e-224D, -1e-223D, 1e-222D, 1e-221D, 1e-220D, 1e-219D, 1e-218D, 1e-217D, 1e-216D, 1e-215D, 1e-214D, -1e-213D, 1e-212D, 1e-211D, 1e-210D, 1e-209D, 1e-208D, 1e-207D, 1e-206D, 1e-205D, 1e-204D, -1e-203D, 1e-202D, 1e-201D, 1e-200D, 1e-199D, 1e-198D, 1e-197D, 1e-196D, 1e-195D, 1e-194D, -1e-193D, 1e-192D, 1e-191D, 1e-190D, 1e-189D, 1e-188D, 1e-187D, 1e-186D, 1e-185D, 1e-184D, -1e-183D, 1e-182D, 1e-181D, 1e-180D, 1e-179D, 1e-178D, 1e-177D, 1e-176D, 1e-175D, 1e-174D, -1e-173D, 1e-172D, 1e-171D, 1e-170D, 1e-169D, 1e-168D, 1e-167D, 1e-166D, 1e-165D, 1e-164D, -1e-163D, 1e-162D, 1e-161D, 1e-160D, 1e-159D, 1e-158D, 1e-157D, 1e-156D, 1e-155D, 1e-154D, -1e-153D, 1e-152D, 1e-151D, 1e-150D, 1e-149D, 1e-148D, 1e-147D, 1e-146D, 1e-145D, 1e-144D, -1e-143D, 1e-142D, 1e-141D, 1e-140D, 1e-139D, 1e-138D, 1e-137D, 1e-136D, 1e-135D, 1e-134D, -1e-133D, 1e-132D, 1e-131D, 1e-130D, 1e-129D, 1e-128D, 1e-127D, 1e-126D, 1e-125D, 1e-124D, -1e-123D, 1e-122D, 1e-121D, 1e-120D, 1e-119D, 1e-118D, 1e-117D, 1e-116D, 1e-115D, 1e-114D, -1e-113D, 1e-112D, 1e-111D, 1e-110D, 1e-109D, 1e-108D, 1e-107D, 1e-106D, 1e-105D, 1e-104D, -1e-103D, 1e-102D, 1e-101D, 1e-100D, 1e-99D, 1e-98D, 1e-97D, 1e-96D, 1e-95D, 1e-94D, -1e-93D, 1e-92D, 1e-91D, 1e-90D, 1e-89D, 1e-88D, 1e-87D, 1e-86D, 1e-85D, 1e-84D, -1e-83D, 1e-82D, 1e-81D, 1e-80D, 1e-79D, 1e-78D, 1e-77D, 1e-76D, 1e-75D, 1e-74D, -1e-73D, 1e-72D, 1e-71D, 1e-70D, 1e-69D, 1e-68D, 1e-67D, 1e-66D, 1e-65D, 1e-64D, -1e-63D, 1e-62D, 1e-61D, 1e-60D, 1e-59D, 1e-58D, 1e-57D, 1e-56D, 1e-55D, 1e-54D, -1e-53D, 1e-52D, 1e-51D, 1e-50D, 1e-49D, 1e-48D, 1e-47D, 1e-46D, 1e-45D, 1e-44D, -1e-43D, 1e-42D, 1e-41D, 1e-40D, 1e-39D, 1e-38D, 1e-37D, 1e-36D, 1e-35D, 1e-34D, -1e-33D, 1e-32D, 1e-31D, 1e-30D, 1e-29D, 1e-28D, 1e-27D, 1e-26D, 1e-25D, 1e-24D, -1e-23D, 1e-22D, 1e-21D, 1e-20D, 1e-19D, 1e-18D, 1e-17D, 1e-16D, 1e-15D, 1e-14D, -1e-13D, 1e-12D, 1e-11D, 1e-10D, 1e-9D, 1e-8D, 1e-7D, 1e-6D, 1e-5D, 1e-4D, -1e-3D, 1e-2D, 1e-1D, 1e0D, 1e1D, 1e2D, 1e3D, 1e4D, -1e5D, 1e6D, 1e7D, 1e8D, 1e9D, 1e10D, 1e11D, 1e12D, 1e13D, 1e14D, -1e15D, 1e16D, 1e17D, 1e18D, 1e19D, 1e20D, 1e21D, 1e22D, 1e23D, 1e24D, -1e25D, 1e26D, 1e27D, 1e28D, 1e29D, 1e30D, 1e31D, 1e32D, 1e33D, 1e34D, -1e35D, 1e36D, 1e37D, 1e38D, 1e39D, 1e40D, 1e41D, 1e42D, 1e43D, 1e44D, -1e45D, 1e46D, 1e47D, 1e48D, 1e49D, 1e50D, 1e51D, 1e52D, 1e53D, 1e54D, -1e55D, 1e56D, 1e57D, 1e58D, 1e59D, 1e60D, 1e61D, 1e62D, 1e63D, 1e64D, -1e65D, 1e66D, 1e67D, 1e68D, 1e69D, 1e70D, 1e71D, 1e72D, 1e73D, 1e74D, -1e75D, 1e76D, 1e77D, 1e78D, 1e79D, 1e80D, 1e81D, 1e82D, 1e83D, 1e84D, -1e85D, 1e86D, 1e87D, 1e88D, 1e89D, 1e90D, 1e91D, 1e92D, 1e93D, 1e94D, -1e95D, 1e96D, 1e97D, 1e98D, 1e99D, 1e100D, 1e101D, 1e102D, 1e103D, 1e104D, -1e105D, 1e106D, 1e107D, 1e108D, 1e109D, 1e110D, 1e111D, 1e112D, 1e113D, 1e114D, -1e115D, 1e116D, 1e117D, 1e118D, 1e119D, 1e120D, 1e121D, 1e122D, 1e123D, 1e124D, -1e125D, 1e126D, 1e127D, 1e128D, 1e129D, 1e130D, 1e131D, 1e132D, 1e133D, 1e134D, -1e135D, 1e136D, 1e137D, 1e138D, 1e139D, 1e140D, 1e141D, 1e142D, 1e143D, 1e144D, -1e145D, 1e146D, 1e147D, 1e148D, 1e149D, 1e150D, 1e151D, 1e152D, 1e153D, 1e154D, -1e155D, 1e156D, 1e157D, 1e158D, 1e159D, 1e160D, 1e161D, 1e162D, 1e163D, 1e164D, -1e165D, 1e166D, 1e167D, 1e168D, 1e169D, 1e170D, 1e171D, 1e172D, 1e173D, 1e174D, -1e175D, 1e176D, 1e177D, 1e178D, 1e179D, 1e180D, 1e181D, 1e182D, 1e183D, 1e184D, -1e185D, 1e186D, 1e187D, 1e188D, 1e189D, 1e190D, 1e191D, 1e192D, 1e193D, 1e194D, -1e195D, 1e196D, 1e197D, 1e198D, 1e199D, 1e200D, 1e201D, 1e202D, 1e203D, 1e204D, -1e205D, 1e206D, 1e207D, 1e208D, 1e209D, 1e210D, 1e211D, 1e212D, 1e213D, 1e214D, -1e215D, 1e216D, 1e217D, 1e218D, 1e219D, 1e220D, 1e221D, 1e222D, 1e223D, 1e224D, -1e225D, 1e226D, 1e227D, 1e228D, 1e229D, 1e230D, 1e231D, 1e232D, 1e233D, 1e234D, -1e235D, 1e236D, 1e237D, 1e238D, 1e239D, 1e240D, 1e241D, 1e242D, 1e243D, 1e244D, -1e245D, 1e246D, 1e247D, 1e248D, 1e249D, 1e250D, 1e251D, 1e252D, 1e253D, 1e254D, -1e255D, 1e256D, 1e257D, 1e258D, 1e259D, 1e260D, 1e261D, 1e262D, 1e263D, 1e264D, -1e265D, 1e266D, 1e267D, 1e268D, 1e269D, 1e270D, 1e271D, 1e272D, 1e273D, 1e274D, -1e275D, 1e276D, 1e277D, 1e278D, 1e279D, 1e280D, 1e281D, 1e282D, 1e283D, 1e284D, -1e285D, 1e286D, 1e287D, 1e288D, 1e289D, 1e290D, 1e291D, 1e292D, 1e293D, 1e294D, -1e295D, 1e296D, 1e297D, 1e298D, 1e299D, 1e300D, 1e301D, 1e302D, 1e303D, 1e304D, -1e305D, 1e306D, 1e307D, 1e308D - }; - - - public void appendFormatted(StringBuffer s, double d, int numFractDigits, - char decimalPoint, char thousandsSeparator, int numDigitsSeparated, - char negativePrefix, char negativeSuffix) - { - //First check for the special cases, +/-infinity, Not-a-number and -0.0 - if (d == Double.NEGATIVE_INFINITY) - { - //d == -Infinity - if (negativePrefix != '\uFFFF') - s.append(negativePrefix); - s.append(INFINITY); - if (negativeSuffix != '\uFFFF') - s.append(negativeSuffix); - } - else if (d == Double.POSITIVE_INFINITY) - //d == Infinity - s.append(INFINITY); - else if (d != d) - //d == NaN - s.append(NaN); - else if (d == 0.0) - { - if ( (Double.doubleToLongBits(d) & DoubleSignMask) != 0) - { - //d == -0.0 - if (negativePrefix != '\uFFFF') - s.append(negativePrefix); - s.append('0').append(decimalPoint).append(ZEROS[numFractDigits]); - if (negativeSuffix != '\uFFFF') - s.append(negativeSuffix); - } - else - //d == 0.0 - s.append('0').append(decimalPoint).append(ZEROS[numFractDigits]); - } - else - { - //convert to a positive format, and record whether we have a negative - //number so that we know later whether to add the negativeSuffix - boolean negative = false; - if (d < 0) - { - //Even if the number is negative, we only need to set the - //negative flag if there is a printable negativeSuffix - if (negativeSuffix != '\uFFFF') - negative = true; - if (negativePrefix != '\uFFFF') - s.append(negativePrefix); - d = -d; - } - - //Find the magnitude. This is basically the exponent in normal form. - int magnitude = magnitude(d); - - //First off, if the number is too small for the given format, we - //only print 0.0..., which makes this real quick - if ( (magnitude + numFractDigits) < 0) - { - appendNearlyZeroNumber(s, d, magnitude, numFractDigits, decimalPoint); - if (negative) - s.append(negativeSuffix); - return; - } - - long l; - //Now scale the double to the biggest long value we need - //We need to handle the smallest magnitudes differently because of rounding errors - - //This test is unlikely to ever be true. It would require numFractDigits - //to be 305 or more, which is pretty unlikely. - if (magnitude < -305) - l = (long) ((d*1E18) / d_tenthPowers[magnitude + 324]); - else - l = (long) (d / d_tenthPowers[magnitude + 323 - 17]); - - //And round up if necessary. Add one to the numFractDigits digit if the - //numFractDigits+1 digit is 5 or greater. It is useful to know that - //given a long, l, the nth digit is obtained using the formula - // nthDigit = (l/(tenthPower(l)/l_tenthPowers[n-1]))%10; - - long l_tenthPower = tenthPower(l); - //The numFractDigits+1 digit of the double is the - //numFractDigits+1+magnitude digit of the long. - //We only need worry about digits within the long. Very large numbers are - //not rounded because all the digits after the decimal points are 0 anyway - if (numFractDigits+magnitude+1 < l_tenthPowers.length) - { - long digit = (l/(l_tenthPower/l_tenthPowers[numFractDigits+magnitude+1]))%10; - if (digit >= 5) - { - l += l_tenthPower/l_tenthPowers[numFractDigits+magnitude]; - } - } - - //And now we just print out our long, with the decimal point character - //inserted in the right place, using as many places as we wanted. - appendAsDouble(s, l, l_tenthPower, magnitude, numFractDigits, decimalPoint, thousandsSeparator, - numDigitsSeparated, negativePrefix, negativeSuffix); - - //Finally, append the negativeSuffix if necessary - if (negative) - s.append(negativeSuffix); - } - } - - public void appendAsDouble(StringBuffer s, long l, long l_mag, int d_magnitude, - int numFractDigits, char decimalPoint, char thousandsSeparator, - int numDigitsSeparated, char negativePrefix, char negativeSuffix) - { - //If the magnitude is negative, we have a 0.xxx number - if (d_magnitude < 0) - { - s.append('0').append(decimalPoint).append(ZEROS[-d_magnitude-1]); - //And just print successive digits until we have reached numFractDigits - //First decrement numFractDigits by the number of digits already printed - numFractDigits += d_magnitude; - - //get the magnitude of l - long c; - while(numFractDigits-- >= 0) - { - //Get the leading character (e.g. '62345/10000 = 6' using integer-divide) - c = l/l_mag; - //Append the digit character for this digit (e.g. number is 6, so character is '6') - s.append(charForDigit[(int) c]); - //Multiply by the leading digit by the magnitude so that we can eliminate the leading digit - //(e.g. 6 * 10000 = 60000) - c *= l_mag; - //and eliminate the leading digit (e.g. 62345-60000 = 2345) - if ( c <= l) - l -= c; - //Decrease the magnitude by 10, and repeat the loop. - l_mag = l_mag/10; - } - } - else - { - //Just keep printing until magnitude is 0 - long c; - while(d_magnitude-- >= 0) - { - if (l_mag == 0) {s.append('0');continue;} - //Get the leading character (e.g. '62345/10000 = 6' using integer-divide) - c = l/l_mag; - //Append the digit character for this digit (e.g. number is 6, so character is '6') - s.append(charForDigit[(int) c]); - - //Don't forget about the thousands separator - if (d_magnitude % numDigitsSeparated == (numDigitsSeparated-1)) - s.append(thousandsSeparator); - - //Multiply by the leading digit by the magnitude so that we can eliminate the leading digit - //(e.g. 6 * 10000 = 60000) - c *= l_mag; - //and eliminate the leading digit (e.g. 62345-60000 = 2345) - if ( c <= l) - l -= c; - //Decrease the magnitude by 10, and repeat the loop. - l_mag = l_mag/10; - } - s.append(decimalPoint); - if (l_mag == 0) - s.append(ZEROS[numFractDigits]); - else - { - while(numFractDigits-- > 0) - { - if (l_mag == 0) {s.append('0');continue;} - //Get the leading character (e.g. '62345/10000 = 6' using integer-divide) - c = l/l_mag; - //Append the digit character for this digit (e.g. number is 6, so character is '6') - s.append(charForDigit[(int) c]); - //Multiply by the leading digit by the magnitude so that we can eliminate the leading digit - //(e.g. 6 * 10000 = 60000) - c *= l_mag; - //and eliminate the leading digit (e.g. 62345-60000 = 2345) - if ( c <= l) - l -= c; - //Decrease the magnitude by 10, and repeat the loop. - l_mag = l_mag/10; - } - } - } - } - - - private void appendNearlyZeroNumber(StringBuffer s, double d, int d_magnitude, - int numFractDigits, char decimalPoint) - { - if (d_magnitude + numFractDigits == -1) - { - //Possibly too small, depends on whether the top digit is 5 or greater - //So we have to scale to get the leading digit - int i; - if (d_magnitude < -305) - //Probably not necessary. Who is going to print 305 places? - i = (int) ((d*1E19) / d_tenthPowers[d_magnitude + 324 + 18]); - else - i = (int) (d / d_tenthPowers[d_magnitude + 323]); - - if (i >= 5) - { - //Not too small, we get to round up - s.append('0').append(decimalPoint).append(ZEROS[numFractDigits-1]); - s.append('1'); - } - else - { - //Definitely too small. Just print zeros - s.append('0').append(decimalPoint).append(ZEROS[numFractDigits]); - } - } - else - { - //Definitely too small - s.append('0').append(decimalPoint).append(ZEROS[numFractDigits]); - } - } - - /** - * Assumes i is positive. Returns the magnitude of i in base 10. - */ - private static long tenthPower(long i) - { - if (i < 10L) return 1; - else if (i < 100L) return 10L; - else if (i < 1000L) return 100L; - else if (i < 10000L) return 1000L; - else if (i < 100000L) return 10000L; - else if (i < 1000000L) return 100000L; - else if (i < 10000000L) return 1000000L; - else if (i < 100000000L) return 10000000L; - else if (i < 1000000000L) return 100000000L; - else if (i < 10000000000L) return 1000000000L; - else if (i < 100000000000L) return 10000000000L; - else if (i < 1000000000000L) return 100000000000L; - else if (i < 10000000000000L) return 1000000000000L; - else if (i < 100000000000000L) return 10000000000000L; - else if (i < 1000000000000000L) return 100000000000000L; - else if (i < 10000000000000000L) return 1000000000000000L; - else if (i < 100000000000000000L) return 10000000000000000L; - else if (i < 1000000000000000000L) return 100000000000000000L; - else return 1000000000000000000L; - } - - - private static int magnitude(double d) - { - //It works. What else can I say. - long doubleToLongBits = Double.doubleToLongBits(d); - int magnitude = - (int) ((((doubleToLongBits & DoubleExpMask) >> DoubleExpShift) - DoubleExpBias) * 0.301029995663981); - - if (magnitude < -323) - magnitude = -323; - else if (magnitude > 308) - magnitude = 308; - - if (d >= d_tenthPowers[magnitude+323]) - { - while(magnitude < 309 && d >= d_tenthPowers[magnitude+323]) - magnitude++; - magnitude--; - return magnitude; - } - else - { - while(magnitude > -324 && d < d_tenthPowers[magnitude+323]) - magnitude--; - return magnitude; - } - } - - static long[] l_tenthPowers = { - 1, - 10L, - 100L, - 1000L, - 10000L, - 100000L, - 1000000L, - 10000000L, - 100000000L, - 1000000000L, - 10000000000L, - 100000000000L, - 1000000000000L, - 10000000000000L, - 100000000000000L, - 1000000000000000L, - 10000000000000000L, - 100000000000000000L, - 1000000000000000000L, - }; - -public static void append(StringBuffer s, double d) -{ - if (d == Double.NEGATIVE_INFINITY) - s.append(NEGATIVE_INFINITY); - else if (d == Double.POSITIVE_INFINITY) - s.append(POSITIVE_INFINITY); - else if (d != d) - s.append(NaN); - else if (d == 0.0) - { - if ( (Double.doubleToLongBits(d) & DoubleSignMask) != 0) - s.append('-'); - s.append(DOUBLE_ZERO); - } - else - { - if (d < 0) - { - s.append('-'); - d = -d; - } - - if (d >= 0.001 && d < 0.01) - { - long i = (long) (d * 1E20); - i = i%100 >= 50 ? (i/100) + 1 : i/100; - s.append(DOUBLE_ZERO2); - appendFractDigits(s, i,-1); - } - else if (d >= 0.01 && d < 0.1) - { - long i = (long) (d * 1E19); - i = i%100 >= 50 ? (i/100) + 1 : i/100; - s.append(DOUBLE_ZERO); - appendFractDigits(s, i,-1); - } - else if (d >= 0.1 && d < 1) - { - long i = (long) (d * 1E18); - i = i%100 >= 50 ? (i/100) + 1 : i/100; - s.append(DOUBLE_ZERO0); - appendFractDigits(s, i,-1); - } - else if (d >= 1 && d < 10) - { - long i = (long) (d * 1E17); - i = i%100 >= 50 ? (i/100) + 1 : i/100; - appendFractDigits(s, i,1); - } - else if (d >= 10 && d < 100) - { - long i = (long) (d * 1E16); - i = i%100 >= 50 ? (i/100) + 1 : i/100; - appendFractDigits(s, i,2); - } - else if (d >= 100 && d < 1000) - { - long i = (long) (d * 1E15); - i = i%100 >= 50 ? (i/100) + 1 : i/100; - appendFractDigits(s, i,3); - } - else if (d >= 1000 && d < 10000) - { - long i = (long) (d * 1E14); - i = i%100 >= 50 ? (i/100) + 1 : i/100; - appendFractDigits(s, i,4); - } - else if (d >= 10000 && d < 100000) - { - long i = (long) (d * 1E13); - i = i%100 >= 50 ? (i/100) + 1 : i/100; - appendFractDigits(s, i,5); - } - else if (d >= 100000 && d < 1000000) - { - long i = (long) (d * 1E12); - i = i%100 >= 50 ? (i/100) + 1 : i/100; - appendFractDigits(s, i,6); - } - else if (d >= 1000000 && d < 10000000) - { - long i = (long) (d * 1E11); - i = i%100 >= 50 ? (i/100) + 1 : i/100; - appendFractDigits(s, i,7); - } - else if (d == Double.MIN_VALUE) - s.append("4.9E-324"); - else if (d == Double.MAX_VALUE) - s.append("1.7976931348623157E308"); - else - { - int magnitude = magnitude(d); - long i; - if (magnitude < -305) - i = (long) (d*1E18 / d_tenthPowers[magnitude + 324]); - else - i = (long) (d / d_tenthPowers[magnitude + 323 - 17]); - i = i%100 >= 50 ? (i/100) + 1 : i/100; - appendFractDigits(s, i, 1); - s.append('E'); - append(s,magnitude); - } - } -} -public static void append(StringBuffer s, int i) -{ - if (i < 0) - { - if (i == Integer.MIN_VALUE) - { - //cannot make this positive due to integer overflow - s.append("-2147483648"); - } - s.append('-'); - i = -i; - } - int mag; - int c; - if (i < 10) - { - //one digit - s.append(charForDigit[i]); - } - else if (i < 100) - { - //two digits - s.append(charForDigit[i/10]); - s.append(charForDigit[i%10]); - } - else if (i < 1000) - { - //three digits - s.append(charForDigit[i/100]); - s.append(charForDigit[(c=i%100)/10]); - s.append(charForDigit[c%10]); - } - else if (i < 10000) - { - //four digits - s.append(charForDigit[i/1000]); - s.append(charForDigit[(c=i%1000)/100]); - s.append(charForDigit[(c%=100)/10]); - s.append(charForDigit[c%10]); - } - else if (i < 100000) - { - //five digits - s.append(charForDigit[i/10000]); - s.append(charForDigit[(c=i%10000)/1000]); - s.append(charForDigit[(c%=1000)/100]); - s.append(charForDigit[(c%=100)/10]); - s.append(charForDigit[c%10]); - } - else if (i < 1000000) - { - //six digits - s.append(charForDigit[i/100000]); - s.append(charForDigit[(c=i%100000)/10000]); - s.append(charForDigit[(c%=10000)/1000]); - s.append(charForDigit[(c%=1000)/100]); - s.append(charForDigit[(c%=100)/10]); - s.append(charForDigit[c%10]); - } - else if (i < 10000000) - { - //seven digits - s.append(charForDigit[i/1000000]); - s.append(charForDigit[(c=i%1000000)/100000]); - s.append(charForDigit[(c%=100000)/10000]); - s.append(charForDigit[(c%=10000)/1000]); - s.append(charForDigit[(c%=1000)/100]); - s.append(charForDigit[(c%=100)/10]); - s.append(charForDigit[c%10]); - } - else if (i < 100000000) - { - //eight digits - s.append(charForDigit[i/10000000]); - s.append(charForDigit[(c=i%10000000)/1000000]); - s.append(charForDigit[(c%=1000000)/100000]); - s.append(charForDigit[(c%=100000)/10000]); - s.append(charForDigit[(c%=10000)/1000]); - s.append(charForDigit[(c%=1000)/100]); - s.append(charForDigit[(c%=100)/10]); - s.append(charForDigit[c%10]); - } - else if (i < 1000000000) - { - //nine digits - s.append(charForDigit[i/100000000]); - s.append(charForDigit[(c=i%100000000)/10000000]); - s.append(charForDigit[(c%=10000000)/1000000]); - s.append(charForDigit[(c%=1000000)/100000]); - s.append(charForDigit[(c%=100000)/10000]); - s.append(charForDigit[(c%=10000)/1000]); - s.append(charForDigit[(c%=1000)/100]); - s.append(charForDigit[(c%=100)/10]); - s.append(charForDigit[c%10]); - } - else - { - //ten digits - s.append(charForDigit[i/1000000000]); - s.append(charForDigit[(c=i%1000000000)/100000000]); - s.append(charForDigit[(c%=100000000)/10000000]); - s.append(charForDigit[(c%=10000000)/1000000]); - s.append(charForDigit[(c%=1000000)/100000]); - s.append(charForDigit[(c%=100000)/10000]); - s.append(charForDigit[(c%=10000)/1000]); - s.append(charForDigit[(c%=1000)/100]); - s.append(charForDigit[(c%=100)/10]); - s.append(charForDigit[c%10]); - } -} -private static void appendFractDigits(StringBuffer s, long i, int decimalOffset) -{ - long mag = tenthPower(i); - long c; - while ( i > 0 ) - { - c = i/mag; - s.append(charForDigit[(int) c]); - decimalOffset--; - if (decimalOffset == 0) - s.append('.'); - c *= mag; - if ( c <= i) - i -= c; - mag = mag/10; - } - if (i != 0) - s.append(charForDigit[(int) i]); - else if (decimalOffset > 0) - { - s.append(ZEROS[decimalOffset]); - decimalOffset = 1; - } - - decimalOffset--; - if (decimalOffset == 0) - s.append(DOT_ZERO); - else if (decimalOffset == -1) - s.append('0'); -} - - public static final char[] NEGATIVE_INFINITY = {'-','I','n','f','i','n','i','t','y'}; - public static final char[] POSITIVE_INFINITY = {'I','n','f','i','n','i','t','y'}; - public static final char[] DOUBLE_ZERO = {'0','.','0'}; - public static final char[] DOUBLE_ZERO2 = {'0','.','0','0'}; - public static final char[] DOUBLE_ZERO0 = {'0','.'}; - public static final char[] DOT_ZERO = {'.','0'}; - -} \ No newline at end of file diff --git a/classpath/java/lang/VMDouble.java b/classpath/java/lang/VMDouble.java index 1f0643c3..e87ac842 100644 --- a/classpath/java/lang/VMDouble.java +++ b/classpath/java/lang/VMDouble.java @@ -1,5 +1,5 @@ /* - Copyright (C) 2003, 2005, 2006 Jeroen Frijters + Copyright (C) 2003, 2005, 2006, 2007 Jeroen Frijters This software is provided 'as-is', without any express or implied warranty. In no event will the authors be held liable for any damages @@ -25,8 +25,7 @@ package java.lang; import cli.System.BitConverter; -import cli.System.Globalization.CultureInfo; -import ikvm.lang.CIL; +import sun.misc.FloatingDecimal; final class VMDouble { @@ -51,55 +50,12 @@ final class VMDouble static String toString(double d, boolean isFloat) { - if(isFloat) - { - float f = (float)d; - if(Float.isNaN(f)) - { - return "NaN"; - } - if(Float.isInfinite(f)) - { - return f < 0f ? "-Infinity" : "Infinity"; - } - if(f == 0f) - { - return BitConverter.DoubleToInt64Bits(d) < 0 ? "-0.0" : "0.0"; - } - // TODO this is not correct, we need to use the Java algorithm of converting a float to string - // HACK really lame hack to approximate the Java behavior a little bit - String s = CIL.box_float(f).ToString(CultureInfo.get_InvariantCulture()); - if(s.indexOf('.') == -1) - { - int e = s.indexOf('E'); - if(e == -1) - { - s += ".0"; - } - else - { - int plus = s.charAt(e + 1) == '+' ? 1 : 0; - s = s.substring(0, e) + ".0E" + Integer.parseInt(s.substring(e + 1 + plus)); - } - } - else - { - int e = s.indexOf('E'); - if(e != -1) - { - int plus = s.charAt(e + 1) == '+' ? 1 : 0; - s = s.substring(0, e) + "E" + Integer.parseInt(s.substring(e + 1 + plus)); - } - } - return s; - } - else - { - StringBuffer sb = new StringBuffer(); - DoubleToString.append(sb, d); - return sb.toString(); - } + assert !isFloat; + return new FloatingDecimal(d).toJavaFormatString(); } - static native double parseDouble(String s); + static double parseDouble(String s) + { + return FloatingDecimal.readJavaFormatString(s).doubleValue(); + } } diff --git a/classpath/java/lang/VMFloat.java b/classpath/java/lang/VMFloat.java index d1f47d53..2fc868d2 100644 --- a/classpath/java/lang/VMFloat.java +++ b/classpath/java/lang/VMFloat.java @@ -1,5 +1,5 @@ /* - Copyright (C) 2003 Jeroen Frijters + Copyright (C) 2003, 2007 Jeroen Frijters This software is provided 'as-is', without any express or implied warranty. In no event will the authors be held liable for any damages @@ -24,6 +24,8 @@ package java.lang; +import sun.misc.FloatingDecimal; + final class VMFloat { static float intBitsToFloat(int v) @@ -44,4 +46,14 @@ final class VMFloat { return cli.System.BitConverter.ToInt32(cli.System.BitConverter.GetBytes(v), 0); } + + static String toString(float f) + { + return new FloatingDecimal(f).toJavaFormatString(); + } + + static float parseFloat(String str) + { + return FloatingDecimal.readJavaFormatString(str).floatValue(); + } } diff --git a/classpath/sun/misc/DoubleConsts.java b/classpath/sun/misc/DoubleConsts.java new file mode 100644 index 00000000..1970e933 --- /dev/null +++ b/classpath/sun/misc/DoubleConsts.java @@ -0,0 +1,118 @@ +/* + * Copyright 2003 Sun Microsystems, Inc. All Rights Reserved. + * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER. + * + * This code is free software; you can redistribute it and/or modify it + * under the terms of the GNU General Public License version 2 only, as + * published by the Free Software Foundation. Sun designates this + * particular file as subject to the "Classpath" exception as provided + * by Sun in the LICENSE file that accompanied this code. + * + * This code is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * version 2 for more details (a copy is included in the LICENSE file that + * accompanied this code). + * + * You should have received a copy of the GNU General Public License version + * 2 along with this work; if not, write to the Free Software Foundation, + * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. + * + * Please contact Sun Microsystems, Inc., 4150 Network Circle, Santa Clara, + * CA 95054 USA or visit www.sun.com if you need additional information or + * have any questions. + */ + +package sun.misc; + +/** + * This class contains additional constants documenting limits of the + * double type. + * + * @author Joseph D. Darcy + * @version 1.7, 05/05/07 + */ + +public class DoubleConsts { + /** + * Don't let anyone instantiate this class. + */ + private DoubleConsts() {} + + public static final double POSITIVE_INFINITY = java.lang.Double.POSITIVE_INFINITY; + public static final double NEGATIVE_INFINITY = java.lang.Double.NEGATIVE_INFINITY; + public static final double NaN = java.lang.Double.NaN; + public static final double MAX_VALUE = java.lang.Double.MAX_VALUE; + public static final double MIN_VALUE = java.lang.Double.MIN_VALUE; + + /** + * A constant holding the smallest positive normal value of type + * double, 2-1022. It is equal to the + * value returned by + * Double.longBitsToDouble(0x0010000000000000L). + * + * @since 1.5 + */ + public static final double MIN_NORMAL = 2.2250738585072014E-308; + + + /** + * The number of logical bits in the significand of a + * double number, including the implicit bit. + */ + public static final int SIGNIFICAND_WIDTH = 53; + + /** + * Maximum exponent a finite double number may have. + * It is equal to the value returned by + * Math.ilogb(Double.MAX_VALUE). + */ + public static final int MAX_EXPONENT = 1023; + + /** + * Minimum exponent a normalized double number may + * have. It is equal to the value returned by + * Math.ilogb(Double.MIN_NORMAL). + */ + public static final int MIN_EXPONENT = -1022; + + /** + * The exponent the smallest positive double + * subnormal value would have if it could be normalized. It is + * equal to the value returned by + * FpUtils.ilogb(Double.MIN_VALUE). + */ + public static final int MIN_SUB_EXPONENT = MIN_EXPONENT - + (SIGNIFICAND_WIDTH - 1); + + /** + * Bias used in representing a double exponent. + */ + public static final int EXP_BIAS = 1023; + + /** + * Bit mask to isolate the sign bit of a double. + */ + public static final long SIGN_BIT_MASK = 0x8000000000000000L; + + /** + * Bit mask to isolate the exponent field of a + * double. + */ + public static final long EXP_BIT_MASK = 0x7FF0000000000000L; + + /** + * Bit mask to isolate the significand field of a + * double. + */ + public static final long SIGNIF_BIT_MASK = 0x000FFFFFFFFFFFFFL; + + static { + // verify bit masks cover all bit positions and that the bit + // masks are non-overlapping + assert(((SIGN_BIT_MASK | EXP_BIT_MASK | SIGNIF_BIT_MASK) == ~0L) && + (((SIGN_BIT_MASK & EXP_BIT_MASK) == 0L) && + ((SIGN_BIT_MASK & SIGNIF_BIT_MASK) == 0L) && + ((EXP_BIT_MASK & SIGNIF_BIT_MASK) == 0L))); + } +} diff --git a/classpath/sun/misc/FloatConsts.java b/classpath/sun/misc/FloatConsts.java new file mode 100644 index 00000000..f97f0962 --- /dev/null +++ b/classpath/sun/misc/FloatConsts.java @@ -0,0 +1,113 @@ +/* + * Copyright 2003 Sun Microsystems, Inc. All Rights Reserved. + * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER. + * + * This code is free software; you can redistribute it and/or modify it + * under the terms of the GNU General Public License version 2 only, as + * published by the Free Software Foundation. Sun designates this + * particular file as subject to the "Classpath" exception as provided + * by Sun in the LICENSE file that accompanied this code. + * + * This code is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * version 2 for more details (a copy is included in the LICENSE file that + * accompanied this code). + * + * You should have received a copy of the GNU General Public License version + * 2 along with this work; if not, write to the Free Software Foundation, + * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. + * + * Please contact Sun Microsystems, Inc., 4150 Network Circle, Santa Clara, + * CA 95054 USA or visit www.sun.com if you need additional information or + * have any questions. + */ + +package sun.misc; + +/** + * This class contains additional constants documenting limits of the + * float type. + * + * @author Joseph D. Darcy + * @version 1.7, 05/05/07 + */ + +public class FloatConsts { + /** + * Don't let anyone instantiate this class. + */ + private FloatConsts() {} + + public static final float POSITIVE_INFINITY = java.lang.Float.POSITIVE_INFINITY; + public static final float NEGATIVE_INFINITY = java.lang.Float.NEGATIVE_INFINITY; + public static final float NaN = java.lang.Float.NaN; + public static final float MAX_VALUE = java.lang.Float.MAX_VALUE; + public static final float MIN_VALUE = java.lang.Float.MIN_VALUE; + + /** + * A constant holding the smallest positive normal value of type + * float, 2-126. It is equal to the value + * returned by Float.intBitsToFloat(0x00800000). + */ + public static final float MIN_NORMAL = 1.17549435E-38f; + + /** + * The number of logical bits in the significand of a + * float number, including the implicit bit. + */ + public static final int SIGNIFICAND_WIDTH = 24; + + /** + * Maximum exponent a finite float number may have. + * It is equal to the value returned by + * Math.ilogb(Float.MAX_VALUE). + */ + public static final int MAX_EXPONENT = 127; + + /** + * Minimum exponent a normalized float number may + * have. It is equal to the value returned by + * Math.ilogb(Float.MIN_NORMAL). + */ + public static final int MIN_EXPONENT = -126; + + /** + * The exponent the smallest positive float subnormal + * value would have if it could be normalized. It is equal to the + * value returned by FpUtils.ilogb(Float.MIN_VALUE). + */ + public static final int MIN_SUB_EXPONENT = MIN_EXPONENT - + (SIGNIFICAND_WIDTH - 1); + + /** + * Bias used in representing a float exponent. + */ + public static final int EXP_BIAS = 127; + + /** + * Bit mask to isolate the sign bit of a float. + */ + public static final int SIGN_BIT_MASK = 0x80000000; + + /** + * Bit mask to isolate the exponent field of a + * float. + */ + public static final int EXP_BIT_MASK = 0x7F800000; + + /** + * Bit mask to isolate the significand field of a + * float. + */ + public static final int SIGNIF_BIT_MASK = 0x007FFFFF; + + static { + // verify bit masks cover all bit positions and that the bit + // masks are non-overlapping + assert(((SIGN_BIT_MASK | EXP_BIT_MASK | SIGNIF_BIT_MASK) == ~0) && + (((SIGN_BIT_MASK & EXP_BIT_MASK) == 0) && + ((SIGN_BIT_MASK & SIGNIF_BIT_MASK) == 0) && + ((EXP_BIT_MASK & SIGNIF_BIT_MASK) == 0))); + } +} diff --git a/classpath/sun/misc/FloatingDecimal.java b/classpath/sun/misc/FloatingDecimal.java new file mode 100644 index 00000000..e9d92b89 --- /dev/null +++ b/classpath/sun/misc/FloatingDecimal.java @@ -0,0 +1,2885 @@ +/* + * Copyright 1996-2004 Sun Microsystems, Inc. All Rights Reserved. + * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER. + * + * This code is free software; you can redistribute it and/or modify it + * under the terms of the GNU General Public License version 2 only, as + * published by the Free Software Foundation. Sun designates this + * particular file as subject to the "Classpath" exception as provided + * by Sun in the LICENSE file that accompanied this code. + * + * This code is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * version 2 for more details (a copy is included in the LICENSE file that + * accompanied this code). + * + * You should have received a copy of the GNU General Public License version + * 2 along with this work; if not, write to the Free Software Foundation, + * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. + * + * Please contact Sun Microsystems, Inc., 4150 Network Circle, Santa Clara, + * CA 95054 USA or visit www.sun.com if you need additional information or + * have any questions. + */ + +/*IKVM*/ +// This file was modified for IKVM by Jeroen Frijters. + +package sun.misc; + +import sun.misc.FpUtils; +import sun.misc.DoubleConsts; +import sun.misc.FloatConsts; +import java.util.regex.*; + +public class FloatingDecimal{ + boolean isExceptional; + boolean isNegative; + int decExponent; + char digits[]; + int nDigits; + int bigIntExp; + int bigIntNBits; + boolean mustSetRoundDir = false; + boolean fromHex = false; + int roundDir = 0; // set by doubleValue + + private FloatingDecimal( boolean negSign, int decExponent, char []digits, int n, boolean e ) + { + isNegative = negSign; + isExceptional = e; + this.decExponent = decExponent; + this.digits = digits; + this.nDigits = n; + } + + /* + * Constants of the implementation + * Most are IEEE-754 related. + * (There are more really boring constants at the end.) + */ + static final long signMask = 0x8000000000000000L; + static final long expMask = 0x7ff0000000000000L; + static final long fractMask= ~(signMask|expMask); + static final int expShift = 52; + static final int expBias = 1023; + static final long fractHOB = ( 1L< 0L ) { // i.e. while ((v&highbit) == 0L ) + v <<= 1; + } + + int n = 0; + while (( v & lowbytes ) != 0L ){ + v <<= 8; + n += 8; + } + while ( v != 0L ){ + v <<= 1; + n += 1; + } + return n; + } + + /* + * Keep big powers of 5 handy for future reference. + */ + private static FDBigInt b5p[]; + + private static synchronized FDBigInt + big5pow( int p ){ + assert p >= 0 : p; // negative power of 5 + if ( b5p == null ){ + b5p = new FDBigInt[ p+1 ]; + }else if (b5p.length <= p ){ + FDBigInt t[] = new FDBigInt[ p+1 ]; + System.arraycopy( b5p, 0, t, 0, b5p.length ); + b5p = t; + } + if ( b5p[p] != null ) + return b5p[p]; + else if ( p < small5pow.length ) + return b5p[p] = new FDBigInt( small5pow[p] ); + else if ( p < long5pow.length ) + return b5p[p] = new FDBigInt( long5pow[p] ); + else { + // construct the value. + // recursively. + int q, r; + // in order to compute 5^p, + // compute its square root, 5^(p/2) and square. + // or, let q = p / 2, r = p -q, then + // 5^p = 5^(q+r) = 5^q * 5^r + q = p >> 1; + r = p - q; + FDBigInt bigq = b5p[q]; + if ( bigq == null ) + bigq = big5pow ( q ); + if ( r < small5pow.length ){ + return (b5p[p] = bigq.mult( small5pow[r] ) ); + }else{ + FDBigInt bigr = b5p[ r ]; + if ( bigr == null ) + bigr = big5pow( r ); + return (b5p[p] = bigq.mult( bigr ) ); + } + } + } + + // + // a common operation + // + private static FDBigInt + multPow52( FDBigInt v, int p5, int p2 ){ + if ( p5 != 0 ){ + if ( p5 < small5pow.length ){ + v = v.mult( small5pow[p5] ); + } else { + v = v.mult( big5pow( p5 ) ); + } + } + if ( p2 != 0 ){ + v.lshiftMe( p2 ); + } + return v; + } + + // + // another common operation + // + private static FDBigInt + constructPow52( int p5, int p2 ){ + FDBigInt v = new FDBigInt( big5pow( p5 ) ); + if ( p2 != 0 ){ + v.lshiftMe( p2 ); + } + return v; + } + + /* + * Make a floating double into a FDBigInt. + * This could also be structured as a FDBigInt + * constructor, but we'd have to build a lot of knowledge + * about floating-point representation into it, and we don't want to. + * + * AS A SIDE EFFECT, THIS METHOD WILL SET THE INSTANCE VARIABLES + * bigIntExp and bigIntNBits + * + */ + private FDBigInt + doubleToBigInt( double dval ){ + long lbits = Double.doubleToLongBits( dval ) & ~signMask; + int binexp = (int)(lbits >>> expShift); + lbits &= fractMask; + if ( binexp > 0 ){ + lbits |= fractHOB; + } else { + assert lbits != 0L : lbits; // doubleToBigInt(0.0) + binexp +=1; + while ( (lbits & fractHOB ) == 0L){ + lbits <<= 1; + binexp -= 1; + } + } + binexp -= expBias; + int nbits = countBits( lbits ); + /* + * We now know where the high-order 1 bit is, + * and we know how many there are. + */ + int lowOrderZeros = expShift+1-nbits; + lbits >>>= lowOrderZeros; + + bigIntExp = binexp+1-nbits; + bigIntNBits = nbits; + return new FDBigInt( lbits ); + } + + /* + * Compute a number that is the ULP of the given value, + * for purposes of addition/subtraction. Generally easy. + * More difficult if subtracting and the argument + * is a normalized a power of 2, as the ULP changes at these points. + */ + private static double + ulp( double dval, boolean subtracting ){ + long lbits = Double.doubleToLongBits( dval ) & ~signMask; + int binexp = (int)(lbits >>> expShift); + double ulpval; + if ( subtracting && ( binexp >= expShift ) && ((lbits&fractMask) == 0L) ){ + // for subtraction from normalized, powers of 2, + // use next-smaller exponent + binexp -= 1; + } + if ( binexp > expShift ){ + ulpval = Double.longBitsToDouble( ((long)(binexp-expShift))< 0L (not zero, nor negative). + * + * The only reason that we develop the digits here, rather than + * calling on Long.toString() is that we can do it a little faster, + * and besides want to treat trailing 0s specially. If Long.toString + * changes, we should re-evaluate this strategy! + */ + private void + developLongDigits( int decExponent, long lvalue, long insignificant ){ + char digits[]; + int ndigits; + int digitno; + int c; + // + // Discard non-significant low-order bits, while rounding, + // up to insignificant value. + int i; + for ( i = 0; insignificant >= 10L; i++ ) + insignificant /= 10L; + if ( i != 0 ){ + long pow10 = long5pow[i] << i; // 10^i == 5^i * 2^i; + long residue = lvalue % pow10; + lvalue /= pow10; + decExponent += i; + if ( residue >= (pow10>>1) ){ + // round up based on the low-order bits we're discarding + lvalue++; + } + } + if ( lvalue <= Integer.MAX_VALUE ){ + assert lvalue > 0L : lvalue; // lvalue <= 0 + // even easier subcase! + // can do int arithmetic rather than long! + int ivalue = (int)lvalue; + ndigits = 10; + digits = (char[])(perThreadBuffer.get()); + digitno = ndigits-1; + c = ivalue%10; + ivalue /= 10; + while ( c == 0 ){ + decExponent++; + c = ivalue%10; + ivalue /= 10; + } + while ( ivalue != 0){ + digits[digitno--] = (char)(c+'0'); + decExponent++; + c = ivalue%10; + ivalue /= 10; + } + digits[digitno] = (char)(c+'0'); + } else { + // same algorithm as above (same bugs, too ) + // but using long arithmetic. + ndigits = 20; + digits = (char[])(perThreadBuffer.get()); + digitno = ndigits-1; + c = (int)(lvalue%10L); + lvalue /= 10L; + while ( c == 0 ){ + decExponent++; + c = (int)(lvalue%10L); + lvalue /= 10L; + } + while ( lvalue != 0L ){ + digits[digitno--] = (char)(c+'0'); + decExponent++; + c = (int)(lvalue%10L); + lvalue /= 10; + } + digits[digitno] = (char)(c+'0'); + } + char result []; + ndigits -= digitno; + result = new char[ ndigits ]; + System.arraycopy( digits, digitno, result, 0, ndigits ); + this.digits = result; + this.decExponent = decExponent+1; + this.nDigits = ndigits; + } + + // + // add one to the least significant digit. + // in the unlikely event there is a carry out, + // deal with it. + // assert that this will only happen where there + // is only one digit, e.g. (float)1e-44 seems to do it. + // + private void + roundup(){ + int i; + int q = digits[ i = (nDigits-1)]; + if ( q == '9' ){ + while ( q == '9' && i > 0 ){ + digits[i] = '0'; + q = digits[--i]; + } + if ( q == '9' ){ + // carryout! High-order 1, rest 0s, larger exp. + decExponent += 1; + digits[0] = '1'; + return; + } + // else fall through. + } + digits[i] = (char)(q+1); + } + + /* + * FIRST IMPORTANT CONSTRUCTOR: DOUBLE + */ + public FloatingDecimal( double d ) + { + long dBits = Double.doubleToLongBits( d ); + long fractBits; + int binExp; + int nSignificantBits; + + // discover and delete sign + if ( (dBits&signMask) != 0 ){ + isNegative = true; + dBits ^= signMask; + } else { + isNegative = false; + } + // Begin to unpack + // Discover obvious special cases of NaN and Infinity. + binExp = (int)( (dBits&expMask) >> expShift ); + fractBits = dBits&fractMask; + if ( binExp == (int)(expMask>>expShift) ) { + isExceptional = true; + if ( fractBits == 0L ){ + digits = infinity; + } else { + digits = notANumber; + isNegative = false; // NaN has no sign! + } + nDigits = digits.length; + return; + } + isExceptional = false; + // Finish unpacking + // Normalize denormalized numbers. + // Insert assumed high-order bit for normalized numbers. + // Subtract exponent bias. + if ( binExp == 0 ){ + if ( fractBits == 0L ){ + // not a denorm, just a 0! + decExponent = 0; + digits = zero; + nDigits = 1; + return; + } + while ( (fractBits&fractHOB) == 0L ){ + fractBits <<= 1; + binExp -= 1; + } + nSignificantBits = expShift + binExp +1; // recall binExp is - shift count. + binExp += 1; + } else { + fractBits |= fractHOB; + nSignificantBits = expShift+1; + } + binExp -= expBias; + // call the routine that actually does all the hard work. + dtoa( binExp, fractBits, nSignificantBits ); + } + + /* + * SECOND IMPORTANT CONSTRUCTOR: SINGLE + */ + public FloatingDecimal( float f ) + { + int fBits = Float.floatToIntBits( f ); + int fractBits; + int binExp; + int nSignificantBits; + + // discover and delete sign + if ( (fBits&singleSignMask) != 0 ){ + isNegative = true; + fBits ^= singleSignMask; + } else { + isNegative = false; + } + // Begin to unpack + // Discover obvious special cases of NaN and Infinity. + binExp = (int)( (fBits&singleExpMask) >> singleExpShift ); + fractBits = fBits&singleFractMask; + if ( binExp == (int)(singleExpMask>>singleExpShift) ) { + isExceptional = true; + if ( fractBits == 0L ){ + digits = infinity; + } else { + digits = notANumber; + isNegative = false; // NaN has no sign! + } + nDigits = digits.length; + return; + } + isExceptional = false; + // Finish unpacking + // Normalize denormalized numbers. + // Insert assumed high-order bit for normalized numbers. + // Subtract exponent bias. + if ( binExp == 0 ){ + if ( fractBits == 0 ){ + // not a denorm, just a 0! + decExponent = 0; + digits = zero; + nDigits = 1; + return; + } + while ( (fractBits&singleFractHOB) == 0 ){ + fractBits <<= 1; + binExp -= 1; + } + nSignificantBits = singleExpShift + binExp +1; // recall binExp is - shift count. + binExp += 1; + } else { + fractBits |= singleFractHOB; + nSignificantBits = singleExpShift+1; + } + binExp -= singleExpBias; + // call the routine that actually does all the hard work. + dtoa( binExp, ((long)fractBits)<<(expShift-singleExpShift), nSignificantBits ); + } + + private void + dtoa( int binExp, long fractBits, int nSignificantBits ) + { + int nFractBits; // number of significant bits of fractBits; + int nTinyBits; // number of these to the right of the point. + int decExp; + + // Examine number. Determine if it is an easy case, + // which we can do pretty trivially using float/long conversion, + // or whether we must do real work. + nFractBits = countBits( fractBits ); + nTinyBits = Math.max( 0, nFractBits - binExp - 1 ); + if ( binExp <= maxSmallBinExp && binExp >= minSmallBinExp ){ + // Look more closely at the number to decide if, + // with scaling by 10^nTinyBits, the result will fit in + // a long. + if ( (nTinyBits < long5pow.length) && ((nFractBits + n5bits[nTinyBits]) < 64 ) ){ + /* + * We can do this: + * take the fraction bits, which are normalized. + * (a) nTinyBits == 0: Shift left or right appropriately + * to align the binary point at the extreme right, i.e. + * where a long int point is expected to be. The integer + * result is easily converted to a string. + * (b) nTinyBits > 0: Shift right by expShift-nFractBits, + * which effectively converts to long and scales by + * 2^nTinyBits. Then multiply by 5^nTinyBits to + * complete the scaling. We know this won't overflow + * because we just counted the number of bits necessary + * in the result. The integer you get from this can + * then be converted to a string pretty easily. + */ + long halfULP; + if ( nTinyBits == 0 ) { + if ( binExp > nSignificantBits ){ + halfULP = 1L << ( binExp-nSignificantBits-1); + } else { + halfULP = 0L; + } + if ( binExp >= expShift ){ + fractBits <<= (binExp-expShift); + } else { + fractBits >>>= (expShift-binExp) ; + } + developLongDigits( 0, fractBits, halfULP ); + return; + } + /* + * The following causes excess digits to be printed + * out in the single-float case. Our manipulation of + * halfULP here is apparently not correct. If we + * better understand how this works, perhaps we can + * use this special case again. But for the time being, + * we do not. + * else { + * fractBits >>>= expShift+1-nFractBits; + * fractBits *= long5pow[ nTinyBits ]; + * halfULP = long5pow[ nTinyBits ] >> (1+nSignificantBits-nFractBits); + * developLongDigits( -nTinyBits, fractBits, halfULP ); + * return; + * } + */ + } + } + /* + * This is the hard case. We are going to compute large positive + * integers B and S and integer decExp, s.t. + * d = ( B / S ) * 10^decExp + * 1 <= B / S < 10 + * Obvious choices are: + * decExp = floor( log10(d) ) + * B = d * 2^nTinyBits * 10^max( 0, -decExp ) + * S = 10^max( 0, decExp) * 2^nTinyBits + * (noting that nTinyBits has already been forced to non-negative) + * I am also going to compute a large positive integer + * M = (1/2^nSignificantBits) * 2^nTinyBits * 10^max( 0, -decExp ) + * i.e. M is (1/2) of the ULP of d, scaled like B. + * When we iterate through dividing B/S and picking off the + * quotient bits, we will know when to stop when the remainder + * is <= M. + * + * We keep track of powers of 2 and powers of 5. + */ + + /* + * Estimate decimal exponent. (If it is small-ish, + * we could double-check.) + * + * First, scale the mantissa bits such that 1 <= d2 < 2. + * We are then going to estimate + * log10(d2) ~=~ (d2-1.5)/1.5 + log(1.5) + * and so we can estimate + * log10(d) ~=~ log10(d2) + binExp * log10(2) + * take the floor and call it decExp. + * FIXME -- use more precise constants here. It costs no more. + */ + double d2 = Double.longBitsToDouble( + expOne | ( fractBits &~ fractHOB ) ); + decExp = (int)Math.floor( + (d2-1.5D)*0.289529654D + 0.176091259 + (double)binExp * 0.301029995663981 ); + int B2, B5; // powers of 2 and powers of 5, respectively, in B + int S2, S5; // powers of 2 and powers of 5, respectively, in S + int M2, M5; // powers of 2 and powers of 5, respectively, in M + int Bbits; // binary digits needed to represent B, approx. + int tenSbits; // binary digits needed to represent 10*S, approx. + FDBigInt Sval, Bval, Mval; + + B5 = Math.max( 0, -decExp ); + B2 = B5 + nTinyBits + binExp; + + S5 = Math.max( 0, decExp ); + S2 = S5 + nTinyBits; + + M5 = B5; + M2 = B2 - nSignificantBits; + + /* + * the long integer fractBits contains the (nFractBits) interesting + * bits from the mantissa of d ( hidden 1 added if necessary) followed + * by (expShift+1-nFractBits) zeros. In the interest of compactness, + * I will shift out those zeros before turning fractBits into a + * FDBigInt. The resulting whole number will be + * d * 2^(nFractBits-1-binExp). + */ + fractBits >>>= (expShift+1-nFractBits); + B2 -= nFractBits-1; + int common2factor = Math.min( B2, S2 ); + B2 -= common2factor; + S2 -= common2factor; + M2 -= common2factor; + + /* + * HACK!! For exact powers of two, the next smallest number + * is only half as far away as we think (because the meaning of + * ULP changes at power-of-two bounds) for this reason, we + * hack M2. Hope this works. + */ + if ( nFractBits == 1 ) + M2 -= 1; + + if ( M2 < 0 ){ + // oops. + // since we cannot scale M down far enough, + // we must scale the other values up. + B2 -= M2; + S2 -= M2; + M2 = 0; + } + /* + * Construct, Scale, iterate. + * Some day, we'll write a stopping test that takes + * account of the asymmetry of the spacing of floating-point + * numbers below perfect powers of 2 + * 26 Sept 96 is not that day. + * So we use a symmetric test. + */ + char digits[] = this.digits = new char[18]; + int ndigit = 0; + boolean low, high; + long lowDigitDifference; + int q; + + /* + * Detect the special cases where all the numbers we are about + * to compute will fit in int or long integers. + * In these cases, we will avoid doing FDBigInt arithmetic. + * We use the same algorithms, except that we "normalize" + * our FDBigInts before iterating. This is to make division easier, + * as it makes our fist guess (quotient of high-order words) + * more accurate! + * + * Some day, we'll write a stopping test that takes + * account of the asymmetry of the spacing of floating-point + * numbers below perfect powers of 2 + * 26 Sept 96 is not that day. + * So we use a symmetric test. + */ + Bbits = nFractBits + B2 + (( B5 < n5bits.length )? n5bits[B5] : ( B5*3 )); + tenSbits = S2+1 + (( (S5+1) < n5bits.length )? n5bits[(S5+1)] : ( (S5+1)*3 )); + if ( Bbits < 64 && tenSbits < 64){ + if ( Bbits < 32 && tenSbits < 32){ + // wa-hoo! They're all ints! + int b = ((int)fractBits * small5pow[B5] ) << B2; + int s = small5pow[S5] << S2; + int m = small5pow[M5] << M2; + int tens = s * 10; + /* + * Unroll the first iteration. If our decExp estimate + * was too high, our first quotient will be zero. In this + * case, we discard it and decrement decExp. + */ + ndigit = 0; + q = (int) ( b / s ); + b = 10 * ( b % s ); + m *= 10; + low = (b < m ); + high = (b+m > tens ); + assert q < 10 : q; // excessively large digit + if ( (q == 0) && ! high ){ + // oops. Usually ignore leading zero. + decExp--; + } else { + digits[ndigit++] = (char)('0' + q); + } + /* + * HACK! Java spec sez that we always have at least + * one digit after the . in either F- or E-form output. + * Thus we will need more than one digit if we're using + * E-form + */ + if ( decExp <= -3 || decExp >= 8 ){ + high = low = false; + } + while( ! low && ! high ){ + q = (int) ( b / s ); + b = 10 * ( b % s ); + m *= 10; + assert q < 10 : q; // excessively large digit + if ( m > 0L ){ + low = (b < m ); + high = (b+m > tens ); + } else { + // hack -- m might overflow! + // in this case, it is certainly > b, + // which won't + // and b+m > tens, too, since that has overflowed + // either! + low = true; + high = true; + } + digits[ndigit++] = (char)('0' + q); + } + lowDigitDifference = (b<<1) - tens; + } else { + // still good! they're all longs! + long b = (fractBits * long5pow[B5] ) << B2; + long s = long5pow[S5] << S2; + long m = long5pow[M5] << M2; + long tens = s * 10L; + /* + * Unroll the first iteration. If our decExp estimate + * was too high, our first quotient will be zero. In this + * case, we discard it and decrement decExp. + */ + ndigit = 0; + q = (int) ( b / s ); + b = 10L * ( b % s ); + m *= 10L; + low = (b < m ); + high = (b+m > tens ); + assert q < 10 : q; // excessively large digit + if ( (q == 0) && ! high ){ + // oops. Usually ignore leading zero. + decExp--; + } else { + digits[ndigit++] = (char)('0' + q); + } + /* + * HACK! Java spec sez that we always have at least + * one digit after the . in either F- or E-form output. + * Thus we will need more than one digit if we're using + * E-form + */ + if ( decExp <= -3 || decExp >= 8 ){ + high = low = false; + } + while( ! low && ! high ){ + q = (int) ( b / s ); + b = 10 * ( b % s ); + m *= 10; + assert q < 10 : q; // excessively large digit + if ( m > 0L ){ + low = (b < m ); + high = (b+m > tens ); + } else { + // hack -- m might overflow! + // in this case, it is certainly > b, + // which won't + // and b+m > tens, too, since that has overflowed + // either! + low = true; + high = true; + } + digits[ndigit++] = (char)('0' + q); + } + lowDigitDifference = (b<<1) - tens; + } + } else { + FDBigInt tenSval; + int shiftBias; + + /* + * We really must do FDBigInt arithmetic. + * Fist, construct our FDBigInt initial values. + */ + Bval = multPow52( new FDBigInt( fractBits ), B5, B2 ); + Sval = constructPow52( S5, S2 ); + Mval = constructPow52( M5, M2 ); + + + // normalize so that division works better + Bval.lshiftMe( shiftBias = Sval.normalizeMe() ); + Mval.lshiftMe( shiftBias ); + tenSval = Sval.mult( 10 ); + /* + * Unroll the first iteration. If our decExp estimate + * was too high, our first quotient will be zero. In this + * case, we discard it and decrement decExp. + */ + ndigit = 0; + q = Bval.quoRemIteration( Sval ); + Mval = Mval.mult( 10 ); + low = (Bval.cmp( Mval ) < 0); + high = (Bval.add( Mval ).cmp( tenSval ) > 0 ); + assert q < 10 : q; // excessively large digit + if ( (q == 0) && ! high ){ + // oops. Usually ignore leading zero. + decExp--; + } else { + digits[ndigit++] = (char)('0' + q); + } + /* + * HACK! Java spec sez that we always have at least + * one digit after the . in either F- or E-form output. + * Thus we will need more than one digit if we're using + * E-form + */ + if ( decExp <= -3 || decExp >= 8 ){ + high = low = false; + } + while( ! low && ! high ){ + q = Bval.quoRemIteration( Sval ); + Mval = Mval.mult( 10 ); + assert q < 10 : q; // excessively large digit + low = (Bval.cmp( Mval ) < 0); + high = (Bval.add( Mval ).cmp( tenSval ) > 0 ); + digits[ndigit++] = (char)('0' + q); + } + if ( high && low ){ + Bval.lshiftMe(1); + lowDigitDifference = Bval.cmp(tenSval); + } else + lowDigitDifference = 0L; // this here only for flow analysis! + } + this.decExponent = decExp+1; + this.digits = digits; + this.nDigits = ndigit; + /* + * Last digit gets rounded based on stopping condition. + */ + if ( high ){ + if ( low ){ + if ( lowDigitDifference == 0L ){ + // it's a tie! + // choose based on which digits we like. + if ( (digits[nDigits-1]&1) != 0 ) roundup(); + } else if ( lowDigitDifference > 0 ){ + roundup(); + } + } else { + roundup(); + } + } + } + + public String + toString(){ + // most brain-dead version + StringBuffer result = new StringBuffer( nDigits+8 ); + if ( isNegative ){ result.append( '-' ); } + if ( isExceptional ){ + result.append( digits, 0, nDigits ); + } else { + result.append( "0."); + result.append( digits, 0, nDigits ); + result.append('e'); + result.append( decExponent ); + } + return new String(result); + } + + public String toJavaFormatString() { + char result[] = (char[])(perThreadBuffer.get()); + int i = getChars(result); + return new String(result, 0, i); + } + + private int getChars(char[] result) { + assert nDigits <= 19 : nDigits; // generous bound on size of nDigits + int i = 0; + if (isNegative) { result[0] = '-'; i = 1; } + if (isExceptional) { + System.arraycopy(digits, 0, result, i, nDigits); + i += nDigits; + } else { + if (decExponent > 0 && decExponent < 8) { + // print digits.digits. + int charLength = Math.min(nDigits, decExponent); + System.arraycopy(digits, 0, result, i, charLength); + i += charLength; + if (charLength < decExponent) { + charLength = decExponent-charLength; + System.arraycopy(zero, 0, result, i, charLength); + i += charLength; + result[i++] = '.'; + result[i++] = '0'; + } else { + result[i++] = '.'; + if (charLength < nDigits) { + int t = nDigits - charLength; + System.arraycopy(digits, charLength, result, i, t); + i += t; + } else { + result[i++] = '0'; + } + } + } else if (decExponent <=0 && decExponent > -3) { + result[i++] = '0'; + result[i++] = '.'; + if (decExponent != 0) { + System.arraycopy(zero, 0, result, i, -decExponent); + i -= decExponent; + } + System.arraycopy(digits, 0, result, i, nDigits); + i += nDigits; + } else { + result[i++] = digits[0]; + result[i++] = '.'; + if (nDigits > 1) { + System.arraycopy(digits, 1, result, i, nDigits-1); + i += nDigits-1; + } else { + result[i++] = '0'; + } + result[i++] = 'E'; + int e; + if (decExponent <= 0) { + result[i++] = '-'; + e = -decExponent+1; + } else { + e = decExponent-1; + } + // decExponent has 1, 2, or 3, digits + if (e <= 9) { + result[i++] = (char)(e+'0'); + } else if (e <= 99) { + result[i++] = (char)(e/10 +'0'); + result[i++] = (char)(e%10 + '0'); + } else { + result[i++] = (char)(e/100+'0'); + e %= 100; + result[i++] = (char)(e/10+'0'); + result[i++] = (char)(e%10 + '0'); + } + } + } + return i; + } + + // Per-thread buffer for string/stringbuffer conversion + private static ThreadLocal perThreadBuffer = new ThreadLocal() { + protected synchronized Object initialValue() { + return new char[26]; + } + }; + + public void appendTo(Appendable buf) { + char result[] = (char[])(perThreadBuffer.get()); + int i = getChars(result); + if (buf instanceof StringBuilder) + ((StringBuilder) buf).append(result, 0, i); + else if (buf instanceof StringBuffer) + ((StringBuffer) buf).append(result, 0, i); + else + assert false; + } + + public static FloatingDecimal + readJavaFormatString( String in ) throws NumberFormatException { + boolean isNegative = false; + boolean signSeen = false; + int decExp; + char c; + + parseNumber: + try{ + in = in.trim(); // don't fool around with white space. + // throws NullPointerException if null + int l = in.length(); + if ( l == 0 ) throw new NumberFormatException("empty String"); + int i = 0; + switch ( c = in.charAt( i ) ){ + case '-': + isNegative = true; + //FALLTHROUGH + case '+': + i++; + signSeen = true; + } + + // Check for NaN and Infinity strings + c = in.charAt(i); + if(c == 'N' || c == 'I') { // possible NaN or infinity + boolean potentialNaN = false; + char targetChars[] = null; // char array of "NaN" or "Infinity" + + if(c == 'N') { + targetChars = notANumber; + potentialNaN = true; + } else { + targetChars = infinity; + } + + // compare Input string to "NaN" or "Infinity" + int j = 0; + while(i < l && j < targetChars.length) { + if(in.charAt(i) == targetChars[j]) { + i++; j++; + } + else // something is amiss, throw exception + break parseNumber; + } + + // For the candidate string to be a NaN or infinity, + // all characters in input string and target char[] + // must be matched ==> j must equal targetChars.length + // and i must equal l + if( (j == targetChars.length) && (i == l) ) { // return NaN or infinity + return (potentialNaN ? new FloatingDecimal(Double.NaN) // NaN has no sign + : new FloatingDecimal(isNegative? + Double.NEGATIVE_INFINITY: + Double.POSITIVE_INFINITY)) ; + } + else { // something went wrong, throw exception + break parseNumber; + } + + } else if (c == '0') { // check for hexadecimal floating-point number + if (l > i+1 ) { + char ch = in.charAt(i+1); + if (ch == 'x' || ch == 'X' ) // possible hex string + return parseHexString(in); + } + } // look for and process decimal floating-point string + + char[] digits = new char[ l ]; + int nDigits= 0; + boolean decSeen = false; + int decPt = 0; + int nLeadZero = 0; + int nTrailZero= 0; + digitLoop: + while ( i < l ){ + switch ( c = in.charAt( i ) ){ + case '0': + if ( nDigits > 0 ){ + nTrailZero += 1; + } else { + nLeadZero += 1; + } + break; // out of switch. + case '1': + case '2': + case '3': + case '4': + case '5': + case '6': + case '7': + case '8': + case '9': + while ( nTrailZero > 0 ){ + digits[nDigits++] = '0'; + nTrailZero -= 1; + } + digits[nDigits++] = c; + break; // out of switch. + case '.': + if ( decSeen ){ + // already saw one ., this is the 2nd. + throw new NumberFormatException("multiple points"); + } + decPt = i; + if ( signSeen ){ + decPt -= 1; + } + decSeen = true; + break; // out of switch. + default: + break digitLoop; + } + i++; + } + /* + * At this point, we've scanned all the digits and decimal + * point we're going to see. Trim off leading and trailing + * zeros, which will just confuse us later, and adjust + * our initial decimal exponent accordingly. + * To review: + * we have seen i total characters. + * nLeadZero of them were zeros before any other digits. + * nTrailZero of them were zeros after any other digits. + * if ( decSeen ), then a . was seen after decPt characters + * ( including leading zeros which have been discarded ) + * nDigits characters were neither lead nor trailing + * zeros, nor point + */ + /* + * special hack: if we saw no non-zero digits, then the + * answer is zero! + * Unfortunately, we feel honor-bound to keep parsing! + */ + if ( nDigits == 0 ){ + digits = zero; + nDigits = 1; + if ( nLeadZero == 0 ){ + // we saw NO DIGITS AT ALL, + // not even a crummy 0! + // this is not allowed. + break parseNumber; // go throw exception + } + + } + + /* Our initial exponent is decPt, adjusted by the number of + * discarded zeros. Or, if there was no decPt, + * then its just nDigits adjusted by discarded trailing zeros. + */ + if ( decSeen ){ + decExp = decPt - nLeadZero; + } else { + decExp = nDigits+nTrailZero; + } + + /* + * Look for 'e' or 'E' and an optionally signed integer. + */ + if ( (i < l) && (((c = in.charAt(i) )=='e') || (c == 'E') ) ){ + int expSign = 1; + int expVal = 0; + int reallyBig = Integer.MAX_VALUE / 10; + boolean expOverflow = false; + switch( in.charAt(++i) ){ + case '-': + expSign = -1; + //FALLTHROUGH + case '+': + i++; + } + int expAt = i; + expLoop: + while ( i < l ){ + if ( expVal >= reallyBig ){ + // the next character will cause integer + // overflow. + expOverflow = true; + } + switch ( c = in.charAt(i++) ){ + case '0': + case '1': + case '2': + case '3': + case '4': + case '5': + case '6': + case '7': + case '8': + case '9': + expVal = expVal*10 + ( (int)c - (int)'0' ); + continue; + default: + i--; // back up. + break expLoop; // stop parsing exponent. + } + } + int expLimit = bigDecimalExponent+nDigits+nTrailZero; + if ( expOverflow || ( expVal > expLimit ) ){ + // + // The intent here is to end up with + // infinity or zero, as appropriate. + // The reason for yielding such a small decExponent, + // rather than something intuitive such as + // expSign*Integer.MAX_VALUE, is that this value + // is subject to further manipulation in + // doubleValue() and floatValue(), and I don't want + // it to be able to cause overflow there! + // (The only way we can get into trouble here is for + // really outrageous nDigits+nTrailZero, such as 2 billion. ) + // + decExp = expSign*expLimit; + } else { + // this should not overflow, since we tested + // for expVal > (MAX+N), where N >= abs(decExp) + decExp = decExp + expSign*expVal; + } + + // if we saw something not a digit ( or end of string ) + // after the [Ee][+-], without seeing any digits at all + // this is certainly an error. If we saw some digits, + // but then some trailing garbage, that might be ok. + // so we just fall through in that case. + // HUMBUG + if ( i == expAt ) + break parseNumber; // certainly bad + } + /* + * We parsed everything we could. + * If there are leftovers, then this is not good input! + */ + if ( i < l && + ((i != l - 1) || + (in.charAt(i) != 'f' && + in.charAt(i) != 'F' && + in.charAt(i) != 'd' && + in.charAt(i) != 'D'))) { + break parseNumber; // go throw exception + } + + return new FloatingDecimal( isNegative, decExp, digits, nDigits, false ); + } catch ( StringIndexOutOfBoundsException e ){ } + throw new NumberFormatException("For input string: \"" + in + "\""); + } + + /* + * Take a FloatingDecimal, which we presumably just scanned in, + * and find out what its value is, as a double. + * + * AS A SIDE EFFECT, SET roundDir TO INDICATE PREFERRED + * ROUNDING DIRECTION in case the result is really destined + * for a single-precision float. + */ + + public double + doubleValue(){ + int kDigits = Math.min( nDigits, maxDecimalDigits+1 ); + long lValue; + double dValue; + double rValue, tValue; + + // First, check for NaN and Infinity values + if(digits == infinity || digits == notANumber) { + if(digits == notANumber) + return Double.NaN; + else + return (isNegative?Double.NEGATIVE_INFINITY:Double.POSITIVE_INFINITY); + } + else { + if (mustSetRoundDir) { + roundDir = 0; + } + /* + * convert the lead kDigits to a long integer. + */ + // (special performance hack: start to do it using int) + int iValue = (int)digits[0]-(int)'0'; + int iDigits = Math.min( kDigits, intDecimalDigits ); + for ( int i=1; i < iDigits; i++ ){ + iValue = iValue*10 + (int)digits[i]-(int)'0'; + } + lValue = (long)iValue; + for ( int i=iDigits; i < kDigits; i++ ){ + lValue = lValue*10L + (long)((int)digits[i]-(int)'0'); + } + dValue = (double)lValue; + int exp = decExponent-kDigits; + /* + * lValue now contains a long integer with the value of + * the first kDigits digits of the number. + * dValue contains the (double) of the same. + */ + + if ( nDigits <= maxDecimalDigits ){ + /* + * possibly an easy case. + * We know that the digits can be represented + * exactly. And if the exponent isn't too outrageous, + * the whole thing can be done with one operation, + * thus one rounding error. + * Note that all our constructors trim all leading and + * trailing zeros, so simple values (including zero) + * will always end up here + */ + if (exp == 0 || dValue == 0.0) + return (isNegative)? -dValue : dValue; // small floating integer + else if ( exp >= 0 ){ + if ( exp <= maxSmallTen ){ + /* + * Can get the answer with one operation, + * thus one roundoff. + */ + rValue = dValue * small10pow[exp]; + if ( mustSetRoundDir ){ + tValue = rValue / small10pow[exp]; + roundDir = ( tValue == dValue ) ? 0 + :( tValue < dValue ) ? 1 + : -1; + } + return (isNegative)? -rValue : rValue; + } + int slop = maxDecimalDigits - kDigits; + if ( exp <= maxSmallTen+slop ){ + /* + * We can multiply dValue by 10^(slop) + * and it is still "small" and exact. + * Then we can multiply by 10^(exp-slop) + * with one rounding. + */ + dValue *= small10pow[slop]; + rValue = dValue * small10pow[exp-slop]; + + if ( mustSetRoundDir ){ + tValue = rValue / small10pow[exp-slop]; + roundDir = ( tValue == dValue ) ? 0 + :( tValue < dValue ) ? 1 + : -1; + } + return (isNegative)? -rValue : rValue; + } + /* + * Else we have a hard case with a positive exp. + */ + } else { + if ( exp >= -maxSmallTen ){ + /* + * Can get the answer in one division. + */ + rValue = dValue / small10pow[-exp]; + tValue = rValue * small10pow[-exp]; + if ( mustSetRoundDir ){ + roundDir = ( tValue == dValue ) ? 0 + :( tValue < dValue ) ? 1 + : -1; + } + return (isNegative)? -rValue : rValue; + } + /* + * Else we have a hard case with a negative exp. + */ + } + } + + /* + * Harder cases: + * The sum of digits plus exponent is greater than + * what we think we can do with one error. + * + * Start by approximating the right answer by, + * naively, scaling by powers of 10. + */ + if ( exp > 0 ){ + if ( decExponent > maxDecimalExponent+1 ){ + /* + * Lets face it. This is going to be + * Infinity. Cut to the chase. + */ + return (isNegative)? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY; + } + if ( (exp&15) != 0 ){ + dValue *= small10pow[exp&15]; + } + if ( (exp>>=4) != 0 ){ + int j; + for( j = 0; exp > 1; j++, exp>>=1 ){ + if ( (exp&1)!=0) + dValue *= big10pow[j]; + } + /* + * The reason for the weird exp > 1 condition + * in the above loop was so that the last multiply + * would get unrolled. We handle it here. + * It could overflow. + */ + double t = dValue * big10pow[j]; + if ( Double.isInfinite( t ) ){ + /* + * It did overflow. + * Look more closely at the result. + * If the exponent is just one too large, + * then use the maximum finite as our estimate + * value. Else call the result infinity + * and punt it. + * ( I presume this could happen because + * rounding forces the result here to be + * an ULP or two larger than + * Double.MAX_VALUE ). + */ + t = dValue / 2.0; + t *= big10pow[j]; + if ( Double.isInfinite( t ) ){ + return (isNegative)? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY; + } + t = Double.MAX_VALUE; + } + dValue = t; + } + } else if ( exp < 0 ){ + exp = -exp; + if ( decExponent < minDecimalExponent-1 ){ + /* + * Lets face it. This is going to be + * zero. Cut to the chase. + */ + /*IKVM*/ + // HACK we use Double.longBitsToDouble(0) to work around .NET 2.0 x64 CLR JIT bug that optimized the expression to -0.0 + return (isNegative) ? -0.0 : Double.longBitsToDouble(0); + } + if ( (exp&15) != 0 ){ + dValue /= small10pow[exp&15]; + } + if ( (exp>>=4) != 0 ){ + int j; + for( j = 0; exp > 1; j++, exp>>=1 ){ + if ( (exp&1)!=0) + dValue *= tiny10pow[j]; + } + /* + * The reason for the weird exp > 1 condition + * in the above loop was so that the last multiply + * would get unrolled. We handle it here. + * It could underflow. + */ + double t = dValue * tiny10pow[j]; + if ( t == 0.0 ){ + /* + * It did underflow. + * Look more closely at the result. + * If the exponent is just one too small, + * then use the minimum finite as our estimate + * value. Else call the result 0.0 + * and punt it. + * ( I presume this could happen because + * rounding forces the result here to be + * an ULP or two less than + * Double.MIN_VALUE ). + */ + t = dValue * 2.0; + t *= tiny10pow[j]; + if ( t == 0.0 ){ + /*IKVM*/ + // HACK we use Double.longBitsToDouble(0) to work around .NET 2.0 x64 CLR JIT bug that optimized the expression to -0.0 + return (isNegative)? -0.0 : Double.longBitsToDouble(0); + } + t = Double.MIN_VALUE; + } + dValue = t; + } + } + + /* + * dValue is now approximately the result. + * The hard part is adjusting it, by comparison + * with FDBigInt arithmetic. + * Formulate the EXACT big-number result as + * bigD0 * 10^exp + */ + FDBigInt bigD0 = new FDBigInt( lValue, digits, kDigits, nDigits ); + exp = decExponent - nDigits; + + correctionLoop: + while(true){ + /* AS A SIDE EFFECT, THIS METHOD WILL SET THE INSTANCE VARIABLES + * bigIntExp and bigIntNBits + */ + FDBigInt bigB = doubleToBigInt( dValue ); + + /* + * Scale bigD, bigB appropriately for + * big-integer operations. + * Naively, we multiply by powers of ten + * and powers of two. What we actually do + * is keep track of the powers of 5 and + * powers of 2 we would use, then factor out + * common divisors before doing the work. + */ + int B2, B5; // powers of 2, 5 in bigB + int D2, D5; // powers of 2, 5 in bigD + int Ulp2; // powers of 2 in halfUlp. + if ( exp >= 0 ){ + B2 = B5 = 0; + D2 = D5 = exp; + } else { + B2 = B5 = -exp; + D2 = D5 = 0; + } + if ( bigIntExp >= 0 ){ + B2 += bigIntExp; + } else { + D2 -= bigIntExp; + } + Ulp2 = B2; + // shift bigB and bigD left by a number s. t. + // halfUlp is still an integer. + int hulpbias; + if ( bigIntExp+bigIntNBits <= -expBias+1 ){ + // This is going to be a denormalized number + // (if not actually zero). + // half an ULP is at 2^-(expBias+expShift+1) + hulpbias = bigIntExp+ expBias + expShift; + } else { + hulpbias = expShift + 2 - bigIntNBits; + } + B2 += hulpbias; + D2 += hulpbias; + // if there are common factors of 2, we might just as well + // factor them out, as they add nothing useful. + int common2 = Math.min( B2, Math.min( D2, Ulp2 ) ); + B2 -= common2; + D2 -= common2; + Ulp2 -= common2; + // do multiplications by powers of 5 and 2 + bigB = multPow52( bigB, B5, B2 ); + FDBigInt bigD = multPow52( new FDBigInt( bigD0 ), D5, D2 ); + // + // to recap: + // bigB is the scaled-big-int version of our floating-point + // candidate. + // bigD is the scaled-big-int version of the exact value + // as we understand it. + // halfUlp is 1/2 an ulp of bigB, except for special cases + // of exact powers of 2 + // + // the plan is to compare bigB with bigD, and if the difference + // is less than halfUlp, then we're satisfied. Otherwise, + // use the ratio of difference to halfUlp to calculate a fudge + // factor to add to the floating value, then go 'round again. + // + FDBigInt diff; + int cmpResult; + boolean overvalue; + if ( (cmpResult = bigB.cmp( bigD ) ) > 0 ){ + overvalue = true; // our candidate is too big. + diff = bigB.sub( bigD ); + if ( (bigIntNBits == 1) && (bigIntExp > -expBias) ){ + // candidate is a normalized exact power of 2 and + // is too big. We will be subtracting. + // For our purposes, ulp is the ulp of the + // next smaller range. + Ulp2 -= 1; + if ( Ulp2 < 0 ){ + // rats. Cannot de-scale ulp this far. + // must scale diff in other direction. + Ulp2 = 0; + diff.lshiftMe( 1 ); + } + } + } else if ( cmpResult < 0 ){ + overvalue = false; // our candidate is too small. + diff = bigD.sub( bigB ); + } else { + // the candidate is exactly right! + // this happens with surprising frequency + break correctionLoop; + } + FDBigInt halfUlp = constructPow52( B5, Ulp2 ); + if ( (cmpResult = diff.cmp( halfUlp ) ) < 0 ){ + // difference is small. + // this is close enough + if (mustSetRoundDir) { + roundDir = overvalue ? -1 : 1; + } + break correctionLoop; + } else if ( cmpResult == 0 ){ + // difference is exactly half an ULP + // round to some other value maybe, then finish + dValue += 0.5*ulp( dValue, overvalue ); + // should check for bigIntNBits == 1 here?? + if (mustSetRoundDir) { + roundDir = overvalue ? -1 : 1; + } + break correctionLoop; + } else { + // difference is non-trivial. + // could scale addend by ratio of difference to + // halfUlp here, if we bothered to compute that difference. + // Most of the time ( I hope ) it is about 1 anyway. + dValue += ulp( dValue, overvalue ); + if ( dValue == 0.0 || dValue == Double.POSITIVE_INFINITY ) + break correctionLoop; // oops. Fell off end of range. + continue; // try again. + } + + } + return (isNegative)? -dValue : dValue; + } + } + + /* + * Take a FloatingDecimal, which we presumably just scanned in, + * and find out what its value is, as a float. + * This is distinct from doubleValue() to avoid the extremely + * unlikely case of a double rounding error, wherein the conversion + * to double has one rounding error, and the conversion of that double + * to a float has another rounding error, IN THE WRONG DIRECTION, + * ( because of the preference to a zero low-order bit ). + */ + + public float + floatValue(){ + int kDigits = Math.min( nDigits, singleMaxDecimalDigits+1 ); + int iValue; + float fValue; + + // First, check for NaN and Infinity values + if(digits == infinity || digits == notANumber) { + if(digits == notANumber) + return Float.NaN; + else + return (isNegative?Float.NEGATIVE_INFINITY:Float.POSITIVE_INFINITY); + } + else { + /* + * convert the lead kDigits to an integer. + */ + iValue = (int)digits[0]-(int)'0'; + for ( int i=1; i < kDigits; i++ ){ + iValue = iValue*10 + (int)digits[i]-(int)'0'; + } + fValue = (float)iValue; + int exp = decExponent-kDigits; + /* + * iValue now contains an integer with the value of + * the first kDigits digits of the number. + * fValue contains the (float) of the same. + */ + + if ( nDigits <= singleMaxDecimalDigits ){ + /* + * possibly an easy case. + * We know that the digits can be represented + * exactly. And if the exponent isn't too outrageous, + * the whole thing can be done with one operation, + * thus one rounding error. + * Note that all our constructors trim all leading and + * trailing zeros, so simple values (including zero) + * will always end up here. + */ + if (exp == 0 || fValue == 0.0f) + return (isNegative)? -fValue : fValue; // small floating integer + else if ( exp >= 0 ){ + if ( exp <= singleMaxSmallTen ){ + /* + * Can get the answer with one operation, + * thus one roundoff. + */ + fValue *= singleSmall10pow[exp]; + return (isNegative)? -fValue : fValue; + } + int slop = singleMaxDecimalDigits - kDigits; + if ( exp <= singleMaxSmallTen+slop ){ + /* + * We can multiply dValue by 10^(slop) + * and it is still "small" and exact. + * Then we can multiply by 10^(exp-slop) + * with one rounding. + */ + fValue *= singleSmall10pow[slop]; + fValue *= singleSmall10pow[exp-slop]; + return (isNegative)? -fValue : fValue; + } + /* + * Else we have a hard case with a positive exp. + */ + } else { + if ( exp >= -singleMaxSmallTen ){ + /* + * Can get the answer in one division. + */ + fValue /= singleSmall10pow[-exp]; + return (isNegative)? -fValue : fValue; + } + /* + * Else we have a hard case with a negative exp. + */ + } + } else if ( (decExponent >= nDigits) && (nDigits+decExponent <= maxDecimalDigits) ){ + /* + * In double-precision, this is an exact floating integer. + * So we can compute to double, then shorten to float + * with one round, and get the right answer. + * + * First, finish accumulating digits. + * Then convert that integer to a double, multiply + * by the appropriate power of ten, and convert to float. + */ + long lValue = (long)iValue; + for ( int i=kDigits; i < nDigits; i++ ){ + lValue = lValue*10L + (long)((int)digits[i]-(int)'0'); + } + double dValue = (double)lValue; + exp = decExponent-nDigits; + dValue *= small10pow[exp]; + fValue = (float)dValue; + return (isNegative)? -fValue : fValue; + + } + /* + * Harder cases: + * The sum of digits plus exponent is greater than + * what we think we can do with one error. + * + * Start by weeding out obviously out-of-range + * results, then convert to double and go to + * common hard-case code. + */ + if ( decExponent > singleMaxDecimalExponent+1 ){ + /* + * Lets face it. This is going to be + * Infinity. Cut to the chase. + */ + return (isNegative)? Float.NEGATIVE_INFINITY : Float.POSITIVE_INFINITY; + } else if ( decExponent < singleMinDecimalExponent-1 ){ + /* + * Lets face it. This is going to be + * zero. Cut to the chase. + */ + /*IKVM*/ + // HACK we use Float.intBitsToFloat(0) to work around .NET 2.0 x64 CLR JIT bug that optimized the expression to -0.0 + return (isNegative) ? -0.0f : Float.intBitsToFloat(0); + } + + /* + * Here, we do 'way too much work, but throwing away + * our partial results, and going and doing the whole + * thing as double, then throwing away half the bits that computes + * when we convert back to float. + * + * The alternative is to reproduce the whole multiple-precision + * algorithm for float precision, or to try to parameterize it + * for common usage. The former will take about 400 lines of code, + * and the latter I tried without success. Thus the semi-hack + * answer here. + */ + mustSetRoundDir = !fromHex; + double dValue = doubleValue(); + return stickyRound( dValue ); + } + } + + + /* + * All the positive powers of 10 that can be + * represented exactly in double/float. + */ + private static final double small10pow[] = { + 1.0e0, + 1.0e1, 1.0e2, 1.0e3, 1.0e4, 1.0e5, + 1.0e6, 1.0e7, 1.0e8, 1.0e9, 1.0e10, + 1.0e11, 1.0e12, 1.0e13, 1.0e14, 1.0e15, + 1.0e16, 1.0e17, 1.0e18, 1.0e19, 1.0e20, + 1.0e21, 1.0e22 + }; + + private static final float singleSmall10pow[] = { + 1.0e0f, + 1.0e1f, 1.0e2f, 1.0e3f, 1.0e4f, 1.0e5f, + 1.0e6f, 1.0e7f, 1.0e8f, 1.0e9f, 1.0e10f + }; + + private static final double big10pow[] = { + 1e16, 1e32, 1e64, 1e128, 1e256 }; + private static final double tiny10pow[] = { + 1e-16, 1e-32, 1e-64, 1e-128, 1e-256 }; + + private static final int maxSmallTen = small10pow.length-1; + private static final int singleMaxSmallTen = singleSmall10pow.length-1; + + private static final int small5pow[] = { + 1, + 5, + 5*5, + 5*5*5, + 5*5*5*5, + 5*5*5*5*5, + 5*5*5*5*5*5, + 5*5*5*5*5*5*5, + 5*5*5*5*5*5*5*5, + 5*5*5*5*5*5*5*5*5, + 5*5*5*5*5*5*5*5*5*5, + 5*5*5*5*5*5*5*5*5*5*5, + 5*5*5*5*5*5*5*5*5*5*5*5, + 5*5*5*5*5*5*5*5*5*5*5*5*5 + }; + + + private static final long long5pow[] = { + 1L, + 5L, + 5L*5, + 5L*5*5, + 5L*5*5*5, + 5L*5*5*5*5, + 5L*5*5*5*5*5, + 5L*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5, + 5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5, + }; + + // approximately ceil( log2( long5pow[i] ) ) + private static final int n5bits[] = { + 0, + 3, + 5, + 7, + 10, + 12, + 14, + 17, + 19, + 21, + 24, + 26, + 28, + 31, + 33, + 35, + 38, + 40, + 42, + 45, + 47, + 49, + 52, + 54, + 56, + 59, + 61, + }; + + private static final char infinity[] = { 'I', 'n', 'f', 'i', 'n', 'i', 't', 'y' }; + private static final char notANumber[] = { 'N', 'a', 'N' }; + private static final char zero[] = { '0', '0', '0', '0', '0', '0', '0', '0' }; + + + /* + * Grammar is compatible with hexadecimal floating-point constants + * described in section 6.4.4.2 of the C99 specification. + */ + private static Pattern hexFloatPattern = Pattern.compile( + //1 234 56 7 8 9 + "([-+])?0[xX](((\\p{XDigit}+)\\.?)|((\\p{XDigit}*)\\.(\\p{XDigit}+)))[pP]([-+])?(\\p{Digit}+)[fFdD]?" + ); + + /* + * Convert string s to a suitable floating decimal; uses the + * double constructor and set the roundDir variable appropriately + * in case the value is later converted to a float. + */ + static FloatingDecimal parseHexString(String s) { + // Verify string is a member of the hexadecimal floating-point + // string language. + Matcher m = hexFloatPattern.matcher(s); + boolean validInput = m.matches(); + + if (!validInput) { + // Input does not match pattern + throw new NumberFormatException("For input string: \"" + s + "\""); + } else { // validInput + /* + * We must isolate the sign, significand, and exponent + * fields. The sign value is straightforward. Since + * floating-point numbers are stored with a normalized + * representation, the significand and exponent are + * interrelated. + * + * After extracting the sign, we normalized the + * significand as a hexadecimal value, calculating an + * exponent adjust for any shifts made during + * normalization. If the significand is zero, the + * exponent doesn't need to be examined since the output + * will be zero. + * + * Next the exponent in the input string is extracted. + * Afterwards, the significand is normalized as a *binary* + * value and the input value's normalized exponent can be + * computed. The significand bits are copied into a + * double significand; if the string has more logical bits + * than can fit in a double, the extra bits affect the + * round and sticky bits which are used to round the final + * value. + */ + + // Extract significand sign + String group1 = m.group(1); + double sign = (( group1 == null ) || group1.equals("+"))? 1.0 : -1.0; + + + // Extract Significand magnitude + /* + * Based on the form of the significand, calculate how the + * binary exponent needs to be adjusted to create a + * normalized *hexadecimal* floating-point number; that + * is, a number where there is one nonzero hex digit to + * the left of the (hexa)decimal point. Since we are + * adjusting a binary, not hexadecimal exponent, the + * exponent is adjusted by a multiple of 4. + * + * There are a number of significand scenarios to consider; + * letters are used in indicate nonzero digits: + * + * 1. 000xxxx => x.xxx normalized + * increase exponent by (number of x's - 1)*4 + * + * 2. 000xxx.yyyy => x.xxyyyy normalized + * increase exponent by (number of x's - 1)*4 + * + * 3. .000yyy => y.yy normalized + * decrease exponent by (number of zeros + 1)*4 + * + * 4. 000.00000yyy => y.yy normalized + * decrease exponent by (number of zeros to right of point + 1)*4 + * + * If the significand is exactly zero, return a properly + * signed zero. + */ + + String significandString =null; + int signifLength = 0; + int exponentAdjust = 0; + { + int leftDigits = 0; // number of meaningful digits to + // left of "decimal" point + // (leading zeros stripped) + int rightDigits = 0; // number of digits to right of + // "decimal" point; leading zeros + // must always be accounted for + /* + * The significand is made up of either + * + * 1. group 4 entirely (integer portion only) + * + * OR + * + * 2. the fractional portion from group 7 plus any + * (optional) integer portions from group 6. + */ + String group4; + if( (group4 = m.group(4)) != null) { // Integer-only significand + // Leading zeros never matter on the integer portion + significandString = stripLeadingZeros(group4); + leftDigits = significandString.length(); + } + else { + // Group 6 is the optional integer; leading zeros + // never matter on the integer portion + String group6 = stripLeadingZeros(m.group(6)); + leftDigits = group6.length(); + + // fraction + String group7 = m.group(7); + rightDigits = group7.length(); + + // Turn "integer.fraction" into "integer"+"fraction" + significandString = + ((group6 == null)?"":group6) + // is the null + // check necessary? + group7; + } + + significandString = stripLeadingZeros(significandString); + signifLength = significandString.length(); + + /* + * Adjust exponent as described above + */ + if (leftDigits >= 1) { // Cases 1 and 2 + exponentAdjust = 4*(leftDigits - 1); + } else { // Cases 3 and 4 + exponentAdjust = -4*( rightDigits - signifLength + 1); + } + + // If the significand is zero, the exponent doesn't + // matter; return a properly signed zero. + + if (signifLength == 0) { // Only zeros in input + return new FloatingDecimal(sign * 0.0); + } + } + + // Extract Exponent + /* + * Use an int to read in the exponent value; this should + * provide more than sufficient range for non-contrived + * inputs. If reading the exponent in as an int does + * overflow, examine the sign of the exponent and + * significand to determine what to do. + */ + String group8 = m.group(8); + boolean positiveExponent = ( group8 == null ) || group8.equals("+"); + long unsignedRawExponent; + try { + unsignedRawExponent = Integer.parseInt(m.group(9)); + } + catch (NumberFormatException e) { + // At this point, we know the exponent is + // syntactically well-formed as a sequence of + // digits. Therefore, if an NumberFormatException + // is thrown, it must be due to overflowing int's + // range. Also, at this point, we have already + // checked for a zero significand. Thus the signs + // of the exponent and significand determine the + // final result: + // + // significand + // + - + // exponent + +infinity -infinity + // - +0.0 -0.0 + return new FloatingDecimal(sign * (positiveExponent ? + Double.POSITIVE_INFINITY : 0.0)); + } + + long rawExponent = + (positiveExponent ? 1L : -1L) * // exponent sign + unsignedRawExponent; // exponent magnitude + + // Calculate partially adjusted exponent + long exponent = rawExponent + exponentAdjust ; + + // Starting copying non-zero bits into proper position in + // a long; copy explicit bit too; this will be masked + // later for normal values. + + boolean round = false; + boolean sticky = false; + int bitsCopied=0; + int nextShift=0; + long significand=0L; + // First iteration is different, since we only copy + // from the leading significand bit; one more exponent + // adjust will be needed... + + // IMPORTANT: make leadingDigit a long to avoid + // surprising shift semantics! + long leadingDigit = getHexDigit(significandString, 0); + + /* + * Left shift the leading digit (53 - (bit position of + * leading 1 in digit)); this sets the top bit of the + * significand to 1. The nextShift value is adjusted + * to take into account the number of bit positions of + * the leadingDigit actually used. Finally, the + * exponent is adjusted to normalize the significand + * as a binary value, not just a hex value. + */ + if (leadingDigit == 1) { + significand |= leadingDigit << 52; + nextShift = 52 - 4; + /* exponent += 0 */ } + else if (leadingDigit <= 3) { // [2, 3] + significand |= leadingDigit << 51; + nextShift = 52 - 5; + exponent += 1; + } + else if (leadingDigit <= 7) { // [4, 7] + significand |= leadingDigit << 50; + nextShift = 52 - 6; + exponent += 2; + } + else if (leadingDigit <= 15) { // [8, f] + significand |= leadingDigit << 49; + nextShift = 52 - 7; + exponent += 3; + } else { + throw new AssertionError("Result from digit conversion too large!"); + } + // The preceding if-else could be replaced by a single + // code block based on the high-order bit set in + // leadingDigit. Given leadingOnePosition, + + // significand |= leadingDigit << (SIGNIFICAND_WIDTH - leadingOnePosition); + // nextShift = 52 - (3 + leadingOnePosition); + // exponent += (leadingOnePosition-1); + + + /* + * Now the exponent variable is equal to the normalized + * binary exponent. Code below will make representation + * adjustments if the exponent is incremented after + * rounding (includes overflows to infinity) or if the + * result is subnormal. + */ + + // Copy digit into significand until the significand can't + // hold another full hex digit or there are no more input + // hex digits. + int i = 0; + for(i = 1; + i < signifLength && nextShift >= 0; + i++) { + long currentDigit = getHexDigit(significandString, i); + significand |= (currentDigit << nextShift); + nextShift-=4; + } + + // After the above loop, the bulk of the string is copied. + // Now, we must copy any partial hex digits into the + // significand AND compute the round bit and start computing + // sticky bit. + + if ( i < signifLength ) { // at least one hex input digit exists + long currentDigit = getHexDigit(significandString, i); + + // from nextShift, figure out how many bits need + // to be copied, if any + switch(nextShift) { // must be negative + case -1: + // three bits need to be copied in; can + // set round bit + significand |= ((currentDigit & 0xEL) >> 1); + round = (currentDigit & 0x1L) != 0L; + break; + + case -2: + // two bits need to be copied in; can + // set round and start sticky + significand |= ((currentDigit & 0xCL) >> 2); + round = (currentDigit &0x2L) != 0L; + sticky = (currentDigit & 0x1L) != 0; + break; + + case -3: + // one bit needs to be copied in + significand |= ((currentDigit & 0x8L)>>3); + // Now set round and start sticky, if possible + round = (currentDigit &0x4L) != 0L; + sticky = (currentDigit & 0x3L) != 0; + break; + + case -4: + // all bits copied into significand; set + // round and start sticky + round = ((currentDigit & 0x8L) != 0); // is top bit set? + // nonzeros in three low order bits? + sticky = (currentDigit & 0x7L) != 0; + break; + + default: + throw new AssertionError("Unexpected shift distance remainder."); + // break; + } + + // Round is set; sticky might be set. + + // For the sticky bit, it suffices to check the + // current digit and test for any nonzero digits in + // the remaining unprocessed input. + i++; + while(i < signifLength && !sticky) { + currentDigit = getHexDigit(significandString,i); + sticky = sticky || (currentDigit != 0); + i++; + } + + } + // else all of string was seen, round and sticky are + // correct as false. + + + // Check for overflow and update exponent accordingly. + + if (exponent > DoubleConsts.MAX_EXPONENT) { // Infinite result + // overflow to properly signed infinity + return new FloatingDecimal(sign * Double.POSITIVE_INFINITY); + } else { // Finite return value + if (exponent <= DoubleConsts.MAX_EXPONENT && // (Usually) normal result + exponent >= DoubleConsts.MIN_EXPONENT) { + + // The result returned in this block cannot be a + // zero or subnormal; however after the + // significand is adjusted from rounding, we could + // still overflow in infinity. + + // AND exponent bits into significand; if the + // significand is incremented and overflows from + // rounding, this combination will update the + // exponent correctly, even in the case of + // Double.MAX_VALUE overflowing to infinity. + + significand = (( ((long)exponent + + (long)DoubleConsts.EXP_BIAS) << + (DoubleConsts.SIGNIFICAND_WIDTH-1)) + & DoubleConsts.EXP_BIT_MASK) | + (DoubleConsts.SIGNIF_BIT_MASK & significand); + + } else { // Subnormal or zero + // (exponent < DoubleConsts.MIN_EXPONENT) + + if (exponent < (DoubleConsts.MIN_SUB_EXPONENT -1 )) { + // No way to round back to nonzero value + // regardless of significand if the exponent is + // less than -1075. + return new FloatingDecimal(sign * 0.0); + } else { // -1075 <= exponent <= MIN_EXPONENT -1 = -1023 + /* + * Find bit position to round to; recompute + * round and sticky bits, and shift + * significand right appropriately. + */ + + sticky = sticky || round; + round = false; + + // Number of bits of significand to preserve is + // exponent - abs_min_exp +1 + // check: + // -1075 +1074 + 1 = 0 + // -1023 +1074 + 1 = 52 + + int bitsDiscarded = 53 - + ((int)exponent - DoubleConsts.MIN_SUB_EXPONENT + 1); + assert bitsDiscarded >= 1 && bitsDiscarded <= 53; + + // What to do here: + // First, isolate the new round bit + round = (significand & (1L << (bitsDiscarded -1))) != 0L; + if (bitsDiscarded > 1) { + // create mask to update sticky bits; low + // order bitsDiscarded bits should be 1 + long mask = ~((~0L) << (bitsDiscarded -1)); + sticky = sticky || ((significand & mask) != 0L ) ; + } + + // Now, discard the bits + significand = significand >> bitsDiscarded; + + significand = (( ((long)(DoubleConsts.MIN_EXPONENT -1) + // subnorm exp. + (long)DoubleConsts.EXP_BIAS) << + (DoubleConsts.SIGNIFICAND_WIDTH-1)) + & DoubleConsts.EXP_BIT_MASK) | + (DoubleConsts.SIGNIF_BIT_MASK & significand); + } + } + + // The significand variable now contains the currently + // appropriate exponent bits too. + + /* + * Determine if significand should be incremented; + * making this determination depends on the least + * significant bit and the round and sticky bits. + * + * Round to nearest even rounding table, adapted from + * table 4.7 in "Computer Arithmetic" by IsraelKoren. + * The digit to the left of the "decimal" point is the + * least significant bit, the digits to the right of + * the point are the round and sticky bits + * + * Number Round(x) + * x0.00 x0. + * x0.01 x0. + * x0.10 x0. + * x0.11 x1. = x0. +1 + * x1.00 x1. + * x1.01 x1. + * x1.10 x1. + 1 + * x1.11 x1. + 1 + */ + boolean incremented = false; + boolean leastZero = ((significand & 1L) == 0L); + if( ( leastZero && round && sticky ) || + ((!leastZero) && round )) { + incremented = true; + significand++; + } + + FloatingDecimal fd = new FloatingDecimal(FpUtils.rawCopySign( + Double.longBitsToDouble(significand), + sign)); + + /* + * Set roundingDir variable field of fd properly so + * that the input string can be properly rounded to a + * float value. There are two cases to consider: + * + * 1. rounding to double discards sticky bit + * information that would change the result of a float + * rounding (near halfway case between two floats) + * + * 2. rounding to double rounds up when rounding up + * would not occur when rounding to float. + * + * For former case only needs to be considered when + * the bits rounded away when casting to float are all + * zero; otherwise, float round bit is properly set + * and sticky will already be true. + * + * The lower exponent bound for the code below is the + * minimum (normalized) subnormal exponent - 1 since a + * value with that exponent can round up to the + * minimum subnormal value and the sticky bit + * information must be preserved (i.e. case 1). + */ + if ((exponent >= FloatConsts.MIN_SUB_EXPONENT-1) && + (exponent <= FloatConsts.MAX_EXPONENT ) ){ + // Outside above exponent range, the float value + // will be zero or infinity. + + /* + * If the low-order 28 bits of a rounded double + * significand are 0, the double could be a + * half-way case for a rounding to float. If the + * double value is a half-way case, the double + * significand may have to be modified to round + * the the right float value (see the stickyRound + * method). If the rounding to double has lost + * what would be float sticky bit information, the + * double significand must be incremented. If the + * double value's significand was itself + * incremented, the float value may end up too + * large so the increment should be undone. + */ + if ((significand & 0xfffffffL) == 0x0L) { + // For negative values, the sign of the + // roundDir is the same as for positive values + // since adding 1 increasing the significand's + // magnitude and subtracting 1 decreases the + // significand's magnitude. If neither round + // nor sticky is true, the double value is + // exact and no adjustment is required for a + // proper float rounding. + if( round || sticky) { + if (leastZero) { // prerounding lsb is 0 + // If round and sticky were both true, + // and the least significant + // significand bit were 0, the rounded + // significand would not have its + // low-order bits be zero. Therefore, + // we only need to adjust the + // significand if round XOR sticky is + // true. + if (round ^ sticky) { + fd.roundDir = 1; + } + } + else { // prerounding lsb is 1 + // If the prerounding lsb is 1 and the + // resulting significand has its + // low-order bits zero, the significand + // was incremented. Here, we undo the + // increment, which will ensure the + // right guard and sticky bits for the + // float rounding. + if (round) + fd.roundDir = -1; + } + } + } + } + + fd.fromHex = true; + return fd; + } + } + } + + /** + * Return s with any leading zeros removed. + */ + static String stripLeadingZeros(String s) { + return s.replaceFirst("^0+", ""); + } + + /** + * Extract a hexadecimal digit from position position + * of string s. + */ + static int getHexDigit(String s, int position) { + int value = Character.digit(s.charAt(position), 16); + if (value <= -1 || value >= 16) { + throw new AssertionError("Unexpected failure of digit conversion of " + + s.charAt(position)); + } + return value; + } + + +} + +/* + * A really, really simple bigint package + * tailored to the needs of floating base conversion. + */ +class FDBigInt { + int nWords; // number of words used + int data[]; // value: data[0] is least significant + + + public FDBigInt( int v ){ + nWords = 1; + data = new int[1]; + data[0] = v; + } + + public FDBigInt( long v ){ + data = new int[2]; + data[0] = (int)v; + data[1] = (int)(v>>>32); + nWords = (data[1]==0) ? 1 : 2; + } + + public FDBigInt( FDBigInt other ){ + data = new int[nWords = other.nWords]; + System.arraycopy( other.data, 0, data, 0, nWords ); + } + + private FDBigInt( int [] d, int n ){ + data = d; + nWords = n; + } + + public FDBigInt( long seed, char digit[], int nd0, int nd ){ + int n= (nd+8)/9; // estimate size needed. + if ( n < 2 ) n = 2; + data = new int[n]; // allocate enough space + data[0] = (int)seed; // starting value + data[1] = (int)(seed>>>32); + nWords = (data[1]==0) ? 1 : 2; + int i = nd0; + int limit = nd-5; // slurp digits 5 at a time. + int v; + while ( i < limit ){ + int ilim = i+5; + v = (int)digit[i++]-(int)'0'; + while( i >5; + int bitcount = c & 0x1f; + int anticount = 32-bitcount; + int t[] = data; + int s[] = data; + if ( nWords+wordcount+1 > t.length ){ + // reallocate. + t = new int[ nWords+wordcount+1 ]; + } + int target = nWords+wordcount; + int src = nWords-1; + if ( bitcount == 0 ){ + // special hack, since an anticount of 32 won't go! + System.arraycopy( s, 0, t, wordcount, nWords ); + target = wordcount-1; + } else { + t[target--] = s[src]>>>anticount; + while ( src >= 1 ){ + t[target--] = (s[src]<>>anticount); + } + t[target--] = s[src]<= 0 ){ + t[target--] = 0; + } + data = t; + nWords += wordcount + 1; + // may have constructed high-order word of 0. + // if so, trim it + while ( nWords > 1 && data[nWords-1] == 0 ) + nWords--; + } + + /* + * normalize this number by shifting until + * the MSB of the number is at 0x08000000. + * This is in preparation for quoRemIteration, below. + * The idea is that, to make division easier, we want the + * divisor to be "normalized" -- usually this means shifting + * the MSB into the high words sign bit. But because we know that + * the quotient will be 0 < q < 10, we would like to arrange that + * the dividend not span up into another word of precision. + * (This needs to be explained more clearly!) + */ + public int + normalizeMe() throws IllegalArgumentException { + int src; + int wordcount = 0; + int bitcount = 0; + int v = 0; + for ( src= nWords-1 ; src >= 0 && (v=data[src]) == 0 ; src--){ + wordcount += 1; + } + if ( src < 0 ){ + // oops. Value is zero. Cannot normalize it! + throw new IllegalArgumentException("zero value"); + } + /* + * In most cases, we assume that wordcount is zero. This only + * makes sense, as we try not to maintain any high-order + * words full of zeros. In fact, if there are zeros, we will + * simply SHORTEN our number at this point. Watch closely... + */ + nWords -= wordcount; + /* + * Compute how far left we have to shift v s.t. its highest- + * order bit is in the right place. Then call lshiftMe to + * do the work. + */ + if ( (v & 0xf0000000) != 0 ){ + // will have to shift up into the next word. + // too bad. + for( bitcount = 32 ; (v & 0xf0000000) != 0 ; bitcount-- ) + v >>>= 1; + } else { + while ( v <= 0x000fffff ){ + // hack: byte-at-a-time shifting + v <<= 8; + bitcount += 8; + } + while ( v <= 0x07ffffff ){ + v <<= 1; + bitcount += 1; + } + } + if ( bitcount != 0 ) + lshiftMe( bitcount ); + return bitcount; + } + + /* + * Multiply a FDBigInt by an int. + * Result is a new FDBigInt. + */ + public FDBigInt + mult( int iv ) { + long v = iv; + int r[]; + long p; + + // guess adequate size of r. + r = new int[ ( v * ((long)data[nWords-1]&0xffffffffL) > 0xfffffffL ) ? nWords+1 : nWords ]; + p = 0L; + for( int i=0; i < nWords; i++ ) { + p += v * ((long)data[i]&0xffffffffL); + r[i] = (int)p; + p >>>= 32; + } + if ( p == 0L){ + return new FDBigInt( r, nWords ); + } else { + r[nWords] = (int)p; + return new FDBigInt( r, nWords+1 ); + } + } + + /* + * Multiply a FDBigInt by an int and add another int. + * Result is computed in place. + * Hope it fits! + */ + public void + multaddMe( int iv, int addend ) { + long v = iv; + long p; + + // unroll 0th iteration, doing addition. + p = v * ((long)data[0]&0xffffffffL) + ((long)addend&0xffffffffL); + data[0] = (int)p; + p >>>= 32; + for( int i=1; i < nWords; i++ ) { + p += v * ((long)data[i]&0xffffffffL); + data[i] = (int)p; + p >>>= 32; + } + if ( p != 0L){ + data[nWords] = (int)p; // will fail noisily if illegal! + nWords++; + } + } + + /* + * Multiply a FDBigInt by another FDBigInt. + * Result is a new FDBigInt. + */ + public FDBigInt + mult( FDBigInt other ){ + // crudely guess adequate size for r + int r[] = new int[ nWords + other.nWords ]; + int i; + // I think I am promised zeros... + + for( i = 0; i < this.nWords; i++ ){ + long v = (long)this.data[i] & 0xffffffffL; // UNSIGNED CONVERSION + long p = 0L; + int j; + for( j = 0; j < other.nWords; j++ ){ + p += ((long)r[i+j]&0xffffffffL) + v*((long)other.data[j]&0xffffffffL); // UNSIGNED CONVERSIONS ALL 'ROUND. + r[i+j] = (int)p; + p >>>= 32; + } + r[i+j] = (int)p; + } + // compute how much of r we actually needed for all that. + for ( i = r.length-1; i> 0; i--) + if ( r[i] != 0 ) + break; + return new FDBigInt( r, i+1 ); + } + + /* + * Add one FDBigInt to another. Return a FDBigInt + */ + public FDBigInt + add( FDBigInt other ){ + int i; + int a[], b[]; + int n, m; + long c = 0L; + // arrange such that a.nWords >= b.nWords; + // n = a.nWords, m = b.nWords + if ( this.nWords >= other.nWords ){ + a = this.data; + n = this.nWords; + b = other.data; + m = other.nWords; + } else { + a = other.data; + n = other.nWords; + b = this.data; + m = this.nWords; + } + int r[] = new int[ n ]; + for ( i = 0; i < n; i++ ){ + c += (long)a[i] & 0xffffffffL; + if ( i < m ){ + c += (long)b[i] & 0xffffffffL; + } + r[i] = (int) c; + c >>= 32; // signed shift. + } + if ( c != 0L ){ + // oops -- carry out -- need longer result. + int s[] = new int[ r.length+1 ]; + System.arraycopy( r, 0, s, 0, r.length ); + s[i++] = (int)c; + return new FDBigInt( s, i ); + } + return new FDBigInt( r, i ); + } + + /* + * Subtract one FDBigInt from another. Return a FDBigInt + * Assert that the result is positive. + */ + public FDBigInt + sub( FDBigInt other ){ + int r[] = new int[ this.nWords ]; + int i; + int n = this.nWords; + int m = other.nWords; + int nzeros = 0; + long c = 0L; + for ( i = 0; i < n; i++ ){ + c += (long)this.data[i] & 0xffffffffL; + if ( i < m ){ + c -= (long)other.data[i] & 0xffffffffL; + } + if ( ( r[i] = (int) c ) == 0 ) + nzeros++; + else + nzeros = 0; + c >>= 32; // signed shift + } + assert c == 0L : c; // borrow out of subtract + assert dataInRangeIsZero(i, m, other); // negative result of subtract + return new FDBigInt( r, n-nzeros ); + } + + private static boolean dataInRangeIsZero(int i, int m, FDBigInt other) { + while ( i < m ) + if (other.data[i++] != 0) + return false; + return true; + } + + /* + * Compare FDBigInt with another FDBigInt. Return an integer + * >0: this > other + * 0: this == other + * <0: this < other + */ + public int + cmp( FDBigInt other ){ + int i; + if ( this.nWords > other.nWords ){ + // if any of my high-order words is non-zero, + // then the answer is evident + int j = other.nWords-1; + for ( i = this.nWords-1; i > j ; i-- ) + if ( this.data[i] != 0 ) return 1; + }else if ( this.nWords < other.nWords ){ + // if any of other's high-order words is non-zero, + // then the answer is evident + int j = this.nWords-1; + for ( i = other.nWords-1; i > j ; i-- ) + if ( other.data[i] != 0 ) return -1; + } else{ + i = this.nWords-1; + } + for ( ; i > 0 ; i-- ) + if ( this.data[i] != other.data[i] ) + break; + // careful! want unsigned compare! + // use brute force here. + int a = this.data[i]; + int b = other.data[i]; + if ( a < 0 ){ + // a is really big, unsigned + if ( b < 0 ){ + return a-b; // both big, negative + } else { + return 1; // b not big, answer is obvious; + } + } else { + // a is not really big + if ( b < 0 ) { + // but b is really big + return -1; + } else { + return a - b; + } + } + } + + /* + * Compute + * q = (int)( this / S ) + * this = 10 * ( this mod S ) + * Return q. + * This is the iteration step of digit development for output. + * We assume that S has been normalized, as above, and that + * "this" has been lshift'ed accordingly. + * Also assume, of course, that the result, q, can be expressed + * as an integer, 0 <= q < 10. + */ + public int + quoRemIteration( FDBigInt S )throws IllegalArgumentException { + // ensure that this and S have the same number of + // digits. If S is properly normalized and q < 10 then + // this must be so. + if ( nWords != S.nWords ){ + throw new IllegalArgumentException("disparate values"); + } + // estimate q the obvious way. We will usually be + // right. If not, then we're only off by a little and + // will re-add. + int n = nWords-1; + long q = ((long)data[n]&0xffffffffL) / (long)S.data[n]; + long diff = 0L; + for ( int i = 0; i <= n ; i++ ){ + diff += ((long)data[i]&0xffffffffL) - q*((long)S.data[i]&0xffffffffL); + data[i] = (int)diff; + diff >>= 32; // N.B. SIGNED shift. + } + if ( diff != 0L ) { + // damn, damn, damn. q is too big. + // add S back in until this turns +. This should + // not be very many times! + long sum = 0L; + while ( sum == 0L ){ + sum = 0L; + for ( int i = 0; i <= n; i++ ){ + sum += ((long)data[i]&0xffffffffL) + ((long)S.data[i]&0xffffffffL); + data[i] = (int) sum; + sum >>= 32; // Signed or unsigned, answer is 0 or 1 + } + /* + * Originally the following line read + * "if ( sum !=0 && sum != -1 )" + * but that would be wrong, because of the + * treatment of the two values as entirely unsigned, + * it would be impossible for a carry-out to be interpreted + * as -1 -- it would have to be a single-bit carry-out, or + * +1. + */ + assert sum == 0 || sum == 1 : sum; // carry out of division correction + q -= 1; + } + } + // finally, we can multiply this by 10. + // it cannot overflow, right, as the high-order word has + // at least 4 high-order zeros! + long p = 0L; + for ( int i = 0; i <= n; i++ ){ + p += 10*((long)data[i]&0xffffffffL); + data[i] = (int)p; + p >>= 32; // SIGNED shift. + } + assert p == 0L : p; // Carry out of *10 + return (int)q; + } + + public long + longValue(){ + // if this can be represented as a long, return the value + assert this.nWords > 0 : this.nWords; // longValue confused + + if (this.nWords == 1) + return ((long)data[0]&0xffffffffL); + + assert dataInRangeIsZero(2, this.nWords, this); // value too big + assert data[1] >= 0; // value too big + return ((long)(data[1]) << 32) | ((long)data[0]&0xffffffffL); + } + + public String + toString() { + StringBuffer r = new StringBuffer(30); + r.append('['); + int i = Math.min( nWords-1, data.length-1) ; + if ( nWords > data.length ){ + r.append( "("+data.length+"<"+nWords+"!)" ); + } + for( ; i> 0 ; i-- ){ + r.append( Integer.toHexString( data[i] ) ); + r.append( (char) ' ' ); + } + r.append( Integer.toHexString( data[0] ) ); + r.append( (char) ']' ); + return new String( r ); + } +} + diff --git a/classpath/sun/misc/FpUtils.java b/classpath/sun/misc/FpUtils.java new file mode 100644 index 00000000..7ab2410c --- /dev/null +++ b/classpath/sun/misc/FpUtils.java @@ -0,0 +1,1230 @@ +/* + * Copyright 2003 Sun Microsystems, Inc. All Rights Reserved. + * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER. + * + * This code is free software; you can redistribute it and/or modify it + * under the terms of the GNU General Public License version 2 only, as + * published by the Free Software Foundation. Sun designates this + * particular file as subject to the "Classpath" exception as provided + * by Sun in the LICENSE file that accompanied this code. + * + * This code is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * version 2 for more details (a copy is included in the LICENSE file that + * accompanied this code). + * + * You should have received a copy of the GNU General Public License version + * 2 along with this work; if not, write to the Free Software Foundation, + * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. + * + * Please contact Sun Microsystems, Inc., 4150 Network Circle, Santa Clara, + * CA 95054 USA or visit www.sun.com if you need additional information or + * have any questions. + */ + +package sun.misc; + +import sun.misc.FloatConsts; +import sun.misc.DoubleConsts; + +/** + * The class FpUtils contains static utility methods for + * manipulating and inspecting float and + * double floating-point numbers. These methods include + * functionality recommended or required by the IEEE 754 + * floating-point standard. + * + * @author Joseph D. Darcy + * @version 1.9, 05/05/07 + */ + +public class FpUtils { + /* + * The methods in this class are reasonably implemented using + * direct or indirect bit-level manipulation of floating-point + * values. However, having access to the IEEE 754 recommended + * functions would obviate the need for most programmers to engage + * in floating-point bit-twiddling. + * + * An IEEE 754 number has three fields, from most significant bit + * to to least significant, sign, exponent, and significand. + * + * msb lsb + * [sign|exponent| fractional_significand] + * + * Using some encoding cleverness, explained below, the high order + * bit of the logical significand does not need to be explicitly + * stored, thus "fractional_significand" instead of simply + * "significand" in the figure above. + * + * For finite normal numbers, the numerical value encoded is + * + * (-1)^sign * 2^(exponent)*(1.fractional_significand) + * + * Most finite floating-point numbers are normalized; the exponent + * value is reduced until the leading significand bit is 1. + * Therefore, the leading 1 is redundant and is not explicitly + * stored. If a numerical value is so small it cannot be + * normalized, it has a subnormal representation. Subnormal + * numbers don't have a leading 1 in their significand; subnormals + * are encoding using a special exponent value. In other words, + * the high-order bit of the logical significand can be elided in + * from the representation in either case since the bit's value is + * implicit from the exponent value. + * + * The exponent field uses a biased representation; if the bits of + * the exponent are interpreted as a unsigned integer E, the + * exponent represented is E - E_bias where E_bias depends on the + * floating-point format. E can range between E_min and E_max, + * constants which depend on the floating-point format. E_min and + * E_max are -126 and +127 for float, -1022 and +1023 for double. + * + * The 32-bit float format has 1 sign bit, 8 exponent bits, and 23 + * bits for the significand (which is logically 24 bits wide + * because of the implicit bit). The 64-bit double format has 1 + * sign bit, 11 exponent bits, and 52 bits for the significand + * (logically 53 bits). + * + * Subnormal numbers and zero have the special exponent value + * E_min -1; the numerical value represented by a subnormal is: + * + * (-1)^sign * 2^(E_min)*(0.fractional_significand) + * + * Zero is represented by all zero bits in the exponent and all + * zero bits in the significand; zero can have either sign. + * + * Infinity and NaN are encoded using the exponent value E_max + + * 1. Signed infinities have all significand bits zero; NaNs have + * at least one non-zero significand bit. + * + * The details of IEEE 754 floating-point encoding will be used in + * the methods below without further comment. For further + * exposition on IEEE 754 numbers, see "IEEE Standard for Binary + * Floating-Point Arithmetic" ANSI/IEEE Std 754-1985 or William + * Kahan's "Lecture Notes on the Status of IEEE Standard 754 for + * Binary Floating-Point Arithmetic", + * http://www.cs.berkeley.edu/~wkahan/ieee754status/ieee754.ps. + * + * Many of this class's methods are members of the set of IEEE 754 + * recommended functions or similar functions recommended or + * required by IEEE 754R. Discussion of various implementation + * techniques for these functions have occurred in: + * + * W.J. Cody and Jerome T. Coonen, "Algorithm 772 Functions to + * Support the IEEE Standard for Binary Floating-Point + * Arithmetic," ACM Transactions on Mathematical Software, + * vol. 19, no. 4, December 1993, pp. 443-451. + * + * Joseph D. Darcy, "Writing robust IEEE recommended functions in + * ``100% Pure Java''(TM)," University of California, Berkeley + * technical report UCB//CSD-98-1009. + */ + + /** + * Don't let anyone instantiate this class. + */ + private FpUtils() {} + + // Constants used in scalb + static double twoToTheDoubleScaleUp = powerOfTwoD(512); + static double twoToTheDoubleScaleDown = powerOfTwoD(-512); + + // Helper Methods + + // The following helper methods are used in the implementation of + // the public recommended functions; they generally omit certain + // tests for exception cases. + + /** + * Returns unbiased exponent of a double. + */ + public static int getExponent(double d){ + /* + * Bitwise convert d to long, mask out exponent bits, shift + * to the right and then subtract out double's bias adjust to + * get true exponent value. + */ + return (int)(((Double.doubleToRawLongBits(d) & DoubleConsts.EXP_BIT_MASK) >> + (DoubleConsts.SIGNIFICAND_WIDTH - 1)) - DoubleConsts.EXP_BIAS); + } + + /** + * Returns unbiased exponent of a float. + */ + public static int getExponent(float f){ + /* + * Bitwise convert f to integer, mask out exponent bits, shift + * to the right and then subtract out float's bias adjust to + * get true exponent value + */ + return ((Float.floatToRawIntBits(f) & FloatConsts.EXP_BIT_MASK) >> + (FloatConsts.SIGNIFICAND_WIDTH - 1)) - FloatConsts.EXP_BIAS; + } + + /** + * Returns a floating-point power of two in the normal range. + */ + static double powerOfTwoD(int n) { + assert(n >= DoubleConsts.MIN_EXPONENT && n <= DoubleConsts.MAX_EXPONENT); + return Double.longBitsToDouble((((long)n + (long)DoubleConsts.EXP_BIAS) << + (DoubleConsts.SIGNIFICAND_WIDTH-1)) + & DoubleConsts.EXP_BIT_MASK); + } + + /** + * Returns a floating-point power of two in the normal range. + */ + static float powerOfTwoF(int n) { + assert(n >= FloatConsts.MIN_EXPONENT && n <= FloatConsts.MAX_EXPONENT); + return Float.intBitsToFloat(((n + FloatConsts.EXP_BIAS) << + (FloatConsts.SIGNIFICAND_WIDTH-1)) + & FloatConsts.EXP_BIT_MASK); + } + + /** + * Returns the first floating-point argument with the sign of the + * second floating-point argument. Note that unlike the {@link + * FpUtils#copySign(double, double) copySign} method, this method + * does not require NaN sign arguments to be treated + * as positive values; implementations are permitted to treat some + * NaN arguments as positive and other NaN arguments as negative + * to allow greater performance. + * + * @param magnitude the parameter providing the magnitude of the result + * @param sign the parameter providing the sign of the result + * @return a value with the magnitude of magnitude + * and the sign of sign. + * @author Joseph D. Darcy + */ + public static double rawCopySign(double magnitude, double sign) { + return Double.longBitsToDouble((Double.doubleToRawLongBits(sign) & + (DoubleConsts.SIGN_BIT_MASK)) | + (Double.doubleToRawLongBits(magnitude) & + (DoubleConsts.EXP_BIT_MASK | + DoubleConsts.SIGNIF_BIT_MASK))); + } + + /** + * Returns the first floating-point argument with the sign of the + * second floating-point argument. Note that unlike the {@link + * FpUtils#copySign(float, float) copySign} method, this method + * does not require NaN sign arguments to be treated + * as positive values; implementations are permitted to treat some + * NaN arguments as positive and other NaN arguments as negative + * to allow greater performance. + * + * @param magnitude the parameter providing the magnitude of the result + * @param sign the parameter providing the sign of the result + * @return a value with the magnitude of magnitude + * and the sign of sign. + * @author Joseph D. Darcy + */ + public static float rawCopySign(float magnitude, float sign) { + return Float.intBitsToFloat((Float.floatToRawIntBits(sign) & + (FloatConsts.SIGN_BIT_MASK)) | + (Float.floatToRawIntBits(magnitude) & + (FloatConsts.EXP_BIT_MASK | + FloatConsts.SIGNIF_BIT_MASK))); + } + + /* ***************************************************************** */ + + /** + * Returns true if the argument is a finite + * floating-point value; returns false otherwise (for + * NaN and infinity arguments). + * + * @param d the double value to be tested + * @return true if the argument is a finite + * floating-point value, false otherwise. + */ + public static boolean isFinite(double d) { + return Math.abs(d) <= DoubleConsts.MAX_VALUE; + } + + /** + * Returns true if the argument is a finite + * floating-point value; returns false otherwise (for + * NaN and infinity arguments). + * + * @param f the float value to be tested + * @return true if the argument is a finite + * floating-point value, false otherwise. + */ + public static boolean isFinite(float f) { + return Math.abs(f) <= FloatConsts.MAX_VALUE; + } + + /** + * Returns true if the specified number is infinitely + * large in magnitude, false otherwise. + * + *

Note that this method is equivalent to the {@link + * Double#isInfinite(double) Double.isInfinite} method; the + * functionality is included in this class for convenience. + * + * @param d the value to be tested. + * @return true if the value of the argument is positive + * infinity or negative infinity; false otherwise. + */ + public static boolean isInfinite(double d) { + return Double.isInfinite(d); + } + + /** + * Returns true if the specified number is infinitely + * large in magnitude, false otherwise. + * + *

Note that this method is equivalent to the {@link + * Float#isInfinite(float) Float.isInfinite} method; the + * functionality is included in this class for convenience. + * + * @param f the value to be tested. + * @return true if the argument is positive infinity or + * negative infinity; false otherwise. + */ + public static boolean isInfinite(float f) { + return Float.isInfinite(f); + } + + /** + * Returns true if the specified number is a + * Not-a-Number (NaN) value, false otherwise. + * + *

Note that this method is equivalent to the {@link + * Double#isNaN(double) Double.isNaN} method; the functionality is + * included in this class for convenience. + * + * @param d the value to be tested. + * @return true if the value of the argument is NaN; + * false otherwise. + */ + public static boolean isNaN(double d) { + return Double.isNaN(d); + } + + /** + * Returns true if the specified number is a + * Not-a-Number (NaN) value, false otherwise. + * + *

Note that this method is equivalent to the {@link + * Float#isNaN(float) Float.isNaN} method; the functionality is + * included in this class for convenience. + * + * @param f the value to be tested. + * @return true if the argument is NaN; + * false otherwise. + */ + public static boolean isNaN(float f) { + return Float.isNaN(f); + } + + /** + * Returns true if the unordered relation holds + * between the two arguments. When two floating-point values are + * unordered, one value is neither less than, equal to, nor + * greater than the other. For the unordered relation to be true, + * at least one argument must be a NaN. + * + * @param arg1 the first argument + * @param arg2 the second argument + * @return true if at least one argument is a NaN, + * false otherwise. + */ + public static boolean isUnordered(double arg1, double arg2) { + return isNaN(arg1) || isNaN(arg2); + } + + /** + * Returns true if the unordered relation holds + * between the two arguments. When two floating-point values are + * unordered, one value is neither less than, equal to, nor + * greater than the other. For the unordered relation to be true, + * at least one argument must be a NaN. + * + * @param arg1 the first argument + * @param arg2 the second argument + * @return true if at least one argument is a NaN, + * false otherwise. + */ + public static boolean isUnordered(float arg1, float arg2) { + return isNaN(arg1) || isNaN(arg2); + } + + /** + * Returns unbiased exponent of a double; for + * subnormal values, the number is treated as if it were + * normalized. That is for all finite, non-zero, positive numbers + * x, scalb(x, -ilogb(x)) is + * always in the range [1, 2). + *

+ * Special cases: + *

    + *
  • If the argument is NaN, then the result is 230. + *
  • If the argument is infinite, then the result is 228. + *
  • If the argument is zero, then the result is -(228). + *
+ * + * @param d floating-point number whose exponent is to be extracted + * @return unbiased exponent of the argument. + * @author Joseph D. Darcy + */ + public static int ilogb(double d) { + int exponent = getExponent(d); + + switch (exponent) { + case DoubleConsts.MAX_EXPONENT+1: // NaN or infinity + if( isNaN(d) ) + return (1<<30); // 2^30 + else // infinite value + return (1<<28); // 2^28 + // break; + + case DoubleConsts.MIN_EXPONENT-1: // zero or subnormal + if(d == 0.0) { + return -(1<<28); // -(2^28) + } + else { + long transducer = Double.doubleToRawLongBits(d); + + /* + * To avoid causing slow arithmetic on subnormals, + * the scaling to determine when d's significand + * is normalized is done in integer arithmetic. + * (there must be at least one "1" bit in the + * significand since zero has been screened out. + */ + + // isolate significand bits + transducer &= DoubleConsts.SIGNIF_BIT_MASK; + assert(transducer != 0L); + + // This loop is simple and functional. We might be + // able to do something more clever that was faster; + // e.g. number of leading zero detection on + // (transducer << (# exponent and sign bits). + while (transducer < + (1L << (DoubleConsts.SIGNIFICAND_WIDTH - 1))) { + transducer *= 2; + exponent--; + } + exponent++; + assert( exponent >= + DoubleConsts.MIN_EXPONENT - (DoubleConsts.SIGNIFICAND_WIDTH-1) && + exponent < DoubleConsts.MIN_EXPONENT); + return exponent; + } + // break; + + default: + assert( exponent >= DoubleConsts.MIN_EXPONENT && + exponent <= DoubleConsts.MAX_EXPONENT); + return exponent; + // break; + } + } + + /** + * Returns unbiased exponent of a float; for + * subnormal values, the number is treated as if it were + * normalized. That is for all finite, non-zero, positive numbers + * x, scalb(x, -ilogb(x)) is + * always in the range [1, 2). + *

+ * Special cases: + *

    + *
  • If the argument is NaN, then the result is 230. + *
  • If the argument is infinite, then the result is 228. + *
  • If the argument is zero, then the result is -(228). + *
+ * + * @param f floating-point number whose exponent is to be extracted + * @return unbiased exponent of the argument. + * @author Joseph D. Darcy + */ + public static int ilogb(float f) { + int exponent = getExponent(f); + + switch (exponent) { + case FloatConsts.MAX_EXPONENT+1: // NaN or infinity + if( isNaN(f) ) + return (1<<30); // 2^30 + else // infinite value + return (1<<28); // 2^28 + // break; + + case FloatConsts.MIN_EXPONENT-1: // zero or subnormal + if(f == 0.0f) { + return -(1<<28); // -(2^28) + } + else { + int transducer = Float.floatToRawIntBits(f); + + /* + * To avoid causing slow arithmetic on subnormals, + * the scaling to determine when f's significand + * is normalized is done in integer arithmetic. + * (there must be at least one "1" bit in the + * significand since zero has been screened out. + */ + + // isolate significand bits + transducer &= FloatConsts.SIGNIF_BIT_MASK; + assert(transducer != 0); + + // This loop is simple and functional. We might be + // able to do something more clever that was faster; + // e.g. number of leading zero detection on + // (transducer << (# exponent and sign bits). + while (transducer < + (1 << (FloatConsts.SIGNIFICAND_WIDTH - 1))) { + transducer *= 2; + exponent--; + } + exponent++; + assert( exponent >= + FloatConsts.MIN_EXPONENT - (FloatConsts.SIGNIFICAND_WIDTH-1) && + exponent < FloatConsts.MIN_EXPONENT); + return exponent; + } + // break; + + default: + assert( exponent >= FloatConsts.MIN_EXPONENT && + exponent <= FloatConsts.MAX_EXPONENT); + return exponent; + // break; + } + } + + + /* + * The scalb operation should be reasonably fast; however, there + * are tradeoffs in writing a method to minimize the worst case + * performance and writing a method to minimize the time for + * expected common inputs. Some processors operate very slowly on + * subnormal operands, taking hundreds or thousands of cycles for + * one floating-point add or multiply as opposed to, say, four + * cycles for normal operands. For processors with very slow + * subnormal execution, scalb would be fastest if written entirely + * with integer operations; in other words, scalb would need to + * include the logic of performing correct rounding of subnormal + * values. This could be reasonably done in at most a few hundred + * cycles. However, this approach may penalize normal operations + * since at least the exponent of the floating-point argument must + * be examined. + * + * The approach taken in this implementation is a compromise. + * Floating-point multiplication is used to do most of the work; + * but knowingly multiplying by a subnormal scaling factor is + * avoided. However, the floating-point argument is not examined + * to see whether or not it is subnormal since subnormal inputs + * are assumed to be rare. At most three multiplies are needed to + * scale from the largest to smallest exponent ranges (scaling + * down, at most two multiplies are needed if subnormal scaling + * factors are allowed). However, in this implementation an + * expensive integer remainder operation is avoided at the cost of + * requiring five floating-point multiplies in the worst case, + * which should still be a performance win. + * + * If scaling of entire arrays is a concern, it would probably be + * more efficient to provide a double[] scalb(double[], int) + * version of scalb to avoid having to recompute the needed + * scaling factors for each floating-point value. + */ + + /** + * Return d × + * 2scale_factor rounded as if performed + * by a single correctly rounded floating-point multiply to a + * member of the double value set. See §4.2.3 + * of the Java + * Language Specification for a discussion of floating-point + * value sets. If the exponent of the result is between the + * double's minimum exponent and maximum exponent, + * the answer is calculated exactly. If the exponent of the + * result would be larger than doubles's maximum + * exponent, an infinity is returned. Note that if the result is + * subnormal, precision may be lost; that is, when scalb(x, + * n) is subnormal, scalb(scalb(x, n), -n) may + * not equal x. When the result is non-NaN, the result has + * the same sign as d. + * + *

+ * Special cases: + *

    + *
  • If the first argument is NaN, NaN is returned. + *
  • If the first argument is infinite, then an infinity of the + * same sign is returned. + *
  • If the first argument is zero, then a zero of the same + * sign is returned. + *
+ * + * @param d number to be scaled by a power of two. + * @param scale_factor power of 2 used to scale d + * @return d * 2scale_factor + * @author Joseph D. Darcy + */ + public static double scalb(double d, int scale_factor) { + /* + * This method does not need to be declared strictfp to + * compute the same correct result on all platforms. When + * scaling up, it does not matter what order the + * multiply-store operations are done; the result will be + * finite or overflow regardless of the operation ordering. + * However, to get the correct result when scaling down, a + * particular ordering must be used. + * + * When scaling down, the multiply-store operations are + * sequenced so that it is not possible for two consecutive + * multiply-stores to return subnormal results. If one + * multiply-store result is subnormal, the next multiply will + * round it away to zero. This is done by first multiplying + * by 2 ^ (scale_factor % n) and then multiplying several + * times by by 2^n as needed where n is the exponent of number + * that is a covenient power of two. In this way, at most one + * real rounding error occurs. If the double value set is + * being used exclusively, the rounding will occur on a + * multiply. If the double-extended-exponent value set is + * being used, the products will (perhaps) be exact but the + * stores to d are guaranteed to round to the double value + * set. + * + * It is _not_ a valid implementation to first multiply d by + * 2^MIN_EXPONENT and then by 2 ^ (scale_factor % + * MIN_EXPONENT) since even in a strictfp program double + * rounding on underflow could occur; e.g. if the scale_factor + * argument was (MIN_EXPONENT - n) and the exponent of d was a + * little less than -(MIN_EXPONENT - n), meaning the final + * result would be subnormal. + * + * Since exact reproducibility of this method can be achieved + * without any undue performance burden, there is no + * compelling reason to allow double rounding on underflow in + * scalb. + */ + + // magnitude of a power of two so large that scaling a finite + // nonzero value by it would be guaranteed to over or + // underflow; due to rounding, scaling down takes takes an + // additional power of two which is reflected here + final int MAX_SCALE = DoubleConsts.MAX_EXPONENT + -DoubleConsts.MIN_EXPONENT + + DoubleConsts.SIGNIFICAND_WIDTH + 1; + int exp_adjust = 0; + int scale_increment = 0; + double exp_delta = Double.NaN; + + // Make sure scaling factor is in a reasonable range + + if(scale_factor < 0) { + scale_factor = Math.max(scale_factor, -MAX_SCALE); + scale_increment = -512; + exp_delta = twoToTheDoubleScaleDown; + } + else { + scale_factor = Math.min(scale_factor, MAX_SCALE); + scale_increment = 512; + exp_delta = twoToTheDoubleScaleUp; + } + + // Calculate (scale_factor % +/-512), 512 = 2^9, using + // technique from "Hacker's Delight" section 10-2. + int t = (scale_factor >> 9-1) >>> 32 - 9; + exp_adjust = ((scale_factor + t) & (512 -1)) - t; + + d *= powerOfTwoD(exp_adjust); + scale_factor -= exp_adjust; + + while(scale_factor != 0) { + d *= exp_delta; + scale_factor -= scale_increment; + } + return d; + } + + /** + * Return f × + * 2scale_factor rounded as if performed + * by a single correctly rounded floating-point multiply to a + * member of the float value set. See §4.2.3 + * of the Java + * Language Specification for a discussion of floating-point + * value set. If the exponent of the result is between the + * float's minimum exponent and maximum exponent, the + * answer is calculated exactly. If the exponent of the result + * would be larger than float's maximum exponent, an + * infinity is returned. Note that if the result is subnormal, + * precision may be lost; that is, when scalb(x, n) + * is subnormal, scalb(scalb(x, n), -n) may not equal + * x. When the result is non-NaN, the result has the same + * sign as f. + * + *

+ * Special cases: + *

    + *
  • If the first argument is NaN, NaN is returned. + *
  • If the first argument is infinite, then an infinity of the + * same sign is returned. + *
  • If the first argument is zero, then a zero of the same + * sign is returned. + *
+ * + * @param f number to be scaled by a power of two. + * @param scale_factor power of 2 used to scale f + * @return f * 2scale_factor + * @author Joseph D. Darcy + */ + public static float scalb(float f, int scale_factor) { + // magnitude of a power of two so large that scaling a finite + // nonzero value by it would be guaranteed to over or + // underflow; due to rounding, scaling down takes takes an + // additional power of two which is reflected here + final int MAX_SCALE = FloatConsts.MAX_EXPONENT + -FloatConsts.MIN_EXPONENT + + FloatConsts.SIGNIFICAND_WIDTH + 1; + + // Make sure scaling factor is in a reasonable range + scale_factor = Math.max(Math.min(scale_factor, MAX_SCALE), -MAX_SCALE); + + /* + * Since + MAX_SCALE for float fits well within the double + * exponent range and + float -> double conversion is exact + * the multiplication below will be exact. Therefore, the + * rounding that occurs when the double product is cast to + * float will be the correctly rounded float result. Since + * all operations other than the final multiply will be exact, + * it is not necessary to declare this method strictfp. + */ + return (float)((double)f*powerOfTwoD(scale_factor)); + } + + /** + * Returns the floating-point number adjacent to the first + * argument in the direction of the second argument. If both + * arguments compare as equal the second argument is returned. + * + *

+ * Special cases: + *

    + *
  • If either argument is a NaN, then NaN is returned. + * + *
  • If both arguments are signed zeros, direction + * is returned unchanged (as implied by the requirement of + * returning the second argument if the arguments compare as + * equal). + * + *
  • If start is + * ±Double.MIN_VALUE and direction + * has a value such that the result should have a smaller + * magnitude, then a zero with the same sign as start + * is returned. + * + *
  • If start is infinite and + * direction has a value such that the result should + * have a smaller magnitude, Double.MAX_VALUE with the + * same sign as start is returned. + * + *
  • If start is equal to ± + * Double.MAX_VALUE and direction has a + * value such that the result should have a larger magnitude, an + * infinity with same sign as start is returned. + *
+ * + * @param start starting floating-point value + * @param direction value indicating which of + * start's neighbors or start should + * be returned + * @return The floating-point number adjacent to start in the + * direction of direction. + * @author Joseph D. Darcy + */ + public static double nextAfter(double start, double direction) { + /* + * The cases: + * + * nextAfter(+infinity, 0) == MAX_VALUE + * nextAfter(+infinity, +infinity) == +infinity + * nextAfter(-infinity, 0) == -MAX_VALUE + * nextAfter(-infinity, -infinity) == -infinity + * + * are naturally handled without any additional testing + */ + + // First check for NaN values + if (isNaN(start) || isNaN(direction)) { + // return a NaN derived from the input NaN(s) + return start + direction; + } else if (start == direction) { + return direction; + } else { // start > direction or start < direction + // Add +0.0 to get rid of a -0.0 (+0.0 + -0.0 => +0.0) + // then bitwise convert start to integer. + long transducer = Double.doubleToRawLongBits(start + 0.0d); + + /* + * IEEE 754 floating-point numbers are lexicographically + * ordered if treated as signed- magnitude integers . + * Since Java's integers are two's complement, + * incrementing" the two's complement representation of a + * logically negative floating-point value *decrements* + * the signed-magnitude representation. Therefore, when + * the integer representation of a floating-point values + * is less than zero, the adjustment to the representation + * is in the opposite direction than would be expected at + * first . + */ + if (direction > start) { // Calculate next greater value + transducer = transducer + (transducer >= 0L ? 1L:-1L); + } else { // Calculate next lesser value + assert direction < start; + if (transducer > 0L) + --transducer; + else + if (transducer < 0L ) + ++transducer; + /* + * transducer==0, the result is -MIN_VALUE + * + * The transition from zero (implicitly + * positive) to the smallest negative + * signed magnitude value must be done + * explicitly. + */ + else + transducer = DoubleConsts.SIGN_BIT_MASK | 1L; + } + + return Double.longBitsToDouble(transducer); + } + } + + /** + * Returns the floating-point number adjacent to the first + * argument in the direction of the second argument. If both + * arguments compare as equal, the second argument is returned. + * + *

+ * Special cases: + *

    + *
  • If either argument is a NaN, then NaN is returned. + * + *
  • If both arguments are signed zeros, a float + * zero with the same sign as direction is returned + * (as implied by the requirement of returning the second argument + * if the arguments compare as equal). + * + *
  • If start is + * ±Float.MIN_VALUE and direction + * has a value such that the result should have a smaller + * magnitude, then a zero with the same sign as start + * is returned. + * + *
  • If start is infinite and + * direction has a value such that the result should + * have a smaller magnitude, Float.MAX_VALUE with the + * same sign as start is returned. + * + *
  • If start is equal to ± + * Float.MAX_VALUE and direction has a + * value such that the result should have a larger magnitude, an + * infinity with same sign as start is returned. + *
+ * + * @param start starting floating-point value + * @param direction value indicating which of + * start's neighbors or start should + * be returned + * @return The floating-point number adjacent to start in the + * direction of direction. + * @author Joseph D. Darcy + */ + public static float nextAfter(float start, double direction) { + /* + * The cases: + * + * nextAfter(+infinity, 0) == MAX_VALUE + * nextAfter(+infinity, +infinity) == +infinity + * nextAfter(-infinity, 0) == -MAX_VALUE + * nextAfter(-infinity, -infinity) == -infinity + * + * are naturally handled without any additional testing + */ + + // First check for NaN values + if (isNaN(start) || isNaN(direction)) { + // return a NaN derived from the input NaN(s) + return start + (float)direction; + } else if (start == direction) { + return (float)direction; + } else { // start > direction or start < direction + // Add +0.0 to get rid of a -0.0 (+0.0 + -0.0 => +0.0) + // then bitwise convert start to integer. + int transducer = Float.floatToRawIntBits(start + 0.0f); + + /* + * IEEE 754 floating-point numbers are lexicographically + * ordered if treated as signed- magnitude integers . + * Since Java's integers are two's complement, + * incrementing" the two's complement representation of a + * logically negative floating-point value *decrements* + * the signed-magnitude representation. Therefore, when + * the integer representation of a floating-point values + * is less than zero, the adjustment to the representation + * is in the opposite direction than would be expected at + * first. + */ + if (direction > start) {// Calculate next greater value + transducer = transducer + (transducer >= 0 ? 1:-1); + } else { // Calculate next lesser value + assert direction < start; + if (transducer > 0) + --transducer; + else + if (transducer < 0 ) + ++transducer; + /* + * transducer==0, the result is -MIN_VALUE + * + * The transition from zero (implicitly + * positive) to the smallest negative + * signed magnitude value must be done + * explicitly. + */ + else + transducer = FloatConsts.SIGN_BIT_MASK | 1; + } + + return Float.intBitsToFloat(transducer); + } + } + + /** + * Returns the floating-point value adjacent to d in + * the direction of positive infinity. This method is + * semantically equivalent to nextAfter(d, + * Double.POSITIVE_INFINITY); however, a nextUp + * implementation may run faster than its equivalent + * nextAfter call. + * + *

Special Cases: + *

    + *
  • If the argument is NaN, the result is NaN. + * + *
  • If the argument is positive infinity, the result is + * positive infinity. + * + *
  • If the argument is zero, the result is + * Double.MIN_VALUE + * + *
+ * + * @param d starting floating-point value + * @return The adjacent floating-point value closer to positive + * infinity. + * @author Joseph D. Darcy + */ + public static double nextUp(double d) { + if( isNaN(d) || d == Double.POSITIVE_INFINITY) + return d; + else { + d += 0.0d; + return Double.longBitsToDouble(Double.doubleToRawLongBits(d) + + ((d >= 0.0d)?+1L:-1L)); + } + } + + /** + * Returns the floating-point value adjacent to f in + * the direction of positive infinity. This method is + * semantically equivalent to nextAfter(f, + * Double.POSITIVE_INFINITY); however, a nextUp + * implementation may run faster than its equivalent + * nextAfter call. + * + *

Special Cases: + *

    + *
  • If the argument is NaN, the result is NaN. + * + *
  • If the argument is positive infinity, the result is + * positive infinity. + * + *
  • If the argument is zero, the result is + * Float.MIN_VALUE + * + *
+ * + * @param f starting floating-point value + * @return The adjacent floating-point value closer to positive + * infinity. + * @author Joseph D. Darcy + */ + public static float nextUp(float f) { + if( isNaN(f) || f == FloatConsts.POSITIVE_INFINITY) + return f; + else { + f += 0.0f; + return Float.intBitsToFloat(Float.floatToRawIntBits(f) + + ((f >= 0.0f)?+1:-1)); + } + } + + /** + * Returns the floating-point value adjacent to d in + * the direction of negative infinity. This method is + * semantically equivalent to nextAfter(d, + * Double.NEGATIVE_INFINITY); however, a + * nextDown implementation may run faster than its + * equivalent nextAfter call. + * + *

Special Cases: + *

    + *
  • If the argument is NaN, the result is NaN. + * + *
  • If the argument is negative infinity, the result is + * negative infinity. + * + *
  • If the argument is zero, the result is + * -Double.MIN_VALUE + * + *
+ * + * @param d starting floating-point value + * @return The adjacent floating-point value closer to negative + * infinity. + * @author Joseph D. Darcy + */ + public static double nextDown(double d) { + if( isNaN(d) || d == Double.NEGATIVE_INFINITY) + return d; + else { + if (d == 0.0) + return -Double.MIN_VALUE; + else + return Double.longBitsToDouble(Double.doubleToRawLongBits(d) + + ((d > 0.0d)?-1L:+1L)); + } + } + + /** + * Returns the floating-point value adjacent to f in + * the direction of negative infinity. This method is + * semantically equivalent to nextAfter(f, + * Float.NEGATIVE_INFINITY); however, a + * nextDown implementation may run faster than its + * equivalent nextAfter call. + * + *

Special Cases: + *

    + *
  • If the argument is NaN, the result is NaN. + * + *
  • If the argument is negative infinity, the result is + * negative infinity. + * + *
  • If the argument is zero, the result is + * -Float.MIN_VALUE + * + *
+ * + * @param f starting floating-point value + * @return The adjacent floating-point value closer to negative + * infinity. + * @author Joseph D. Darcy + */ + public static double nextDown(float f) { + if( isNaN(f) || f == Float.NEGATIVE_INFINITY) + return f; + else { + if (f == 0.0f) + return -Float.MIN_VALUE; + else + return Float.intBitsToFloat(Float.floatToRawIntBits(f) + + ((f > 0.0f)?-1:+1)); + } + } + + /** + * Returns the first floating-point argument with the sign of the + * second floating-point argument. For this method, a NaN + * sign argument is always treated as if it were + * positive. + * + * @param magnitude the parameter providing the magnitude of the result + * @param sign the parameter providing the sign of the result + * @return a value with the magnitude of magnitude + * and the sign of sign. + * @author Joseph D. Darcy + * @since 1.5 + */ + public static double copySign(double magnitude, double sign) { + return rawCopySign(magnitude, (isNaN(sign)?1.0d:sign)); + } + + /** + * Returns the first floating-point argument with the sign of the + * second floating-point argument. For this method, a NaN + * sign argument is always treated as if it were + * positive. + * + * @param magnitude the parameter providing the magnitude of the result + * @param sign the parameter providing the sign of the result + * @return a value with the magnitude of magnitude + * and the sign of sign. + * @author Joseph D. Darcy + */ + public static float copySign(float magnitude, float sign) { + return rawCopySign(magnitude, (isNaN(sign)?1.0f:sign)); + } + + /** + * Returns the size of an ulp of the argument. An ulp of a + * double value is the positive distance between this + * floating-point value and the double value next + * larger in magnitude. Note that for non-NaN x, + * ulp(-x) == ulp(x). + * + *

Special Cases: + *

    + *
  • If the argument is NaN, then the result is NaN. + *
  • If the argument is positive or negative infinity, then the + * result is positive infinity. + *
  • If the argument is positive or negative zero, then the result is + * Double.MIN_VALUE. + *
  • If the argument is ±Double.MAX_VALUE, then + * the result is equal to 2971. + *
+ * + * @param d the floating-point value whose ulp is to be returned + * @return the size of an ulp of the argument + * @author Joseph D. Darcy + * @since 1.5 + */ + public static double ulp(double d) { + int exp = getExponent(d); + + switch(exp) { + case DoubleConsts.MAX_EXPONENT+1: // NaN or infinity + return Math.abs(d); + // break; + + case DoubleConsts.MIN_EXPONENT-1: // zero or subnormal + return Double.MIN_VALUE; + // break + + default: + assert exp <= DoubleConsts.MAX_EXPONENT && exp >= DoubleConsts.MIN_EXPONENT; + + // ulp(x) is usually 2^(SIGNIFICAND_WIDTH-1)*(2^ilogb(x)) + exp = exp - (DoubleConsts.SIGNIFICAND_WIDTH-1); + if (exp >= DoubleConsts.MIN_EXPONENT) { + return powerOfTwoD(exp); + } + else { + // return a subnormal result; left shift integer + // representation of Double.MIN_VALUE appropriate + // number of positions + return Double.longBitsToDouble(1L << + (exp - (DoubleConsts.MIN_EXPONENT - (DoubleConsts.SIGNIFICAND_WIDTH-1)) )); + } + // break + } + } + + /** + * Returns the size of an ulp of the argument. An ulp of a + * float value is the positive distance between this + * floating-point value and the float value next + * larger in magnitude. Note that for non-NaN x, + * ulp(-x) == ulp(x). + * + *

Special Cases: + *

    + *
  • If the argument is NaN, then the result is NaN. + *
  • If the argument is positive or negative infinity, then the + * result is positive infinity. + *
  • If the argument is positive or negative zero, then the result is + * Float.MIN_VALUE. + *
  • If the argument is ±Float.MAX_VALUE, then + * the result is equal to 2104. + *
+ * + * @param f the floating-point value whose ulp is to be returned + * @return the size of an ulp of the argument + * @author Joseph D. Darcy + * @since 1.5 + */ + public static float ulp(float f) { + int exp = getExponent(f); + + switch(exp) { + case FloatConsts.MAX_EXPONENT+1: // NaN or infinity + return Math.abs(f); + // break; + + case FloatConsts.MIN_EXPONENT-1: // zero or subnormal + return FloatConsts.MIN_VALUE; + // break + + default: + assert exp <= FloatConsts.MAX_EXPONENT && exp >= FloatConsts.MIN_EXPONENT; + + // ulp(x) is usually 2^(SIGNIFICAND_WIDTH-1)*(2^ilogb(x)) + exp = exp - (FloatConsts.SIGNIFICAND_WIDTH-1); + if (exp >= FloatConsts.MIN_EXPONENT) { + return powerOfTwoF(exp); + } + else { + // return a subnormal result; left shift integer + // representation of FloatConsts.MIN_VALUE appropriate + // number of positions + return Float.intBitsToFloat(1 << + (exp - (FloatConsts.MIN_EXPONENT - (FloatConsts.SIGNIFICAND_WIDTH-1)) )); + } + // break + } + } + + /** + * Returns the signum function of the argument; zero if the argument + * is zero, 1.0 if the argument is greater than zero, -1.0 if the + * argument is less than zero. + * + *

Special Cases: + *

    + *
  • If the argument is NaN, then the result is NaN. + *
  • If the argument is positive zero or negative zero, then the + * result is the same as the argument. + *
+ * + * @param d the floating-point value whose signum is to be returned + * @return the signum function of the argument + * @author Joseph D. Darcy + * @since 1.5 + */ + public static double signum(double d) { + return (d == 0.0 || isNaN(d))?d:copySign(1.0, d); + } + + /** + * Returns the signum function of the argument; zero if the argument + * is zero, 1.0f if the argument is greater than zero, -1.0f if the + * argument is less than zero. + * + *

Special Cases: + *

    + *
  • If the argument is NaN, then the result is NaN. + *
  • If the argument is positive zero or negative zero, then the + * result is the same as the argument. + *
+ * + * @param f the floating-point value whose signum is to be returned + * @return the signum function of the argument + * @author Joseph D. Darcy + * @since 1.5 + */ + public static float signum(float f) { + return (f == 0.0f || isNaN(f))?f:copySign(1.0f, f); + } + +} diff --git a/runtime/classpath.cs b/runtime/classpath.cs index 084c8727..a7c4bca3 100644 --- a/runtime/classpath.cs +++ b/runtime/classpath.cs @@ -369,109 +369,6 @@ namespace IKVM.NativeCode.java } } - public class VMDouble - { - public static double parseDouble(string s) - { - if(s.Length == 0) goto error; - int first = 0; - while(s[first] <= ' ') - { - first++; - if(first == s.Length) goto error; - } - int last = s.Length - 1; - while(s[last] <= ' ') - { - last--; - if(first > last) goto error; - } - bool sign = false; - if(s[first] == '-') - { - sign = true; - first++; - if(first > last) goto error; - } - else if(s[first] == '+') - { - first++; - if(first > last) goto error; - } - if(last - first == 7 && string.CompareOrdinal(s, first, "Infinity", 0, 8) == 0) - { - return sign ? double.NegativeInfinity : double.PositiveInfinity; - } - if(last - first == 2 && string.CompareOrdinal(s, first, "NaN", 0, 3) == 0) - { - return double.NaN; - } - // Java allows 'f' or 'd' at the end - if("dfDF".IndexOf(s[last]) >= 0) - { - last--; - if(first > last) goto error; - } - bool dot = false; - bool exp = false; - for(int i = first; i <= last; i++) - { - char c = s[i]; - if(c >= '0' && c <= '9') - { - // ok - } - else if(c == '.') - { - if(dot || exp) goto error; - dot = true; - } - else if(c == 'e' || c == 'E') - { - if(i == first || i == last || exp) goto error; - if(s[i + 1] == '-' || s[i + 1] == '+') - { - i++; - if(i == last) goto error; - } - exp = true; - } - else goto error; - } - if(first != 0 || last != s.Length - 1) - { - s = s.Substring(first, last - first + 1); - } - try - { - // I doubt that this is fully correct, but since we don't implement FP properly, - // we can probably get away with this as well. - double d = double.Parse(s, System.Globalization.CultureInfo.InvariantCulture); - if(d == 0.0 && (s == "4.9E-324" || s == "4.9e-324")) - { - // HACK special support for Double.MIN_VALUE - return sign ? -double.Epsilon : double.Epsilon; - } - return sign ? -d : d; - } - catch(OverflowException) - { - return sign ? double.NegativeInfinity : double.PositiveInfinity; - } - catch(FormatException x) - { - // this can't happen, since we already validated the format of the string - System.Diagnostics.Debug.Fail(x.ToString()); - } - error: -#if !FIRST_PASS - throw new NumberFormatException(string.Format("For input string: \"{0}\"", s)); -#else - return 0; -#endif - } - } - public class VMSystem { public static void arraycopy(object src, int srcStart, object dest, int destStart, int len)