~drizzle-trunk/drizzle/development : revision 1089

42

43

using namespace std;

44

45

/**

46

Appears to suffice to not call malloc() in most cases.

47

@todo

48

see if it is possible to get rid of malloc().

49

this constant is sufficient to avoid malloc() on all inputs I have tried.

50

*/

51

#define DTOA_BUFF_SIZE (420 * sizeof(void *))

52

53

45

/* Magic value returned by dtoa() to indicate overflow */

54

46

#define DTOA_OVERFLOW 9999

55

47

56

static double my_strtod_int(const char *, char **, int *, char *, size_t);

57

static char *dtoa(double, int, int, int *, int *, char **, char *, size_t);

58

static void dtoa_free(char *, char *, size_t);

48

static double my_strtod_int(const char *, char **, int *);

49

static char *dtoa(double, int, int, int *, int *, char **);

59

50

60

51

/**

61

52

@brief

92

83

{

93

84

int decpt, sign, len, i;

94

85

char *res, *src, *end, *dst= to;

95

double alignedbuf[(DTOA_BUFF_SIZE+1)/sizeof(double)];

96

char *buf= (char*)alignedbuf;

97

86

assert(precision >= 0 && precision < NOT_FIXED_DEC && to != NULL);

98

87

99

res= dtoa(x, 5, precision, &decpt, &sign, &end, buf, sizeof(buf));

88

res= dtoa(x, 5, precision, &decpt, &sign, &end);

100

89

101

90

if (decpt == DTOA_OVERFLOW)

102

91

{

103

dtoa_free(res, buf, sizeof(buf));

92

free(res);

104

93

*to++= '0';

105

94

*to= '\0';

106

95

if (error != NULL)

144

133

if (error != NULL)

145

134

*error= false;

146

135

147

dtoa_free(res, buf, sizeof(buf));

136

free(res);

148

137

149

138

return dst - to;

150

139

}

217

206

{

218

207

int decpt, sign, len, exp_len;

219

208

char *res, *src, *end, *dst= to, *dend= dst + width;

220

double alignedbuf[(DTOA_BUFF_SIZE+1)/sizeof(double)];

221

char *buf= (char*)alignedbuf;

222

209

bool have_space, force_e_format;

223

210

assert(width > 0 && to != NULL);

224

211

227

214

width--;

228

215

229

216

res= dtoa(x, 4, type == MY_GCVT_ARG_DOUBLE ? width : min(width, FLT_DIG),

230

&decpt, &sign, &end, buf, sizeof(buf));

217

&decpt, &sign, &end);

218

231

219

if (decpt == DTOA_OVERFLOW)

232

220

{

233

dtoa_free(res, buf, sizeof(buf));

221

free(res);

234

222

*to++= '0';

235

223

*to= '\0';

236

224

if (error != NULL)

331

319

decimal point. For this we are calling dtoa with mode=5, passing the

332

320

number of significant digits = (len-decpt) - (len-width) = width-decpt

333

321

*/

334

dtoa_free(res, buf, sizeof(buf));

335

res= dtoa(x, 5, width - decpt, &decpt, &sign, &end, buf, sizeof(buf));

322

free(res);

323

res= dtoa(x, 5, width - decpt, &decpt, &sign, &end);

336

324

src= res;

337

325

len= end - res;

338

326

}

397

385

if (width < len)

398

386

{

399

387

/* Yes, re-convert with a smaller width */

400

dtoa_free(res, buf, sizeof(buf));

401

res= dtoa(x, 4, width, &decpt, &sign, &end, buf, sizeof(buf));

388

free(res);

389

res= dtoa(x, 4, width, &decpt, &sign, &end);

402

390

src= res;

403

391

len= end - res;

404

392

if (--decpt < 0)

438

426

}

439

427

440

428

end:

441

dtoa_free(res, buf, sizeof(buf));

429

free(res);

442

430

*dst= '\0';

443

431

444

432

return dst - to;

464

452

465

453

double my_strtod(const char *str, char **end, int *error)

466

454

