/*
Written 2025-11-07 by Brendan G Bohannon.
Do whatever with this, just assume Public Domain or MIT-0 or CC0 or whatever.
I take no responsibility for anything that may result from using this code.
*/
/*
Decimal128 DPD test stuff...
*/
#include
#include
#include
#include
#include
#include
#ifndef BTIC2F_BYTE
#define BTIC2F_BYTE
typedef unsigned char byte;
typedef unsigned short u16;
typedef unsigned int u32;
typedef unsigned long long u64;
typedef signed char sbyte;
typedef signed short s16;
typedef signed int s32;
typedef signed long long s64;
#endif
#define TKD128_MASK30 0x3FFFFFFFU
#define TKD128_MASK20 0x000FFFFFU
#define TKD128_MASK46 0x00003FFFFFFFFFFFULL
// #define TKD128_BIAS 6176
#define TKD128_BIAS 6143
typedef struct TKD128_Packed128_s TKD128_Packed128;
typedef struct TKD128_Unpacked4x30_s TKD128_Unpacked4x30;
typedef struct TKD128_Unpacked8x30_s TKD128_Unpacked8x30;
struct TKD128_Packed128_s {
u64 lo;
u64 hi;
};
struct TKD128_Unpacked4x30_s {
u32 m[4];
s16 exp;
byte sgn;
};
struct TKD128_Unpacked8x30_s {
u32 m[8];
s16 exp;
byte sgn;
};
char *TKD128_OprDumpAsString(
TKD128_Unpacked4x30 val);
u16 tkd128_unpackdpd[1024];
u16 tkd128_packdpd[1024];
int TKD128_InitDPD()
{
int d0, d1, d2;
int i, j, k;
if(tkd128_packdpd[999])
return(0);
for(i=0; i<1024; i++)
{
if(i&8)
{
if(((i>>1)&3)==0)
{
d0=(i>>7)&7;
d1=(i>>4)&7;
d2=8+(i&1);
}else
if(((i>>1)&3)==1)
{
d0=(i>>7)&7;
d1=8+((i>>4)&1);
d2=((i>>4)&6)|(i&1);
}else
if(((i>>1)&3)==2)
{
d0=8+((i>>7)&1);
d1=(i>>4)&7;
d2=((i>>7)&6)|(i&1);
}else
{
if(((i>>5)&3)==0)
{
d0=8+((i>>7)&1);
d1=8+((i>>4)&1);
d2=((i>>7)&6)|(i&1);
}else
if(((i>>5)&3)==1)
{
d0=8+((i>>7)&1);
d1=((i>>7)&6)+((i>>4)&1);
d2=8+((i>>0)&1);
}else
if(((i>>5)&3)==2)
{
d0=((i>>7)&6)+((i>>7)&1);
d1=8+((i>>4)&1);
d2=8+((i>>0)&1);
}else
{
d0=8+((i>>7)&1);
d1=8+((i>>4)&1);
d2=8+((i>>0)&1);
}
}
}else
{
d0=(i>>7)&7;
d1=(i>>4)&7;
d2=(i>>0)&7;
}
j=d0*100+d1*10+d2;
tkd128_unpackdpd[i]=j;
if(!tkd128_packdpd[j])
tkd128_packdpd[j]=i;
}
return(1);
}
TKD128_Unpacked4x30 TKD128_UnpackValueX30(TKD128_Packed128 val)
{
TKD128_Unpacked4x30 tmp;
int mc, mcd, ex;
tmp.m[0]=(val.lo>> 0)&TKD128_MASK30;
tmp.m[1]=(val.lo>>30)&TKD128_MASK30;
tmp.m[2]=((val.hi<< 4)|(val.lo>>60))&TKD128_MASK30;
tmp.m[3]=(val.hi>>26)&TKD128_MASK20;
mc=(val.hi>>46)&0x1FFFF;
if((mc>>15)==3)
{ mcd=8+(mc>>12)&1; ex=(mc&0xFFF)|(((mc>>13)&3)<<12); }
else
{ mcd=(mc>>12)&7; ex=(mc&0xFFF)|(((mc>>15)&3)<<12); }
tmp.m[3]+=mcd*1000000;
tmp.exp=ex;
tmp.sgn=val.hi>>63;
return(tmp);
}
TKD128_Packed128 TKD128_PackValueX30(TKD128_Unpacked4x30 val)
{
TKD128_Packed128 tmp;
int m3n, mcd, mc, ex;
m3n=val.m[3]%1000000000;
mcd=val.m[3]/1000000000;
ex=val.exp;
if(mcd>=8)
{ mc=0x18000|(ex&0xFFF)|((mcd&1)<<12)|(((ex>>12)&3)<<13); }
else
{ mc=(ex&0xFFF)|((mcd&7)<<12)|(((ex>>12)&3)<<15); }
tmp.lo=val.m[0]|(((u64)(val.m[1]))<<30)|(((u64)(val.m[2]))<<60);
tmp.hi=(val.m[2]>>4)|(((u64)m3n)<<26)|
(((u64)mcd)<<46)|(((u64)val.sgn)<<63);
return(tmp);
}
#if 1
void TKD128_UnpackValueX30R(
TKD128_Unpacked4x30 *rtmp, TKD128_Packed128 *rval)
{
// TKD128_Unpacked4x30 tmp;
u64 lo, hi;
int mc, mcd, ex;
lo=rval->lo;
hi=rval->hi;
rtmp->m[0]=(lo>> 0)&TKD128_MASK30;
rtmp->m[1]=(lo>>30)&TKD128_MASK30;
rtmp->m[2]=((hi<< 4)|(lo>>60))&TKD128_MASK30;
rtmp->m[3]=(hi>>26)&TKD128_MASK20;
mc=(hi>>46)&0x1FFFF;
if((mc>>15)==3)
{
mcd=8+(mc>>12)&1;
ex=(mc&0xFFF)|(((mc>>13)&3)<<12);
}else
{
mcd=(mc>>12)&7;
ex=(mc&0xFFF)|(((mc>>15)&3)<<12);
}
rtmp->m[3]+=mcd*1000000;
rtmp->exp=ex;
rtmp->sgn=hi>>63;
}
void TKD128_PackValueX30R(
TKD128_Unpacked4x30 *rval, TKD128_Packed128 *rtmp)
{
u32 m0, m1, m2, m3;
int m3n, mcd, mc, ex, sgn;
m0=rval->m[0];
m1=rval->m[1];
m2=rval->m[2];
m3=rval->m[3];
m3n=m3%1000000000;
mcd=m3/1000000000;
sgn=rval->sgn;
ex=rval->exp;
if(mcd>=8)
{ mc=0x18000|(ex&0xFFF)|((mcd&1)<<12)|(((ex>>12)&3)<<13); }
else
{ mc=(ex&0xFFF)|((mcd&7)<<12)|(((ex>>12)&3)<<15); }
rtmp->lo=m0|(((u64)(m1))<<30)|(((u64)(m2))<<60);
rtmp->hi=(m2>>4)|(((u64)m3n)<<26)|
(((u64)mcd)<<46)|(((u64)sgn)<<63);
}
#endif
u32 TKD128_UnpackBitsDPD30(u32 val)
{
u32 v0, v1, v2, v;
v0=tkd128_unpackdpd[(val>> 0)&1023];
v1=tkd128_unpackdpd[(val>>10)&1023];
v2=tkd128_unpackdpd[(val>>20)&1023];
v=v0+(v1*1000)+(v2*1000000);
return(v);
}
u32 TKD128_PackBitsDPD30(u32 val)
{
u32 v0, v1, v2, v;
v=val;
v1=v/1000;
v0=v-(v1*1000);
v2=v1/1000;
v1=v1-(v2*1000);
v0=tkd128_packdpd[v0];
v1=tkd128_packdpd[v1];
v2=tkd128_packdpd[v2];
v=v0|(v1<<10)|(v2<<20);
return(v);
}
TKD128_Unpacked4x30 TKD128_UnpackValueDPD(TKD128_Packed128 val)
{
TKD128_Unpacked4x30 tmp;
int mc, mcd, ex;
tmp.m[0]=TKD128_UnpackBitsDPD30((val.lo>> 0)&TKD128_MASK30);
tmp.m[1]=TKD128_UnpackBitsDPD30((val.lo>>30)&TKD128_MASK30);
tmp.m[2]=TKD128_UnpackBitsDPD30(((val.hi<< 4)|(val.lo>>60))&TKD128_MASK30);
tmp.m[3]=TKD128_UnpackBitsDPD30((val.hi>>26)&TKD128_MASK20);
mc=(val.hi>>46)&0x1FFFF;
if((mc>>15)==3)
{
mcd=8+(mc>>12)&1;
ex=(mc&0xFFF)|(((mc>>13)&3)<<12);
}else
{
mcd=(mc>>12)&7;
ex=(mc&0xFFF)|(((mc>>15)&3)<<12);
}
tmp.m[3]+=mcd*1000000;
tmp.exp=ex;
tmp.sgn=val.hi>>63;
return(tmp);
}
TKD128_Packed128 TKD128_PackValueDPD(TKD128_Unpacked4x30 val)
{
TKD128_Packed128 tmp;
u32 m0, m1, m2, m3;
int m3n, mcd, mc, ex;
mcd=val.m[3]/1000000;
m3n=val.m[3]%1000000;
ex=val.exp;
if(mcd>=8)
{ mc=0x18000|(ex&0xFFF)|((mcd&1)<<12)|(((ex>>12)&3)<<13); }
else
{ mc=(ex&0xFFF)|((mcd&7)<<12)|(((ex>>12)&3)<<15); }
m0=TKD128_PackBitsDPD30(val.m[0]);
m1=TKD128_PackBitsDPD30(val.m[1]);
m2=TKD128_PackBitsDPD30(val.m[2]);
m3=TKD128_PackBitsDPD30(m3n);
tmp.lo=m0|(((u64)m1)<<30)|(((u64)m2)<<60);
tmp.hi=(m2>>4)|(((u64)m3)<<26)|
(((u64)mc)<<46)|(((u64)val.sgn)<<63);
return(tmp);
}
int TKD128_ScaleLeftArray8(u32 *arr, int shl)
{
int sh;
sh=shl;
while(sh>=9)
{
sh-=9;
arr[7]=arr[6];
arr[6]=arr[5];
arr[5]=arr[4];
arr[4]=arr[3];
arr[3]=arr[2];
arr[2]=arr[1];
arr[1]=arr[0];
arr[0]=0;
}
while(sh>=3)
{
sh-=3;
arr[7]=((arr[7]%1000000)*1000)+(arr[6]/1000000);
arr[6]=((arr[6]%1000000)*1000)+(arr[5]/1000000);
arr[5]=((arr[5]%1000000)*1000)+(arr[4]/1000000);
arr[4]=((arr[4]%1000000)*1000)+(arr[3]/1000000);
arr[3]=((arr[3]%1000000)*1000)+(arr[2]/1000000);
arr[2]=((arr[2]%1000000)*1000)+(arr[1]/1000000);
arr[1]=((arr[1]%1000000)*1000)+(arr[0]/1000000);
arr[0]=((arr[0]%1000000)*1000);
}
while(sh>=1)
{
sh--;
arr[7]=((arr[7]%100000000)*10)+(arr[6]/100000000);
arr[6]=((arr[6]%100000000)*10)+(arr[5]/100000000);
arr[5]=((arr[5]%100000000)*10)+(arr[4]/100000000);
arr[4]=((arr[4]%100000000)*10)+(arr[3]/100000000);
arr[3]=((arr[3]%100000000)*10)+(arr[2]/100000000);
arr[2]=((arr[2]%100000000)*10)+(arr[1]/100000000);
arr[1]=((arr[1]%100000000)*10)+(arr[0]/100000000);
arr[0]=((arr[0]%100000000)*10);
}
return(0);
}
int TKD128_ScaleRightArray8(u32 *arr, int shl)
{
int sh;
sh=shl;
while(sh>=9)
{
sh-=9;
arr[0]=arr[1];
arr[1]=arr[2];
arr[2]=arr[3];
arr[3]=arr[4];
arr[4]=arr[5];
arr[5]=arr[6];
arr[6]=arr[7];
arr[7]=0;
}
while(sh>=3)
{
sh-=3;
arr[0]=((arr[1]%1000)*1000000)+(arr[0]/1000);
arr[1]=((arr[2]%1000)*1000000)+(arr[1]/1000);
arr[2]=((arr[3]%1000)*1000000)+(arr[2]/1000);
arr[3]=((arr[4]%1000)*1000000)+(arr[3]/1000);
arr[4]=((arr[5]%1000)*1000000)+(arr[4]/1000);
arr[5]=((arr[6]%1000)*1000000)+(arr[5]/1000);
arr[6]=((arr[7]%1000)*1000000)+(arr[6]/1000);
arr[7]=(arr[7]/1000);
}
while(sh>=1)
{
sh--;
arr[0]=((arr[1]%10)*100000000)+(arr[0]/10);
arr[1]=((arr[2]%10)*100000000)+(arr[1]/10);
arr[2]=((arr[3]%10)*100000000)+(arr[2]/10);
arr[3]=((arr[4]%10)*100000000)+(arr[3]/10);
arr[4]=((arr[5]%10)*100000000)+(arr[4]/10);
arr[5]=((arr[6]%10)*100000000)+(arr[5]/10);
arr[6]=((arr[7]%10)*100000000)+(arr[6]/10);
arr[7]=(arr[7]/10);
}
return(0);
}
int TKD128_AddArray8(u32 *arr_a, u32 *arr_b, u32 *arr_c)
{
u32 v, c;
int i;
c=0;
for(i=0; i<8; i++)
{
v=arr_a[i]+arr_b[i]+c;
c=0;
if(v>=1000000000U)
{ c=1; v-=1000000000U; }
arr_c[i]=v;
}
return(0);
}
int TKD128_SubArray8(u32 *arr_a, u32 *arr_b, u32 *arr_c)
{
u32 v, c;
int i;
c=0;
for(i=0; i<8; i++)
{
v=arr_a[i]-arr_b[i]-c;
c=0;
if(v>=1000000000U)
{ c=1; v+=1000000000U; }
arr_c[i]=v;
}
return(0);
}
TKD128_Unpacked8x30 TKD128_Conv4to8x30(
TKD128_Unpacked4x30 val)
{
TKD128_Unpacked8x30 tmp;
tmp.sgn=val.sgn;
tmp.exp=val.exp;
tmp.m[7]=val.m[3];
tmp.m[6]=val.m[2];
tmp.m[5]=val.m[1];
tmp.m[4]=val.m[0];
tmp.m[3]=0;
tmp.m[2]=0;
tmp.m[1]=0;
tmp.m[0]=0;
return(tmp);
}
TKD128_Unpacked4x30 TKD128_OprNormalize8to4x30(TKD128_Unpacked8x30 valc0)
{
TKD128_Unpacked4x30 valc;
u32 m0, m1, m2, m3;
if(valc0.m[7]>=10000000)
{ valc0.exp++; TKD128_ScaleRightArray8(valc0.m, 1); }
while(!valc0.m[7] && !valc0.m[6] && (valc0.exp>=9))
{ valc0.exp-=9; TKD128_ScaleLeftArray8(valc0.m, 9); }
while(!valc0.m[7] && (valc0.exp>=6))
{ valc0.exp-=6; TKD128_ScaleLeftArray8(valc0.m, 6); }
while((valc0.m[7]<=1000) && (valc0.exp>=3))
{ valc0.exp-=3; TKD128_ScaleLeftArray8(valc0.m, 3); }
while((valc0.m[7]<1000000) && (valc0.exp>=1))
{ valc0.exp--; TKD128_ScaleLeftArray8(valc0.m, 1); }
m0=valc0.m[4];
m1=valc0.m[5];
m2=valc0.m[6];
m3=valc0.m[7];
m0+=(valc0.m[3]>500000000);
if(m0>=1000000000)
{ m1++; m0-=1000000000; }
if(m1>=1000000000)
{ m2++; m1-=1000000000; }
if(m2>=1000000000)
{ m3++; m2-=1000000000; }
valc.sgn=valc0.sgn;
valc.exp=valc0.exp;
valc.m[0]=m0; valc.m[1]=m1;
valc.m[2]=m2; valc.m[3]=m3;
return(valc);
}
TKD128_Unpacked4x30 TKD128_OprAdd(
TKD128_Unpacked4x30 vala,
TKD128_Unpacked4x30 valb)
{
TKD128_Unpacked8x30 vala0, valb0, valc0;
TKD128_Unpacked4x30 valc;
int maxexp, flip;
vala0=TKD128_Conv4to8x30(vala);
valb0=TKD128_Conv4to8x30(valb);
maxexp=vala0.exp;
flip=0;
if(valb0.exp>maxexp)
{
maxexp=valb0.exp;
flip=1;
}else
if(vala0.exp==valb0.exp)
{
if(valb0.m[7]>vala0.m[7])
{ flip=1; }
else if(valb0.m[7]==vala0.m[7])
{ if(valb0.m[6]>vala0.m[6])
{ flip=1; }
else if(valb0.m[6]==vala0.m[6])
{ if(valb0.m[5]>vala0.m[5])
{ flip=1; }
else if(valb0.m[5]==vala0.m[5])
{ if(valb0.m[4]>vala0.m[4])
{ flip=1; }
}
}
}
}
if(flip)
{
valc0=vala0;
vala0=valb0;
valb0=valc0;
}
TKD128_ScaleRightArray8(vala0.m, maxexp-vala0.exp);
TKD128_ScaleRightArray8(valb0.m, maxexp-valb0.exp);
vala0.exp=maxexp;
valb0.exp=maxexp;
valc0.exp=maxexp;
valc0.sgn=vala0.sgn;
if(vala0.sgn==valb0.sgn)
{ TKD128_AddArray8(vala0.m, valb0.m, valc0.m); }
else
{ TKD128_SubArray8(vala0.m, valb0.m, valc0.m); }
valc=TKD128_OprNormalize8to4x30(valc0);
return(valc);
}
TKD128_Unpacked4x30 TKD128_OprSub(
TKD128_Unpacked4x30 vala,
TKD128_Unpacked4x30 valb)
{
TKD128_Unpacked4x30 valb1;
valb1=valb;
valb1.sgn=!valb.sgn;
return(TKD128_OprAdd(vala, valb1));
}
TKD128_Unpacked4x30 TKD128_OprMul(
TKD128_Unpacked4x30 vala,
TKD128_Unpacked4x30 valb)
{
TKD128_Unpacked8x30 valc0;
TKD128_Unpacked4x30 valc;
u64 m0_aa, m0_ab, m0_ac, m0_ad;
u64 m0_ba, m0_bb, m0_bc, m0_bd;
u64 m0_ca, m0_cb, m0_cc, m0_cd;
u64 m0_da, m0_db, m0_dc, m0_dd;
u64 m1a, m1b, m1c, m1d;
valc0.exp=vala.exp+valb.exp-TKD128_BIAS;
valc0.sgn=vala.sgn^valb.sgn;
m0_aa=((u64)vala.m[3])*((u64)valb.m[3]);
m0_ab=((u64)vala.m[3])*((u64)valb.m[2]);
m0_ac=((u64)vala.m[3])*((u64)valb.m[1]);
m0_ad=((u64)vala.m[3])*((u64)valb.m[0]);
m0_ba=((u64)vala.m[2])*((u64)valb.m[3]);
m0_bb=((u64)vala.m[2])*((u64)valb.m[2]);
m0_bc=((u64)vala.m[2])*((u64)valb.m[1]);
m0_bd=((u64)vala.m[2])*((u64)valb.m[0]);
m0_ca=((u64)vala.m[1])*((u64)valb.m[3]);
m0_cb=((u64)vala.m[1])*((u64)valb.m[2]);
m0_cc=((u64)vala.m[1])*((u64)valb.m[1]);
m0_cd=((u64)vala.m[1])*((u64)valb.m[0]);
m0_da=((u64)vala.m[0])*((u64)valb.m[3]);
m0_db=((u64)vala.m[0])*((u64)valb.m[2]);
m0_dc=((u64)vala.m[0])*((u64)valb.m[1]);
m0_dd=((u64)vala.m[0])*((u64)valb.m[0]);
/*
* AA
* AB
* BA
* BB
* AC
* CA
* AD
* DA
* BC
* CB
* CC
* BD
* DB
* CD
* DC
* DD
*/
m1a=m0_aa+
(m0_ba/1000000000ULL) +
(m0_ab/1000000000ULL) ;
m1b=m0_bb+
m0_ac+
m0_ca+
(m0_cb/1000000000ULL) +
(m0_bc/1000000000ULL) +
(m0_ad/1000000000ULL) +
(m0_da/1000000000ULL) +
((m0_ba%1000000000ULL)*1000000000ULL) +
((m0_ab%1000000000ULL)*1000000000ULL) ;
m1c=m0_cc+
m0_bd+
m0_db+
(m0_dc/1000000000ULL) +
(m0_cd/1000000000ULL) +
((m0_ad%1000000000ULL)*1000000000ULL) +
((m0_da%1000000000ULL)*1000000000ULL) +
((m0_cb%1000000000ULL)*1000000000ULL) +
((m0_bc%1000000000ULL)*1000000000ULL) ;
m1d=m0_dd+
((m0_dc%1000000000ULL)*1000000000ULL) +
((m0_cd%1000000000ULL)*1000000000ULL) ;
while(m1d>=1000000000000000000ULL)
{ m1d-=1000000000000000000ULL; m1c++; }
while(m1c>=1000000000000000000ULL)
{ m1c-=1000000000000000000ULL; m1b++; }
while(m1b>=1000000000000000000ULL)
{ m1b-=1000000000000000000ULL; m1a++; }
valc0.m[7]=m1a/1000000000ULL;
valc0.m[6]=m1a%1000000000ULL;
valc0.m[5]=m1b/1000000000ULL;
valc0.m[4]=m1b%1000000000ULL;
valc0.m[3]=m1c/1000000000ULL;
valc0.m[2]=m1c%1000000000ULL;
valc0.m[1]=m1d/1000000000ULL;
valc0.m[0]=m1d%1000000000ULL;
TKD128_ScaleLeftArray8(valc0.m, 3);
valc=TKD128_OprNormalize8to4x30(valc0);
return(valc);
}
TKD128_Unpacked4x30 TKD128_OprRecipApx(
TKD128_Unpacked4x30 val)
{
TKD128_Unpacked4x30 valc;
valc.sgn=val.sgn;
valc.exp=(2*TKD128_BIAS-1)-val.exp;
valc.m[3]=1000000000-val.m[3];
valc.m[2]=1000000000-val.m[2];
valc.m[1]=1000000000-val.m[1];
valc.m[0]=1000000000-val.m[0];
return(valc);
}
int TKD128_OprCompare(
TKD128_Unpacked4x30 vala,
TKD128_Unpacked4x30 valb)
{
if(vala.sgn!=valb.sgn)
{
if(!vala.sgn)
return(1);
return(-1);
}
if(vala.exp!=valb.exp)
{
if(vala.exp>valb.exp)
return(1);
return(-1);
}
if(vala.m[3]!=valb.m[3])
{
if(vala.m[3]>valb.m[3])
return(1);
return(-1);
}
if(vala.m[2]!=valb.m[2])
{
if(vala.m[2]>valb.m[2])
return(1);
return(-1);
}
if(vala.m[1]!=valb.m[1])
{
if(vala.m[1]>valb.m[1])
return(1);
return(-1);
}
if(vala.m[0]!=valb.m[0])
{
if(vala.m[0]>valb.m[0])
return(1);
return(-1);
}
return(0);
}
TKD128_Unpacked4x30 TKD128_OprRecip(
TKD128_Unpacked4x30 val)
{
static TKD128_Unpacked4x30 cnst_0p5 ={{0,0,0, 500000}, TKD128_BIAS, 0};
static TKD128_Unpacked4x30 cnst_2p0 ={{0,0,0,2000000}, TKD128_BIAS, 0};
static TKD128_Unpacked4x30 cnst_1p0 ={{0,0,0,1000000}, TKD128_BIAS, 0};
static TKD128_Unpacked4x30 cnst_3p0 ={{0,0,0,3000000}, TKD128_BIAS, 0};
static TKD128_Unpacked4x30 cnst_0p75={{0,0,0, 750000}, TKD128_BIAS, 0};
static TKD128_Unpacked4x30 cnst_1p25={{0,0,0,1250000}, TKD128_BIAS, 0};
static TKD128_Unpacked4x30 cnst_0p875={{0,0,0, 875000}, TKD128_BIAS, 0};
static TKD128_Unpacked4x30 cnst_1p125={{0,0,0,1125000}, TKD128_BIAS, 0};
TKD128_Unpacked4x30 valc, valb, valt;
int d0, d1, d2, d3, de;
int i;
valc=TKD128_OprRecipApx(val);
// printf("Rcp Apx-0:\n %s\n %s\n",
// TKD128_OprDumpAsString(val),
// TKD128_OprDumpAsString(valc));
for(i=0; i<32; i++)
{
valb=TKD128_OprMul(val, valc);
de=valb.exp-TKD128_BIAS;
if(de>0)
{ valc=TKD128_OprMul(valc, cnst_0p5); continue; }
if(de<0)
{ valc=TKD128_OprMul(valc, cnst_2p0); continue; }
if(valb.m[3]>=2000000)
{ valc=TKD128_OprMul(valc, cnst_0p5); continue; }
if(valb.m[3]<= 500000)
{ valc=TKD128_OprMul(valc, cnst_2p0); continue; }
if(valb.m[3]>=1500000)
{ valc=TKD128_OprMul(valc, cnst_0p75); continue; }
if(valb.m[3]<= 750000)
{ valc=TKD128_OprMul(valc, cnst_1p25); continue; }
if(valb.m[3]>=1250000)
{ valc=TKD128_OprMul(valc, cnst_0p875); continue; }
if(valb.m[3]<= 875000)
{ valc=TKD128_OprMul(valc, cnst_1p125); continue; }
break;
}
// printf("Rcp Apx-1:\n %s\n %s\n",
// TKD128_OprDumpAsString(val),
// TKD128_OprDumpAsString(valc));
for(i=0; i<32; i++)
{
valb=TKD128_OprMul(val, valc);
valt=TKD128_OprSub(cnst_2p0, valb);
de=TKD128_OprCompare(valt, cnst_1p0);
if(!de)
break;
// valt=TKD128_OprSub(cnst_3p0, valb);
// valt=TKD128_OprMul(valt, cnst_0p5);
// printf("Rcp Apx-I:\n %s\n %s\n",
// TKD128_OprDumpAsString(valb),
// TKD128_OprDumpAsString(valt));
valc=TKD128_OprMul(valc, valt);
}
// printf("Rcp Apx-2:\n %s\n %s\n",
// TKD128_OprDumpAsString(val),
// TKD128_OprDumpAsString(valc));
return(valc);
}
TKD128_Unpacked4x30 TKD128_OprDiv(
TKD128_Unpacked4x30 vala,
TKD128_Unpacked4x30 valb)
{
TKD128_Unpacked4x30 valb1;
valb1=TKD128_OprRecip(valb);
return(TKD128_OprMul(vala, valb1));
}
int TKD128_OprDump30ToString(char *ct, u32 v)
{
u32 v0, v1, vt;
int i;
vt=v;
for(i=0; i<9; i++)
{
v1=vt/10;
v0=vt-(v1*10);
vt=v1;
ct[8-i]='0'+v0;
}
return(0);
}
char *TKD128_RAlloc(int sz)
{
static char buf[4096];
static u16 pos;
char *ct;
if((pos+sz)>=4096)
pos=0;
ct=buf+pos;
pos+=sz;
return(ct);
}
char *TKD128_RStrDup(char *str)
{
char *ct;
ct=TKD128_RAlloc(strlen(str)+1);
strcpy(ct, str);
return(ct);
}
char *TKD128_OprDumpAsString(
TKD128_Unpacked4x30 val)
{
static char tbuf[256], tb1[48];
char *ct, *cs;
int i, j, k, ex, exn;
if((val.exp==0) && !val.m[3] && !val.m[2] && !val.m[1] && !val.m[0])
{
if(val.sgn)
return("-0.0");
return("0.0");
}
TKD128_OprDump30ToString(tb1+ 0, val.m[3]);
TKD128_OprDump30ToString(tb1+ 9, val.m[2]);
TKD128_OprDump30ToString(tb1+18, val.m[1]);
TKD128_OprDump30ToString(tb1+27, val.m[0]);
tb1[36]=0;
ex=val.exp-TKD128_BIAS;
exn=0;
if((ex>32) || (ex<(-32)))
exn=1;
ct=tbuf;
if(val.sgn)
*ct++='-';
if(exn)
{
*ct++='0';
*ct++='.';
strcpy(ct, tb1+2);
ct+=strlen(ct);
*ct++='E';
sprintf(ct, "%+d", ex);
}else
if(ex<0)
{
k=(-ex)-1;
*ct++='0';
*ct++='.';
while(k--)
*ct++='0';
strcpy(ct, tb1+2);
}else
{
k=ex+1;
cs=tb1+2;
while(k && *cs=='0')
{ *cs++; k--; }
while(k && *cs)
{ *ct++=*cs++; k--; }
if(*cs)
{
*ct++='.';
strcpy(ct, cs);
}else
{
while(k--)
*ct++='0';
*ct++='.';
*ct++='0';
*ct=0;
}
}
return(TKD128_RStrDup(tbuf));
}
u32 TKD128_OprParse30ToString(char *str)
{
u32 v;
int i;
v=0;
for(i=0; i<9; i++)
{ v=(v*10)+(str[i]-'0'); }
return(v);
}
TKD128_Unpacked4x30 TKD128_OprParseFromString(char *str)
{
TKD128_Unpacked4x30 val;
char tb[64];
char *cs, *ct;
int exh, exn, exsz, ex, exadj, sgn;
cs=str;
ct=tb;
sgn=0;
if(*cs=='+')
{ cs++; sgn=0; }
if(*cs=='-')
{ cs++; sgn=1; }
while(*cs=='0')
cs++;
*ct++='0';
*ct++='0';
exh=0;
while((*cs>='0') && (*cs<='9') && (exh<34))
{ *ct++=*cs++; exh++; }
if(exh>=34)
{
while((*cs>='0') && (*cs<='9'))
{ cs++; exh++; }
}else
{
exn=0;
exsz=0;
if(*cs=='.')
{
cs++;
if(!exh)
{
while(*cs=='0')
{ cs++; exsz++; }
}
while((*cs>='0') && (*cs<='9') && ((exh+exn)<34))
{ *ct++=*cs++; exn++; }
}
while((exh+exn)<34)
{ *ct++='0'; exn++; }
while((*cs>='0') && (*cs<='9'))
{ cs++; }
}
*ct=0;
exadj=0;
if((*cs=='e') || (*cs=='E'))
{
if(cs[1]=='+')
{ exadj=atoi(cs+2); }
else
{ exadj=atoi(cs+1); }
}
ex=TKD128_BIAS+exh+exadj-1-exsz;
val.sgn=sgn;
val.exp=ex;
val.m[3]=TKD128_OprParse30ToString(tb+ 0);
val.m[2]=TKD128_OprParse30ToString(tb+ 9);
val.m[1]=TKD128_OprParse30ToString(tb+18);
val.m[0]=TKD128_OprParse30ToString(tb+27);
return(val);
}
TKD128_Packed128 TKD128_OprAddDPD(
TKD128_Packed128 vala,
TKD128_Packed128 valb)
{
TKD128_Packed128 vpc;
TKD128_Unpacked4x30 vna, vnb, vnc;
vna=TKD128_UnpackValueDPD(vala);
vnb=TKD128_UnpackValueDPD(valb);
vnc=TKD128_OprAdd(vna, vnb);
vpc=TKD128_PackValueDPD(vnc);
return(vpc);
}
TKD128_Packed128 TKD128_OprMulDPD(
TKD128_Packed128 vala,
TKD128_Packed128 valb)
{
TKD128_Packed128 vpc;
TKD128_Unpacked4x30 vna, vnb, vnc;
vna=TKD128_UnpackValueDPD(vala);
vnb=TKD128_UnpackValueDPD(valb);
vnc=TKD128_OprMul(vna, vnb);
vpc=TKD128_PackValueDPD(vnc);
return(vpc);
}
TKD128_Packed128 TKD128_OprSubDPD(
TKD128_Packed128 vala,
TKD128_Packed128 valb)
{
TKD128_Packed128 vpc;
TKD128_Unpacked4x30 vna, vnb, vnc;
vna=TKD128_UnpackValueDPD(vala);
vnb=TKD128_UnpackValueDPD(valb);
vnc=TKD128_OprSub(vna, vnb);
vpc=TKD128_PackValueDPD(vnc);
return(vpc);
}
TKD128_Packed128 TKD128_OprDivDPD(
TKD128_Packed128 vala,
TKD128_Packed128 valb)
{
TKD128_Packed128 vpc;
TKD128_Unpacked4x30 vna, vnb, vnc;
vna=TKD128_UnpackValueDPD(vala);
vnb=TKD128_UnpackValueDPD(valb);
vnc=TKD128_OprDiv(vna, vnb);
vpc=TKD128_PackValueDPD(vnc);
return(vpc);
}
TKD128_Packed128 TKD128_OprAddX30(
TKD128_Packed128 vala,
TKD128_Packed128 valb)
{
TKD128_Packed128 vpc;
TKD128_Unpacked4x30 vna, vnb, vnc;
// vna=TKD128_UnpackValueX30(vala);
// vnb=TKD128_UnpackValueX30(valb);
TKD128_UnpackValueX30R(&vna, &vala);
TKD128_UnpackValueX30R(&vnb, &valb);
vnc=TKD128_OprAdd(vna, vnb);
// vpc=TKD128_PackValueX30(vnc);
TKD128_PackValueX30R(&vnc, &vpc);
return(vpc);
}
TKD128_Packed128 TKD128_OprMulX30(
TKD128_Packed128 vala,
TKD128_Packed128 valb)
{
TKD128_Packed128 vpc;
TKD128_Unpacked4x30 vna, vnb, vnc;
// vna=TKD128_UnpackValueX30(vala);
// vnb=TKD128_UnpackValueX30(valb);
TKD128_UnpackValueX30R(&vna, &vala);
TKD128_UnpackValueX30R(&vnb, &valb);
vnc=TKD128_OprMul(vna, vnb);
// vpc=TKD128_PackValueX30(vnc);
TKD128_PackValueX30R(&vnc, &vpc);
return(vpc);
}
char *TKD128_OprDumpAsStringDPD(TKD128_Packed128 val)
{
return(TKD128_OprDumpAsString(TKD128_UnpackValueDPD(val)));
}
TKD128_Packed128 TKD128_OprParseFromStringDPD(char *str)
{
return(TKD128_PackValueDPD(TKD128_OprParseFromString(str)));
}
// #define TEST_PI_DPD_HI 0x4329F0D84E95AF97
// #define TEST_PI_DPD_LO 0x7B1E09EE2B4D7F06
// 2F3F 1D4C 5D9C 1E8A
// 5E7F 1A5B 2C3A 5E3A
#define TEST_PI_DPD_HI 0x2200431857ACAD71ULL
#define TEST_PI_DPD_LO 0xF6F1CAC9A334F3DBULL
int main()
{
TKD128_Packed128 vp0, vp1, vp2;
TKD128_Unpacked4x30 vn0, vn1, vn2, vn3;
u32 v0, v1, v2, v3;
int t0, t1, t2, t0e;
int i, n;
TKD128_InitDPD();
v0=0x12345678;
v2=v0&1023;
if(v2>=1000) v2-=1000;
v3=v0%1000;
printf("A %d %d\n", v2, v3);
v1=((((u64)v0)<<10)+(((u64)v0)<<4)+(((u64)v0)<<3))>>20;
// v1=v0>>10;
// v1=(v0+(v0>>6)+(v0>>7))>>10;
// v1=v0/1000;
v2=v0-v1*1000;
if(v2>1000)
{
v1++;
v2=v0-v1*1000;
}
printf("B %d (%X) %d %d\n", v2, v2, v1, v0/1000);
#if 1
vp0.hi=TEST_PI_DPD_HI;
vp0.lo=TEST_PI_DPD_LO;
vn0=TKD128_UnpackValueDPD(vp0);
printf("DPD Test: %s\n",
TKD128_OprDumpAsString(vn0));
#endif
vn0.sgn=0;
vn1.sgn=0;
vn0.exp=TKD128_BIAS;
vn1.exp=TKD128_BIAS;
vn0.m[3]= 1234567;
vn0.m[2]=890000000;
vn0.m[1]=0;
vn0.m[0]=0;
vn1.m[3]= 3141592;
vn1.m[2]=653589793;
vn1.m[1]=238462643;
vn1.m[0]=383279500;
vp0=TKD128_PackValueDPD(vn0);
vp1=TKD128_PackValueDPD(vn1);
printf("%016llX %016llX\n", vp1.hi, vp1.lo);
vn1=TKD128_UnpackValueDPD(vp1);
for(i=0; i<16; i++)
{
printf("%02X ", ((byte *)(&vp1))[i]);
}
printf("\n");
vn2=TKD128_OprParseFromString("987.654321");
printf("Parse Test 1: %s\n", TKD128_OprDumpAsString(vn2));
vn2=TKD128_OprParseFromString("0.0707123");
printf("Parse Test 2: %s\n", TKD128_OprDumpAsString(vn2));
t0=clock(); t0e=t0+CLOCKS_PER_SEC;
t1=t0; n=0;
while(t1