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