{

467

double alignedbuf[(DTOA_BUFF_SIZE+1)/sizeof(double)];

468

char *buf= (char*)alignedbuf;

469

455

double res;

470

456

assert(str != NULL && end != NULL && *end != NULL && error != NULL);

471

457

472

res= my_strtod_int(str, end, error, buf, sizeof(buf));

458

res= my_strtod_int(str, end, error);

473

459

return (*error == 0) ? res : (res < 0 ? -DBL_MAX : DBL_MAX);

474

460

}

475

461

598

584

#define Big1 0xffffffff

599

585

#define FFFFFFFF 0xffffffffUL

600

586

601

/* This is tested to be enough for dtoa */

602

603

#define Kmax 15

604

605

#define Bcopy(x,y) memcpy(&x->sign, &y->sign, \

606

2*sizeof(int) + y->wds*sizeof(ULong))

607

608

587

/* Arbitrary-length integer */

609

588

610

589

typedef struct Bigint

619

598

int wds; /* current length in 32-bit words */

620

599

} Bigint;

621

600

622

623

/* A simple stack-memory based allocator for Bigints */

624

625

typedef struct Stack_alloc

601

static Bigint *Bcopy(Bigint* dst, Bigint* src)

626

602

{

627

char *begin;

628

char *free;

629

char *end;

630

/*

631

Having list of free blocks lets us reduce maximum required amount

632

of memory from ~4000 bytes to < 1680 (tested on x86).

633

*/

634

Bigint *freelist[Kmax+1];

635

} Stack_alloc;

636

637

638

/*

639

Try to allocate object on stack, and resort to malloc if all

640

stack memory is used.

641

*/

642

643

static Bigint *Balloc(int k, Stack_alloc *alloc)

603

dst->sign= src->sign;

604

dst->wds= src->wds;

605

606

assert(dst->maxwds >= src->wds);

607

608

memcpy(dst->p.x, src->p.x, src->wds*sizeof(ULong));

609

610

return dst;

611

}

612

613

static Bigint *Balloc(int k)

644

614

{

645

615

Bigint *rv;

646

if (k <= Kmax && alloc->freelist[k])

647

{

648

rv= alloc->freelist[k];

649

alloc->freelist[k]= rv->p.next;

650

}

651

else

652

{

653

int x, len;

654

655

x= 1 << k;

656

len= ALIGN_SIZE(sizeof(Bigint) + x * sizeof(ULong));

657

658

if (alloc->free + len <= alloc->end)

659

{

660

rv= (Bigint*) alloc->free;

661

alloc->free+= len;

662

}

663

else

664

rv= (Bigint*) malloc(len);

665

666

rv->k= k;

667

rv->maxwds= x;

668

}

616

617

/* TODO: some malloc failure checking */

618

619

int x= 1 << k;

620

rv= (Bigint*) malloc(sizeof(Bigint));

621

622

rv->p.x= (ULong*)malloc(x * sizeof(ULong));

623

624

rv->k= k;

625

rv->maxwds= x;

626

669

627

rv->sign= rv->wds= 0;

670

rv->p.x= (ULong*) (rv + 1);

628

671

629

return rv;

672

630

}

673

631

677

635

list. Otherwise call free().

678

636

*/

679

637

680

static void Bfree(Bigint *v, Stack_alloc *alloc)

681

{

682

char *gptr= (char*) v; /* generic pointer */

683

if (gptr < alloc->begin || gptr >= alloc->end)

684

free(gptr);

685

else if (v->k <= Kmax)

686

{

687

/*

688

Maintain free lists only for stack objects: this way we don't

689

have to bother with freeing lists in the end of dtoa;

690

heap should not be used normally anyway.

691

*/

692

v->p.next= alloc->freelist[v->k];

693

alloc->freelist[v->k]= v;

694

}

695

}

696

697

698

/*

699

This is to place return value of dtoa in: tries to use stack

700

as well, but passes by free lists management and just aligns len by

701

sizeof(ULong).

702

*/

703

704

static char *dtoa_alloc(int i, Stack_alloc *alloc)

705

