模拟浮点运算器进行双精度浮点运算

D3D设备初始化以后,双精度浮点运算的结果出现了不同。

在FPU中,存在着三种运算精度:single precision(24bits),double precision(53bits),double extended precision(64bits)。而默认精度是53bits的double precision,也就是双精度浮点。D3D出于性能考虑,会将fpu的计算精度改为单精度。因为fpu线程相关的特性,渲染线程中所有的浮点运算都会保持与D3D一致。这种转变体现在fpu的控制寄存器(CTRL)的变化上,CTRL寄存器的值从007F变成027F。
RC字段,这个字段控制浮点转整型的转换方式,
00 = 朝最接近或者偶数舍入 01 = 朝负无穷大方向舍入
10 = 朝正无穷大方向舍入 11 = 超0方向截断
PC 字段精度控制
00 = 单精度 01 = 保留 10 = 双精度 11 = 扩展精度

D3D针对这种情况也有相应的处理办法,IDirect3D9::CreateDevice()接口的第4个参数,DWORD BehaviorFlags,可以选择D3DCREATE_FPU_PRESERVE,D3DSDK上原文描述:Indicates that the application needs either double-precision floating-point unit (FPU) or FPU exceptions enabled. Microsoft® Direct3D® sets the FPU state each time it is called.
By default, the pipeline uses single precision. Be sure to use this flag to get double precision. Setting the flag will reduce Direct3D performance.
引入这个flag,则D3D将不会修改fpu的状态。 

以下是几种解决方法:
1、 创建D3D设备的时候,使用D3DCREATE_FPU_PRESERVE,这个可以根本解决问题,但是会影响D3D的性能;而且目前底层是封装了,这个方法对于客户端编程的帮助不大。
CreateDevice(D3DADAPTER_DEFAULT, D3DDEVTYPE_HAL, m_hWnd, D3DCREATE_HARDWARE_VERTEXPROCESSING|D3DCREATE_FPU_PRESERVE, &d3dpp, &pD3DDevice);

2、 如果没有精度要求,尽量使用单精度的浮点运算。
3、 采用多线程,因为fpu的状态仅在一个线程中起作用,所以,可以创建一个独立的运算线程进行双精度的浮点运算。

4、如果只是单纯除于整百整千的数值进行输出的话,可以将大数转换成带小数点的字符串。

5、模拟浮点运算器进行双精度浮点运算,不过效率比较低,有待优化:

#include <cfloat>
#include <windows.h>

const UINT32 _FLT_EXP_INIT = 127; // float指数起始值
const UINT32 _FLT_MARK_SGN = 0x80000000; // float符号位掩码
const UINT32 _FLT_MARK_EXP = 0x7F800000; // float指数位掩码
const UINT32 _FLT_MARK_MAN = 0x007FFFFF; // float尾数位掩码
const UINT32 _FLT_OFFSET_SGN = 31; // float符号位偏移量
const UINT32 _FLT_OFFSET_EXP = 23; // float指数位偏移量
const UINT32 _FLT_OFFSET_MAN = 0; // float尾数位偏移量
const UINT64 _DBL_EXP_INIT = 1023; // double指数起始值
const UINT64 _DBL_MARK_SGN = 0x8000000000000000; // double符号位掩码
const UINT64 _DBL_MARK_EXP = 0x7FF0000000000000; // double指数位掩码
const UINT64 _DBL_MARK_MAN = 0x000FFFFFFFFFFFFF; // double尾数位掩码
const UINT64 _DBL_OFFSET_SGN = 63; // double符号位偏移量
const UINT64 _DBL_OFFSET_EXP = 52; // double指数位偏移量
const UINT64 _DBL_OFFSET_MAN = 0; // double尾数位偏移量
const UINT32 _UINT32_MARK_SGN = 0x80000000; // UINT32符号位掩码
const UINT64 _UINT64_MARK_SGN = 0x8000000000000000; // UINT64符号位掩码

//
// 单精度浮点数结构
//
struct FLT_PNT
{
UINT32 sgn : 1; // 符号
UINT32 exp : 8; // 指数
UINT32 man : 23; // 尾数
};

//
// 双精度浮点数结构
//
struct DBL_PNT
{
UINT64 sgn : 1; // 符号
UINT64 exp : 11; // 指数
UINT64 man : 52; // 尾数
};

