C-XSC - A C++ Class Library for Extended Scientific Computing
2.5.4
|
00001 /* 00002 ** CXSC is a C++ library for eXtended Scientific Computing (V 2.5.4) 00003 ** 00004 ** Copyright (C) 1990-2000 Institut fuer Angewandte Mathematik, 00005 ** Universitaet Karlsruhe, Germany 00006 ** (C) 2000-2014 Wiss. Rechnen/Softwaretechnologie 00007 ** Universitaet Wuppertal, Germany 00008 ** 00009 ** This library is free software; you can redistribute it and/or 00010 ** modify it under the terms of the GNU Library General Public 00011 ** License as published by the Free Software Foundation; either 00012 ** version 2 of the License, or (at your option) any later version. 00013 ** 00014 ** This library is distributed in the hope that it will be useful, 00015 ** but WITHOUT ANY WARRANTY; without even the implied warranty of 00016 ** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00017 ** Library General Public License for more details. 00018 ** 00019 ** You should have received a copy of the GNU Library General Public 00020 ** License along with this library; if not, write to the Free 00021 ** Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA 00022 */ 00023 00024 /* CVS $Id: realio.cpp,v 1.31 2014/01/30 17:23:48 cxsc Exp $ */ 00025 00026 #include <cstring> 00027 #include "real.hpp" 00028 #include "ioflags.hpp" 00029 #include "RtsFunc.h" 00030 #include "dot.hpp" 00031 00032 namespace cxsc { 00033 00034 int d_init_dm (void); 00035 void d_outp(char *buffer, Dotprecision c, 00036 int FormatFlag, int FracDigits, int rnd, 00037 int *length); 00038 00039 #if _WIN32 00040 extern __declspec(thread) char *dm; 00041 #elif __APPLE__ && !CXSC_FORCE_TLS 00042 extern char *dm; 00043 #else 00044 extern __thread char *dm; 00045 #endif 00046 00047 00048 char* GetHexDigits (char* s, a_btyp& x, int count); 00049 00050 int r_outpx (char *buffer, a_real s,a_intg FormatFlag, a_intg FracDigits, a_intg rnd, a_intg *length); 00051 00052 #define WASGOOD 0 00053 #define ISINFINITYPLUS -1 00054 #define ISINFINITYMINUS -2 00055 #define ISQUIETNAN -3 00056 #define ISSIGNALINGNAN -4 00057 00058 std::string realToHex(const real& a) 00059 { 00060 a_btyp* b = (a_btyp*) &a; 00061 int i; 00062 char cs[256]; 00063 00064 for (i=0; i < realwidth-19; i++) cs[i] = ' '; 00065 cs[i] = 0; 00066 sprintf (&cs[strlen(cs)], "%s", ((b[HIGHREAL] & 0x80000000L) ? "-" : "+")); 00067 sprintf (&cs[strlen(cs)],"%c", '1'); 00068 sprintf (&cs[strlen(cs)],"%05lX", (b[HIGHREAL] & 0xFFFFFL)); 00069 sprintf (&cs[strlen(cs)],"%08lX",(long unsigned int)(b[LOWREAL])); 00070 sprintf (&cs[strlen(cs)],"e%03X",(unsigned int)((b[HIGHREAL] >> 20) & 0x7FF)); 00071 return string(cs); 00072 } 00073 00074 std::string & operator<<(std::string &s, const real& a) throw() 00075 { 00076 #if _WIN32 00077 static __declspec(thread) char cs[256]; 00078 #elif __APPLE__ && !CXSC_FORCE_TLS 00079 static char cs[256]; 00080 #else 00081 static __thread char cs[256]; 00082 #endif 00083 00084 if (ioflags.isset(IOFlags::hex)) 00085 { 00086 a_btyp* b = (a_btyp*) &a; 00087 int i; 00088 char cs[256]; 00089 00090 for (i=0; i < realwidth-19; i++) cs[i] = ' '; 00091 cs[i] = 0; 00092 sprintf (&cs[strlen(cs)], "%s", ((b[HIGHREAL] & 0x80000000L) ? "-" : "+")); 00093 sprintf (&cs[strlen(cs)],"%c", '1'); 00094 sprintf (&cs[strlen(cs)],"%05lX", (b[HIGHREAL] & 0xFFFFFL)); 00095 sprintf (&cs[strlen(cs)],"%08lX",(long unsigned int)(b[LOWREAL])); 00096 sprintf (&cs[strlen(cs)],"e%03X",(unsigned int)((b[HIGHREAL] >> 20) & 0x7FF)); 00097 s+=cs; 00098 } else if (ioflags.isset(IOFlags::rndnone)) 00099 { 00100 if (IsSignalingNaN(a)) s+="<SignallingNaN>"; 00101 else if (IsQuietNaN(a)) s+="<QuietNaN>"; 00102 else if (IsInfinity(a)) s+="<Infinity>"; 00103 else 00104 { 00105 if (realdigits && realwidth) 00106 sprintf (cs, "%*.*g", realwidth, realdigits, a.w); // no need for "lg" 00107 else if (realwidth) 00108 sprintf (cs, "%*g", realwidth, a.w); 00109 else 00110 sprintf (cs, "%g", a.w); 00111 s+=cs; 00112 } 00113 } else 00114 { 00115 rndtype rnd; 00116 a_intg length, formatflag, addblanks; 00117 a_intg digits = realdigits; 00118 char *str; 00119 00120 if (d_init_dm () == -1) 00121 { 00122 // throw 00123 // errmon (ERR_ALL(NOMOREMEMORY)); 00124 // errmon (ERR_ALL(NOCONTINUEPOSSIBLE)); 00125 } 00126 00127 if (ioflags.isset(IOFlags::rndup)) 00128 rnd = RND_UP; 00129 else if (ioflags.isset(IOFlags::rnddown)) 00130 rnd = RND_DOWN; 00131 else 00132 rnd = RND_NEXT; 00133 00134 if (ioflags.isset(IOFlags::variable)) 00135 formatflag = realwidth; 00136 else if (ioflags.isset(IOFlags::varfixwidth)) 00137 formatflag = realwidth, digits = -digits; 00138 else 00139 formatflag = (ioflags.isset(IOFlags::fixed)) ? 0 : -1; 00140 00141 switch (r_outpx (dm, a.w, formatflag, digits, rnd, &length)) 00142 { 00143 case WASGOOD: 00144 dm[length] = 0; 00145 str = dm; 00146 if (*str == '+') 00147 { 00148 if (ioflags.isset(IOFlags::blank)) 00149 *str = ' '; 00150 else if (ioflags.isset(IOFlags::noblank)) 00151 str++; 00152 } 00153 break; 00154 case ISINFINITYPLUS: 00155 str = (char*)"<+Infinity>"; 00156 break; 00157 case ISINFINITYMINUS: 00158 str = (char*)"<-Infinity>"; 00159 break; 00160 case ISQUIETNAN: 00161 str = (char*)"<QuietNaN>"; 00162 break; 00163 case ISSIGNALINGNAN: 00164 str = (char*)"<SignalingNaN>"; 00165 break; 00166 default: 00167 str = (char*)"<ERROR>"; 00168 break; 00169 } 00170 length = strlen(str); 00171 addblanks = (length < realwidth) ? realwidth - length : 0; 00172 00173 if (ioflags.isset(IOFlags::rightjust)) 00174 for (;addblanks; addblanks--) 00175 s+= ' '; 00176 00177 s+=str; 00178 00179 for (;addblanks; addblanks--) 00180 s+= ' '; 00181 } 00182 return s; 00183 } 00184 00185 std::ostream & operator <<(std::ostream &o,const real &a) throw() 00186 { 00187 std::string s=""; 00188 s << a; 00189 o << s; 00190 return o; 00191 } 00192 00193 std::string & operator>> (std::string & str, real& a) throw() 00194 { 00195 char *s=new char[str.size()+1]; 00196 char *orgs=s; 00197 strcpy(s,str.c_str()); 00198 00199 if (ioflags.isset(IOFlags::hex)) 00200 { 00201 a_btyp* b = (a_btyp*) &a; 00202 a_btyp x; 00203 00204 b[0] = b[1] = 0; 00205 s = cskipwhitespaces (s); 00206 if (*s == '-') 00207 { 00208 b[HIGHREAL] |= 0x80000000L; 00209 s++; 00210 } else if (*s == '+') 00211 s++; 00212 00213 // ---------------------------------------------------- 00214 // es wird ohne Pruefung folgendes Format vorausgestzt : 00215 // 1xxxxxxxxxxxxxeXXX 00216 // wobei x = HexDigits der Mantisse 00217 // X = HexDigits des Exponents 00218 // saemmtliche HexDigits muessen gross geschrieben sein ! 00219 00220 if (*s) 00221 s++; // skip '1' 00222 s = GetHexDigits (s, x, 5); 00223 b[HIGHREAL] |= x; // get mantissa 00224 s = GetHexDigits (s, x, 8); 00225 b[LOWREAL] = x; 00226 if (*s) 00227 s++; // skip 'e' 00228 s = GetHexDigits (s, x, 3); 00229 b[HIGHREAL] |= x << 20; // get exponent 00230 if (*s) 00231 s++; // skip at least one more char 00232 } else 00233 { 00234 rndtype rndfl; 00235 00236 if (ioflags.isset(IOFlags::rndup)) 00237 rndfl = RND_UP; 00238 else if (ioflags.isset(IOFlags::rnddown)) 00239 rndfl = RND_DOWN; 00240 else 00241 rndfl = RND_NEXT; 00242 00243 str=s; 00244 dotprecision dot; 00245 str >> dot; 00246 strcpy(s,str.c_str()); // Ooooooooohhhh... :(( 00247 a = rnd (dot, rndfl); 00248 } 00249 str=s; 00250 delete [] orgs; 00251 return str; 00252 } 00253 void operator >>(const char *a,real &b) throw() 00254 { 00255 std::string c(a); 00256 c>>b; 00257 } 00258 void operator >>(const string &a,real &b) throw() 00259 { 00260 std::string c(a); 00261 c>>b; 00262 } 00263 std::istream& operator>> (std::istream& s, real& a) throw() 00264 { 00265 if (ioflags.isset(IOFlags::hex)) 00266 { 00267 char inp[20]; 00268 int i; 00269 char c; 00270 00271 // !! There are no checks about the right input-format, it is assumed 00272 // the fixed format : S1xxxxxxxxxxxxxeXXX 00273 // whereby S means the sign 00274 // '1' is the leading implicit binary one of the mantissa 00275 // x are the mantissa hexdigits 00276 // 'e' signals the beginning of the exponent 00277 // X are the exponent hexdigits 00278 00279 // skip white spaces and 00280 // get the sign, the 14 Mantissadigits and the exponent 00281 // (total 19 chars) 00282 00283 c = skipwhitespaces (s); 00284 for (i=0; i < 19; i++) 00285 { 00286 inp[i] = c; 00287 if (s.good()) s.get(c); else c = '\0'; 00288 } 00289 inp[i] = '\0'; 00290 00291 inp >> a; 00292 } else 00293 { 00294 rndtype rndfl; 00295 00296 if (ioflags.isset(IOFlags::rndup)) 00297 rndfl = RND_UP; 00298 else if (ioflags.isset(IOFlags::rnddown)) 00299 rndfl = RND_DOWN; 00300 else 00301 rndfl = RND_NEXT; 00302 00303 dotprecision dot; 00304 s >> dot; 00305 a = rnd (dot, rndfl); 00306 } 00307 return s; 00308 } 00309 00310 //---------------------------------------------------------------------------- 00311 // GetHexDigits 00312 // 00313 // Interpretiert die naechsten "count" Zeichen aus dem String "s" 00314 // als HexDigits, der binaere Wert wird dann in "x" zurueckgegeben. 00315 // Als Returnwert wird ein Zeiger auf die Stringposition nach den 00316 // HexDigits zurueckgegeben. 00317 00318 char* GetHexDigits (char* s, a_btyp& x, int count) 00319 { 00320 int i, c; 00321 00322 for (x=0,i=0; i < count && *s; s++,i++) { 00323 if ((c = *s) >= 'A') c -= 'A' - 10; else c -= '0'; 00324 if (c < 0 || c > 0xF) c = 0; 00325 x = (x << 4) | c; 00326 } 00327 return s; 00328 } 00329 00330 } // namespace cxsc 00331 00332 /****************************************************************/ 00333 /* */ 00334 /* Filename : r_outpx.c */ 00335 /* */ 00336 /* Entries : void r_outp */ 00337 /* (buffer,s,FormatFlag, */ 00338 /* FracDigits,rnd,length) */ 00339 /* char *buffer; */ 00340 /* a_real s; */ 00341 /* a_intg rnd,FormatFlag,FracDigits; */ 00342 /* a_intg *length; */ 00343 /* */ 00344 /* Arguments : buffer - output buffer holding string */ 00345 /* s - IEEE value */ 00346 /* FormatFlag - format selection */ 00347 /* -1 = scientific format */ 00348 /* 0 = fixed format */ 00349 /* > 0 = variable format */ 00350 /* value is the total field */ 00351 /* width */ 00352 /* FracDigits - number of fraction digits*/ 00353 /* rnd - rounding mode */ 00354 /* -1 = round downwards */ 00355 /* 0 = round to nearest */ 00356 /* 1 = round upwards */ 00357 /* length - size of buffer string */ 00358 /* */ 00359 /* Description : Convert an IEEE double format number */ 00360 /* to a character string. */ 00361 /* */ 00362 /****************************************************************/ 00363 00364 /*#ifndef ALL_IN_ONE 00365 #ifdef AIX 00366 #include "/u/p88c/runtime/o_defs.h" 00367 #else 00368 #include "o_defs.h" 00369 #endif 00370 #define local 00371 #endif*/ 00372 00373 namespace cxsc { 00374 00375 #define WASGOOD 0 00376 #define ISINFINITYPLUS -1 00377 #define ISINFINITYMINUS -2 00378 #define ISQUIETNAN -3 00379 #define ISSIGNALINGNAN -4 00380 00381 #define MANT_INFINITY(a) ((a)[0]==HIDDEN_BIT && (a)[1]==ZERO) 00382 #define SIGNALING(a) ((a) & SIGNAL_BIT) 00383 #define SIGNAL_BIT ((a_btyp)0x00080000L) 00384 00385 int r_outpx(char *buffer, a_real s, a_intg FormatFlag, 00386 a_intg FracDigits, a_intg rnd, a_intg *length) 00387 { 00388 a_intg ActWidth,DecPlaces,expo,IntDigits,MinNumChars; 00389 a_intg dexpo,digits,bdp,k,l,addpoint; 00390 a_intg HoldWidth = (FracDigits < 0); 00391 a_bool vz,zero; 00392 a_btyp mant[BSIZE]; 00393 00394 if (HoldWidth) 00395 FracDigits = -FracDigits; 00396 *length = 0; 00397 00398 zero = b_deko(s,&expo,mant,&vz); 00399 00400 /* */ 00401 /* infinity or NaN */ 00402 /* */ 00403 if (expo>EXPO_MAX) 00404 { 00405 if (MANT_INFINITY(mant)) 00406 { /* Infinity */ 00407 k = (vz) ? ISINFINITYMINUS : ISINFINITYPLUS; 00408 } else 00409 { /* NaN */ 00410 k = (SIGNALING(mant[0])) ? ISSIGNALINGNAN : ISQUIETNAN; 00411 } 00412 return k; 00413 } 00414 00415 for (k=D_U_RATIO;k<BSIZE;k++) 00416 mant[k] = ZERO; 00417 00418 if (vz) 00419 rnd = -rnd; 00420 00421 if (FormatFlag > 0 && FormatFlag-FracDigits <= 2) 00422 { 00423 FormatFlag = FracDigits+3; 00424 } 00425 if (!zero && FormatFlag > 0) 00426 { 00427 dexpo = (a_intg) ((expo*30103L)/100000L); 00428 if (dexpo < -((FracDigits+1)/2) || 00429 (dexpo > 0 && dexpo >= FormatFlag-FracDigits-2)) 00430 { 00431 FormatFlag = -1; 00432 if (HoldWidth) 00433 { 00434 FracDigits = (FracDigits < 5) ? 0 : FracDigits-5; 00435 } 00436 } 00437 } 00438 00439 if (expo>800) 00440 bdp = (BUFFERSIZE-8)-(FormatFlag == -1 ? 0 : FracDigits); 00441 else if (expo<-800) 00442 bdp = 8; 00443 else 00444 bdp = B_D_P; 00445 00446 addpoint = (FracDigits == 0) ? 0 : 1; 00447 /* */ 00448 /* floating-point representation */ 00449 /* */ 00450 if (FormatFlag == -1) 00451 { 00452 DecPlaces = FracDigits; 00453 ActWidth = DecPlaces+ExpDigits+4+addpoint; 00454 /* */ 00455 /* number is zero */ 00456 /* */ 00457 if (zero) 00458 { 00459 *buffer++ = (vz) ? '-' : '+'; 00460 *buffer++ = '0'; 00461 if (FracDigits) 00462 { 00463 *buffer++ = XSC_DECIMAL_POINT; 00464 for (k=0;k<DecPlaces;k++) *buffer++ = '0'; 00465 } 00466 *buffer++ = EXPONENT_E; 00467 *buffer++ = PLUS_SIGN; 00468 for (k=0;k<ExpDigits;k++) 00469 *buffer++ = '0'; 00470 00471 *length = ActWidth; 00472 00473 return WASGOOD; 00474 } 00475 /* */ 00476 /* determine output string */ 00477 /* */ 00478 digits = DecPlaces+(1+2); 00479 dexpo = -1; 00480 b_out(mant,expo,digits,buffer,&bdp,&dexpo); 00481 00482 if (dexpo>0 && dexpo>DecPlaces+2) 00483 digits = dexpo+1; 00484 00485 b_rnd(rnd,buffer,digits,DecPlaces+1,&bdp,&dexpo); 00486 00487 if (FracDigits) 00488 { 00489 buffer[bdp-dexpo-1] = buffer[bdp-dexpo]; 00490 buffer[bdp-dexpo] = XSC_DECIMAL_POINT; 00491 } 00492 buffer[bdp-dexpo+DecPlaces+1] = EXPONENT_E; 00493 buffer[bdp-dexpo+DecPlaces+2] = (dexpo<0) ? MINUS_SIGN : PLUS_SIGN; 00494 expo = (dexpo<0) ? -dexpo : dexpo; 00495 for (k=ExpDigits;k>0;k--) 00496 { 00497 buffer[bdp-dexpo+DecPlaces+2+k] = expo%10+'0'; 00498 expo /= 10; 00499 } 00500 *length = 2+addpoint+DecPlaces+2+ExpDigits; 00501 00502 l = bdp - 1 - addpoint - dexpo; 00503 buffer[l] = (vz) ? '-' : '+'; 00504 for (k=0;k<*length;k++,l++) 00505 buffer[k] = buffer[l]; 00506 } else 00507 /* */ 00508 /* fixed-point representation */ 00509 /* */ 00510 { 00511 /* */ 00512 /* number is zero */ 00513 /* */ 00514 if (zero) 00515 { 00516 *length = FracDigits+2+addpoint; 00517 *buffer++ = (vz) ? '-' : '+'; 00518 *buffer++ = '0'; 00519 if (FracDigits) 00520 { 00521 *buffer++ = XSC_DECIMAL_POINT; 00522 for (k=0;k<FracDigits;k++) 00523 *buffer++ = '0'; 00524 } 00525 return WASGOOD; 00526 } 00527 00528 /* estimate number of decimal digits */ 00529 if (expo>=0) 00530 { 00531 IntDigits = ((expo+1)*61)/200+1; 00532 } else 00533 { 00534 IntDigits = 0; 00535 dexpo = 0; 00536 } 00537 00538 /* fill fractional part with zeros */ 00539 for (k=0;k<=FracDigits+2;k++) 00540 buffer[bdp+k] = '0'; 00541 00542 digits = IntDigits+FracDigits+2; 00543 b_out(mant,expo,digits,buffer,&bdp,&dexpo); 00544 00545 /* correct setting of IntDigits */ 00546 if (expo>=0) 00547 { 00548 IntDigits = dexpo+1; 00549 } else 00550 { 00551 IntDigits = 1; 00552 digits++; 00553 } 00554 00555 b_rnd(rnd,buffer,digits,IntDigits+FracDigits,&bdp,&dexpo); 00556 00557 /* correct setting of IntDigits after rounding */ 00558 IntDigits = dexpo+1; 00559 00560 /* value is zero after rounding; */ 00561 if (vz) 00562 { 00563 for (k=(bdp+1)-IntDigits;k<(bdp+1)+FracDigits;k++) 00564 { 00565 if (buffer[k]!='0') 00566 break; 00567 } 00568 if (k==(bdp+1)+FracDigits) 00569 vz = FALSE; 00570 } 00571 00572 MinNumChars = IntDigits+FracDigits+1+addpoint; 00573 00574 if (FracDigits) 00575 { 00576 /* generate decimal point */ 00577 for (k=bdp-IntDigits;k<bdp;k++) 00578 buffer[k] = buffer[k+1]; 00579 buffer[bdp] = XSC_DECIMAL_POINT; 00580 } 00581 00582 *length = MinNumChars; 00583 /* */ 00584 /* sign of number */ 00585 /* */ 00586 l = bdp - 1 - addpoint - dexpo; 00587 buffer[l] = (vz) ? '-' : '+'; 00588 for (k=0;k<*length;k++,l++) 00589 buffer[k] = buffer[l]; 00590 } 00591 00592 return WASGOOD; 00593 } 00594 00595 } // namespace cxsc 00596