{

706

char *rv;

707

int aligned_size= (i + sizeof(double) - 1) / sizeof(double) * sizeof(double);

708

if (alloc->free + aligned_size <= alloc->end)

709

{

710

rv= alloc->free;

711

alloc->free+= aligned_size;

712

}

713

else

714

rv= (char *)malloc(i);

715

return rv;

716

}

717

718

719

/*

720

dtoa_free() must be used to free values s returned by dtoa()

721

This is the counterpart of dtoa_alloc()

722

*/

723

724

static void dtoa_free(char *gptr, char *buf, size_t buf_size)

725

{

726

if (gptr < buf || gptr >= buf + buf_size)

727

free(gptr);

728

}

729

638

static void Bfree(Bigint *v)

639

{

640

if(!v)

641

return;

642

free(v->p.x);

643

free(v);

644

}

730

645

731

646

/* Bigint arithmetic functions */

732

647

733

648

/* Multiply by m and add a */

734

649

735

static Bigint *multadd(Bigint *b, int m, int a, Stack_alloc *alloc)

650

static Bigint *multadd(Bigint *b, int m, int a)

736

651

{

737

652

int i, wds;

738

653

ULong *x;

754

669

{

755

670

if (wds >= b->maxwds)

756

671

{

757

b1= Balloc(b->k+1, alloc);

672

b1= Balloc(b->k+1);

758

673

Bcopy(b1, b);

759

Bfree(b, alloc);

674

Bfree(b);

760

675

b= b1;

761

676

}

762

677

b->p.x[wds++]= (ULong) carry;

766

681

}

767

682

768

683

769

static Bigint *s2b(const char *s, int nd0, int nd, ULong y9, Stack_alloc *alloc)

684

static Bigint *s2b(const char *s, int nd0, int nd, ULong y9)

770

685

{

771

686

Bigint *b;

772

687

int i, k;

774

689

775

690

x= (nd + 8) / 9;

776

691

for (k= 0, y= 1; x > y; y <<= 1, k++) ;

777

b= Balloc(k, alloc);

692

b= Balloc(k);

778

693

b->p.x[0]= y9;

779

694

b->wds= 1;

780

695

783

698

{

784

699

s+= 9;

785

700

do

786

b= multadd(b, 10, *s++ - '0', alloc);

701

b= multadd(b, 10, *s++ - '0');

787

702

while (++i < nd0);

788

703

s++;

789

704

}

790

705

else

791

706

s+= 10;

792

707

for(; i < nd; i++)

793

b= multadd(b, 10, *s++ - '0', alloc);

708

b= multadd(b, 10, *s++ - '0');

794

709

return b;

795

710

}

796

711

881

796

882

797

/* Convert integer to Bigint number */

883

798

884

static Bigint *i2b(int i, Stack_alloc *alloc)

799

static Bigint *i2b(int i)

885

800

{

886

801

Bigint *b;

887

802

888

b= Balloc(1, alloc);

803

b= Balloc(1);

889

804

b->p.x[0]= i;

890

805

b->wds= 1;

891

806

return b;

894

809

895

810

/* Multiply two Bigint numbers */

896

811

897

static Bigint *mult(Bigint *a, Bigint *b, Stack_alloc *alloc)

812

static Bigint *mult(Bigint *a, Bigint *b)

898

813

{

899

814

Bigint *c;

900

815

int k, wa, wb, wc;

914

829

wc= wa + wb;

915

830

if (wc > a->maxwds)

916

831

k++;

917

c= Balloc(k, alloc);

832

c= Balloc(k);

918

833

for (x= c->p.x, xa= x + wc; x < xa; x++)

919

834

*x= 0;

920

835

xa= a->p.x;

986

901

987

902

#define P5A_MAX (sizeof(p5_a)/sizeof(*p5_a) - 1)

988

903

989

static Bigint *pow5mult(Bigint *b, int k, Stack_alloc *alloc)

904

static Bigint *pow5mult(Bigint *b, int k)

990

905

{

991

906

Bigint *b1, *p5, *p51;

992

907

int i;

993

908

static int p05[3]= { 5, 25, 125 };

994

909

995

910

if ((i= k & 3))

996

b= multadd(b, p05[i-1], 0, alloc);

911

b= multadd(b, p05[i-1], 0);

997

912

998

913

if (!(k>>= 2))

999

914

return b;

1002

917

{

1003

918

if (k & 1)

1004

919

{

1005

b1= mult(b, p5, alloc);

1006

Bfree(b, alloc);

920

b1= mult(b, p5);

921

Bfree(b);

1007

922

b= b1;

1008

923

}

1009

924

if (!(k>>= 1))

1012

927

if (p5 < p5_a + P5A_MAX)

1013

928

++p5;

1014

929

else if (p5 == p5_a + P5A_MAX)

1015

p5= mult(p5, p5, alloc);

930

p5= mult(p5, p5);

1016

931

else

1017

932

{

1018

p51= mult(p5, p5, alloc);

1019

Bfree(p5, alloc);

933

p51= mult(p5, p5);

934

Bfree(p5);

1020

935

p5= p51;

1021

936

}

1022

937

}

1024

939

}

1025

940

1026

941

1027

static Bigint *lshift(Bigint *b, int k, Stack_alloc *alloc)

942

static Bigint *lshift(Bigint *b, int k)

1028

943

{

1029

944

int i, k1, n, n1;

1030

945

Bigint *b1;

1035

950

n1= n + b->wds + 1;

1036

951

for (i= b->maxwds; n1 > i; i<<= 1)

1037

952

k1++;

1038

b1= Balloc(k1, alloc);

953

b1= Balloc(k1);

1039

954

x1= b1->p.x;

1040

955

for (i= 0; i < n; i++)

1041

956

*x1++= 0;

1059

974

*x1++= *x++;

1060

975

while (x < xe);

1061

976

b1->wds= n1 - 1;

1062

Bfree(b, alloc);

977

Bfree(b);

1063

978

return b1;

1064

979

}

1065

980

1088

1003

}

1089

1004

1090

1005

1091

static Bigint *diff(Bigint *a, Bigint *b, Stack_alloc *alloc)

1006

static Bigint *diff(Bigint *a, Bigint *b)

1092

1007

{

1093

1008

Bigint *c;

1094

1009

int i, wa, wb;

1098

1013

i= cmp(a,b);

1099

1014

if (!i)

1100

1015

{

1101

c= Balloc(0, alloc);

1016

c= Balloc(0);

1102

1017

c->wds= 1;

1103

1018

c->p.x[0]= 0;

1104

1019

return c;

1112

1027

}

1113

1028

else

1114

1029

i= 0;

1115

c= Balloc(a->k, alloc);

1030

c= Balloc(a->k);

1116

1031

c->sign= i;

1117

1032

wa= a->wds;

1118

1033

xa= a->p.x;

1193

1108

}

1194

1109

1195

1110

1196

static Bigint *d2b(double d, int *e, int *bits, Stack_alloc *alloc)

1111

static Bigint *d2b(double d, int *e, int *bits)

1197

1112

{

1198

1113

Bigint *b;

1199

1114

int de, k;

1202

1117

#define d0 word0(d)

1203

1118

#define d1 word1(d)

1204

1119

1205

b= Balloc(1, alloc);

1120

b= Balloc(1);

1206

1121

x= b->p.x;

1207

1122

1208

1123

z= d0 & Frac_mask;

1309

1224

for 0 <= k <= 22).

1310

1225

*/

1311

1226

1312

static double my_strtod_int(const char *s00, char **se, int *error, char *buf, size_t buf_size)

1227

static double my_strtod_int(const char *s00, char **se, int *error)

1313

1228

{

1314

1229

int scale;

1315

1230

int bb2, bb5, bbe, bd2, bd5, bbbits, bs2, c, dsign,

1322

1237

#ifdef SET_INEXACT

1323

1238

int inexact, oldinexact;

1324

1239

#endif

1325

Stack_alloc alloc;

1326

1240

1327

1241

c= 0;

1328

1242

*error= 0;

1329

1243

1330

alloc.begin= alloc.free= buf;

1331

alloc.end= buf + buf_size;

1332

memset(alloc.freelist, 0, sizeof(alloc.freelist));

1333

1334

1244

sign= nz0= nz= 0;

1335

1245

dval(rv)= 0.;

1336

1246

for (s= s00; s < end; s++)

1606

1516

1607

1517

/* Put digits into bd: true value = bd * 10^e */

1608

1518

1609

bd0= s2b(s0, nd0, nd, y, &alloc);

1519

bd0= s2b(s0, nd0, nd, y);

1610

1520

1611

1521

for(;;)

1612

1522

{

1613

bd= Balloc(bd0->k, &alloc);

1523

bd= Balloc(bd0->k);

1614

1524

Bcopy(bd, bd0);

1615

bb= d2b(dval(rv), &bbe, &bbbits, &alloc); /* rv = bb * 2^bbe */

1616

bs= i2b(1, &alloc);

1525

bb= d2b(dval(rv), &bbe, &bbbits); /* rv = bb * 2^bbe */

1526

bs= i2b(1);

1617

1527

1618

1528

if (e >= 0)

1619

1529

{

1650

1560

}

1651

1561

if (bb5 > 0)

1652

1562

{

1653

bs= pow5mult(bs, bb5, &alloc);

1654

bb1= mult(bs, bb, &alloc);

1655

Bfree(bb, &alloc);

1563

bs= pow5mult(bs, bb5);

1564

bb1= mult(bs, bb);

1565

Bfree(bb);

1656

1566

bb= bb1;

1657

1567

}

1658

1568

if (bb2 > 0)

1659

bb= lshift(bb, bb2, &alloc);

1569

bb= lshift(bb, bb2);

1660

1570

if (bd5 > 0)

1661

bd= pow5mult(bd, bd5, &alloc);

1571

bd= pow5mult(bd, bd5);

1662

1572

if (bd2 > 0)

1663

bd= lshift(bd, bd2, &alloc);

1573

bd= lshift(bd, bd2);

1664

1574

if (bs2 > 0)

1665

bs= lshift(bs, bs2, &alloc);

1666

delta= diff(bb, bd, &alloc);

1575

bs= lshift(bs, bs2);

1576

delta= diff(bb, bd);

1667

1577

dsign= delta->sign;

1668

1578

delta->sign= 0;

1669

1579

i= cmp(delta, bs);

1691

1601

#endif

1692

1602

break;

1693

1603

}

1694

delta= lshift(delta, Log2P, &alloc);

1604

delta= lshift(delta, Log2P);

1695

1605

if (cmp(delta, bs) > 0)

1696

1606

goto drop_down;

1697

1607

break;

1832

1742

}

1833

1743

#endif

1834

1744

cont:

1835

Bfree(bb, &alloc);

1836

Bfree(bd, &alloc);

1837

Bfree(bs, &alloc);

1838

Bfree(delta, &alloc);

1745

Bfree(bb);

1746

Bfree(bd);

1747

Bfree(bs);

1748

Bfree(delta);

1839

1749

}

1840

1750

#ifdef SET_INEXACT

1841

1751

if (inexact)

1865

1775

}