//
// 从float类型数据中解析出符号位、指数位和尾数位
//
inline void _ParseFloat(OUT FLT_PNT& fltPnt, FLOAT flt)
{
UINT32* pun32Flt = reinterpret_cast<UINT32*>(&flt);

fltPnt.sgn = (*pun32Flt & _FLT_MARK_SGN) >> _FLT_OFFSET_SGN;
fltPnt.exp = (*pun32Flt & _FLT_MARK_EXP) >> _FLT_OFFSET_EXP;
fltPnt.man = (*pun32Flt & _FLT_MARK_MAN) >> _FLT_OFFSET_MAN;
}

void ParseFloat(OUT FLT_PNT& fltPnt, INT32 n32)
{
if (0 == n32)
{
fltPnt.sgn = 0;
fltPnt.exp = 0;
fltPnt.man = 0;
}
else
{
fltPnt.sgn = n32 & _FLT_MARK_SGN;
n32 &= ~_FLT_MARK_SGN;

for (fltPnt.exp = _FLT_EXP_INIT;
(n32 & 1) != 1;
++fltPnt.exp, n32 >>= 1)
{;}

UINT32 un32CountMan = 0;
for (INT32 n32Temp = n32;
1 < n32Temp;
++un32CountMan, n32Temp >>= 1)
{;}

fltPnt.exp += un32CountMan;
if (_FLT_OFFSET_EXP < un32CountMan)
{
fltPnt.man = n32 >> (un32CountMan - _FLT_OFFSET_EXP);
}
else
{
fltPnt.man = n32 << (_FLT_OFFSET_EXP - un32CountMan);
}
}
}

//
// 从double类型数据中解析出符号位、指数位和尾数位
//
inline void _ParseDouble(OUT DBL_PNT& dblPnt, DOUBLE dbl)
{
UINT64* pun64Dbl = reinterpret_cast<UINT64*>(&dbl);

dblPnt.sgn = (*pun64Dbl & _DBL_MARK_SGN) >> _DBL_OFFSET_SGN;
dblPnt.exp = (*pun64Dbl & _DBL_MARK_EXP) >> _DBL_OFFSET_EXP;
dblPnt.man = (*pun64Dbl & _DBL_MARK_MAN) >> _DBL_OFFSET_MAN;
}

//
// 从double类型数据中解析出符号位、指数位和尾数位
//
void ParseDouble(OUT DBL_PNT& dblPnt, INT64 n64)
{
if (0 == n64)
{
dblPnt.sgn = 0;
dblPnt.exp = 0;
dblPnt.man = 0;
}
else
{
dblPnt.sgn = n64 & _DBL_MARK_SGN;
n64 &= ~_DBL_MARK_SGN;

for (dblPnt.exp = _DBL_EXP_INIT;
(n64 & 1) != 1;
++dblPnt.exp, n64 >>= 1)
{;}

UINT64 un64CountMan = 0;
for (INT64 n64Temp = n64;
1 < n64Temp;
++un64CountMan, n64Temp >>= 1)
{;}

dblPnt.exp += un64CountMan;
if (_DBL_OFFSET_EXP < un64CountMan)
{
dblPnt.man = n64 >> (un64CountMan - _DBL_OFFSET_EXP);
}
else
{
dblPnt.man = n64 << (_DBL_OFFSET_EXP - un64CountMan);
}
}
}

//
// float类型数据转换为int32类型数据
//
INT32 FloatToInt32(FLOAT fOrgn)
{
if (0 == fOrgn)
{
return 0;
}

FLT_PNT fltPnt = {0};
_ParseFloat(fltPnt, fOrgn);

const UINT32 un32RealExp = fltPnt.exp - _FLT_EXP_INIT;
UINT32 un32Mantissa = fltPnt.man;
if (_FLT_OFFSET_EXP < un32RealExp)
{
un32Mantissa <<= un32RealExp - _FLT_OFFSET_EXP;
}
else
{
un32Mantissa >>= _FLT_OFFSET_EXP - un32RealExp;
}
un32Mantissa |= 1 << un32RealExp;

return (fltPnt.sgn ? -1 : 1) * un32Mantissa;
}

