1 /**
2  * 80-bit floating point value implementation if the C/D compiler does not support them natively.
3  *
4  * Copyright (C) 1999-2021 by The D Language Foundation, All Rights Reserved
5  * All Rights Reserved, written by Rainer Schuetze
6  * http://www.digitalmars.com
7  * Distributed under the Boost Software License, Version 1.0.
8  * (See accompanying file LICENSE or copy at http://www.boost.org/LICENSE_1_0.txt)
9  * https://github.com/dlang/dmd/blob/master/src/root/longdouble.d
10  */
11 
12 module dmd.root.longdouble;
13 
14 static if (real.sizeof > 8)
15     alias longdouble = real;
16 else
17     alias longdouble = longdouble_soft;
18 
19 // longdouble_soft needed when building the backend with
20 // Visual C or the frontend with LDC on Windows
21 version(CRuntime_Microsoft):
22 extern(C++):
23 nothrow:
24 @nogc:
25 
26 version(D_InlineAsm_X86_64)
27     version = AsmX86;
28 else version(D_InlineAsm_X86)
29     version = AsmX86;
30 else
31     static assert(false, "longdouble_soft not supported on this platform");
32 
33 bool initFPU()
34 {
35     version(D_InlineAsm_X86_64)
36     {
37         // set precision to 64-bit mantissa and rounding control to nearest
38         asm nothrow @nogc @trusted
39         {
40             push    RAX;                 // add space on stack
41             fstcw   word ptr [RSP];
42             movzx   EAX,word ptr [RSP];  // also return old CW in EAX
43             and     EAX, ~0xF00;         // mask for PC and RC
44             or      EAX, 0x300;
45             mov     dword ptr [RSP],EAX;
46             fldcw   word ptr [RSP];
47             pop     RAX;
48         }
49     }
50     else version(D_InlineAsm_X86)
51     {
52         // set precision to 64-bit mantissa and rounding control to nearest
53         asm nothrow @nogc @trusted
54         {
55             push    EAX;                 // add space on stack
56             fstcw   word ptr [ESP];
57             movzx   EAX,word ptr [ESP];  // also return old CW in EAX
58             and     EAX, ~0xF00;         // mask for PC and RC
59             or      EAX, 0x300;
60             mov     dword ptr [ESP],EAX;
61             fldcw   word ptr [ESP];
62             pop     EAX;
63         }
64     }
65 
66     return true;
67 }
68 
69 version(unittest) version(CRuntime_Microsoft)
70 extern(D) shared static this()
71 {
72     initFPU(); // otherwise not guaranteed to be run before pure unittest below
73 }
74 
75 void ld_clearfpu()
76 {
77     version(AsmX86)
78     {
79         asm nothrow @nogc @trusted
80         {
81             fclex;
82         }
83     }
84 }
85 
86 pure:
87 @trusted: // LDC: LLVM __asm is @system AND requires taking the address of variables
88 
89 struct longdouble_soft
90 {
91 nothrow @nogc pure:
92     // DMD's x87 `real` on Windows is packed (alignof = 2 -> sizeof = 10).
93     align(2) ulong mantissa = 0xC000000000000000UL; // default to qnan
94     ushort exp_sign = 0x7fff; // sign is highest bit
95 
96     this(ulong m, ushort es) { mantissa = m; exp_sign = es; }
97     this(longdouble_soft ld) { mantissa = ld.mantissa; exp_sign = ld.exp_sign; }
98     this(int i) { ld_set(&this, i); }
99     this(uint i) { ld_set(&this, i); }
100     this(long i) { ld_setll(&this, i); }
101     this(ulong i) { ld_setull(&this, i); }
102     this(float f) { ld_set(&this, f); }
103     this(double d)
104     {
105         // allow zero initialization at compile time
106         if (__ctfe && d == 0)
107         {
108             mantissa = 0;
109             exp_sign = 0;
110         }
111         else
112             ld_set(&this, d);
113     }
114     this(real r)
115     {
116         static if (real.sizeof > 8)
117             *cast(real*)&this = r;
118         else
119             this(cast(double)r);
120     }
121 
122     ushort exponent() const { return exp_sign & 0x7fff; }
123     bool sign() const { return (exp_sign & 0x8000) != 0; }
124 
125     extern(D)
126     {
127         ref longdouble_soft opAssign(longdouble_soft ld) return { mantissa = ld.mantissa; exp_sign = ld.exp_sign; return this; }
128         ref longdouble_soft opAssign(T)(T rhs) { this = longdouble_soft(rhs); return this; }
129 
130         longdouble_soft opUnary(string op)() const
131         {
132             static if (op == "-") return longdouble_soft(mantissa, exp_sign ^ 0x8000);
133             else static assert(false, "Operator `"~op~"` is not implemented");
134         }
135 
136         bool opEquals(T)(T rhs) const { return this.ld_cmpe(longdouble_soft(rhs)); }
137         int  opCmp(T)(T rhs) const { return this.ld_cmp(longdouble_soft(rhs)); }
138 
139         longdouble_soft opBinary(string op, T)(T rhs) const
140         {
141             static if      (op == "+") return this.ld_add(longdouble_soft(rhs));
142             else static if (op == "-") return this.ld_sub(longdouble_soft(rhs));
143             else static if (op == "*") return this.ld_mul(longdouble_soft(rhs));
144             else static if (op == "/") return this.ld_div(longdouble_soft(rhs));
145             else static if (op == "%") return this.ld_mod(longdouble_soft(rhs));
146             else static assert(false, "Operator `"~op~"` is not implemented");
147         }
148 
149         longdouble_soft opBinaryRight(string op, T)(T rhs) const
150         {
151             static if      (op == "+") return longdouble_soft(rhs).ld_add(this);
152             else static if (op == "-") return longdouble_soft(rhs).ld_sub(this);
153             else static if (op == "*") return longdouble_soft(rhs).ld_mul(this);
154             else static if (op == "%") return longdouble_soft(rhs).ld_mod(this);
155             else static assert(false, "Operator `"~op~"` is not implemented");
156         }
157 
158         ref longdouble_soft opOpAssign(string op)(longdouble_soft rhs)
159         {
160             mixin("this = this " ~ op ~ " rhs;");
161             return this;
162         }
163 
164         T opCast(T)() const @trusted
165         {
166             static      if (is(T == bool))   return mantissa != 0 || (exp_sign & 0x7fff) != 0;
167             else static if (is(T == byte))   return cast(T)ld_read(&this);
168             else static if (is(T == ubyte))  return cast(T)ld_read(&this);
169             else static if (is(T == short))  return cast(T)ld_read(&this);
170             else static if (is(T == ushort)) return cast(T)ld_read(&this);
171             else static if (is(T == int))    return cast(T)ld_read(&this);
172             else static if (is(T == uint))   return cast(T)ld_read(&this);
173             else static if (is(T == float))  return cast(T)ld_read(&this);
174             else static if (is(T == double)) return cast(T)ld_read(&this);
175             else static if (is(T == long))   return ld_readll(&this);
176             else static if (is(T == ulong))  return ld_readull(&this);
177             else static if (is(T == real))
178             {
179                 // convert to front end real if built with dmd
180                 if (real.sizeof > 8)
181                     return *cast(real*)&this;
182                 else
183                     return ld_read(&this);
184             }
185             else static assert(false, "usupported type");
186         }
187     }
188 
189     // a qnan
190     static longdouble_soft nan() { return longdouble_soft(0xC000000000000000UL, 0x7fff); }
191     static longdouble_soft infinity() { return longdouble_soft(0x8000000000000000UL, 0x7fff); }
192     static longdouble_soft zero() { return longdouble_soft(0, 0); }
193     static longdouble_soft max() { return longdouble_soft(0xffffffffffffffffUL, 0x7ffe); }
194     static longdouble_soft min_normal() { return longdouble_soft(0x8000000000000000UL, 1); }
195     static longdouble_soft epsilon() { return longdouble_soft(0x8000000000000000UL, 0x3fff - 63); }
196 
197     static uint dig() { return 18; }
198     static uint mant_dig() { return 64; }
199     static uint max_exp() { return 16_384; }
200     static uint min_exp() { return -16_381; }
201     static uint max_10_exp() { return 4932; }
202     static uint min_10_exp() { return -4932; }
203 }
204 
205 static assert(longdouble_soft.alignof == longdouble.alignof);
206 static assert(longdouble_soft.sizeof == longdouble.sizeof);
207 
208 version(LDC)
209 {
210     import ldc.llvmasm;
211 
212     extern(D):
213     private:
214     string fld_arg  (string arg)() { return `__asm("fldt $0",  "*m,~{st}",  &` ~ arg ~ `);`; }
215     string fstp_arg (string arg)() { return `__asm("fstpt $0", "=*m,~{st}", &` ~ arg ~ `);`; }
216     string fld_parg (string arg)() { return `__asm("fldt $0",  "*m,~{st}",   ` ~ arg ~ `);`; }
217     string fstp_parg(string arg)() { return `__asm("fstpt $0", "=*m,~{st}",  ` ~ arg ~ `);`; }
218 }
219 else version(D_InlineAsm_X86_64)
220 {
221     // longdouble_soft passed by reference
222     extern(D):
223     private:
224     string fld_arg(string arg)()
225     {
226         return "asm nothrow @nogc pure @trusted { mov RAX, " ~ arg ~ "; fld real ptr [RAX]; }";
227     }
228     string fstp_arg(string arg)()
229     {
230         return "asm nothrow @nogc pure @trusted { mov RAX, " ~ arg ~ "; fstp real ptr [RAX]; }";
231     }
232     alias fld_parg = fld_arg;
233     alias fstp_parg = fstp_arg;
234 }
235 else version(D_InlineAsm_X86)
236 {
237     // longdouble_soft passed by value
238     extern(D):
239     private:
240     string fld_arg(string arg)()
241     {
242         return "asm nothrow @nogc pure @trusted { lea EAX, " ~ arg ~ "; fld real ptr [EAX]; }";
243     }
244     string fstp_arg(string arg)()
245     {
246         return "asm nothrow @nogc pure @trusted { lea EAX, " ~ arg ~ "; fstp real ptr [EAX]; }";
247     }
248     string fld_parg(string arg)()
249     {
250         return "asm nothrow @nogc pure @trusted { mov EAX, " ~ arg ~ "; fld real ptr [EAX]; }";
251     }
252     string fstp_parg(string arg)()
253     {
254         return "asm nothrow @nogc pure @trusted { mov EAX, " ~ arg ~ "; fstp real ptr [EAX]; }";
255     }
256 }
257 
258 double ld_read(const longdouble_soft* pthis)
259 {
260     double res;
261     version(AsmX86)
262     {
263         mixin(fld_parg!("pthis"));
264         asm nothrow @nogc pure @trusted
265         {
266             fstp res;
267         }
268     }
269     return res;
270 }
271 
272 long ld_readll(const longdouble_soft* pthis)
273 {
274     return ld_readull(pthis);
275 }
276 
277 ulong ld_readull(const longdouble_soft* pthis)
278 {
279     // somehow the FPU does not respect the CHOP mode of the rounding control
280     // in 64-bit mode
281     // so we roll our own conversion (it also allows the usual C wrap-around
282     // instead of the "invalid value" created by the FPU)
283     int expo = pthis.exponent - 0x3fff;
284     ulong u;
285     if(expo < 0 || expo > 127)
286         return 0;
287     if(expo < 64)
288         u = pthis.mantissa >> (63 - expo);
289     else
290         u = pthis.mantissa << (expo - 63);
291     if(pthis.sign)
292         u = ~u + 1;
293     return u;
294 }
295 
296 int ld_statusfpu()
297 {
298     int res = 0;
299     version(AsmX86)
300     {
301         asm nothrow @nogc pure @trusted
302         {
303             fstsw word ptr [res];
304         }
305     }
306     return res;
307 }
308 
309 void ld_set(longdouble_soft* pthis, double d)
310 {
311     version(AsmX86)
312     {
313         asm nothrow @nogc pure @trusted
314         {
315             fld d;
316         }
317         mixin(fstp_parg!("pthis"));
318     }
319 }
320 
321 void ld_setll(longdouble_soft* pthis, long d)
322 {
323     version(AsmX86)
324     {
325         asm nothrow @nogc pure @trusted
326         {
327             fild qword ptr d;
328         }
329         mixin(fstp_parg!("pthis"));
330     }
331 }
332 
333 void ld_setull(longdouble_soft* pthis, ulong d)
334 {
335     d ^= (1L << 63);
336     version(AsmX86)
337     {
338         auto pTwoPow63 = &twoPow63;
339         mixin(fld_parg!("pTwoPow63"));
340         asm nothrow @nogc pure @trusted
341         {
342             fild qword ptr d;
343             faddp;
344         }
345         mixin(fstp_parg!("pthis"));
346     }
347 }
348 
349 // using an argument as result to avoid RVO, see https://issues.dlang.org/show_bug.cgi?id=18758
350 longdouble_soft ldexpl(longdouble_soft ld, int exp)
351 {
352     version(AsmX86)
353     {
354         asm nothrow @nogc pure @trusted
355         {
356             fild    dword ptr exp;
357         }
358         mixin(fld_arg!("ld"));
359         asm nothrow @nogc pure @trusted
360         {
361             fscale;                 // ST(0) = ST(0) * (2**ST(1))
362             fstp    ST(1);
363         }
364         mixin(fstp_arg!("ld"));
365     }
366     return ld;
367 }
368 
369 ///////////////////////////////////////////////////////////////////////
370 longdouble_soft ld_add(longdouble_soft ld1, longdouble_soft ld2)
371 {
372     version(AsmX86)
373     {
374         mixin(fld_arg!("ld1"));
375         mixin(fld_arg!("ld2"));
376         asm nothrow @nogc pure @trusted
377         {
378             fadd;
379         }
380         mixin(fstp_arg!("ld1"));
381     }
382     return ld1;
383 }
384 
385 longdouble_soft ld_sub(longdouble_soft ld1, longdouble_soft ld2)
386 {
387     version(AsmX86)
388     {
389         mixin(fld_arg!("ld1"));
390         mixin(fld_arg!("ld2"));
391         asm nothrow @nogc pure @trusted
392         {
393             fsub;
394         }
395         mixin(fstp_arg!("ld1"));
396     }
397     return ld1;
398 }
399 
400 longdouble_soft ld_mul(longdouble_soft ld1, longdouble_soft ld2)
401 {
402     version(AsmX86)
403     {
404         mixin(fld_arg!("ld1"));
405         mixin(fld_arg!("ld2"));
406         asm nothrow @nogc pure @trusted
407         {
408             fmul;
409         }
410         mixin(fstp_arg!("ld1"));
411     }
412     return ld1;
413 }
414 
415 longdouble_soft ld_div(longdouble_soft ld1, longdouble_soft ld2)
416 {
417     version(AsmX86)
418     {
419         mixin(fld_arg!("ld1"));
420         mixin(fld_arg!("ld2"));
421         asm nothrow @nogc pure @trusted
422         {
423             fdiv;
424         }
425         mixin(fstp_arg!("ld1"));
426     }
427     return ld1;
428 }
429 
430 bool ld_cmpb(longdouble_soft x, longdouble_soft y)
431 {
432     short sw;
433     bool res;
434     version(AsmX86)
435     {
436         mixin(fld_arg!("y"));
437         mixin(fld_arg!("x"));
438         asm nothrow @nogc pure @trusted
439         {
440             fucomip ST(1);
441             setb    AL;
442             setnp   AH;
443             and     AL,AH;
444             mov     res,AL;
445             fstp    ST(0);
446         }
447     }
448     return res;
449 }
450 
451 bool ld_cmpbe(longdouble_soft x, longdouble_soft y)
452 {
453     short sw;
454     bool res;
455     version(AsmX86)
456     {
457         mixin(fld_arg!("y"));
458         mixin(fld_arg!("x"));
459         asm nothrow @nogc pure @trusted
460         {
461             fucomip ST(1);
462             setbe   AL;
463             setnp   AH;
464             and     AL,AH;
465             mov     res,AL;
466             fstp    ST(0);
467         }
468     }
469     return res;
470 }
471 
472 bool ld_cmpa(longdouble_soft x, longdouble_soft y)
473 {
474     short sw;
475     bool res;
476     version(AsmX86)
477     {
478         mixin(fld_arg!("y"));
479         mixin(fld_arg!("x"));
480         asm nothrow @nogc pure @trusted
481         {
482             fucomip ST(1);
483             seta    AL;
484             setnp   AH;
485             and     AL,AH;
486             mov     res,AL;
487             fstp    ST(0);
488         }
489     }
490     return res;
491 }
492 
493 bool ld_cmpae(longdouble_soft x, longdouble_soft y)
494 {
495     short sw;
496     bool res;
497     version(AsmX86)
498     {
499         mixin(fld_arg!("y"));
500         mixin(fld_arg!("x"));
501         asm nothrow @nogc pure @trusted
502         {
503             fucomip ST(1);
504             setae   AL;
505             setnp   AH;
506             and     AL,AH;
507             mov     res,AL;
508             fstp    ST(0);
509         }
510     }
511     return res;
512 }
513 
514 bool ld_cmpe(longdouble_soft x, longdouble_soft y)
515 {
516     short sw;
517     bool res;
518     version(AsmX86)
519     {
520         mixin(fld_arg!("y"));
521         mixin(fld_arg!("x"));
522         asm nothrow @nogc pure @trusted
523         {
524             fucomip ST(1);
525             sete    AL;
526             setnp   AH;
527             and     AL,AH;
528             mov     res,AL;
529             fstp    ST(0);
530         }
531     }
532     return res;
533 }
534 
535 bool ld_cmpne(longdouble_soft x, longdouble_soft y)
536 {
537     short sw;
538     bool res;
539     version(AsmX86)
540     {
541         mixin(fld_arg!("y"));
542         mixin(fld_arg!("x"));
543         asm nothrow @nogc pure @trusted
544         {
545             fucomip ST(1);
546             setne   AL;
547             setp    AH;
548             or      AL,AH;
549             mov     res,AL;
550             fstp    ST(0);
551         }
552     }
553     return res;
554 }
555 
556 int ld_cmp(longdouble_soft x, longdouble_soft y)
557 {
558     // return -1 if x < y, 0 if x == y or unordered, 1 if x > y
559     short sw;
560     int res;
561     version(AsmX86)
562     {
563         mixin(fld_arg!("y"));
564         mixin(fld_arg!("x"));
565         asm nothrow @nogc pure @trusted
566         {
567             fucomip ST(1);
568             seta    AL;
569             setb    AH;
570             setp    DL;
571             or      AL, DL;
572             or      AH, DL;
573             sub     AL, AH;
574             movsx   EAX, AL;
575             fstp    ST(0);
576             mov     res, EAX;
577         }
578     }
579 }
580 
581 
582 int _isnan(longdouble_soft ld)
583 {
584     return (ld.exponent == 0x7fff && ld.mantissa != 0 && ld.mantissa != (1L << 63)); // exclude pseudo-infinity and infinity, but not FP Indefinite
585 }
586 
587 longdouble_soft fabsl(longdouble_soft ld)
588 {
589     ld.exp_sign = ld.exponent;
590     return ld;
591 }
592 
593 longdouble_soft sqrtl(longdouble_soft ld)
594 {
595     version(AsmX86)
596     {
597         mixin(fld_arg!("ld"));
598         asm nothrow @nogc pure @trusted
599         {
600             fsqrt;
601         }
602         mixin(fstp_arg!("ld"));
603     }
604     return ld;
605 }
606 
607 longdouble_soft sqrt(longdouble_soft ld) { return sqrtl(ld); }
608 
609 longdouble_soft sinl (longdouble_soft ld)
610 {
611     version(AsmX86)
612     {
613         mixin(fld_arg!("ld"));
614         asm nothrow @nogc pure @trusted
615         {
616             fsin; // exact for |x|<=PI/4
617         }
618         mixin(fstp_arg!("ld"));
619     }
620     return ld;
621 }
622 longdouble_soft cosl (longdouble_soft ld)
623 {
624     version(AsmX86)
625     {
626         mixin(fld_arg!("ld"));
627         asm nothrow @nogc pure @trusted
628         {
629             fcos; // exact for |x|<=PI/4
630         }
631         mixin(fstp_arg!("ld"));
632     }
633     return ld;
634 }
635 longdouble_soft tanl (longdouble_soft ld)
636 {
637     version(AsmX86)
638     {
639         mixin(fld_arg!("ld"));
640         asm nothrow @nogc pure @trusted
641         {
642             fptan;
643             fstp ST(0); // always 1
644         }
645         mixin(fstp_arg!("ld"));
646     }
647     return ld;
648 }
649 
650 longdouble_soft fmodl(longdouble_soft x, longdouble_soft y)
651 {
652     return ld_mod(x, y);
653 }
654 
655 longdouble_soft ld_mod(longdouble_soft x, longdouble_soft y)
656 {
657     short sw;
658     version(AsmX86)
659     {
660         mixin(fld_arg!("y"));
661         mixin(fld_arg!("x"));
662         asm nothrow @nogc pure @trusted
663         {
664         FM1:    // We don't use fprem1 because for some inexplicable
665                 // reason we get -5 when we do _modulo(15, 10)
666             fprem;                          // ST = ST % ST1
667             fstsw   word ptr sw;
668             fwait;
669             mov     AH,byte ptr sw+1;       // get msb of status word in AH
670             sahf;                           // transfer to flags
671             jp      FM1;                    // continue till ST < ST1
672             fstp    ST(1);                  // leave remainder on stack
673         }
674         mixin(fstp_arg!("x"));
675     }
676     return x;
677 }
678 
679 //////////////////////////////////////////////////////////////
680 
681 @safe:
682 
683 __gshared const
684 {
685     longdouble_soft ld_qnan = longdouble_soft(0xC000000000000000UL, 0x7fff);
686     longdouble_soft ld_inf  = longdouble_soft(0x8000000000000000UL, 0x7fff);
687 
688     longdouble_soft ld_zero  = longdouble_soft(0, 0);
689     longdouble_soft ld_one   = longdouble_soft(0x8000000000000000UL, 0x3fff);
690     longdouble_soft ld_pi    = longdouble_soft(0xc90fdaa22168c235UL, 0x4000);
691     longdouble_soft ld_log2t = longdouble_soft(0xd49a784bcd1b8afeUL, 0x4000);
692     longdouble_soft ld_log2e = longdouble_soft(0xb8aa3b295c17f0bcUL, 0x3fff);
693     longdouble_soft ld_log2  = longdouble_soft(0x9a209a84fbcff799UL, 0x3ffd);
694     longdouble_soft ld_ln2   = longdouble_soft(0xb17217f7d1cf79acUL, 0x3ffe);
695 
696     longdouble_soft ld_pi2     = longdouble_soft(0xc90fdaa22168c235UL, 0x4001);
697     longdouble_soft ld_piOver2 = longdouble_soft(0xc90fdaa22168c235UL, 0x3fff);
698     longdouble_soft ld_piOver4 = longdouble_soft(0xc90fdaa22168c235UL, 0x3ffe);
699 
700     longdouble_soft twoPow63 = longdouble_soft(1UL << 63, 0x3fff + 63);
701 }
702 
703 //////////////////////////////////////////////////////////////
704 
705 enum LD_TYPE_OTHER    = 0;
706 enum LD_TYPE_ZERO     = 1;
707 enum LD_TYPE_INFINITE = 2;
708 enum LD_TYPE_SNAN     = 3;
709 enum LD_TYPE_QNAN     = 4;
710 
711 int ld_type(longdouble_soft x)
712 {
713     // see https://en.wikipedia.org/wiki/Extended_precision
714     if(x.exponent == 0)
715         return x.mantissa == 0 ? LD_TYPE_ZERO : LD_TYPE_OTHER; // dnormal if not zero
716     if(x.exponent != 0x7fff)
717         return LD_TYPE_OTHER;    // normal or denormal
718     uint  upper2  = x.mantissa >> 62;
719     ulong lower62 = x.mantissa & ((1L << 62) - 1);
720     if(upper2 == 0 && lower62 == 0)
721         return LD_TYPE_INFINITE; // pseudo-infinity
722     if(upper2 == 2 && lower62 == 0)
723         return LD_TYPE_INFINITE; // infinity
724     if(upper2 == 2 && lower62 != 0)
725         return LD_TYPE_SNAN;
726     return LD_TYPE_QNAN;         // qnan, indefinite, pseudo-nan
727 }
728 
729 // consider sprintf pure
730 private extern(C) int sprintf(scope char* s, scope const char* format, ...) pure @nogc nothrow;
731 
732 size_t ld_sprint(char* str, int fmt, longdouble_soft x) @system
733 {
734     // ensure dmc compatible strings for nan and inf
735     switch(ld_type(x))
736     {
737         case LD_TYPE_QNAN:
738         case LD_TYPE_SNAN:
739             return sprintf(str, "nan");
740         case LD_TYPE_INFINITE:
741             return sprintf(str, x.sign ? "-inf" : "inf");
742         default:
743             break;
744     }
745 
746     // fmt is 'a','A','f' or 'g'
747     if(fmt != 'a' && fmt != 'A')
748     {
749         char[3] format = ['%', cast(char)fmt, 0];
750         return sprintf(str, format.ptr, ld_read(&x));
751     }
752 
753     ushort exp = x.exponent;
754     ulong mantissa = x.mantissa;
755 
756     if(ld_type(x) == LD_TYPE_ZERO)
757         return sprintf(str, fmt == 'a' ? "0x0.0L" : "0X0.0L");
758 
759     size_t len = 0;
760     if(x.sign)
761         str[len++] = '-';
762     str[len++] = '0';
763     str[len++] = cast(char)('X' + fmt - 'A');
764     str[len++] = mantissa & (1L << 63) ? '1' : '0';
765     str[len++] = '.';
766     mantissa = mantissa << 1;
767     while(mantissa)
768     {
769         int dig = (mantissa >> 60) & 0xf;
770         dig += dig < 10 ? '0' : fmt - 10;
771         str[len++] = cast(char)dig;
772         mantissa = mantissa << 4;
773     }
774     str[len++] = cast(char)('P' + fmt - 'A');
775     if(exp < 0x3fff)
776     {
777         str[len++] = '-';
778         exp = cast(ushort)(0x3fff - exp);
779     }
780     else
781     {
782         str[len++] = '+';
783         exp = cast(ushort)(exp - 0x3fff);
784     }
785     size_t exppos = len;
786     for(int i = 12; i >= 0; i -= 4)
787     {
788         int dig = (exp >> i) & 0xf;
789         if(dig != 0 || len > exppos || i == 0)
790             str[len++] = cast(char)(dig + (dig < 10 ? '0' : fmt - 10));
791     }
792     str[len] = 0;
793     return len;
794 }
795 
796 //////////////////////////////////////////////////////////////
797 
798 @system unittest
799 {
800     import core.stdc.string;
801     import core.stdc.stdio;
802 
803     char[32] buffer;
804     ld_sprint(buffer.ptr, 'a', ld_pi);
805     assert(strcmp(buffer.ptr, "0x1.921fb54442d1846ap+1") == 0);
806 
807     auto len = ld_sprint(buffer.ptr, 'g', longdouble_soft(2.0));
808     assert(buffer[0 .. len] == "2.00000" || buffer[0 .. len] == "2"); // Win10 - 64bit
809 
810     ld_sprint(buffer.ptr, 'g', longdouble_soft(1_234_567.89));
811     assert(strcmp(buffer.ptr, "1.23457e+06") == 0);
812 
813     ld_sprint(buffer.ptr, 'g', ld_inf);
814     assert(strcmp(buffer.ptr, "inf") == 0);
815 
816     ld_sprint(buffer.ptr, 'g', ld_qnan);
817     assert(strcmp(buffer.ptr, "nan") == 0);
818 
819     longdouble_soft ldb = longdouble_soft(0.4);
820     long b = cast(long)ldb;
821     assert(b == 0);
822 
823     b = cast(long)longdouble_soft(0.9);
824     assert(b == 0);
825 
826     long x = 0x12345678abcdef78L;
827     longdouble_soft ldx = longdouble_soft(x);
828     assert(ldx > ld_zero);
829     long y = cast(long)ldx;
830     assert(x == y);
831 
832     x = -0x12345678abcdef78L;
833     ldx = longdouble_soft(x);
834     assert(ldx < ld_zero);
835     y = cast(long)ldx;
836     assert(x == y);
837 
838     ulong u = 0x12345678abcdef78L;
839     longdouble_soft ldu = longdouble_soft(u);
840     assert(ldu > ld_zero);
841     ulong v = cast(ulong)ldu;
842     assert(u == v);
843 
844     u = 0xf234567812345678UL;
845     ldu = longdouble_soft(u);
846     assert(ldu > ld_zero);
847     v = cast(ulong)ldu;
848     assert(u == v);
849 
850     u = 0xf2345678;
851     ldu = longdouble_soft(u);
852     ldu = ldu * ldu;
853     ldu = sqrt(ldu);
854     v = cast(ulong)ldu;
855     assert(u == v);
856 
857     u = 0x123456789A;
858     ldu = longdouble_soft(u);
859     ldu = ldu * longdouble_soft(1L << 23);
860     v = cast(ulong)ldu;
861     u = u * (1L << 23);
862     assert(u == v);
863 }