1866

1776

#endif

1867

1777

retfree:

1868

Bfree(bb, &alloc);

1869

Bfree(bd, &alloc);

1870

Bfree(bs, &alloc);

1871

Bfree(bd0, &alloc);

1872

Bfree(delta, &alloc);

1778

Bfree(bb);

1779

Bfree(bd);

1780

Bfree(bs);

1781

Bfree(bd0);

1782

Bfree(delta);

1873

1783

ret:

1874

1784

*se= (char *)s;

1875

1785

return sign ? -dval(rv) : dval(rv);

1976

1886

*/

1977

1887

1978

1888

static char *dtoa(double d, int mode, int ndigits, int *decpt, int *sign,

1979

char **rve, char *buf, size_t buf_size)

1889

char **rve)

1980

1890

{

1981

1891

/*

1982

1892

Arguments ndigits, decpt, sign are similar to those

1988

1898

mode:

1989

1899

0 ==> shortest string that yields d when read in

1990

1900

and rounded to nearest.

1901

1991

1902

1 ==> like 0, but with Steele & White stopping rule;

1992

1903

e.g. with IEEE P754 arithmetic , mode 0 gives

1993

1904

1e23 whereas mode 1 gives 9.999999999999999e22.

2019

1930

Bigint *b, *b1, *delta, *mlo = NULL, *mhi, *S;

2020

1931

double d2, ds, eps;

2021

1932

char *s, *s0;

2022

Stack_alloc alloc;

2023

2024

alloc.begin= alloc.free= buf;

2025

alloc.end= buf + buf_size;

2026

memset(alloc.freelist, 0, sizeof(alloc.freelist));

2027

1933

2028

1934

if (word0(d) & Sign_bit)

2029

1935

{

2039

1945

(!dval(d) && (*decpt= 1)))

2040

1946

{

2041

1947

/* Infinity, NaN, 0 */

2042

char *res= (char*) dtoa_alloc(2, &alloc);

1948

char *res= (char*) malloc(2);

2043

1949

res[0]= '0';

2044

1950

res[1]= '\0';

2045

1951

if (rve)

2048

1954

}