//
// double类型数据转换为int64类型数据
//
INT64 DoubleToInt64(DOUBLE dOrgn)
{
if (0 == dOrgn)
{
return 0;
}

DBL_PNT dblPnt;
_ParseDouble(dblPnt, dOrgn);

const UINT64 un64RealExp = dblPnt.exp - _DBL_EXP_INIT;
UINT64 un64Mantissa = dblPnt.man;
if (_DBL_OFFSET_EXP < un64RealExp)
{
un64Mantissa <<= un64RealExp - _DBL_OFFSET_EXP;
}
else
{
un64Mantissa >>= _DBL_OFFSET_EXP - un64RealExp;
}
un64Mantissa |= static_cast<UINT64>(1) << un64RealExp;

return (dblPnt.sgn ? -1 : 1) * un64Mantissa;
}

//
// Int32类型数据转换为float类型数据
//
FLOAT Int32ToFloat(INT32 n32Orgn)
{
if (0 == n32Orgn)
{
return 0;
}

FLT_PNT fltPnt = {0};
ParseFloat(fltPnt, n32Orgn);

UINT32 sgn = fltPnt.sgn;
UINT32 exp = fltPnt.exp;
UINT32 man = fltPnt.man;

UINT32 un32Result = (sgn << _FLT_OFFSET_SGN) & _FLT_MARK_SGN
| (exp << _FLT_OFFSET_EXP) & _FLT_MARK_EXP
| (man << _FLT_OFFSET_MAN) & _FLT_MARK_MAN;

return *reinterpret_cast<FLOAT*>(&un32Result);
}

//
// Int32类型数据转换为float类型数据
//
DOUBLE Int64ToDouble(INT64 n64Orgn)
{
if (0 == n64Orgn)
{
return 0;
}

DBL_PNT dblPnt = {0};
ParseDouble(dblPnt, n64Orgn);

UINT64 sgn = dblPnt.sgn;
UINT64 exp = dblPnt.exp;
UINT64 man = dblPnt.man;

UINT64 un64Result = (sgn << _DBL_OFFSET_SGN) & _DBL_MARK_SGN
| (exp << _DBL_OFFSET_EXP) & _DBL_MARK_EXP
| (man << _DBL_OFFSET_MAN) & _DBL_MARK_MAN;
return *reinterpret_cast<DOUBLE*>(&un64Result);
}

//
// double类型数据相乘运算
//
DOUBLE DoubleDivision(INT64 n64Dividend, INT64 n64Divisor)
{
if (0 == n64Dividend)
{
return 0;
}

if (0 == n64Divisor)
{
return DBL_MAX;
}

UINT64 sgn = (n64Dividend & _DBL_MARK_SGN) ^ (n64Divisor & _DBL_MARK_SGN);
n64Dividend &= ~_DBL_MARK_SGN;
n64Divisor &= ~_DBL_MARK_SGN;

UINT64 exp = _DBL_EXP_INIT;
UINT64 exp1 = 0;
UINT64 exp2 = 0;

for (UINT64 un64Mark = ~0 << 1;
0 != (n64Dividend & un64Mark);
++exp1, un64Mark <<= 1)
{;}

for (UINT64 un64Mark = ~0 << 1;
0 != (n64Divisor & un64Mark);
++exp2, un64Mark <<= 1)
{;}

if (n64Dividend < n64Divisor && n64Divisor > (n64Dividend << (exp2 - exp1))
|| n64Dividend > n64Divisor && n64Dividend < (n64Divisor << (exp1 - exp2)))
{
--exp;
}
exp += exp1 - exp2;

for (; (n64Dividend & _DBL_MARK_SGN) == 0; n64Dividend <<= 1)
{;}

for (; (n64Divisor & 1) == 0; n64Divisor >>= 1)
{;}

UINT64* pun64Dividend = reinterpret_cast<UINT64*>(&n64Dividend);
UINT64* pun64Divisor = reinterpret_cast<UINT64*>(&n64Divisor);
UINT64 man = *pun64Dividend / *pun64Divisor;

for (; 0 != (man & ~_DBL_MARK_MAN); man >>= 1)
{;}

for (; 0 == (man & ~_DBL_MARK_MAN); man <<= 1)
{;}

UINT64 un64Result =
sgn & _DBL_MARK_SGN
| (exp << _DBL_OFFSET_EXP) & _DBL_MARK_EXP
| man & _DBL_MARK_MAN;
return *reinterpret_cast<DOUBLE*>(&un64Result);
}
原文地址:https://www.cnblogs.com/wxxweb/p/2285565.html