2049

1955

2050

1956

2051

b= d2b(dval(d), &be, &bbits, &alloc);

1957

b= d2b(dval(d), &be, &bbits);

2052

1958

if ((i= (int)(word0(d) >> Exp_shift1 & (Exp_mask>>Exp_shift1))))

2053

1959

{

2054

1960

dval(d2)= dval(d);

2163

2069

if (i <= 0)

2164

2070

i= 1;

2165

2071

}

2166

s= s0= dtoa_alloc(i, &alloc);

2072

s= s0= (char*)malloc(i+1); /* +1 for trailing '\0' appended

2073

at end of function */

2167

2074

2168

2075

if (ilim >= 0 && ilim <= Quick_max && try_quick)

2169

2076

{

2328

2235

i = denorm ? be + (Bias + (P-1) - 1 + 1) : 1 + P - bbits;

2329

2236

b2+= i;

2330

2237

s2+= i;

2331

mhi= i2b(1, &alloc);

2238

mhi= i2b(1);

2332

2239

}

2333

2240

if (m2 > 0 && s2 > 0)

2334

2241

{

2343

2250

{

2344

2251

if (m5 > 0)

2345

2252

{

2346

mhi= pow5mult(mhi, m5, &alloc);

2347

b1= mult(mhi, b, &alloc);

2348

Bfree(b, &alloc);

2253

mhi= pow5mult(mhi, m5);

2254

b1= mult(mhi, b);

2255

Bfree(b);

2349

2256

b= b1;

2350

2257

}

2351

2258

if ((j= b5 - m5))

2352

b= pow5mult(b, j, &alloc);

2259

b= pow5mult(b, j);

2353

2260

}

2354

2261

else

2355

b= pow5mult(b, b5, &alloc);

2262

b= pow5mult(b, b5);

2356

2263

}

2357

S= i2b(1, &alloc);

2264

S= i2b(1);

2358

2265

if (s5 > 0)

2359

S= pow5mult(S, s5, &alloc);

2266

S= pow5mult(S, s5);

2360

2267

2361

2268

/* Check for special case that d is a normalized power of 2. */

2362

2269

2400

2307

s2+= i;

2401

2308

}

2402

2309

if (b2 > 0)

2403

b= lshift(b, b2, &alloc);

2310

b= lshift(b, b2);

2404

2311

if (s2 > 0)

2405

S= lshift(S, s2, &alloc);

2312

S= lshift(S, s2);

2406

2313

if (k_check)

2407

2314

{

2408

2315

if (cmp(b,S) < 0)

2409

2316

{

2410

2317

k--;

2411

2318

/* we botched the k estimate */

2412

b= multadd(b, 10, 0, &alloc);

2319

b= multadd(b, 10, 0);

2413

2320

if (leftright)

2414

mhi= multadd(mhi, 10, 0, &alloc);

2321

mhi= multadd(mhi, 10, 0);

2415

2322

ilim= ilim1;

2416

2323

}

2417

2324

}

2418

2325

if (ilim <= 0 && (mode == 3 || mode == 5))

2419

2326

{

2420

if (ilim < 0 || cmp(b,S= multadd(S,5,0, &alloc)) <= 0)

2327

if (ilim < 0 || cmp(b,S= multadd(S,5,0)) <= 0)

2421

2328

{

2422

2329

/* no digits, fcvt style */

2423

2330

no_digits:

2432

2339

if (leftright)

2433

2340

{

2434

2341

if (m2 > 0)

2435

mhi= lshift(mhi, m2, &alloc);

2342

mhi= lshift(mhi, m2);

2436

2343

2437

2344

/*

2438

2345

Compute mlo -- check for special case that d is a normalized power of 2.

2441

2348

mlo= mhi;

2442

2349

if (spec_case)

2443

2350

{

2444

mhi= Balloc(mhi->k, &alloc);

2351

mhi= Balloc(mhi->k);

2445

2352

Bcopy(mhi, mlo);

2446

mhi= lshift(mhi, Log2P, &alloc);

2353

mhi= lshift(mhi, Log2P);

2447

2354

}

2448

2355

2449

2356

for (i= 1;;i++)

2451

2358

dig= quorem(b,S) + '0';

2452

2359

/* Do we yet have the shortest decimal string that will round to d? */

2453

2360

j= cmp(b, mlo);

2454

delta= diff(S, mhi, &alloc);

2361

delta= diff(S, mhi);

2455

2362

j1= delta->sign ? 1 : cmp(b, delta);

2456

Bfree(delta, &alloc);

2363

Bfree(delta);

2457

2364

if (j1 == 0 && mode != 1 && !(word1(d) & 1)

2458

2365

)

2459

2366

{

2472

2379

}

2473

2380

if (j1 > 0)

2474

2381

{

2475

b= lshift(b, 1, &alloc);

2382

b= lshift(b, 1);

2476

2383

j1= cmp(b, S);

2477

2384

if ((j1 > 0 || (j1 == 0 && dig & 1))

2478

2385

&& dig++ == '9')

2496

2403

*s++= dig;

2497

2404

if (i == ilim)

2498

2405

break;

2499

b= multadd(b, 10, 0, &alloc);

2406

b= multadd(b, 10, 0);

2500

2407

if (mlo == mhi)

2501

mlo= mhi= multadd(mhi, 10, 0, &alloc);

2408

mlo= mhi= multadd(mhi, 10, 0);

2502

2409

else

2503

2410

{

2504

mlo= multadd(mlo, 10, 0, &alloc);

2505

mhi= multadd(mhi, 10, 0, &alloc);

2411

mlo= multadd(mlo, 10, 0);

2412

mhi= multadd(mhi, 10, 0);

2506

2413

}

2507

2414

}

2508

2415

}

2516

2423

}

2517

2424

if (i >= ilim)

2518

2425

break;

2519

b= multadd(b, 10, 0, &alloc);

2426

b= multadd(b, 10, 0);

2520

2427

}

2521

2428

2522

2429

/* Round off last digit */

2523

2430

2524

b= lshift(b, 1, &alloc);

2431

b= lshift(b, 1);

2525

2432

j= cmp(b, S);

2526

2433

if (j > 0 || (j == 0 && dig & 1))

2527

2434

{

2541

2448

s++;

2542

2449

}

2543

2450

ret:

2544

Bfree(S, &alloc);

2451

Bfree(S);

2545

2452

if (mhi)

2546

2453

{

2547

2454

if (mlo && mlo != mhi)

2548

Bfree(mlo, &alloc);

2549

Bfree(mhi, &alloc);

2455

Bfree(mlo);

2456

Bfree(mhi);

2550

2457

}

2551

2458

ret1:

2552

Bfree(b, &alloc);

2459

Bfree(b);

2553

2460

*s= 0;

2554

2461

*decpt= k + 1;

2555

2462

if (rve)